LAPACK 3.12.1
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
ssygv.f
Go to the documentation of this file.
1*> \brief \b SSYGV
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8*> Download SSYGV + dependencies
9*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/ssygv.f">
10*> [TGZ]</a>
11*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/ssygv.f">
12*> [ZIP]</a>
13*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/ssygv.f">
14*> [TXT]</a>
15*
16* Definition:
17* ===========
18*
19* SUBROUTINE SSYGV( ITYPE, JOBZ, UPLO, N, A, LDA, B, LDB, W, WORK,
20* LWORK, INFO )
21*
22* .. Scalar Arguments ..
23* CHARACTER JOBZ, UPLO
24* INTEGER INFO, ITYPE, LDA, LDB, LWORK, N
25* ..
26* .. Array Arguments ..
27* REAL A( LDA, * ), B( LDB, * ), W( * ), WORK( * )
28* ..
29*
30*
31*> \par Purpose:
32* =============
33*>
34*> \verbatim
35*>
36*> SSYGV computes all the eigenvalues, and optionally, the eigenvectors
37*> of a real generalized symmetric-definite eigenproblem, of the form
38*> A*x=(lambda)*B*x, A*Bx=(lambda)*x, or B*A*x=(lambda)*x.
39*> Here A and B are assumed to be symmetric and B is also
40*> positive definite.
41*> \endverbatim
42*
43* Arguments:
44* ==========
45*
46*> \param[in] ITYPE
47*> \verbatim
48*> ITYPE is INTEGER
49*> Specifies the problem type to be solved:
50*> = 1: A*x = (lambda)*B*x
51*> = 2: A*B*x = (lambda)*x
52*> = 3: B*A*x = (lambda)*x
53*> \endverbatim
54*>
55*> \param[in] JOBZ
56*> \verbatim
57*> JOBZ is CHARACTER*1
58*> = 'N': Compute eigenvalues only;
59*> = 'V': Compute eigenvalues and eigenvectors.
60*> \endverbatim
61*>
62*> \param[in] UPLO
63*> \verbatim
64*> UPLO is CHARACTER*1
65*> = 'U': Upper triangles of A and B are stored;
66*> = 'L': Lower triangles of A and B are stored.
67*> \endverbatim
68*>
69*> \param[in] N
70*> \verbatim
71*> N is INTEGER
72*> The order of the matrices A and B. N >= 0.
73*> \endverbatim
74*>
75*> \param[in,out] A
76*> \verbatim
77*> A is REAL array, dimension (LDA, N)
78*> On entry, the symmetric matrix A. If UPLO = 'U', the
79*> leading N-by-N upper triangular part of A contains the
80*> upper triangular part of the matrix A. If UPLO = 'L',
81*> the leading N-by-N lower triangular part of A contains
82*> the lower triangular part of the matrix A.
83*>
84*> On exit, if JOBZ = 'V', then if INFO = 0, A contains the
85*> matrix Z of eigenvectors. The eigenvectors are normalized
86*> as follows:
87*> if ITYPE = 1 or 2, Z**T*B*Z = I;
88*> if ITYPE = 3, Z**T*inv(B)*Z = I.
89*> If JOBZ = 'N', then on exit the upper triangle (if UPLO='U')
90*> or the lower triangle (if UPLO='L') of A, including the
91*> diagonal, is destroyed.
92*> \endverbatim
93*>
94*> \param[in] LDA
95*> \verbatim
96*> LDA is INTEGER
97*> The leading dimension of the array A. LDA >= max(1,N).
98*> \endverbatim
99*>
100*> \param[in,out] B
101*> \verbatim
102*> B is REAL array, dimension (LDB, N)
103*> On entry, the symmetric positive definite matrix B.
104*> If UPLO = 'U', the leading N-by-N upper triangular part of B
105*> contains the upper triangular part of the matrix B.
106*> If UPLO = 'L', the leading N-by-N lower triangular part of B
107*> contains the lower triangular part of the matrix B.
108*>
109*> On exit, if INFO <= N, the part of B containing the matrix is
110*> overwritten by the triangular factor U or L from the Cholesky
111*> factorization B = U**T*U or B = L*L**T.
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] W
121*> \verbatim
122*> W is REAL array, dimension (N)
123*> If INFO = 0, the eigenvalues in ascending order.
124*> \endverbatim
125*>
126*> \param[out] WORK
127*> \verbatim
128*> WORK is REAL array, dimension (MAX(1,LWORK))
129*> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
130*> \endverbatim
131*>
132*> \param[in] LWORK
133*> \verbatim
134*> LWORK is INTEGER
135*> The length of the array WORK. LWORK >= max(1,3*N-1).
136*> For optimal efficiency, LWORK >= (NB+2)*N,
137*> where NB is the blocksize for SSYTRD returned by ILAENV.
138*>
139*> If LWORK = -1, then a workspace query is assumed; the routine
140*> only calculates the optimal size of the WORK array, returns
141*> this value as the first entry of the WORK array, and no error
142*> message related to LWORK is issued by XERBLA.
143*> \endverbatim
144*>
145*> \param[out] INFO
146*> \verbatim
147*> INFO is INTEGER
148*> = 0: successful exit
149*> < 0: if INFO = -i, the i-th argument had an illegal value
150*> > 0: SPOTRF or SSYEV returned an error code:
151*> <= N: if INFO = i, SSYEV failed to converge;
152*> i off-diagonal elements of an intermediate
153*> tridiagonal form did not converge to zero;
154*> > N: if INFO = N + i, for 1 <= i <= N, then the leading
155*> principal minor of order i of B is not positive.
156*> The factorization of B could not be completed and
157*> no eigenvalues or eigenvectors were computed.
158*> \endverbatim
159*
160* Authors:
161* ========
162*
163*> \author Univ. of Tennessee
164*> \author Univ. of California Berkeley
165*> \author Univ. of Colorado Denver
166*> \author NAG Ltd.
167*
168*> \ingroup hegv
169*
170* =====================================================================
171 SUBROUTINE ssygv( ITYPE, JOBZ, UPLO, N, A, LDA, B, LDB, W,
172 $ WORK,
173 $ LWORK, INFO )
174*
175* -- LAPACK driver routine --
176* -- LAPACK is a software package provided by Univ. of Tennessee, --
177* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
178*
179* .. Scalar Arguments ..
180 CHARACTER JOBZ, UPLO
181 INTEGER INFO, ITYPE, LDA, LDB, LWORK, N
182* ..
183* .. Array Arguments ..
184 REAL A( LDA, * ), B( LDB, * ), W( * ), WORK( * )
185* ..
186*
187* =====================================================================
188*
189* .. Parameters ..
190 REAL ONE
191 PARAMETER ( ONE = 1.0e+0 )
192* ..
193* .. Local Scalars ..
194 LOGICAL LQUERY, UPPER, WANTZ
195 CHARACTER TRANS
196 INTEGER LWKMIN, LWKOPT, NB, NEIG
197* ..
198* .. External Functions ..
199 LOGICAL LSAME
200 INTEGER ILAENV
201 REAL SROUNDUP_LWORK
202 EXTERNAL ilaenv, lsame, sroundup_lwork
203* ..
204* .. External Subroutines ..
205 EXTERNAL spotrf, ssyev, ssygst, strmm, strsm,
206 $ xerbla
207* ..
208* .. Intrinsic Functions ..
209 INTRINSIC max
210* ..
211* .. Executable Statements ..
212*
213* Test the input parameters.
214*
215 wantz = lsame( jobz, 'V' )
216 upper = lsame( uplo, 'U' )
217 lquery = ( lwork.EQ.-1 )
218*
219 info = 0
220 IF( itype.LT.1 .OR. itype.GT.3 ) THEN
221 info = -1
222 ELSE IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
223 info = -2
224 ELSE IF( .NOT.( upper .OR. lsame( uplo, 'L' ) ) ) THEN
225 info = -3
226 ELSE IF( n.LT.0 ) THEN
227 info = -4
228 ELSE IF( lda.LT.max( 1, n ) ) THEN
229 info = -6
230 ELSE IF( ldb.LT.max( 1, n ) ) THEN
231 info = -8
232 END IF
233*
234 IF( info.EQ.0 ) THEN
235 lwkmin = max( 1, 3*n - 1 )
236 nb = ilaenv( 1, 'SSYTRD', uplo, n, -1, -1, -1 )
237 lwkopt = max( lwkmin, ( nb + 2 )*n )
238 work( 1 ) = sroundup_lwork(lwkopt)
239*
240 IF( lwork.LT.lwkmin .AND. .NOT.lquery ) THEN
241 info = -11
242 END IF
243 END IF
244*
245 IF( info.NE.0 ) THEN
246 CALL xerbla( 'SSYGV ', -info )
247 RETURN
248 ELSE IF( lquery ) THEN
249 RETURN
250 END IF
251*
252* Quick return if possible
253*
254 IF( n.EQ.0 )
255 $ RETURN
256*
257* Form a Cholesky factorization of B.
258*
259 CALL spotrf( uplo, n, b, ldb, info )
260 IF( info.NE.0 ) THEN
261 info = n + info
262 RETURN
263 END IF
264*
265* Transform problem to standard eigenvalue problem and solve.
266*
267 CALL ssygst( itype, uplo, n, a, lda, b, ldb, info )
268 CALL ssyev( jobz, uplo, n, a, lda, w, work, lwork, info )
269*
270 IF( wantz ) THEN
271*
272* Backtransform eigenvectors to the original problem.
273*
274 neig = n
275 IF( info.GT.0 )
276 $ neig = info - 1
277 IF( itype.EQ.1 .OR. itype.EQ.2 ) THEN
278*
279* For A*x=(lambda)*B*x and A*B*x=(lambda)*x;
280* backtransform eigenvectors: x = inv(L)**T*y or inv(U)*y
281*
282 IF( upper ) THEN
283 trans = 'N'
284 ELSE
285 trans = 'T'
286 END IF
287*
288 CALL strsm( 'Left', uplo, trans, 'Non-unit', n, neig,
289 $ one,
290 $ b, ldb, a, lda )
291*
292 ELSE IF( itype.EQ.3 ) THEN
293*
294* For B*A*x=(lambda)*x;
295* backtransform eigenvectors: x = L*y or U**T*y
296*
297 IF( upper ) THEN
298 trans = 'T'
299 ELSE
300 trans = 'N'
301 END IF
302*
303 CALL strmm( 'Left', uplo, trans, 'Non-unit', n, neig,
304 $ one,
305 $ b, ldb, a, lda )
306 END IF
307 END IF
308*
309 work( 1 ) = sroundup_lwork(lwkopt)
310 RETURN
311*
312* End of SSYGV
313*
314 END
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine ssyev(jobz, uplo, n, a, lda, w, work, lwork, info)
SSYEV computes the eigenvalues and, optionally, the left and/or right eigenvectors for SY matrices
Definition ssyev.f:130
subroutine ssygst(itype, uplo, n, a, lda, b, ldb, info)
SSYGST
Definition ssygst.f:125
subroutine ssygv(itype, jobz, uplo, n, a, lda, b, ldb, w, work, lwork, info)
SSYGV
Definition ssygv.f:174
subroutine spotrf(uplo, n, a, lda, info)
SPOTRF
Definition spotrf.f:105
subroutine strmm(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb)
STRMM
Definition strmm.f:177
subroutine strsm(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb)
STRSM
Definition strsm.f:181