LAPACK  3.10.0
LAPACK: Linear Algebra PACKage
ssyev.f
Go to the documentation of this file.
1 *> \brief <b> SSYEV computes the eigenvalues and, optionally, the left and/or right eigenvectors 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 SSYEV + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/ssyev.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/ssyev.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/ssyev.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE SSYEV( JOBZ, UPLO, N, A, LDA, W, WORK, LWORK, INFO )
22 *
23 * .. Scalar Arguments ..
24 * CHARACTER JOBZ, UPLO
25 * INTEGER INFO, LDA, LWORK, N
26 * ..
27 * .. Array Arguments ..
28 * REAL A( LDA, * ), W( * ), WORK( * )
29 * ..
30 *
31 *
32 *> \par Purpose:
33 * =============
34 *>
35 *> \verbatim
36 *>
37 *> SSYEV computes all eigenvalues and, optionally, eigenvectors of a
38 *> real symmetric matrix A.
39 *> \endverbatim
40 *
41 * Arguments:
42 * ==========
43 *
44 *> \param[in] JOBZ
45 *> \verbatim
46 *> JOBZ is CHARACTER*1
47 *> = 'N': Compute eigenvalues only;
48 *> = 'V': Compute eigenvalues and eigenvectors.
49 *> \endverbatim
50 *>
51 *> \param[in] UPLO
52 *> \verbatim
53 *> UPLO is CHARACTER*1
54 *> = 'U': Upper triangle of A is stored;
55 *> = 'L': Lower triangle of A is stored.
56 *> \endverbatim
57 *>
58 *> \param[in] N
59 *> \verbatim
60 *> N is INTEGER
61 *> The order of the matrix A. N >= 0.
62 *> \endverbatim
63 *>
64 *> \param[in,out] A
65 *> \verbatim
66 *> A is REAL array, dimension (LDA, N)
67 *> On entry, the symmetric matrix A. If UPLO = 'U', the
68 *> leading N-by-N upper triangular part of A contains the
69 *> upper triangular part of the matrix A. If UPLO = 'L',
70 *> the leading N-by-N lower triangular part of A contains
71 *> the lower triangular part of the matrix A.
72 *> On exit, if JOBZ = 'V', then if INFO = 0, A contains the
73 *> orthonormal eigenvectors of the matrix A.
74 *> If JOBZ = 'N', then on exit the lower triangle (if UPLO='L')
75 *> or the upper triangle (if UPLO='U') of A, including the
76 *> diagonal, is destroyed.
77 *> \endverbatim
78 *>
79 *> \param[in] LDA
80 *> \verbatim
81 *> LDA is INTEGER
82 *> The leading dimension of the array A. LDA >= max(1,N).
83 *> \endverbatim
84 *>
85 *> \param[out] W
86 *> \verbatim
87 *> W is REAL array, dimension (N)
88 *> If INFO = 0, the eigenvalues in ascending order.
89 *> \endverbatim
90 *>
91 *> \param[out] WORK
92 *> \verbatim
93 *> WORK is REAL array, dimension (MAX(1,LWORK))
94 *> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
95 *> \endverbatim
96 *>
97 *> \param[in] LWORK
98 *> \verbatim
99 *> LWORK is INTEGER
100 *> The length of the array WORK. LWORK >= max(1,3*N-1).
101 *> For optimal efficiency, LWORK >= (NB+2)*N,
102 *> where NB is the blocksize for SSYTRD returned by ILAENV.
103 *>
104 *> If LWORK = -1, then a workspace query is assumed; the routine
105 *> only calculates the optimal size of the WORK array, returns
106 *> this value as the first entry of the WORK array, and no error
107 *> message related to LWORK is issued by XERBLA.
108 *> \endverbatim
109 *>
110 *> \param[out] INFO
111 *> \verbatim
112 *> INFO is INTEGER
113 *> = 0: successful exit
114 *> < 0: if INFO = -i, the i-th argument had an illegal value
115 *> > 0: if INFO = i, the algorithm failed to converge; i
116 *> off-diagonal elements of an intermediate tridiagonal
117 *> form did not converge to zero.
118 *> \endverbatim
119 *
120 * Authors:
121 * ========
122 *
123 *> \author Univ. of Tennessee
124 *> \author Univ. of California Berkeley
125 *> \author Univ. of Colorado Denver
126 *> \author NAG Ltd.
127 *
128 *> \ingroup realSYeigen
129 *
130 * =====================================================================
131  SUBROUTINE ssyev( JOBZ, UPLO, N, A, LDA, W, WORK, LWORK, INFO )
132 *
133 * -- LAPACK driver routine --
134 * -- LAPACK is a software package provided by Univ. of Tennessee, --
135 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
136 *
137 * .. Scalar Arguments ..
138  CHARACTER JOBZ, UPLO
139  INTEGER INFO, LDA, LWORK, N
140 * ..
141 * .. Array Arguments ..
142  REAL A( LDA, * ), W( * ), WORK( * )
143 * ..
144 *
145 * =====================================================================
146 *
147 * .. Parameters ..
148  REAL ZERO, ONE
149  parameter( zero = 0.0e0, one = 1.0e0 )
150 * ..
151 * .. Local Scalars ..
152  LOGICAL LOWER, LQUERY, WANTZ
153  INTEGER IINFO, IMAX, INDE, INDTAU, INDWRK, ISCALE,
154  $ LLWORK, LWKOPT, NB
155  REAL ANRM, BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA,
156  $ SMLNUM
157 * ..
158 * .. External Functions ..
159  LOGICAL LSAME
160  INTEGER ILAENV
161  REAL SLAMCH, SLANSY
162  EXTERNAL ilaenv, lsame, slamch, slansy
163 * ..
164 * .. External Subroutines ..
165  EXTERNAL slascl, sorgtr, sscal, ssteqr, ssterf, ssytrd,
166  $ xerbla
167 * ..
168 * .. Intrinsic Functions ..
169  INTRINSIC max, sqrt
170 * ..
171 * .. Executable Statements ..
172 *
173 * Test the input parameters.
174 *
175  wantz = lsame( jobz, 'V' )
176  lower = lsame( uplo, 'L' )
177  lquery = ( lwork.EQ.-1 )
178 *
179  info = 0
180  IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
181  info = -1
182  ELSE IF( .NOT.( lower .OR. lsame( uplo, 'U' ) ) ) THEN
183  info = -2
184  ELSE IF( n.LT.0 ) THEN
185  info = -3
186  ELSE IF( lda.LT.max( 1, n ) ) THEN
187  info = -5
188  END IF
189 *
190  IF( info.EQ.0 ) THEN
191  nb = ilaenv( 1, 'SSYTRD', uplo, n, -1, -1, -1 )
192  lwkopt = max( 1, ( nb+2 )*n )
193  work( 1 ) = lwkopt
194 *
195  IF( lwork.LT.max( 1, 3*n-1 ) .AND. .NOT.lquery )
196  $ info = -8
197  END IF
198 *
199  IF( info.NE.0 ) THEN
200  CALL xerbla( 'SSYEV ', -info )
201  RETURN
202  ELSE IF( lquery ) THEN
203  RETURN
204  END IF
205 *
206 * Quick return if possible
207 *
208  IF( n.EQ.0 ) THEN
209  RETURN
210  END IF
211 *
212  IF( n.EQ.1 ) THEN
213  w( 1 ) = a( 1, 1 )
214  work( 1 ) = 2
215  IF( wantz )
216  $ a( 1, 1 ) = one
217  RETURN
218  END IF
219 *
220 * Get machine constants.
221 *
222  safmin = slamch( 'Safe minimum' )
223  eps = slamch( 'Precision' )
224  smlnum = safmin / eps
225  bignum = one / smlnum
226  rmin = sqrt( smlnum )
227  rmax = sqrt( bignum )
228 *
229 * Scale matrix to allowable range, if necessary.
230 *
231  anrm = slansy( 'M', uplo, n, a, lda, work )
232  iscale = 0
233  IF( anrm.GT.zero .AND. anrm.LT.rmin ) THEN
234  iscale = 1
235  sigma = rmin / anrm
236  ELSE IF( anrm.GT.rmax ) THEN
237  iscale = 1
238  sigma = rmax / anrm
239  END IF
240  IF( iscale.EQ.1 )
241  $ CALL slascl( uplo, 0, 0, one, sigma, n, n, a, lda, info )
242 *
243 * Call SSYTRD to reduce symmetric matrix to tridiagonal form.
244 *
245  inde = 1
246  indtau = inde + n
247  indwrk = indtau + n
248  llwork = lwork - indwrk + 1
249  CALL ssytrd( uplo, n, a, lda, w, work( inde ), work( indtau ),
250  $ work( indwrk ), llwork, iinfo )
251 *
252 * For eigenvalues only, call SSTERF. For eigenvectors, first call
253 * SORGTR to generate the orthogonal matrix, then call SSTEQR.
254 *
255  IF( .NOT.wantz ) THEN
256  CALL ssterf( n, w, work( inde ), info )
257  ELSE
258  CALL sorgtr( uplo, n, a, lda, work( indtau ), work( indwrk ),
259  $ llwork, iinfo )
260  CALL ssteqr( jobz, n, w, work( inde ), a, lda, work( indtau ),
261  $ info )
262  END IF
263 *
264 * If matrix was scaled, then rescale eigenvalues appropriately.
265 *
266  IF( iscale.EQ.1 ) THEN
267  IF( info.EQ.0 ) THEN
268  imax = n
269  ELSE
270  imax = info - 1
271  END IF
272  CALL sscal( imax, one / sigma, w, 1 )
273  END IF
274 *
275 * Set WORK(1) to optimal workspace size.
276 *
277  work( 1 ) = lwkopt
278 *
279  RETURN
280 *
281 * End of SSYEV
282 *
283  END
subroutine slascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
SLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: slascl.f:143
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine ssteqr(COMPZ, N, D, E, Z, LDZ, WORK, INFO)
SSTEQR
Definition: ssteqr.f:131
subroutine ssterf(N, D, E, INFO)
SSTERF
Definition: ssterf.f:86
subroutine sorgtr(UPLO, N, A, LDA, TAU, WORK, LWORK, INFO)
SORGTR
Definition: sorgtr.f:123
subroutine ssytrd(UPLO, N, A, LDA, D, E, TAU, WORK, LWORK, INFO)
SSYTRD
Definition: ssytrd.f:192
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:132
subroutine sscal(N, SA, SX, INCX)
SSCAL
Definition: sscal.f:79