naginterfaces.library.lapacklin.zgesvx

naginterfaces.library.lapacklin.zgesvx(fact, trans, nrhs, a, af, ipiv, equed, r, c, b)

zgesvx uses the $LU$ factorization to compute the solution to a complex system of linear equations

$$AX=B\text{or}{A}^{T}X=B\text{or}{A}^{H}X=B\text{,}$$

where $A$ is an $n\times n$ matrix 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 f07ap

https://support.nag.com/numeric/nl/nagdoc_30.3/flhtml/f07/f07apf.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'}$

$af$ and $ipiv$ contain the factorized form of $A$. If $equed\ne \text{'N'}$, the matrix $A$ has been equilibrated with scaling factors given by $r$ and $c$. $a$, $af$ and $ipiv$ are not modified.

$fact=\text{'N'}$

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

$fact=\text{'E'}$

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

transstr, length 1

Specifies the form of the system of equations.

$trans=\text{'N'}$

$AX=B$ (No transpose).

$trans=\text{'T'}$

$A^{T}X=B$ (Transpose).

$trans=\text{'C'}$

$A^{H}X=B$ (Conjugate transpose).

nrhsint

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

acomplex, array-like, shape $(n,n)$

The $n\times n$ matrix $A$.

If $fact=\text{'F'}$ and $equed\ne \text{'N'}$, $a$ must have been equilibrated by the scaling factors in $r$ and/or $c$.

afcomplex, array-like, shape $(n,n)$

If $fact=\text{'F'}$, $af$ contains the factors $L$ and $U$ from the factorization $A=PLU$ as computed by zgetrf(). If $equed\ne \text{'N'}$, $af$ is the factorized form of the equilibrated matrix $A$.

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

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

If $fact=\text{'F'}$, $ipiv$ contains the pivot indices from the factorization $A=PLU$ as computed by zgetrf(); at the $i$th step row $i$ of the matrix was interchanged with row $ipiv[i-1]$. $ipiv[i-1]=i$ indicates a row interchange was not required.

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

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{'R'}$, row equilibration, i.e., $A$ has been premultiplied by $D_R$;

if $equed=\text{'C'}$, column equilibration, i.e., $A$ has been postmultiplied by $D_C$;

if $equed=\text{'B'}$, both row and column equilibration, i.e., $A$ has been replaced by $D_{R}A{D}_C$.

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

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

If $fact=\text{'F'}$ and $equed=\text{'R'}$ or $\text{'B'}$, $r$ must contain the row scale factors for $A$, $D_R$; each element of $r$ must be positive.

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

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

If $fact=\text{'F'}$ and $equed=\text{'C'}$ or $\text{'B'}$, $c$ must contain the column scale factors for $A$, $D_C$; each element of $c$ must be positive.

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

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

Returns

acomplex, ndarray, shape $(n,n)$

If $fact=\text{'F'}$ or $\text{'N'}$, or if $fact=\text{'E'}$ and $equed=\text{'N'}$, $a$ is not modified.

If $fact=\text{'E'}$ or $equed\ne \text{'N'}$, $A$ is scaled as follows:

if $equed=\text{'R'}$, $A=D_{R}A$;

if $equed=\text{'C'}$, $A=A{D}_C$;

if $equed=\text{'B'}$, $A=D_{R}A{D}_C$.

afcomplex, ndarray, shape $(n,n)$

If $fact=\text{'N'}$, $af$ returns the factors $L$ and $U$ from the factorization $A=PLU$ of the original matrix $A$.

If $fact=\text{'E'}$, $af$ returns the factors $L$ and $U$ from the factorization $A=PLU$ of the equilibrated matrix $A$ (see the description of $a$ for the form of the equilibrated matrix).

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

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

If $fact=\text{'N'}$, $ipiv$ contains the pivot indices from the factorization $A=PLU$ of the original matrix $A$.

If $fact=\text{'E'}$, $ipiv$ contains the pivot indices from the factorization $A=PLU$ of the equilibrated matrix $A$.

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

equedstr, length 1

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

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

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

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

Otherwise, if no constraints are violated and $equed=\text{'R'}$ or $\text{'B'}$, $r$ contains the row scale factors for $A$, $D_R$, such that $A$ is multiplied on the left by $D_R$; each element of $r$ is positive.

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

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

Otherwise, if no constraints are violated and $equed=\text{'C'}$ or $\text{'B'}$, $c$ contains the row scale factors for $A$, $D_C$; each element of $c$ is positive.

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

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

