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