1: *> \brief <b> ZHESV_AA_2STAGE computes the solution to system of linear equations A * X = B for HE 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 ZHESV_AA_2STAGE + dependencies
10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zhesv_aa_2stage.f">
11: *> [TGZ]</a>
12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zhesv_aa_2stage.f">
13: *> [ZIP]</a>
14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zhesv_aa_2stage.f">
15: *> [TXT]</a>
16: *> \endhtmlonly
17: *
18: * Definition:
19: * ===========
20: *
21: * SUBROUTINE ZHESV_AA_2STAGE( UPLO, N, NRHS, A, LDA, TB, LTB,
22: * IPIV, IPIV2, B, LDB, WORK, LWORK,
23: * INFO )
24: *
25: * .. Scalar Arguments ..
26: * CHARACTER UPLO
27: * INTEGER N, NRHS, LDA, LTB, LDB, LWORK, INFO
28: * ..
29: * .. Array Arguments ..
30: * INTEGER IPIV( * ), IPIV2( * )
31: * COMPLEX*16 A( LDA, * ), TB( * ), B( LDB, *), WORK( * )
32: * ..
33: *
34: *
35: *> \par Purpose:
36: * =============
37: *>
38: *> \verbatim
39: *>
40: *> ZHESV_AA_2STAGE computes the solution to a complex system of
41: *> linear equations
42: *> A * X = B,
43: *> where A is an N-by-N Hermitian matrix and X and B are N-by-NRHS
44: *> matrices.
45: *>
46: *> Aasen's 2-stage algorithm is used to factor A as
47: *> A = U**H * T * U, if UPLO = 'U', or
48: *> A = L * T * L**H, if UPLO = 'L',
49: *> where U (or L) is a product of permutation and unit upper (lower)
50: *> triangular matrices, and T is Hermitian and band. The matrix T is
51: *> then LU-factored with partial pivoting. The factored form of A
52: *> is then used to solve the system of equations A * X = B.
53: *>
54: *> This is the blocked version of the algorithm, calling Level 3 BLAS.
55: *> \endverbatim
56: *
57: * Arguments:
58: * ==========
59: *
60: *> \param[in] UPLO
61: *> \verbatim
62: *> UPLO is CHARACTER*1
63: *> = 'U': Upper triangle of A is stored;
64: *> = 'L': Lower triangle of A is stored.
65: *> \endverbatim
66: *>
67: *> \param[in] N
68: *> \verbatim
69: *> N is INTEGER
70: *> The order of the matrix A. N >= 0.
71: *> \endverbatim
72: *>
73: *> \param[in] NRHS
74: *> \verbatim
75: *> NRHS is INTEGER
76: *> The number of right hand sides, i.e., the number of columns
77: *> of the matrix B. NRHS >= 0.
78: *> \endverbatim
79: *>
80: *> \param[in,out] A
81: *> \verbatim
82: *> A is COMPLEX*16 array, dimension (LDA,N)
83: *> On entry, the hermitian matrix A. If UPLO = 'U', the leading
84: *> N-by-N upper triangular part of A contains the upper
85: *> triangular part of the matrix A, and the strictly lower
86: *> triangular part of A is not referenced. If UPLO = 'L', the
87: *> leading N-by-N lower triangular part of A contains the lower
88: *> triangular part of the matrix A, and the strictly upper
89: *> triangular part of A is not referenced.
90: *>
91: *> On exit, L is stored below (or above) the subdiaonal blocks,
92: *> when UPLO is 'L' (or 'U').
93: *> \endverbatim
94: *>
95: *> \param[in] LDA
96: *> \verbatim
97: *> LDA is INTEGER
98: *> The leading dimension of the array A. LDA >= max(1,N).
99: *> \endverbatim
100: *>
101: *> \param[out] TB
102: *> \verbatim
103: *> TB is COMPLEX*16 array, dimension (LTB)
104: *> On exit, details of the LU factorization of the band matrix.
105: *> \endverbatim
106: *>
107: *> \param[in] LTB
108: *> \verbatim
109: *> LTB is INTEGER
110: *> The size of the array TB. LTB >= 4*N, internally
111: *> used to select NB such that LTB >= (3*NB+1)*N.
112: *>
113: *> If LTB = -1, then a workspace query is assumed; the
114: *> routine only calculates the optimal size of LTB,
115: *> returns this value as the first entry of TB, and
116: *> no error message related to LTB is issued by XERBLA.
117: *> \endverbatim
118: *>
119: *> \param[out] IPIV
120: *> \verbatim
121: *> IPIV is INTEGER array, dimension (N)
122: *> On exit, it contains the details of the interchanges, i.e.,
123: *> the row and column k of A were interchanged with the
124: *> row and column IPIV(k).
125: *> \endverbatim
126: *>
127: *> \param[out] IPIV2
128: *> \verbatim
129: *> IPIV2 is INTEGER array, dimension (N)
130: *> On exit, it contains the details of the interchanges, i.e.,
131: *> the row and column k of T were interchanged with the
132: *> row and column IPIV(k).
133: *> \endverbatim
134: *>
135: *> \param[in,out] B
136: *> \verbatim
137: *> B is COMPLEX*16 array, dimension (LDB,NRHS)
138: *> On entry, the right hand side matrix B.
139: *> On exit, the solution matrix X.
140: *> \endverbatim
141: *>
142: *> \param[in] LDB
143: *> \verbatim
144: *> LDB is INTEGER
145: *> The leading dimension of the array B. LDB >= max(1,N).
146: *> \endverbatim
147: *>
148: *> \param[out] WORK
149: *> \verbatim
150: *> WORK is COMPLEX*16 workspace of size LWORK
151: *> \endverbatim
152: *>
153: *> \param[in] LWORK
154: *> \verbatim
155: *> LWORK is INTEGER
156: *> The size of WORK. LWORK >= N, internally used to select NB
157: *> such that LWORK >= N*NB.
158: *>
159: *> If LWORK = -1, then a workspace query is assumed; the
160: *> routine only calculates the optimal size of the WORK array,
161: *> returns this value as the first entry of the WORK array, and
162: *> no error message related to LWORK is issued by XERBLA.
163: *> \endverbatim
164: *>
165: *> \param[out] INFO
166: *> \verbatim
167: *> INFO is INTEGER
168: *> = 0: successful exit
169: *> < 0: if INFO = -i, the i-th argument had an illegal value.
170: *> > 0: if INFO = i, band LU factorization failed on i-th column
171: *> \endverbatim
172: *
173: * Authors:
174: * ========
175: *
176: *> \author Univ. of Tennessee
177: *> \author Univ. of California Berkeley
178: *> \author Univ. of Colorado Denver
179: *> \author NAG Ltd.
180: *
181: *> \ingroup complex16HEsolve
182: *
183: * =====================================================================
184: SUBROUTINE ZHESV_AA_2STAGE( UPLO, N, NRHS, A, LDA, TB, LTB,
185: $ IPIV, IPIV2, B, LDB, WORK, LWORK,
186: $ INFO )
187: *
188: * -- LAPACK driver routine --
189: * -- LAPACK is a software package provided by Univ. of Tennessee, --
190: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
191: *
192: IMPLICIT NONE
193: *
194: * .. Scalar Arguments ..
195: CHARACTER UPLO
196: INTEGER N, NRHS, LDA, LDB, LTB, LWORK, INFO
197: * ..
198: * .. Array Arguments ..
199: INTEGER IPIV( * ), IPIV2( * )
200: COMPLEX*16 A( LDA, * ), B( LDB, * ), TB( * ), WORK( * )
201: * ..
202: *
203: * =====================================================================
204: * .. Parameters ..
205: COMPLEX*16 ZERO, ONE
206: PARAMETER ( ZERO = ( 0.0D+0, 0.0D+0 ),
207: $ ONE = ( 1.0D+0, 0.0D+0 ) )
208: *
209: * .. Local Scalars ..
210: LOGICAL UPPER, TQUERY, WQUERY
211: INTEGER LWKOPT
212: * ..
213: * .. External Functions ..
214: LOGICAL LSAME
215: INTEGER ILAENV
216: EXTERNAL LSAME, ILAENV
217: * ..
218: * .. External Subroutines ..
219: EXTERNAL XERBLA, ZHETRF_AA_2STAGE, ZHETRS_AA_2STAGE
220: * ..
221: * .. Intrinsic Functions ..
222: INTRINSIC MAX
223: * ..
224: * .. Executable Statements ..
225: *
226: * Test the input parameters.
227: *
228: INFO = 0
229: UPPER = LSAME( UPLO, 'U' )
230: WQUERY = ( LWORK.EQ.-1 )
231: TQUERY = ( LTB.EQ.-1 )
232: IF( .NOT.UPPER .AND. .NOT.LSAME( UPLO, 'L' ) ) THEN
233: INFO = -1
234: ELSE IF( N.LT.0 ) THEN
235: INFO = -2
236: ELSE IF( NRHS.LT.0 ) THEN
237: INFO = -3
238: ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
239: INFO = -5
240: ELSE IF( LTB.LT.( 4*N ) .AND. .NOT.TQUERY ) THEN
241: INFO = -7
242: ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
243: INFO = -11
244: ELSE IF( LWORK.LT.N .AND. .NOT.WQUERY ) THEN
245: INFO = -13
246: END IF
247: *
248: IF( INFO.EQ.0 ) THEN
249: CALL ZHETRF_AA_2STAGE( UPLO, N, A, LDA, TB, -1, IPIV,
250: $ IPIV2, WORK, -1, INFO )
251: LWKOPT = INT( WORK(1) )
252: END IF
253: *
254: IF( INFO.NE.0 ) THEN
255: CALL XERBLA( 'ZHESV_AA_2STAGE', -INFO )
256: RETURN
257: ELSE IF( WQUERY .OR. TQUERY ) THEN
258: RETURN
259: END IF
260: *
261: * Compute the factorization A = U**H*T*U or A = L*T*L**H.
262: *
263: CALL ZHETRF_AA_2STAGE( UPLO, N, A, LDA, TB, LTB, IPIV, IPIV2,
264: $ WORK, LWORK, INFO )
265: IF( INFO.EQ.0 ) THEN
266: *
267: * Solve the system A*X = B, overwriting B with X.
268: *
269: CALL ZHETRS_AA_2STAGE( UPLO, N, NRHS, A, LDA, TB, LTB, IPIV,
270: $ IPIV2, B, LDB, INFO )
271: *
272: END IF
273: *
274: WORK( 1 ) = LWKOPT
275: *
276: RETURN
277: *
278: * End of ZHESV_AA_2STAGE
279: *
280: END
CVSweb interface <joel.bertrand@systella.fr>