1DSYSVX(1) LAPACK driver routine (version 3.1) DSYSVX(1)
2
6 DSYSVX - the diagonal pivoting factorization to compute the solution to
7 a real system of linear equations A * X = B,
8
10 SUBROUTINE DSYSVX( FACT, UPLO, N, NRHS, A, LDA, AF, LDAF, IPIV, B, LDB,
11 X, LDX, RCOND, FERR, BERR, WORK, LWORK, IWORK, INFO
12 )
13 CHARACTER FACT, UPLO
15 INTEGER INFO, LDA, LDAF, LDB, LDX, LWORK, N, NRHS
17 DOUBLE PRECISION RCOND
19 INTEGER IPIV( * ), IWORK( * )
21 DOUBLE PRECISION A( LDA, * ), AF( LDAF, * ), B( LDB, * ),
23 BERR( * ), FERR( * ), WORK( * ), X( LDX, * )
24
26 DSYSVX uses the diagonal pivoting factorization to compute the solution
27 to a real system of linear equations A * X = B, where A is an N-by-N
28 symmetric matrix and X and B are N-by-NRHS matrices.
29 Error bounds on the solution and a condition estimate are also pro‐
31 vided.
32
35 The following steps are performed:
36 1. If FACT = 'N', the diagonal pivoting method is used to factor A.
38 The form of the factorization is
39 A = U * D * U**T, if UPLO = 'U', or
40 A = L * D * L**T, if UPLO = 'L',
41 where U (or L) is a product of permutation and unit upper (lower)
42 triangular matrices, and D is symmetric and block diagonal with
43 1-by-1 and 2-by-2 diagonal blocks.
44 2. If some D(i,i)=0, so that D is exactly singular, then the routine
46 returns with INFO = i. Otherwise, the factored form of A is used
47 to estimate the condition number of the matrix A. If the
48 reciprocal of the condition number is less than machine precision,
49 INFO = N+1 is returned as a warning, but the routine still goes on
50 to solve for X and compute error bounds as described below.
51 3. The system of equations is solved for X using the factored form
53 of A.
54 4. Iterative refinement is applied to improve the computed solution
56 matrix and calculate error bounds and backward error estimates
57 for it.
58
61 FACT (input) CHARACTER*1
62 Specifies whether or not the factored form of A has been sup‐
63 plied on entry. = 'F': On entry, AF and IPIV contain the fac‐
64 tored form of A. AF and IPIV will not be modified. = 'N':
65 The matrix A will be copied to AF and factored.
66 UPLO (input) CHARACTER*1
68 = 'U': Upper triangle of A is stored;
69 = 'L': Lower triangle of A is stored.
70 N (input) INTEGER
72 The number of linear equations, i.e., the order of the matrix
73 A. N >= 0.
74 NRHS (input) INTEGER
76 The number of right hand sides, i.e., the number of columns of
77 the matrices B and X. NRHS >= 0.
78 A (input) DOUBLE PRECISION array, dimension (LDA,N)
80 The symmetric matrix A. If UPLO = 'U', the leading N-by-N
81 upper triangular part of A contains the upper triangular part
82 of the matrix A, and the strictly lower triangular part of A is
83 not referenced. If UPLO = 'L', the leading N-by-N lower trian‐
84 gular part of A contains the lower triangular part of the
85 matrix A, and the strictly upper triangular part of A is not
86 referenced.
87 LDA (input) INTEGER
89 The leading dimension of the array A. LDA >= max(1,N).
90 AF (input or output) DOUBLE PRECISION array, dimension (LDAF,N)
92 If FACT = 'F', then AF is an input argument and on entry con‐
93 tains the block diagonal matrix D and the multipliers used to
94 obtain the factor U or L from the factorization A = U*D*U**T or
95 A = L*D*L**T as computed by DSYTRF.
96 If FACT = 'N', then AF is an output argument and on exit
98 returns the block diagonal matrix D and the multipliers used to
99 obtain the factor U or L from the factorization A = U*D*U**T or
100 A = L*D*L**T.
101 LDAF (input) INTEGER
103 The leading dimension of the array AF. LDAF >= max(1,N).
104 IPIV (input or output) INTEGER array, dimension (N)
106 If FACT = 'F', then IPIV is an input argument and on entry con‐
107 tains details of the interchanges and the block structure of D,
108 as determined by DSYTRF. If IPIV(k) > 0, then rows and columns
109 k and IPIV(k) were interchanged and D(k,k) is a 1-by-1 diagonal
110 block. If UPLO = 'U' and IPIV(k) = IPIV(k-1) < 0, then rows
111 and columns k-1 and -IPIV(k) were interchanged and
112 D(k-1:k,k-1:k) is a 2-by-2 diagonal block. If UPLO = 'L' and
113 IPIV(k) = IPIV(k+1) < 0, then rows and columns k+1 and -IPIV(k)
114 were interchanged and D(k:k+1,k:k+1) is a 2-by-2 diagonal
115 block.
116 If FACT = 'N', then IPIV is an output argument and on exit con‐
118 tains details of the interchanges and the block structure of D,
119 as determined by DSYTRF.
120 B (input) DOUBLE PRECISION array, dimension (LDB,NRHS)
122 The N-by-NRHS right hand side matrix B.
123 LDB (input) INTEGER
125 The leading dimension of the array B. LDB >= max(1,N).
126 X (output) DOUBLE PRECISION array, dimension (LDX,NRHS)
128 If INFO = 0 or INFO = N+1, the N-by-NRHS solution matrix X.
129 LDX (input) INTEGER
131 The leading dimension of the array X. LDX >= max(1,N).
132 RCOND (output) DOUBLE PRECISION
134 The estimate of the reciprocal condition number of the matrix
135 A. If RCOND is less than the machine precision (in particular,
136 if RCOND = 0), the matrix is singular to working precision.
137 This condition is indicated by a return code of INFO > 0.
138 FERR (output) DOUBLE PRECISION array, dimension (NRHS)
140 The estimated forward error bound for each solution vector X(j)
141 (the j-th column of the solution matrix X). If XTRUE is the
142 true solution corresponding to X(j), FERR(j) is an estimated
143 upper bound for the magnitude of the largest element in (X(j) -
144 XTRUE) divided by the magnitude of the largest element in X(j).
145 The estimate is as reliable as the estimate for RCOND, and is
146 almost always a slight overestimate of the true error.
147 BERR (output) DOUBLE PRECISION array, dimension (NRHS)
149 The componentwise relative backward error of each solution vec‐
150 tor X(j) (i.e., the smallest relative change in any element of
151 A or B that makes X(j) an exact solution).
152 WORK (workspace/output) DOUBLE PRECISION array, dimension
154 (MAX(1,LWORK))
155 On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
156 LWORK (input) INTEGER
158 The length of WORK. LWORK >= max(1,3*N), and for best perfor‐
159 mance, when FACT = 'N', LWORK >= max(1,3*N,N*NB), where NB is
160 the optimal blocksize for DSYTRF.
161 If LWORK = -1, then a workspace query is assumed; the routine
163 only calculates the optimal size of the WORK array, returns
164 this value as the first entry of the WORK array, and no error
165 message related to LWORK is issued by XERBLA.
166 IWORK (workspace) INTEGER array, dimension (N)
168 INFO (output) INTEGER
170 = 0: successful exit
171 < 0: if INFO = -i, the i-th argument had an illegal value
172 > 0: if INFO = i, and i is
173 <= N: D(i,i) is exactly zero. The factorization has been com‐
174 pleted but the factor D is exactly singular, so the solution
175 and error bounds could not be computed. RCOND = 0 is returned.
176 = N+1: D is nonsingular, but RCOND is less than machine preci‐
177 sion, meaning that the matrix is singular to working precision.
178 Nevertheless, the solution and error bounds are computed
179 because there are a number of situations where the computed
180 solution can be more accurate than the value of RCOND would
181 suggest.
182 LAPACK driver routine (version 3.N1o)vember 2006 DSYSVX(1)