naginterfaces.library.lapacklin.dpbsvx

naginterfaces.library.lapacklin.dpbsvx(fact, uplo, n, kd, nrhs, ab, afb, equed, s, b)

dpbsvx uses the Cholesky factorization

$$A=U^{T}U\text{or}A=L{L}^T$$

to compute the solution to a real system of linear equations

$$AX=B\text{,}$$

where $A$ is an $n\times n$ symmetric positive definite band matrix of bandwidth $(2k_d+1)$ and $X$ and $B$ are $n\times r$ matrices. Error bounds on the solution and a condition estimate are also provided.

For full information please refer to the NAG Library document for f07hb

https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/f07/f07hbf.html

Parameters

factstr, length 1

Specifies whether or not the factorized form of the matrix $A$ is supplied on entry, and if not, whether the matrix $A$ should be equilibrated before it is factorized.

$fact=\text{'F'}$

$afb$ contains the factorized form of $A$. If $equed=\text{'Y'}$, the matrix $A$ has been equilibrated with scaling factors given by $s$. $ab$ and $afb$ will not be modified.

$fact=\text{'N'}$

The matrix $A$ will be copied to $afb$ and factorized.

$fact=\text{'E'}$

The matrix $A$ will be equilibrated if necessary, then copied to $afb$ and factorized.

uplostr, length 1

If $uplo=\text{'U'}$, the upper triangle of $A$ is stored.

If $uplo=\text{'L'}$, the lower triangle of $A$ is stored.

nint

$n$, the number of linear equations, i.e., the order of the matrix $A$.

kdint

$k_d$, the number of superdiagonals of the matrix $A$ if $uplo=\text{'U'}$, or the number of subdiagonals if $uplo=\text{'L'}$.

nrhsint

$r$, the number of right-hand sides, i.e., the number of columns of the matrix $B$.

abfloat, array-like, shape $(kd+1,n)$

The upper or lower triangle of the symmetric band matrix $A$, except if $fact=\text{'F'}$ and $equed=\text{'Y'}$, in which case $ab$ must contain the equilibrated matrix $D_{S}A{D}_S$.

afbfloat, array-like, shape $(kd+1,n)$

If $fact=\text{'F'}$, $afb$ contains the triangular factor $U$ or $L$ from the Cholesky factorization $A=U^{T}U$ or $A=L{L}^T$ of the band matrix $A$, in the same storage format as $A$. If $equed=\text{'Y'}$, $afb$ is the factorized form of the equilibrated matrix $A$.

equedstr, length 1

If $fact=\text{'N'}$ or $\text{'E'}$, $equed$ need not be set.

If $fact=\text{'F'}$, $equed$ must specify the form of the equilibration that was performed as follows:

if $equed=\text{'N'}$, no equilibration;

if $equed=\text{'Y'}$, equilibration was performed, i.e., $A$ has been replaced by $D_{S}A{D}_S$.

sfloat, array-like, shape $\left(n\right)$

If $fact=\text{'N'}$ or $\text{'E'}$, $s$ need not be set.

If $fact=\text{'F'}$ and $equed=\text{'Y'}$, $s$ must contain the scale factors, $D_S$, for $A$; each element of $s$ must be positive.

bfloat, array-like, shape $(n,nrhs)$

The $n\times r$ right-hand side matrix $B$.

Returns

abfloat, ndarray, shape $(kd+1,n)$

If $fact=\text{'E'}$ and $equed=\text{'Y'}$, $ab$ is overwritten by $D_{S}A{D}_S$.

afbfloat, ndarray, shape $(kd+1,n)$

If $fact=\text{'N'}$, $afb$ returns the triangular factor $U$ or $L$ from the Cholesky factorization $A=U^{T}U$ or $A=L{L}^T$.

If $fact=\text{'E'}$, $afb$ returns the triangular factor $U$ or $L$ from the Cholesky factorization $A=U^{T}U$ or $A=L{L}^T$ of the equilibrated matrix $A$ (see the description of $ab$ for the form of the equilibrated matrix).

equedstr, length 1

If $fact=\text{'F'}$, $equed$ is unchanged from entry.

Otherwise, if no constraints are violated, $equed$ specifies the form of the equilibration that was performed as specified above.

sfloat, ndarray, shape $\left(n\right)$

If $fact=\text{'F'}$, $s$ is unchanged from entry.

Otherwise, if no constraints are violated and $equed=\text{'Y'}$, $s$ contains the scale factors, $D_S$, for $A$; each element of $s$ is positive.

bfloat, ndarray, shape $(n,nrhs)$

If $equed=\text{'N'}$, $b$ is not modified.

If $equed=\text{'Y'}$, $b$ is overwritten by $D_{S}B$.

xfloat, ndarray, shape $(n,nrhs)$

