1ZGESVX(1) LAPACK driver routine (version 3.2) ZGESVX(1)
2
6 ZGESVX - uses the LU factorization to compute the solution to a complex
7 system of linear equations A * X = B,
8
10 SUBROUTINE ZGESVX( FACT, TRANS, N, NRHS, A, LDA, AF, LDAF, IPIV, EQUED,
11 R, C, B, LDB, X, LDX, RCOND, FERR, BERR, WORK,
12 RWORK, INFO )
13 CHARACTER EQUED, FACT, TRANS
15 INTEGER INFO, LDA, LDAF, LDB, LDX, N, NRHS
17 DOUBLE PRECISION RCOND
19 INTEGER IPIV( * )
21 DOUBLE PRECISION BERR( * ), C( * ), FERR( * ), R( * ),
23 RWORK( * )
24 COMPLEX*16 A( LDA, * ), AF( LDAF, * ), B( LDB, * ), WORK( * ),
26 X( LDX, * )
27
29 ZGESVX uses the LU factorization to compute the solution to a complex
30 system of linear equations
31 A * X = B, where A is an N-by-N matrix and X and B are N-by-NRHS
32 matrices. Error bounds on the solution and a condition estimate are
33 also provided.
34
36 The following steps are performed:
37 1. If FACT = 'E', real scaling factors are computed to equilibrate
38 the system:
39 TRANS = 'N': diag(R)*A*diag(C) *inv(diag(C))*X = diag(R)*B
40 TRANS = 'T': (diag(R)*A*diag(C))**T *inv(diag(R))*X = diag(C)*B
41 TRANS = 'C': (diag(R)*A*diag(C))**H *inv(diag(R))*X = diag(C)*B
42 Whether or not the system will be equilibrated depends on the
43 scaling of the matrix A, but if equilibration is used, A is
44 overwritten by diag(R)*A*diag(C) and B by diag(R)*B (if TRANS='N')
45 or diag(C)*B (if TRANS = 'T' or 'C').
46 2. If FACT = 'N' or 'E', the LU decomposition is used to factor the
47 matrix A (after equilibration if FACT = 'E') as
48 A = P * L * U,
49 where P is a permutation matrix, L is a unit lower triangular
50 matrix, and U is upper triangular.
51 3. If some U(i,i)=0, so that U is exactly singular, then the routine
52 returns with INFO = i. Otherwise, the factored form of A is used
53 to estimate the condition number of the matrix A. If the
54 reciprocal of the condition number is less than machine precision,
55 INFO = N+1 is returned as a warning, but the routine still goes on
56 to solve for X and compute error bounds as described below. 4. The
57 system of equations is solved for X using the factored form
58 of A.
59 5. Iterative refinement is applied to improve the computed solution
60 matrix and calculate error bounds and backward error estimates
61 for it.
62 6. If equilibration was used, the matrix X is premultiplied by
63 diag(C) (if TRANS = 'N') or diag(R) (if TRANS = 'T' or 'C') so
64 that it solves the original system before equilibration.
65
67 FACT (input) CHARACTER*1
68 Specifies whether or not the factored form of the matrix A is
69 supplied on entry, and if not, whether the matrix A should be
70 equilibrated before it is factored. = 'F': On entry, AF and
71 IPIV contain the factored form of A. If EQUED is not 'N', the
72 matrix A has been equilibrated with scaling factors given by R
73 and C. A, AF, and IPIV are not modified. = 'N': The matrix A
74 will be copied to AF and factored.
75 = 'E': The matrix A will be equilibrated if necessary, then
76 copied to AF and factored.
77 TRANS (input) CHARACTER*1
79 Specifies the form of the system of equations:
80 = 'N': A * X = B (No transpose)
81 = 'T': A**T * X = B (Transpose)
82 = 'C': A**H * X = B (Conjugate transpose)
83 N (input) INTEGER
85 The number of linear equations, i.e., the order of the matrix
86 A. N >= 0.
87 NRHS (input) INTEGER
89 The number of right hand sides, i.e., the number of columns of
90 the matrices B and X. NRHS >= 0.
91 A (input/output) COMPLEX*16 array, dimension (LDA,N)
93 On entry, the N-by-N matrix A. If FACT = 'F' and EQUED is not
94 'N', then A must have been equilibrated by the scaling factors
95 in R and/or C. A is not modified if FACT = 'F' or 'N', or if
96 FACT = 'E' and EQUED = 'N' on exit. On exit, if EQUED .ne.
97 'N', A is scaled as follows: EQUED = 'R': A := diag(R) * A
98 EQUED = 'C': A := A * diag(C)
99 EQUED = 'B': A := diag(R) * A * diag(C).
100 LDA (input) INTEGER
102 The leading dimension of the array A. LDA >= max(1,N).
103 AF (input or output) COMPLEX*16 array, dimension (LDAF,N)
105 If FACT = 'F', then AF is an input argument and on entry con‐
106 tains the factors L and U from the factorization A = P*L*U as
107 computed by ZGETRF. If EQUED .ne. 'N', then AF is the factored
108 form of the equilibrated matrix A. If FACT = 'N', then AF is
109 an output argument and on exit returns the factors L and U from
110 the factorization A = P*L*U of the original matrix A. If FACT
111 = 'E', then AF is an output argument and on exit returns the
112 factors L and U from the factorization A = P*L*U of the equili‐
113 brated matrix A (see the description of A for the form of the
114 equilibrated matrix).
115 LDAF (input) INTEGER
117 The leading dimension of the array AF. LDAF >= max(1,N).
118 IPIV (input or output) INTEGER array, dimension (N)
120 If FACT = 'F', then IPIV is an input argument and on entry con‐
121 tains the pivot indices from the factorization A = P*L*U as
122 computed by ZGETRF; row i of the matrix was interchanged with
123 row IPIV(i). If FACT = 'N', then IPIV is an output argument
124 and on exit contains the pivot indices from the factorization A
125 = P*L*U of the original matrix A. If FACT = 'E', then IPIV is
126 an output argument and on exit contains the pivot indices from
127 the factorization A = P*L*U of the equilibrated matrix A.
128 EQUED (input or output) CHARACTER*1
130 Specifies the form of equilibration that was done. = 'N': No
131 equilibration (always true if FACT = 'N').
132 = 'R': Row equilibration, i.e., A has been premultiplied by
133 diag(R). = 'C': Column equilibration, i.e., A has been post‐
134 multiplied by diag(C). = 'B': Both row and column equilibra‐
135 tion, i.e., A has been replaced by diag(R) * A * diag(C).
136 EQUED is an input argument if FACT = 'F'; otherwise, it is an
137 output argument.
138 R (input or output) DOUBLE PRECISION array, dimension (N)
140 The row scale factors for A. If EQUED = 'R' or 'B', A is mul‐
141 tiplied on the left by diag(R); if EQUED = 'N' or 'C', R is not
142 accessed. R is an input argument if FACT = 'F'; otherwise, R
143 is an output argument. If FACT = 'F' and EQUED = 'R' or 'B',
144 each element of R must be positive.
145 C (input or output) DOUBLE PRECISION array, dimension (N)
147 The column scale factors for A. If EQUED = 'C' or 'B', A is
148 multiplied on the right by diag(C); if EQUED = 'N' or 'R', C is
149 not accessed. C is an input argument if FACT = 'F'; otherwise,
150 C is an output argument. If FACT = 'F' and EQUED = 'C' or 'B',
151 each element of C must be positive.
152 B (input/output) COMPLEX*16 array, dimension (LDB,NRHS)
154 On entry, the N-by-NRHS right hand side matrix B. On exit, if
155 EQUED = 'N', B is not modified; if TRANS = 'N' and EQUED = 'R'
156 or 'B', B is overwritten by diag(R)*B; if TRANS = 'T' or 'C'
157 and EQUED = 'C' or 'B', B is overwritten by diag(C)*B.
158 LDB (input) INTEGER
160 The leading dimension of the array B. LDB >= max(1,N).
161 X (output) COMPLEX*16 array, dimension (LDX,NRHS)
163 If INFO = 0 or INFO = N+1, the N-by-NRHS solution matrix X to
164 the original system of equations. Note that A and B are modi‐
165 fied on exit if EQUED .ne. 'N', and the solution to the equili‐
166 brated system is inv(diag(C))*X if TRANS = 'N' and EQUED = 'C'
167 or 'B', or inv(diag(R))*X if TRANS = 'T' or 'C' and EQUED = 'R'
168 or 'B'.
169 LDX (input) INTEGER
171 The leading dimension of the array X. LDX >= max(1,N).
172 RCOND (output) DOUBLE PRECISION
174 The estimate of the reciprocal condition number of the matrix A
175 after equilibration (if done). If RCOND is less than the
176 machine precision (in particular, if RCOND = 0), the matrix is
177 singular to working precision. This condition is indicated by
178 a return code of INFO > 0.
179 FERR (output) DOUBLE PRECISION array, dimension (NRHS)
181 The estimated forward error bound for each solution vector X(j)
182 (the j-th column of the solution matrix X). If XTRUE is the
183 true solution corresponding to X(j), FERR(j) is an estimated
184 upper bound for the magnitude of the largest element in (X(j) -
185 XTRUE) divided by the magnitude of the largest element in X(j).
186 The estimate is as reliable as the estimate for RCOND, and is
187 almost always a slight overestimate of the true error.
188 BERR (output) DOUBLE PRECISION array, dimension (NRHS)
190 The componentwise relative backward error of each solution vec‐
191 tor X(j) (i.e., the smallest relative change in any element of
192 A or B that makes X(j) an exact solution).
193 WORK (workspace) COMPLEX*16 array, dimension (2*N)
195 RWORK (workspace/output) DOUBLE PRECISION array, dimension (2*N)
197 On exit, RWORK(1) contains the reciprocal pivot growth factor
198 norm(A)/norm(U). The "max absolute element" norm is used. If
199 RWORK(1) is much less than 1, then the stability of the LU fac‐
200 torization of the (equilibrated) matrix A could be poor. This
201 also means that the solution X, condition estimator RCOND, and
202 forward error bound FERR could be unreliable. If factorization
203 fails with 0<INFO<=N, then RWORK(1) contains the reciprocal
204 pivot growth factor for the leading INFO columns of A.
205 INFO (output) INTEGER
207 = 0: successful exit
208 < 0: if INFO = -i, the i-th argument had an illegal value
209 > 0: if INFO = i, and i is
210 <= N: U(i,i) is exactly zero. The factorization has been com‐
211 pleted, but the factor U is exactly singular, so the solution
212 and error bounds could not be computed. RCOND = 0 is returned.
213 = N+1: U is nonsingular, but RCOND is less than machine preci‐
214 sion, meaning that the matrix is singular to working precision.
215 Nevertheless, the solution and error bounds are computed
216 because there are a number of situations where the computed
217 solution can be more accurate than the value of RCOND would
218 suggest.
219 LAPACK driver routine (version 3.N2o)vember 2008 ZGESVX(1)