|
|
◆ Zcgesv()
| Sub Zcgesv |
( |
N As |
Long, |
|
|
A() As |
Complex, |
|
|
IPiv() As |
Long, |
|
|
B() As |
Complex, |
|
|
X() As |
Complex, |
|
|
Iter As |
Long, |
|
|
Info As |
Long, |
|
|
Optional Nrhs As |
Long = 1 |
|
) |
| |
(Simple driver) Solution to system of linear equations AX = B for a complex matrix (mixed precision with iterative refinement)
- Purpose
- This routine computes the solution to a complex system of linear equations where A is an n x n matrix and X and B are n x nrhs matrices.
Zcgesv first attempts to factorize the matrix in single precision and use this factorization within an iterative refinement procedure to produce a solution with double precision normwise backward error quality (see below). If the approach fails the method switches to a double precision factorization and solve.
The iterative refinement is not going to be a winning strategy if the ratio single precision performance over double precision performance is too small. A reasonable strategy should take the number of right-hand sides and the size of the matrix into account. Up to now, we always try iterative refinement.
The iterative refinement process is stopped if
iter > itermax
or for all the rhs we have:
rnrm < sqrt(n)*xnrm*anrm*eps*bwdmax
where
iter is the number of the current iteration in the iterative refinement process
rnrm is the infinity-norm of the residual
xnrm is the infinity-norm of the solution
anrm is the infinity-operator-norm of the matrix A
eps is the machine epsilon returned by Dlamch('E')
The value itermax and bwdmax are fixed to 30 and 1.0D+00 respectively.
- Parameters
-
| [in] | N | Number of linear equations, i.e., order of the matrix A. (N >= 0) (If N = 0, returns without computation) |
| [in,out] | A() | Array A(LA1 - 1, LA2 - 1) (LA1 >= N, LA2 >= N)
[in] N x N coefficient matrix A.
[out] If iterative refinement has been successfully used (Info = 0 and Iter >= 0, see description below), then A() is unchanged. If double precision factorization has been used (Info = 0 and Iter < 0, see description below), then the array A() contains the factors L and U from the factorization A = P*L*U; the unit diagonal elements of L are not stored. |
| [out] | IPiv() | Array IPiv(LIPiv - 1) (LIPiv >= N)
Pivot indices that define the permutation matrix P; row i of the matrix was interchanged with row IPiv(i-1). Corresponds either to the single precision factorization (if Info = 0 and Iter >= 0) or the double precision factorization (if Info = 0 and Iter < 0). |
| [in] | B() | Array B(LB1 - 1, LB2 - 1) (LB1 >= max(1, N), LB2 >= Nrhs) (2D array) or B(LB - 1) (LB >= max(1, N), Nrhs = 1) (1D array)
N x Nrhs matrix of right hand side matrix B. |
| [out] | X() | Array X(LX1 - 1, LX2 - 1) (LX1 >= max(1, N), LX2 >= Nrhs) (2D array) or X(LX - 1) (LX >= max(1, N), Nrhs = 1) (1D array)
If Info = 0, the N x Nrhs solution matrix X. |
| [out] | Iter | < 0: Iterative refinement has failed, double precision factorization has been performed.
= -1: The routine fell back to full precision for implementation- or machine-specific reasons.
= -2: Narrowing the precision induced an overflow, the routine fell back to full precision.
= -3: Failure of Sgetrf.
= -31: Stop the iterative refinement after the 30th iterations.
> 0: Iterative refinement has been successfully used. Returns the number of iterations. |
| [out] | Info | = 0: Successful exit.
= -1: The argument N had an illegal value. (N < 0)
= -2: The argument A() is invalid.
= -3: The argument IPiv() is invalid.
= -4: The argument B() is invalid.
= -5: The argument X() is invalid.
= -8: The argument Nrhs had an illegal value. (Nrhs < 0)
= i > 0: The i-th diagonal element of the factor U computed in double precision is exactly zero. The factorization has been completed, but the factor U is exactly singular, so the solution could not be computed. |
| [in] | Nrhs | (Optional)
Number of right hand sides, i.e., number of columns of the matrix B. (Nrhs >= 0) (If Nrhs = 0, returns without computation) (default = 1) |
- Reference
- LAPACK
- Example Program
- Solve the system of linear equations Ax = B, where
( 0.2-0.11i -0.93-0.32i 0.81+0.37i )
A = ( -0.8-0.92i -0.29+0.86i 0.64+0.51i )
( 0.71+0.59i -0.15+0.19i 0.2+0.94i )
( -0.5853-0.9457i )
B = ( -2.1697-1.0006i )
( 0.0116-0.5094i )
Sub Ex_Zcgesv()
Const N = 3
Dim A(N - 1, N - 1) As Complex, B(N - 1) As Complex, X(N - 1) As Complex
Dim IPiv(N - 1) As Long, Iter As Long, Info As Long
A(0, 0) = Cmplx(0.2, -0.11): A(0, 1) = Cmplx(-0.93, -0.32): A(0, 2) = Cmplx(0.81, 0.37)
A(1, 0) = Cmplx(-0.8, -0.92): A(1, 1) = Cmplx(-0.29, 0.86): A(1, 2) = Cmplx(0.64, 0.51)
A(2, 0) = Cmplx(0.71, 0.59): A(2, 1) = Cmplx(-0.15, 0.19): A(2, 2) = Cmplx(0.2, 0.94)
B(0) = Cmplx(-0.5853, -0.9457): B(1) = Cmplx(-2.1697, -1.0006): B(2) = Cmplx(0.0116, -0.5094)
Call Zcgesv(N, A(), IPiv(), B(), X(), Iter, Info)
Debug.Print "Iter =", Iter, "Info =", Info
End Sub
Function Cmplx(R As Double, Optional I As Double=0) As Complex Building complex number
Function Cimag(A As Complex) As Double Imaginary part of complex number
Function Creal(A As Complex) As Double Real part of complex number
Sub Zcgesv(N As Long, A() As Complex, IPiv() As Long, B() As Complex, X() As Complex, Iter As Long, Info As Long, Optional Nrhs As Long=1) (Simple driver) Solution to system of linear equations AX = B for a complex matrix (mixed precision w...
- Example Results
X = 0.79 -0.13 0.13 0.75 -0.91 0.3
Iter = 2 Info = 0
|