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 * U**T, 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*U**T 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: *> \date December 2016
158: *
159: * @precisions fortran d -> z c
160: *
161: *> \ingroup doubleSYsolve
162: *
163: * =====================================================================
164: SUBROUTINE DSYSV_AA( UPLO, N, NRHS, A, LDA, IPIV, B, LDB, WORK,
165: $ LWORK, INFO )
166: *
167: * -- LAPACK driver routine (version 3.7.0) --
168: * -- LAPACK is a software package provided by Univ. of Tennessee, --
169: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
170: * December 2016
171: *
172: * .. Scalar Arguments ..
173: CHARACTER UPLO
174: INTEGER INFO, LDA, LDB, LWORK, N, NRHS
175: * ..
176: * .. Array Arguments ..
177: INTEGER IPIV( * )
178: DOUBLE PRECISION A( LDA, * ), B( LDB, * ), WORK( * )
179: * ..
180: *
181: * =====================================================================
182: *
183: * .. Local Scalars ..
184: LOGICAL LQUERY
185: INTEGER LWKOPT, LWKOPT_SYTRF, LWKOPT_SYTRS
186: * ..
187: * .. External Functions ..
188: LOGICAL LSAME
189: INTEGER ILAENV
190: EXTERNAL ILAENV, LSAME
191: * ..
192: * .. External Subroutines ..
193: EXTERNAL XERBLA, DSYTRF, DSYTRS, DSYTRS2
194: * ..
195: * .. Intrinsic Functions ..
196: INTRINSIC MAX
197: * ..
198: * .. Executable Statements ..
199: *
200: * Test the input parameters.
201: *
202: INFO = 0
203: LQUERY = ( LWORK.EQ.-1 )
204: IF( .NOT.LSAME( UPLO, 'U' ) .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
205: INFO = -1
206: ELSE IF( N.LT.0 ) THEN
207: INFO = -2
208: ELSE IF( NRHS.LT.0 ) THEN
209: INFO = -3
210: ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
211: INFO = -5
212: ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
213: INFO = -8
214: ELSE IF( LWORK.LT.MAX(2*N, 3*N-2) .AND. .NOT.LQUERY ) THEN
215: INFO = -10
216: END IF
217: *
218: IF( INFO.EQ.0 ) THEN
219: CALL DSYTRF_AA( UPLO, N, A, LDA, IPIV, WORK, -1, INFO )
220: LWKOPT_SYTRF = INT( WORK(1) )
221: CALL DSYTRS_AA( UPLO, N, NRHS, A, LDA, IPIV, B, LDB, WORK,
222: $ -1, INFO )
223: LWKOPT_SYTRS = INT( WORK(1) )
224: LWKOPT = MAX( LWKOPT_SYTRF, LWKOPT_SYTRS )
225: WORK( 1 ) = LWKOPT
226: IF( LWORK.LT.LWKOPT .AND. .NOT.LQUERY ) THEN
227: INFO = -10
228: END IF
229: END IF
230: *
231: IF( INFO.NE.0 ) THEN
232: CALL XERBLA( 'DSYSV_AA ', -INFO )
233: RETURN
234: ELSE IF( LQUERY ) THEN
235: RETURN
236: END IF
237: *
238: * Compute the factorization A = U*T*U**T or A = L*T*L**T.
239: *
240: CALL DSYTRF_AA( UPLO, N, A, LDA, IPIV, WORK, LWORK, INFO )
241: IF( INFO.EQ.0 ) THEN
242: *
243: * Solve the system A*X = B, overwriting B with X.
244: *
245: CALL DSYTRS_AA( UPLO, N, NRHS, A, LDA, IPIV, B, LDB, WORK,
246: $ LWORK, INFO )
247: *
248: END IF
249: *
250: WORK( 1 ) = LWKOPT
251: *
252: RETURN
253: *
254: * End of DSYSV_AA
255: *
256: END
CVSweb interface <joel.bertrand@systella.fr>