1: *> \brief <b> DSYSV_AA computes the solution to system of linear equations A * X = B for SY matrices</b>
2: *
3: * =========== DOCUMENTATION ===========
4: *
5: * Online html documentation available at
6: * http://www.netlib.org/lapack/explore-html/
7: *
8: *> \htmlonly
9: *> Download DSYSV_AA + dependencies
10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dsysv_aa.f">
11: *> [TGZ]</a>
12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dsysv_aa.f">
13: *> [ZIP]</a>
14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dsysv_aa.f">
15: *> [TXT]</a>
16: *> \endhtmlonly
17: *
18: * Definition:
19: * ===========
20: *
21: * SUBROUTINE DSYSV_AA( UPLO, N, NRHS, A, LDA, IPIV, B, LDB, WORK,
22: * LWORK, INFO )
23: *
24: * .. Scalar Arguments ..
25: * CHARACTER UPLO
26: * INTEGER N, NRHS, LDA, LDB, LWORK, INFO
27: * ..
28: * .. Array Arguments ..
29: * INTEGER IPIV( * )
30: * DOUBLE PRECISION A( LDA, * ), B( LDB, * ), WORK( * )
31: * ..
32: *
33: *
34: *> \par Purpose:
35: * =============
36: *>
37: *> \verbatim
38: *>
39: *> DSYSV computes the solution to a real system of linear equations
40: *> A * X = B,
41: *> where A is an N-by-N symmetric matrix and X and B are N-by-NRHS
42: *> matrices.
43: *>
44: *> Aasen's algorithm is used to factor A as
45: *> A = U**T * T * U, if UPLO = 'U', or
46: *> A = L * T * L**T, if UPLO = 'L',
47: *> where U (or L) is a product of permutation and unit upper (lower)
48: *> triangular matrices, and T is symmetric tridiagonal. The factored
49: *> form of A is then used to solve the system of equations A * X = B.
50: *> \endverbatim
51: *
52: * Arguments:
53: * ==========
54: *
55: *> \param[in] UPLO
56: *> \verbatim
57: *> UPLO is CHARACTER*1
58: *> = 'U': Upper triangle of A is stored;
59: *> = 'L': Lower triangle of A is stored.
60: *> \endverbatim
61: *>
62: *> \param[in] N
63: *> \verbatim
64: *> N is INTEGER
65: *> The number of linear equations, i.e., the order of the
66: *> matrix A. N >= 0.
67: *> \endverbatim
68: *>
69: *> \param[in] NRHS
70: *> \verbatim
71: *> NRHS is INTEGER
72: *> The number of right hand sides, i.e., the number of columns
73: *> of the matrix B. NRHS >= 0.
74: *> \endverbatim
75: *>
76: *> \param[in,out] A
77: *> \verbatim
78: *> A is DOUBLE PRECISION array, dimension (LDA,N)
79: *> On entry, the symmetric matrix A. If UPLO = 'U', the leading
80: *> N-by-N upper triangular part of A contains the upper
81: *> triangular part of the matrix A, and the strictly lower
82: *> triangular part of A is not referenced. If UPLO = 'L', the
83: *> leading N-by-N lower triangular part of A contains the lower
84: *> triangular part of the matrix A, and the strictly upper
85: *> triangular part of A is not referenced.
86: *>
87: *> On exit, if INFO = 0, the tridiagonal matrix T and the
88: *> multipliers used to obtain the factor U or L from the
89: *> factorization A = U**T*T*U or A = L*T*L**T as computed by
90: *> DSYTRF.
91: *> \endverbatim
92: *>
93: *> \param[in] LDA
94: *> \verbatim
95: *> LDA is INTEGER
96: *> The leading dimension of the array A. LDA >= max(1,N).
97: *> \endverbatim
98: *>
99: *> \param[out] IPIV
100: *> \verbatim
101: *> IPIV is INTEGER array, dimension (N)
102: *> On exit, it contains the details of the interchanges, i.e.,
103: *> the row and column k of A were interchanged with the
104: *> row and column IPIV(k).
105: *> \endverbatim
106: *>
107: *> \param[in,out] B
108: *> \verbatim
109: *> B is DOUBLE PRECISION array, dimension (LDB,NRHS)
110: *> On entry, the N-by-NRHS right hand side matrix B.
111: *> On exit, if INFO = 0, the N-by-NRHS solution matrix X.
112: *> \endverbatim
113: *>
114: *> \param[in] LDB
115: *> \verbatim
116: *> LDB is INTEGER
117: *> The leading dimension of the array B. LDB >= max(1,N).
118: *> \endverbatim
119: *>
120: *> \param[out] WORK
121: *> \verbatim
122: *> WORK is DOUBLE PRECISION array, dimension (MAX(1,LWORK))
123: *> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
124: *> \endverbatim
125: *>
126: *> \param[in] LWORK
127: *> \verbatim
128: *> LWORK is INTEGER
129: *> The length of WORK. LWORK >= MAX(1,2*N,3*N-2), and for
130: *> the best performance, LWORK >= MAX(1,N*NB), where NB is
131: *> the optimal blocksize for DSYTRF_AA.
132: *>
133: *> If LWORK = -1, then a workspace query is assumed; the routine
134: *> only calculates the optimal size of the WORK array, returns
135: *> this value as the first entry of the WORK array, and no error
136: *> message related to LWORK is issued by XERBLA.
137: *> \endverbatim
138: *>
139: *> \param[out] INFO
140: *> \verbatim
141: *> INFO is INTEGER
142: *> = 0: successful exit
143: *> < 0: if INFO = -i, the i-th argument had an illegal value
144: *> > 0: if INFO = i, D(i,i) is exactly zero. The factorization
145: *> has been completed, but the block diagonal matrix D is
146: *> exactly singular, so the solution could not be computed.
147: *> \endverbatim
148: *
149: * Authors:
150: * ========
151: *
152: *> \author Univ. of Tennessee
153: *> \author Univ. of California Berkeley
154: *> \author Univ. of Colorado Denver
155: *> \author NAG Ltd.
156: *
157: *> \ingroup doubleSYsolve
158: *
159: * =====================================================================
160: SUBROUTINE DSYSV_AA( UPLO, N, NRHS, A, LDA, IPIV, B, LDB, WORK,
161: $ LWORK, INFO )
162: *
163: * -- LAPACK driver routine --
164: * -- LAPACK is a software package provided by Univ. of Tennessee, --
165: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
166: *
167: * .. Scalar Arguments ..
168: CHARACTER UPLO
169: INTEGER INFO, LDA, LDB, LWORK, N, NRHS
170: * ..
171: * .. Array Arguments ..
172: INTEGER IPIV( * )
173: DOUBLE PRECISION A( LDA, * ), B( LDB, * ), WORK( * )
174: * ..
175: *
176: * =====================================================================
177: *
178: * .. Local Scalars ..
179: LOGICAL LQUERY
180: INTEGER LWKOPT, LWKOPT_SYTRF, LWKOPT_SYTRS
181: * ..
182: * .. External Functions ..
183: LOGICAL LSAME
184: INTEGER ILAENV
185: EXTERNAL ILAENV, LSAME
186: * ..
187: * .. External Subroutines ..
188: EXTERNAL XERBLA, DSYTRF_AA, DSYTRS_AA
189: * ..
190: * .. Intrinsic Functions ..
191: INTRINSIC MAX
192: * ..
193: * .. Executable Statements ..
194: *
195: * Test the input parameters.
196: *
197: INFO = 0
198: LQUERY = ( LWORK.EQ.-1 )
199: IF( .NOT.LSAME( UPLO, 'U' ) .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
200: INFO = -1
201: ELSE IF( N.LT.0 ) THEN
202: INFO = -2
203: ELSE IF( NRHS.LT.0 ) THEN
204: INFO = -3
205: ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
206: INFO = -5
207: ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
208: INFO = -8
209: ELSE IF( LWORK.LT.MAX(2*N, 3*N-2) .AND. .NOT.LQUERY ) THEN
210: INFO = -10
211: END IF
212: *
213: IF( INFO.EQ.0 ) THEN
214: CALL DSYTRF_AA( UPLO, N, A, LDA, IPIV, WORK, -1, INFO )
215: LWKOPT_SYTRF = INT( WORK(1) )
216: CALL DSYTRS_AA( UPLO, N, NRHS, A, LDA, IPIV, B, LDB, WORK,
217: $ -1, INFO )
218: LWKOPT_SYTRS = INT( WORK(1) )
219: LWKOPT = MAX( LWKOPT_SYTRF, LWKOPT_SYTRS )
220: WORK( 1 ) = LWKOPT
221: END IF
222: *
223: IF( INFO.NE.0 ) THEN
224: CALL XERBLA( 'DSYSV_AA ', -INFO )
225: RETURN
226: ELSE IF( LQUERY ) THEN
227: RETURN
228: END IF
229: *
230: * Compute the factorization A = U**T*T*U or A = L*T*L**T.
231: *
232: CALL DSYTRF_AA( UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO )
233: IF( INFO.EQ.0 ) THEN
234: *
235: * Solve the system A*X = B, overwriting B with X.
236: *
237: CALL DSYTRS_AA( UPLO, N, NRHS, A, LDA, IPIV, B, LDB, WORK,
238: $ LWORK, INFO )
239: *
240: END IF
241: *
242: WORK( 1 ) = LWKOPT
243: *
244: RETURN
245: *
246: * End of DSYSV_AA
247: *
248: END
CVSweb interface <joel.bertrand@systella.fr>