If $trans=\text{'N'}$ and $equed=\text{'R'}$ or $\text{'B'}$, $b$ is overwritten by $D_{R}B$.

If $trans=\text{'T'}$ or $\text{'C'}$ and $equed=\text{'C'}$ or $\text{'B'}$, $b$ is overwritten by $D_{C}B$.

xcomplex, 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\ne \text{'N'}$, and the solution to the equilibrated system is $D_C^{-1}X$ if $trans=\text{'N'}$ and $equed=\text{'C'}$ or $\text{'B'}$, or $D_R^{-1}X$ if $trans=\text{'T'}$ or $\text{'C'}$ and $equed=\text{'R'}$ or $\text{'B'}$.

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).

recip_growth_factorfloat

If no exception or warning is raised, the reciprocal pivot growth factor $\parallel A\parallel /\parallel U\parallel$, where $\parallel .\parallel$ denotes the maximum absolute element norm. If $recip\_growth\_factor\ll 1$, the stability of the $LU$ factorization of (equilibrated) $A$ could be poor. This also means that the solution $x$, condition estimate $rcond$, and forward error bound $ferr$ could be unreliable. If the factorization fails with $errno$ in 1 … $n$, then $recip\_growth\_factor$ contains the reciprocal pivot growth factor for the leading $\text{errno}$ columns of $A$.

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 $trans$.

Constraint: $trans=\text{'N'}$, $\text{'T'}$ or $\text{'C'}$.

(errno $-3$)

On entry, error in parameter $\text{n}$.

Constraint: $n\ge 0$.

(errno $-4$)

On entry, error in parameter $nrhs$.

Constraint: $nrhs\ge 0$.

(errno $-10$)

On entry, error in parameter $equed$.

Constraint: $equed=\text{'N'}$, $\text{'R'}$, $\text{'C'}$ or $\text{'B'}$.

(errno $-11$)

On entry, error in parameter $r$.

(errno $-12$)

On entry, error in parameter $c$.

Warns

NagAlgorithmicWarning

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

Element $⟨value⟩$ of the diagonal is exactly zero. The factorization has been completed, but the factor $U$ is exactly singular, so the solution and error bounds could not be computed. $rcond=0.0$ is returned.

(errno $n+1$)

$U$ 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

zgesvx performs the following steps:

1. Equilibration

   The linear system to be solved may be badly scaled. However, the system can be equilibrated as a first stage by setting $fact=\text{'E'}$. In this case, real scaling factors are computed and these factors then determine whether the system is to be equilibrated. Equilibrated forms of the systems $AX=B$, $A^{T}X=B$ and $A^{H}X=B$ are

   $$\left(D_{R}A{D}_C\right)\left(D_C^{-1}X\right)=D_{R}B\text{,}$$
   $${\left(D_{R}A{D}_C\right)}^T\left(D_R^{-1}X\right)=D_{C}B\text{,}$$

   and

   $${\left(D_{R}A{D}_C\right)}^H\left(D_R^{-1}X\right)=D_{C}B\text{,}$$

   respectively, where $D_R$ and $D_C$ are diagonal matrices, with positive diagonal elements, formed from the computed scaling factors.

   When equilibration is used, $A$ will be overwritten by $D_{R}A{D}_C$ and $B$ will be overwritten by $D_{R}B$ (or $D_{C}B$ when the solution of $A^{T}X=B$ or $A^{H}X=B$ is sought).

2. Factorization

   The matrix $A$, or its scaled form, is copied and factored using the $LU$ decomposition

   $$A=PLU\text{,}$$

   where $P$ is a permutation matrix, $L$ is a unit lower triangular matrix, and $U$ is upper triangular.

   This stage can be by-passed when a factored matrix (with scaled matrices and scaling factors) are supplied; for example, as provided by a previous call to zgesvx with the same matrix $A$.

3. Condition Number Estimation

   The $LU$ factorization of $A$ determines whether a solution to the linear system exists. If some diagonal element of $U$ is zero, then $U$ is exactly singular, no solution exists and the function returns with a failure. Otherwise the factorized 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 then a warning code is returned on final exit.

4. Solution

   The (equilibrated) system is solved for $X$ ($D_C^{-1}X$ or $D_R^{-1}X$) using the factored form of $A$ ($D_{R}A{D}_C$).

5. Iterative Refinement

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

6. Construct Solution Matrix $X$

   If equilibration was used, the matrix $X$ is premultiplied by $D_C$ (if $trans=\text{'N'}$) or $D_R$ (if $trans=\text{'T'}$ or $\text{'C'}$) 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