If the function exits successfully or $errno$ = $n$ + 1, the $n\times r$ solution matrix $X$ to the original system of equations. Note that the arrays $A$ and $B$ are modified on exit if $equed=\text{'Y'}$, and the solution to the equilibrated system is $D_S^{-1}X$.

rcondfloat

If no constraints are violated, an estimate of the reciprocal condition number of the matrix $A$ (after equilibration if that is performed), computed as $rcond=1.0/\left({\parallel A\parallel }_1{\Vert A^{-1}\Vert }_1\right)$.

ferrfloat, ndarray, shape $\left(nrhs\right)$

If the function exits successfully or $errno$ = $n$ + 1, an estimate of the forward error bound for each computed solution vector, such that ${\Vert {\hat{x}}_j-x_j\Vert }_{\infty }/{\Vert x_j\Vert }_{\infty }\le ferr[j-1]$ where ${\hat{x}}_j$ is the $j$th column of the computed solution returned in the array $x$ and $x_j$ is the corresponding column of the exact solution $X$. The estimate is as reliable as the estimate for $rcond$, and is almost always a slight overestimate of the true error.

berrfloat, ndarray, shape $\left(nrhs\right)$

If the function exits successfully or $errno$ = $n$ + 1, an estimate of the component-wise relative backward error of each computed solution vector ${\hat{x}}_j$ (i.e., the smallest relative change in any element of $A$ or $B$ that makes ${\hat{x}}_j$ an exact solution).

Raises

NagValueError

(errno $-1$)

On entry, error in parameter $fact$.

Constraint: $fact=\text{'F'}$, $\text{'N'}$ or $\text{'E'}$.

(errno $-2$)

On entry, error in parameter $uplo$.

Constraint: $uplo=\text{'U'}$ or $\text{'L'}$.

(errno $-3$)

On entry, error in parameter $n$.

Constraint: $n\ge 0$.

(errno $-4$)

On entry, error in parameter $kd$.

Constraint: $kd\ge 0$.

(errno $-5$)

On entry, error in parameter $nrhs$.

Constraint: $nrhs\ge 0$.

(errno $-10$)

On entry, error in parameter $equed$.

Constraint: $equed=\text{'N'}$ or $\text{'Y'}$.

(errno $-11$)

On entry, error in parameter $s$.

(errno $(i>0)\text{and}(i\le n)$)

The leading minor of order $⟨value⟩$ of $A$ is not positive definite, so the factorization could not be completed, and the solution has not been computed. $rcond=0.0$ is returned.

Warns

NagAlgorithmicWarning

(errno $n+1$)

$U$ (or $L$) is nonsingular, but $rcond$ is less than machine precision, meaning that the matrix is singular to working precision. Nevertheless, the solution and error bounds are computed because there are a number of situations where the computed solution can be more accurate than the value of $rcond$ would suggest.

Notes

dpbsvx performs the following steps:

1. If $fact=\text{'E'}$, real diagonal scaling factors, $D_S$, are computed to equilibrate the system:

   $$\left(D_{S}A{D}_S\right)\left(D_S^{-1}X\right)=D_{S}B\text{.}$$

   Whether or not the system will be equilibrated depends on the scaling of the matrix $A$, but if equilibration is used, $A$ is overwritten by $D_{S}A{D}_S$ and $B\times D_{S}B$.

2. If $fact=\text{'N'}$ or $\text{'E'}$, the Cholesky decomposition is used to factor the matrix $A$ (after equilibration if $fact=\text{'E'}$) as $A=U^{T}U$ if $uplo=\text{'U'}$ or $A=L{L}^T$ if $uplo=\text{'L'}$, where $U$ is an upper triangular matrix and $L$ is a lower triangular matrix.

3. If the leading $i\times i$ principal minor of $A$ is not positive definite, then the function returns with $\text{errno}=i$. Otherwise, the factored form of $A$ is used to estimate the condition number of the matrix $A$. If the reciprocal of the condition number is less than machine precision, $errno$ = $n$ + 1 is returned as a warning, but the function still goes on to solve for $X$ and compute error bounds as described below.

4. The system of equations is solved for $X$ using the factored form of $A$.

5. Iterative refinement is applied to improve the computed solution matrix and to calculate error bounds and backward error estimates for it.

6. If equilibration was used, the matrix $X$ is premultiplied by $D_S$ so that it solves the original system before equilibration.

References

Anderson, E, Bai, Z, Bischof, C, Blackford, S, Demmel, J, Dongarra, J J, Du Croz, J J, Greenbaum, A, Hammarling, S, McKenney, A and Sorensen, D, 1999, LAPACK Users’ Guide, (3rd Edition), SIAM, Philadelphia, https://www.netlib.org/lapack/lug

Golub, G H and Van Loan, C F, 1996, Matrix Computations, (3rd Edition), Johns Hopkins University Press, Baltimore

Higham, N J, 2002, Accuracy and Stability of Numerical Algorithms, (2nd Edition), SIAM, Philadelphia