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