LAPACK  3.6.0
LAPACK: Linear Algebra PACKage
Collaboration diagram for real:

Functions

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 More...
 
subroutine ssyevd (JOBZ, UPLO, N, A, LDA, W, WORK, LWORK, IWORK, LIWORK, INFO)
  SSYEVD computes the eigenvalues and, optionally, the left and/or right eigenvectors for SY matrices More...
 
subroutine ssyevr (JOBZ, RANGE, UPLO, N, A, LDA, VL, VU, IL, IU, ABSTOL, M, W, Z, LDZ, ISUPPZ, WORK, LWORK, IWORK, LIWORK, INFO)
  SSYEVR computes the eigenvalues and, optionally, the left and/or right eigenvectors for SY matrices More...
 
subroutine ssyevx (JOBZ, RANGE, UPLO, N, A, LDA, VL, VU, IL, IU, ABSTOL, M, W, Z, LDZ, WORK, LWORK, IWORK, IFAIL, INFO)
  SSYEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for SY matrices More...
 
subroutine ssygv (ITYPE, JOBZ, UPLO, N, A, LDA, B, LDB, W, WORK, LWORK, INFO)
 SSYGV More...
 
subroutine ssygvd (ITYPE, JOBZ, UPLO, N, A, LDA, B, LDB, W, WORK, LWORK, IWORK, LIWORK, INFO)
 SSYGVD More...
 
subroutine ssygvx (ITYPE, JOBZ, RANGE, UPLO, N, A, LDA, B, LDB, VL, VU, IL, IU, ABSTOL, M, W, Z, LDZ, WORK, LWORK, IWORK, IFAIL, INFO)
 SSYGVX More...
 

Detailed Description

This is the group of real eigenvalue driver functions for SY matrices

Function Documentation

subroutine ssyev ( character  JOBZ,
character  UPLO,
integer  N,
real, dimension( lda, * )  A,
integer  LDA,
real, dimension( * )  W,
real, dimension( * )  WORK,
integer  LWORK,
integer  INFO 
)

SSYEV computes the eigenvalues and, optionally, the left and/or right eigenvectors for SY matrices

Download SSYEV + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 SSYEV computes all eigenvalues and, optionally, eigenvectors of a
 real symmetric matrix A.
Parameters
[in]JOBZ
          JOBZ is CHARACTER*1
          = 'N':  Compute eigenvalues only;
          = 'V':  Compute eigenvalues and eigenvectors.
[in]UPLO
          UPLO is CHARACTER*1
          = 'U':  Upper triangle of A is stored;
          = 'L':  Lower triangle of A is stored.
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in,out]A
          A is REAL array, dimension (LDA, N)
          On entry, the symmetric matrix A.  If UPLO = 'U', the
          leading N-by-N upper triangular part of A contains the
          upper triangular part of the matrix A.  If UPLO = 'L',
          the leading N-by-N lower triangular part of A contains
          the lower triangular part of the matrix A.
          On exit, if JOBZ = 'V', then if INFO = 0, A contains the
          orthonormal eigenvectors of the matrix A.
          If JOBZ = 'N', then on exit the lower triangle (if UPLO='L')
          or the upper triangle (if UPLO='U') of A, including the
          diagonal, is destroyed.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[out]W
          W is REAL array, dimension (N)
          If INFO = 0, the eigenvalues in ascending order.
[out]WORK
          WORK is REAL array, dimension (MAX(1,LWORK))
          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
[in]LWORK
          LWORK is INTEGER
          The length of the array WORK.  LWORK >= max(1,3*N-1).
          For optimal efficiency, LWORK >= (NB+2)*N,
          where NB is the blocksize for SSYTRD returned by ILAENV.

          If LWORK = -1, then a workspace query is assumed; the routine
          only calculates the optimal size of the WORK array, returns
          this value as the first entry of the WORK array, and no error
          message related to LWORK is issued by XERBLA.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  if INFO = i, the algorithm failed to converge; i
                off-diagonal elements of an intermediate tridiagonal
                form did not converge to zero.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011

Definition at line 134 of file ssyev.f.

134 *
135 * -- LAPACK driver routine (version 3.4.0) --
136 * -- LAPACK is a software package provided by Univ. of Tennessee, --
137 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
138 * November 2011
139 *
140 * .. Scalar Arguments ..
141  CHARACTER jobz, uplo
142  INTEGER info, lda, lwork, n
143 * ..
144 * .. Array Arguments ..
145  REAL a( lda, * ), w( * ), work( * )
146 * ..
147 *
148 * =====================================================================
149 *
150 * .. Parameters ..
151  REAL zero, one
152  parameter( zero = 0.0e0, one = 1.0e0 )
153 * ..
154 * .. Local Scalars ..
155  LOGICAL lower, lquery, wantz
156  INTEGER iinfo, imax, inde, indtau, indwrk, iscale,
157  $ llwork, lwkopt, nb
158  REAL anrm, bignum, eps, rmax, rmin, safmin, sigma,
159  $ smlnum
160 * ..
161 * .. External Functions ..
162  LOGICAL lsame
163  INTEGER ilaenv
164  REAL slamch, slansy
165  EXTERNAL ilaenv, lsame, slamch, slansy
166 * ..
167 * .. External Subroutines ..
168  EXTERNAL slascl, sorgtr, sscal, ssteqr, ssterf, ssytrd,
169  $ xerbla
170 * ..
171 * .. Intrinsic Functions ..
172  INTRINSIC max, sqrt
173 * ..
174 * .. Executable Statements ..
175 *
176 * Test the input parameters.
177 *
178  wantz = lsame( jobz, 'V' )
179  lower = lsame( uplo, 'L' )
180  lquery = ( lwork.EQ.-1 )
181 *
182  info = 0
183  IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
184  info = -1
185  ELSE IF( .NOT.( lower .OR. lsame( uplo, 'U' ) ) ) THEN
186  info = -2
187  ELSE IF( n.LT.0 ) THEN
188  info = -3
189  ELSE IF( lda.LT.max( 1, n ) ) THEN
190  info = -5
191  END IF
192 *
193  IF( info.EQ.0 ) THEN
194  nb = ilaenv( 1, 'SSYTRD', uplo, n, -1, -1, -1 )
195  lwkopt = max( 1, ( nb+2 )*n )
196  work( 1 ) = lwkopt
197 *
198  IF( lwork.LT.max( 1, 3*n-1 ) .AND. .NOT.lquery )
199  $ info = -8
200  END IF
201 *
202  IF( info.NE.0 ) THEN
203  CALL xerbla( 'SSYEV ', -info )
204  RETURN
205  ELSE IF( lquery ) THEN
206  RETURN
207  END IF
208 *
209 * Quick return if possible
210 *
211  IF( n.EQ.0 ) THEN
212  RETURN
213  END IF
214 *
215  IF( n.EQ.1 ) THEN
216  w( 1 ) = a( 1, 1 )
217  work( 1 ) = 2
218  IF( wantz )
219  $ a( 1, 1 ) = one
220  RETURN
221  END IF
222 *
223 * Get machine constants.
224 *
225  safmin = slamch( 'Safe minimum' )
226  eps = slamch( 'Precision' )
227  smlnum = safmin / eps
228  bignum = one / smlnum
229  rmin = sqrt( smlnum )
230  rmax = sqrt( bignum )
231 *
232 * Scale matrix to allowable range, if necessary.
233 *
234  anrm = slansy( 'M', uplo, n, a, lda, work )
235  iscale = 0
236  IF( anrm.GT.zero .AND. anrm.LT.rmin ) THEN
237  iscale = 1
238  sigma = rmin / anrm
239  ELSE IF( anrm.GT.rmax ) THEN
240  iscale = 1
241  sigma = rmax / anrm
242  END IF
243  IF( iscale.EQ.1 )
244  $ CALL slascl( uplo, 0, 0, one, sigma, n, n, a, lda, info )
245 *
246 * Call SSYTRD to reduce symmetric matrix to tridiagonal form.
247 *
248  inde = 1
249  indtau = inde + n
250  indwrk = indtau + n
251  llwork = lwork - indwrk + 1
252  CALL ssytrd( uplo, n, a, lda, w, work( inde ), work( indtau ),
253  $ work( indwrk ), llwork, iinfo )
254 *
255 * For eigenvalues only, call SSTERF. For eigenvectors, first call
256 * SORGTR to generate the orthogonal matrix, then call SSTEQR.
257 *
258  IF( .NOT.wantz ) THEN
259  CALL ssterf( n, w, work( inde ), info )
260  ELSE
261  CALL sorgtr( uplo, n, a, lda, work( indtau ), work( indwrk ),
262  $ llwork, iinfo )
263  CALL ssteqr( jobz, n, w, work( inde ), a, lda, work( indtau ),
264  $ info )
265  END IF
266 *
267 * If matrix was scaled, then rescale eigenvalues appropriately.
268 *
269  IF( iscale.EQ.1 ) THEN
270  IF( info.EQ.0 ) THEN
271  imax = n
272  ELSE
273  imax = info - 1
274  END IF
275  CALL sscal( imax, one / sigma, w, 1 )
276  END IF
277 *
278 * Set WORK(1) to optimal workspace size.
279 *
280  work( 1 ) = lwkopt
281 *
282  RETURN
283 *
284 * End of SSYEV
285 *
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
subroutine sorgtr(UPLO, N, A, LDA, TAU, WORK, LWORK, INFO)
SORGTR
Definition: sorgtr.f:125
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine ssytrd(UPLO, N, A, LDA, D, E, TAU, WORK, LWORK, INFO)
SSYTRD
Definition: ssytrd.f:194
subroutine sscal(N, SA, SX, INCX)
SSCAL
Definition: sscal.f:55
real function slansy(NORM, UPLO, N, A, LDA, WORK)
SLANSY returns the value of the 1-norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a real symmetric matrix.
Definition: slansy.f:124
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
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:141
subroutine ssterf(N, D, E, INFO)
SSTERF
Definition: ssterf.f:88
subroutine ssteqr(COMPZ, N, D, E, Z, LDZ, WORK, INFO)
SSTEQR
Definition: ssteqr.f:133

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine ssyevd ( character  JOBZ,
character  UPLO,
integer  N,
real, dimension( lda, * )  A,
integer  LDA,
real, dimension( * )  W,
real, dimension( * )  WORK,
integer  LWORK,
integer, dimension( * )  IWORK,
integer  LIWORK,
integer  INFO 
)

SSYEVD computes the eigenvalues and, optionally, the left and/or right eigenvectors for SY matrices

Download SSYEVD + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 SSYEVD computes all eigenvalues and, optionally, eigenvectors of a
 real symmetric matrix A. If eigenvectors are desired, it uses a
 divide and conquer algorithm.

 The divide and conquer algorithm makes very mild assumptions about
 floating point arithmetic. It will work on machines with a guard
 digit in add/subtract, or on those binary machines without guard
 digits which subtract like the Cray X-MP, Cray Y-MP, Cray C-90, or
 Cray-2. It could conceivably fail on hexadecimal or decimal machines
 without guard digits, but we know of none.

 Because of large use of BLAS of level 3, SSYEVD needs N**2 more
 workspace than SSYEVX.
Parameters
[in]JOBZ
          JOBZ is CHARACTER*1
          = 'N':  Compute eigenvalues only;
          = 'V':  Compute eigenvalues and eigenvectors.
[in]UPLO
          UPLO is CHARACTER*1
          = 'U':  Upper triangle of A is stored;
          = 'L':  Lower triangle of A is stored.
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in,out]A
          A is REAL array, dimension (LDA, N)
          On entry, the symmetric matrix A.  If UPLO = 'U', the
          leading N-by-N upper triangular part of A contains the
          upper triangular part of the matrix A.  If UPLO = 'L',
          the leading N-by-N lower triangular part of A contains
          the lower triangular part of the matrix A.
          On exit, if JOBZ = 'V', then if INFO = 0, A contains the
          orthonormal eigenvectors of the matrix A.
          If JOBZ = 'N', then on exit the lower triangle (if UPLO='L')
          or the upper triangle (if UPLO='U') of A, including the
          diagonal, is destroyed.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[out]W
          W is REAL array, dimension (N)
          If INFO = 0, the eigenvalues in ascending order.
[out]WORK
          WORK is REAL array,
                                         dimension (LWORK)
          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
[in]LWORK
          LWORK is INTEGER
          The dimension of the array WORK.
          If N <= 1,               LWORK must be at least 1.
          If JOBZ = 'N' and N > 1, LWORK must be at least 2*N+1.
          If JOBZ = 'V' and N > 1, LWORK must be at least 
                                                1 + 6*N + 2*N**2.

          If LWORK = -1, then a workspace query is assumed; the routine
          only calculates the optimal sizes of the WORK and IWORK
          arrays, returns these values as the first entries of the WORK
          and IWORK arrays, and no error message related to LWORK or
          LIWORK is issued by XERBLA.
[out]IWORK
          IWORK is INTEGER array, dimension (MAX(1,LIWORK))
          On exit, if INFO = 0, IWORK(1) returns the optimal LIWORK.
[in]LIWORK
          LIWORK is INTEGER
          The dimension of the array IWORK.
          If N <= 1,                LIWORK must be at least 1.
          If JOBZ  = 'N' and N > 1, LIWORK must be at least 1.
          If JOBZ  = 'V' and N > 1, LIWORK must be at least 3 + 5*N.

          If LIWORK = -1, then a workspace query is assumed; the
          routine only calculates the optimal sizes of the WORK and
          IWORK arrays, returns these values as the first entries of
          the WORK and IWORK arrays, and no error message related to
          LWORK or LIWORK is issued by XERBLA.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  if INFO = i and JOBZ = 'N', then the algorithm failed
                to converge; i off-diagonal elements of an intermediate
                tridiagonal form did not converge to zero;
                if INFO = i and JOBZ = 'V', then the algorithm failed
                to compute an eigenvalue while working on the submatrix
                lying in rows and columns INFO/(N+1) through
                mod(INFO,N+1).
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011
Contributors:
Jeff Rutter, Computer Science Division, University of California at Berkeley, USA
Modified by Francoise Tisseur, University of Tennessee
Modified description of INFO. Sven, 16 Feb 05.

Definition at line 185 of file ssyevd.f.

185 *
186 * -- LAPACK driver routine (version 3.4.0) --
187 * -- LAPACK is a software package provided by Univ. of Tennessee, --
188 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
189 * November 2011
190 *
191 * .. Scalar Arguments ..
192  CHARACTER jobz, uplo
193  INTEGER info, lda, liwork, lwork, n
194 * ..
195 * .. Array Arguments ..
196  INTEGER iwork( * )
197  REAL a( lda, * ), w( * ), work( * )
198 * ..
199 *
200 * =====================================================================
201 *
202 * .. Parameters ..
203  REAL zero, one
204  parameter( zero = 0.0e+0, one = 1.0e+0 )
205 * ..
206 * .. Local Scalars ..
207 *
208  LOGICAL lower, lquery, wantz
209  INTEGER iinfo, inde, indtau, indwk2, indwrk, iscale,
210  $ liopt, liwmin, llwork, llwrk2, lopt, lwmin
211  REAL anrm, bignum, eps, rmax, rmin, safmin, sigma,
212  $ smlnum
213 * ..
214 * .. External Functions ..
215  LOGICAL lsame
216  INTEGER ilaenv
217  REAL slamch, slansy
218  EXTERNAL ilaenv, lsame, slamch, slansy
219 * ..
220 * .. External Subroutines ..
221  EXTERNAL slacpy, slascl, sormtr, sscal, sstedc, ssterf,
222  $ ssytrd, xerbla
223 * ..
224 * .. Intrinsic Functions ..
225  INTRINSIC max, sqrt
226 * ..
227 * .. Executable Statements ..
228 *
229 * Test the input parameters.
230 *
231  wantz = lsame( jobz, 'V' )
232  lower = lsame( uplo, 'L' )
233  lquery = ( lwork.EQ.-1 .OR. liwork.EQ.-1 )
234 *
235  info = 0
236  IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
237  info = -1
238  ELSE IF( .NOT.( lower .OR. lsame( uplo, 'U' ) ) ) THEN
239  info = -2
240  ELSE IF( n.LT.0 ) THEN
241  info = -3
242  ELSE IF( lda.LT.max( 1, n ) ) THEN
243  info = -5
244  END IF
245 *
246  IF( info.EQ.0 ) THEN
247  IF( n.LE.1 ) THEN
248  liwmin = 1
249  lwmin = 1
250  lopt = lwmin
251  liopt = liwmin
252  ELSE
253  IF( wantz ) THEN
254  liwmin = 3 + 5*n
255  lwmin = 1 + 6*n + 2*n**2
256  ELSE
257  liwmin = 1
258  lwmin = 2*n + 1
259  END IF
260  lopt = max( lwmin, 2*n +
261  $ ilaenv( 1, 'SSYTRD', uplo, n, -1, -1, -1 ) )
262  liopt = liwmin
263  END IF
264  work( 1 ) = lopt
265  iwork( 1 ) = liopt
266 *
267  IF( lwork.LT.lwmin .AND. .NOT.lquery ) THEN
268  info = -8
269  ELSE IF( liwork.LT.liwmin .AND. .NOT.lquery ) THEN
270  info = -10
271  END IF
272  END IF
273 *
274  IF( info.NE.0 ) THEN
275  CALL xerbla( 'SSYEVD', -info )
276  RETURN
277  ELSE IF( lquery ) THEN
278  RETURN
279  END IF
280 *
281 * Quick return if possible
282 *
283  IF( n.EQ.0 )
284  $ RETURN
285 *
286  IF( n.EQ.1 ) THEN
287  w( 1 ) = a( 1, 1 )
288  IF( wantz )
289  $ a( 1, 1 ) = one
290  RETURN
291  END IF
292 *
293 * Get machine constants.
294 *
295  safmin = slamch( 'Safe minimum' )
296  eps = slamch( 'Precision' )
297  smlnum = safmin / eps
298  bignum = one / smlnum
299  rmin = sqrt( smlnum )
300  rmax = sqrt( bignum )
301 *
302 * Scale matrix to allowable range, if necessary.
303 *
304  anrm = slansy( 'M', uplo, n, a, lda, work )
305  iscale = 0
306  IF( anrm.GT.zero .AND. anrm.LT.rmin ) THEN
307  iscale = 1
308  sigma = rmin / anrm
309  ELSE IF( anrm.GT.rmax ) THEN
310  iscale = 1
311  sigma = rmax / anrm
312  END IF
313  IF( iscale.EQ.1 )
314  $ CALL slascl( uplo, 0, 0, one, sigma, n, n, a, lda, info )
315 *
316 * Call SSYTRD to reduce symmetric matrix to tridiagonal form.
317 *
318  inde = 1
319  indtau = inde + n
320  indwrk = indtau + n
321  llwork = lwork - indwrk + 1
322  indwk2 = indwrk + n*n
323  llwrk2 = lwork - indwk2 + 1
324 *
325  CALL ssytrd( uplo, n, a, lda, w, work( inde ), work( indtau ),
326  $ work( indwrk ), llwork, iinfo )
327 *
328 * For eigenvalues only, call SSTERF. For eigenvectors, first call
329 * SSTEDC to generate the eigenvector matrix, WORK(INDWRK), of the
330 * tridiagonal matrix, then call SORMTR to multiply it by the
331 * Householder transformations stored in A.
332 *
333  IF( .NOT.wantz ) THEN
334  CALL ssterf( n, w, work( inde ), info )
335  ELSE
336  CALL sstedc( 'I', n, w, work( inde ), work( indwrk ), n,
337  $ work( indwk2 ), llwrk2, iwork, liwork, info )
338  CALL sormtr( 'L', uplo, 'N', n, n, a, lda, work( indtau ),
339  $ work( indwrk ), n, work( indwk2 ), llwrk2, iinfo )
340  CALL slacpy( 'A', n, n, work( indwrk ), n, a, lda )
341  END IF
342 *
343 * If matrix was scaled, then rescale eigenvalues appropriately.
344 *
345  IF( iscale.EQ.1 )
346  $ CALL sscal( n, one / sigma, w, 1 )
347 *
348  work( 1 ) = lopt
349  iwork( 1 ) = liopt
350 *
351  RETURN
352 *
353 * End of SSYEVD
354 *
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
subroutine sormtr(SIDE, UPLO, TRANS, M, N, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
SORMTR
Definition: sormtr.f:174
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine ssytrd(UPLO, N, A, LDA, D, E, TAU, WORK, LWORK, INFO)
SSYTRD
Definition: ssytrd.f:194
subroutine sstedc(COMPZ, N, D, E, Z, LDZ, WORK, LWORK, IWORK, LIWORK, INFO)
SSTEDC
Definition: sstedc.f:190
subroutine sscal(N, SA, SX, INCX)
SSCAL
Definition: sscal.f:55
subroutine slacpy(UPLO, M, N, A, LDA, B, LDB)
SLACPY copies all or part of one two-dimensional array to another.
Definition: slacpy.f:105
real function slansy(NORM, UPLO, N, A, LDA, WORK)
SLANSY returns the value of the 1-norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a real symmetric matrix.
Definition: slansy.f:124
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
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:141
subroutine ssterf(N, D, E, INFO)
SSTERF
Definition: ssterf.f:88

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine ssyevr ( character  JOBZ,
character  RANGE,
character  UPLO,
integer  N,
real, dimension( lda, * )  A,
integer  LDA,
real  VL,
real  VU,
integer  IL,
integer  IU,
real  ABSTOL,
integer  M,
real, dimension( * )  W,
real, dimension( ldz, * )  Z,
integer  LDZ,
integer, dimension( * )  ISUPPZ,
real, dimension( * )  WORK,
integer  LWORK,
integer, dimension( * )  IWORK,
integer  LIWORK,
integer  INFO 
)

SSYEVR computes the eigenvalues and, optionally, the left and/or right eigenvectors for SY matrices

Download SSYEVR + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 SSYEVR computes selected eigenvalues and, optionally, eigenvectors
 of a real symmetric matrix A.  Eigenvalues and eigenvectors can be
 selected by specifying either a range of values or a range of
 indices for the desired eigenvalues.

 SSYEVR first reduces the matrix A to tridiagonal form T with a call
 to SSYTRD.  Then, whenever possible, SSYEVR calls SSTEMR to compute
 the eigenspectrum using Relatively Robust Representations.  SSTEMR
 computes eigenvalues by the dqds algorithm, while orthogonal
 eigenvectors are computed from various "good" L D L^T representations
 (also known as Relatively Robust Representations). Gram-Schmidt
 orthogonalization is avoided as far as possible. More specifically,
 the various steps of the algorithm are as follows.

 For each unreduced block (submatrix) of T,
    (a) Compute T - sigma I  = L D L^T, so that L and D
        define all the wanted eigenvalues to high relative accuracy.
        This means that small relative changes in the entries of D and L
        cause only small relative changes in the eigenvalues and
        eigenvectors. The standard (unfactored) representation of the
        tridiagonal matrix T does not have this property in general.
    (b) Compute the eigenvalues to suitable accuracy.
        If the eigenvectors are desired, the algorithm attains full
        accuracy of the computed eigenvalues only right before
        the corresponding vectors have to be computed, see steps c) and d).
    (c) For each cluster of close eigenvalues, select a new
        shift close to the cluster, find a new factorization, and refine
        the shifted eigenvalues to suitable accuracy.
    (d) For each eigenvalue with a large enough relative separation compute
        the corresponding eigenvector by forming a rank revealing twisted
        factorization. Go back to (c) for any clusters that remain.

 The desired accuracy of the output can be specified by the input
 parameter ABSTOL.

 For more details, see SSTEMR's documentation and:
 - Inderjit S. Dhillon and Beresford N. Parlett: "Multiple representations
   to compute orthogonal eigenvectors of symmetric tridiagonal matrices,"
   Linear Algebra and its Applications, 387(1), pp. 1-28, August 2004.
 - Inderjit Dhillon and Beresford Parlett: "Orthogonal Eigenvectors and
   Relative Gaps," SIAM Journal on Matrix Analysis and Applications, Vol. 25,
   2004.  Also LAPACK Working Note 154.
 - Inderjit Dhillon: "A new O(n^2) algorithm for the symmetric
   tridiagonal eigenvalue/eigenvector problem",
   Computer Science Division Technical Report No. UCB/CSD-97-971,
   UC Berkeley, May 1997.


 Note 1 : SSYEVR calls SSTEMR when the full spectrum is requested
 on machines which conform to the ieee-754 floating point standard.
 SSYEVR calls SSTEBZ and SSTEIN on non-ieee machines and
 when partial spectrum requests are made.

 Normal execution of SSTEMR may create NaNs and infinities and
 hence may abort due to a floating point exception in environments
 which do not handle NaNs and infinities in the ieee standard default
 manner.
Parameters
[in]JOBZ
          JOBZ is CHARACTER*1
          = 'N':  Compute eigenvalues only;
          = 'V':  Compute eigenvalues and eigenvectors.
[in]RANGE
          RANGE is CHARACTER*1
          = 'A': all eigenvalues will be found.
          = 'V': all eigenvalues in the half-open interval (VL,VU]
                 will be found.
          = 'I': the IL-th through IU-th eigenvalues will be found.
          For RANGE = 'V' or 'I' and IU - IL < N - 1, SSTEBZ and
          SSTEIN are called
[in]UPLO
          UPLO is CHARACTER*1
          = 'U':  Upper triangle of A is stored;
          = 'L':  Lower triangle of A is stored.
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in,out]A
          A is REAL array, dimension (LDA, N)
          On entry, the symmetric matrix A.  If UPLO = 'U', the
          leading N-by-N upper triangular part of A contains the
          upper triangular part of the matrix A.  If UPLO = 'L',
          the leading N-by-N lower triangular part of A contains
          the lower triangular part of the matrix A.
          On exit, the lower triangle (if UPLO='L') or the upper
          triangle (if UPLO='U') of A, including the diagonal, is
          destroyed.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[in]VL
          VL is REAL
[in]VU
          VU is REAL
          If RANGE='V', the lower and upper bounds of the interval to
          be searched for eigenvalues. VL < VU.
          Not referenced if RANGE = 'A' or 'I'.
[in]IL
          IL is INTEGER
[in]IU
          IU is INTEGER
          If RANGE='I', the indices (in ascending order) of the
          smallest and largest eigenvalues to be returned.
          1 <= IL <= IU <= N, if N > 0; IL = 1 and IU = 0 if N = 0.
          Not referenced if RANGE = 'A' or 'V'.
[in]ABSTOL
          ABSTOL is REAL
          The absolute error tolerance for the eigenvalues.
          An approximate eigenvalue is accepted as converged
          when it is determined to lie in an interval [a,b]
          of width less than or equal to

                  ABSTOL + EPS *   max( |a|,|b| ) ,

          where EPS is the machine precision.  If ABSTOL is less than
          or equal to zero, then  EPS*|T|  will be used in its place,
          where |T| is the 1-norm of the tridiagonal matrix obtained
          by reducing A to tridiagonal form.

          See "Computing Small Singular Values of Bidiagonal Matrices
          with Guaranteed High Relative Accuracy," by Demmel and
          Kahan, LAPACK Working Note #3.

          If high relative accuracy is important, set ABSTOL to
          SLAMCH( 'Safe minimum' ).  Doing so will guarantee that
          eigenvalues are computed to high relative accuracy when
          possible in future releases.  The current code does not
          make any guarantees about high relative accuracy, but
          future releases will. See J. Barlow and J. Demmel,
          "Computing Accurate Eigensystems of Scaled Diagonally
          Dominant Matrices", LAPACK Working Note #7, for a discussion
          of which matrices define their eigenvalues to high relative
          accuracy.
[out]M
          M is INTEGER
          The total number of eigenvalues found.  0 <= M <= N.
          If RANGE = 'A', M = N, and if RANGE = 'I', M = IU-IL+1.
[out]W
          W is REAL array, dimension (N)
          The first M elements contain the selected eigenvalues in
          ascending order.
[out]Z
          Z is REAL array, dimension (LDZ, max(1,M))
          If JOBZ = 'V', then if INFO = 0, the first M columns of Z
          contain the orthonormal eigenvectors of the matrix A
          corresponding to the selected eigenvalues, with the i-th
          column of Z holding the eigenvector associated with W(i).
          If JOBZ = 'N', then Z is not referenced.
          Note: the user must ensure that at least max(1,M) columns are
          supplied in the array Z; if RANGE = 'V', the exact value of M
          is not known in advance and an upper bound must be used.
          Supplying N columns is always safe.
[in]LDZ
          LDZ is INTEGER
          The leading dimension of the array Z.  LDZ >= 1, and if
          JOBZ = 'V', LDZ >= max(1,N).
[out]ISUPPZ
          ISUPPZ is INTEGER array, dimension ( 2*max(1,M) )
          The support of the eigenvectors in Z, i.e., the indices
          indicating the nonzero elements in Z. The i-th eigenvector
          is nonzero only in elements ISUPPZ( 2*i-1 ) through
          ISUPPZ( 2*i ).
          Implemented only for RANGE = 'A' or 'I' and IU - IL = N - 1
[out]WORK
          WORK is REAL array, dimension (MAX(1,LWORK))
          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
[in]LWORK
          LWORK is INTEGER
          The dimension of the array WORK.  LWORK >= max(1,26*N).
          For optimal efficiency, LWORK >= (NB+6)*N,
          where NB is the max of the blocksize for SSYTRD and SORMTR
          returned by ILAENV.

          If LWORK = -1, then a workspace query is assumed; the routine
          only calculates the optimal sizes of the WORK and IWORK
          arrays, returns these values as the first entries of the WORK
          and IWORK arrays, and no error message related to LWORK or
          LIWORK is issued by XERBLA.
[out]IWORK
          IWORK is INTEGER array, dimension (MAX(1,LIWORK))
          On exit, if INFO = 0, IWORK(1) returns the optimal LWORK.
[in]LIWORK
          LIWORK is INTEGER
          The dimension of the array IWORK.  LIWORK >= max(1,10*N).

          If LIWORK = -1, then a workspace query is assumed; the
          routine only calculates the optimal sizes of the WORK and
          IWORK arrays, returns these values as the first entries of
          the WORK and IWORK arrays, and no error message related to
          LWORK or LIWORK is issued by XERBLA.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  Internal error
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
September 2012
Contributors:
Inderjit Dhillon, IBM Almaden, USA
Osni Marques, LBNL/NERSC, USA
Ken Stanley, Computer Science Division, University of California at Berkeley, USA
Jason Riedy, Computer Science Division, University of California at Berkeley, USA

Definition at line 329 of file ssyevr.f.

329 *
330 * -- LAPACK driver routine (version 3.4.2) --
331 * -- LAPACK is a software package provided by Univ. of Tennessee, --
332 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
333 * September 2012
334 *
335 * .. Scalar Arguments ..
336  CHARACTER jobz, range, uplo
337  INTEGER il, info, iu, lda, ldz, liwork, lwork, m, n
338  REAL abstol, vl, vu
339 * ..
340 * .. Array Arguments ..
341  INTEGER isuppz( * ), iwork( * )
342  REAL a( lda, * ), w( * ), work( * ), z( ldz, * )
343 * ..
344 *
345 * =====================================================================
346 *
347 * .. Parameters ..
348  REAL zero, one, two
349  parameter( zero = 0.0e+0, one = 1.0e+0, two = 2.0e+0 )
350 * ..
351 * .. Local Scalars ..
352  LOGICAL alleig, indeig, lower, lquery, test, valeig,
353  $ wantz, tryrac
354  CHARACTER order
355  INTEGER i, ieeeok, iinfo, imax, indd, inddd, inde,
356  $ indee, indibl, indifl, indisp, indiwo, indtau,
357  $ indwk, indwkn, iscale, j, jj, liwmin,
358  $ llwork, llwrkn, lwkopt, lwmin, nb, nsplit
359  REAL abstll, anrm, bignum, eps, rmax, rmin, safmin,
360  $ sigma, smlnum, tmp1, vll, vuu
361 * ..
362 * .. External Functions ..
363  LOGICAL lsame
364  INTEGER ilaenv
365  REAL slamch, slansy
366  EXTERNAL lsame, ilaenv, slamch, slansy
367 * ..
368 * .. External Subroutines ..
369  EXTERNAL scopy, sormtr, sscal, sstebz, sstemr, sstein,
371 * ..
372 * .. Intrinsic Functions ..
373  INTRINSIC max, min, sqrt
374 * ..
375 * .. Executable Statements ..
376 *
377 * Test the input parameters.
378 *
379  ieeeok = ilaenv( 10, 'SSYEVR', 'N', 1, 2, 3, 4 )
380 *
381  lower = lsame( uplo, 'L' )
382  wantz = lsame( jobz, 'V' )
383  alleig = lsame( range, 'A' )
384  valeig = lsame( range, 'V' )
385  indeig = lsame( range, 'I' )
386 *
387  lquery = ( ( lwork.EQ.-1 ) .OR. ( liwork.EQ.-1 ) )
388 *
389  lwmin = max( 1, 26*n )
390  liwmin = max( 1, 10*n )
391 *
392  info = 0
393  IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
394  info = -1
395  ELSE IF( .NOT.( alleig .OR. valeig .OR. indeig ) ) THEN
396  info = -2
397  ELSE IF( .NOT.( lower .OR. lsame( uplo, 'U' ) ) ) THEN
398  info = -3
399  ELSE IF( n.LT.0 ) THEN
400  info = -4
401  ELSE IF( lda.LT.max( 1, n ) ) THEN
402  info = -6
403  ELSE
404  IF( valeig ) THEN
405  IF( n.GT.0 .AND. vu.LE.vl )
406  $ info = -8
407  ELSE IF( indeig ) THEN
408  IF( il.LT.1 .OR. il.GT.max( 1, n ) ) THEN
409  info = -9
410  ELSE IF( iu.LT.min( n, il ) .OR. iu.GT.n ) THEN
411  info = -10
412  END IF
413  END IF
414  END IF
415  IF( info.EQ.0 ) THEN
416  IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) ) THEN
417  info = -15
418  END IF
419  END IF
420 *
421  IF( info.EQ.0 ) THEN
422  nb = ilaenv( 1, 'SSYTRD', uplo, n, -1, -1, -1 )
423  nb = max( nb, ilaenv( 1, 'SORMTR', uplo, n, -1, -1, -1 ) )
424  lwkopt = max( ( nb+1 )*n, lwmin )
425  work( 1 ) = lwkopt
426  iwork( 1 ) = liwmin
427 *
428  IF( lwork.LT.lwmin .AND. .NOT.lquery ) THEN
429  info = -18
430  ELSE IF( liwork.LT.liwmin .AND. .NOT.lquery ) THEN
431  info = -20
432  END IF
433  END IF
434 *
435  IF( info.NE.0 ) THEN
436  CALL xerbla( 'SSYEVR', -info )
437  RETURN
438  ELSE IF( lquery ) THEN
439  RETURN
440  END IF
441 *
442 * Quick return if possible
443 *
444  m = 0
445  IF( n.EQ.0 ) THEN
446  work( 1 ) = 1
447  RETURN
448  END IF
449 *
450  IF( n.EQ.1 ) THEN
451  work( 1 ) = 26
452  IF( alleig .OR. indeig ) THEN
453  m = 1
454  w( 1 ) = a( 1, 1 )
455  ELSE
456  IF( vl.LT.a( 1, 1 ) .AND. vu.GE.a( 1, 1 ) ) THEN
457  m = 1
458  w( 1 ) = a( 1, 1 )
459  END IF
460  END IF
461  IF( wantz ) THEN
462  z( 1, 1 ) = one
463  isuppz( 1 ) = 1
464  isuppz( 2 ) = 1
465  END IF
466  RETURN
467  END IF
468 *
469 * Get machine constants.
470 *
471  safmin = slamch( 'Safe minimum' )
472  eps = slamch( 'Precision' )
473  smlnum = safmin / eps
474  bignum = one / smlnum
475  rmin = sqrt( smlnum )
476  rmax = min( sqrt( bignum ), one / sqrt( sqrt( safmin ) ) )
477 *
478 * Scale matrix to allowable range, if necessary.
479 *
480  iscale = 0
481  abstll = abstol
482  IF (valeig) THEN
483  vll = vl
484  vuu = vu
485  END IF
486  anrm = slansy( 'M', uplo, n, a, lda, work )
487  IF( anrm.GT.zero .AND. anrm.LT.rmin ) THEN
488  iscale = 1
489  sigma = rmin / anrm
490  ELSE IF( anrm.GT.rmax ) THEN
491  iscale = 1
492  sigma = rmax / anrm
493  END IF
494  IF( iscale.EQ.1 ) THEN
495  IF( lower ) THEN
496  DO 10 j = 1, n
497  CALL sscal( n-j+1, sigma, a( j, j ), 1 )
498  10 CONTINUE
499  ELSE
500  DO 20 j = 1, n
501  CALL sscal( j, sigma, a( 1, j ), 1 )
502  20 CONTINUE
503  END IF
504  IF( abstol.GT.0 )
505  $ abstll = abstol*sigma
506  IF( valeig ) THEN
507  vll = vl*sigma
508  vuu = vu*sigma
509  END IF
510  END IF
511 
512 * Initialize indices into workspaces. Note: The IWORK indices are
513 * used only if SSTERF or SSTEMR fail.
514 
515 * WORK(INDTAU:INDTAU+N-1) stores the scalar factors of the
516 * elementary reflectors used in SSYTRD.
517  indtau = 1
518 * WORK(INDD:INDD+N-1) stores the tridiagonal's diagonal entries.
519  indd = indtau + n
520 * WORK(INDE:INDE+N-1) stores the off-diagonal entries of the
521 * tridiagonal matrix from SSYTRD.
522  inde = indd + n
523 * WORK(INDDD:INDDD+N-1) is a copy of the diagonal entries over
524 * -written by SSTEMR (the SSTERF path copies the diagonal to W).
525  inddd = inde + n
526 * WORK(INDEE:INDEE+N-1) is a copy of the off-diagonal entries over
527 * -written while computing the eigenvalues in SSTERF and SSTEMR.
528  indee = inddd + n
529 * INDWK is the starting offset of the left-over workspace, and
530 * LLWORK is the remaining workspace size.
531  indwk = indee + n
532  llwork = lwork - indwk + 1
533 
534 * IWORK(INDIBL:INDIBL+M-1) corresponds to IBLOCK in SSTEBZ and
535 * stores the block indices of each of the M<=N eigenvalues.
536  indibl = 1
537 * IWORK(INDISP:INDISP+NSPLIT-1) corresponds to ISPLIT in SSTEBZ and
538 * stores the starting and finishing indices of each block.
539  indisp = indibl + n
540 * IWORK(INDIFL:INDIFL+N-1) stores the indices of eigenvectors
541 * that corresponding to eigenvectors that fail to converge in
542 * SSTEIN. This information is discarded; if any fail, the driver
543 * returns INFO > 0.
544  indifl = indisp + n
545 * INDIWO is the offset of the remaining integer workspace.
546  indiwo = indifl + n
547 
548 *
549 * Call SSYTRD to reduce symmetric matrix to tridiagonal form.
550 *
551  CALL ssytrd( uplo, n, a, lda, work( indd ), work( inde ),
552  $ work( indtau ), work( indwk ), llwork, iinfo )
553 *
554 * If all eigenvalues are desired
555 * then call SSTERF or SSTEMR and SORMTR.
556 *
557  test = .false.
558  IF( indeig ) THEN
559  IF( il.EQ.1 .AND. iu.EQ.n ) THEN
560  test = .true.
561  END IF
562  END IF
563  IF( ( alleig.OR.test ) .AND. ( ieeeok.EQ.1 ) ) THEN
564  IF( .NOT.wantz ) THEN
565  CALL scopy( n, work( indd ), 1, w, 1 )
566  CALL scopy( n-1, work( inde ), 1, work( indee ), 1 )
567  CALL ssterf( n, w, work( indee ), info )
568  ELSE
569  CALL scopy( n-1, work( inde ), 1, work( indee ), 1 )
570  CALL scopy( n, work( indd ), 1, work( inddd ), 1 )
571 *
572  IF (abstol .LE. two*n*eps) THEN
573  tryrac = .true.
574  ELSE
575  tryrac = .false.
576  END IF
577  CALL sstemr( jobz, 'A', n, work( inddd ), work( indee ),
578  $ vl, vu, il, iu, m, w, z, ldz, n, isuppz,
579  $ tryrac, work( indwk ), lwork, iwork, liwork,
580  $ info )
581 *
582 *
583 *
584 * Apply orthogonal matrix used in reduction to tridiagonal
585 * form to eigenvectors returned by SSTEIN.
586 *
587  IF( wantz .AND. info.EQ.0 ) THEN
588  indwkn = inde
589  llwrkn = lwork - indwkn + 1
590  CALL sormtr( 'L', uplo, 'N', n, m, a, lda,
591  $ work( indtau ), z, ldz, work( indwkn ),
592  $ llwrkn, iinfo )
593  END IF
594  END IF
595 *
596 *
597  IF( info.EQ.0 ) THEN
598 * Everything worked. Skip SSTEBZ/SSTEIN. IWORK(:) are
599 * undefined.
600  m = n
601  GO TO 30
602  END IF
603  info = 0
604  END IF
605 *
606 * Otherwise, call SSTEBZ and, if eigenvectors are desired, SSTEIN.
607 * Also call SSTEBZ and SSTEIN if SSTEMR fails.
608 *
609  IF( wantz ) THEN
610  order = 'B'
611  ELSE
612  order = 'E'
613  END IF
614 
615  CALL sstebz( range, order, n, vll, vuu, il, iu, abstll,
616  $ work( indd ), work( inde ), m, nsplit, w,
617  $ iwork( indibl ), iwork( indisp ), work( indwk ),
618  $ iwork( indiwo ), info )
619 *
620  IF( wantz ) THEN
621  CALL sstein( n, work( indd ), work( inde ), m, w,
622  $ iwork( indibl ), iwork( indisp ), z, ldz,
623  $ work( indwk ), iwork( indiwo ), iwork( indifl ),
624  $ info )
625 *
626 * Apply orthogonal matrix used in reduction to tridiagonal
627 * form to eigenvectors returned by SSTEIN.
628 *
629  indwkn = inde
630  llwrkn = lwork - indwkn + 1
631  CALL sormtr( 'L', uplo, 'N', n, m, a, lda, work( indtau ), z,
632  $ ldz, work( indwkn ), llwrkn, iinfo )
633  END IF
634 *
635 * If matrix was scaled, then rescale eigenvalues appropriately.
636 *
637 * Jump here if SSTEMR/SSTEIN succeeded.
638  30 CONTINUE
639  IF( iscale.EQ.1 ) THEN
640  IF( info.EQ.0 ) THEN
641  imax = m
642  ELSE
643  imax = info - 1
644  END IF
645  CALL sscal( imax, one / sigma, w, 1 )
646  END IF
647 *
648 * If eigenvalues are not in order, then sort them, along with
649 * eigenvectors. Note: We do not sort the IFAIL portion of IWORK.
650 * It may not be initialized (if SSTEMR/SSTEIN succeeded), and we do
651 * not return this detailed information to the user.
652 *
653  IF( wantz ) THEN
654  DO 50 j = 1, m - 1
655  i = 0
656  tmp1 = w( j )
657  DO 40 jj = j + 1, m
658  IF( w( jj ).LT.tmp1 ) THEN
659  i = jj
660  tmp1 = w( jj )
661  END IF
662  40 CONTINUE
663 *
664  IF( i.NE.0 ) THEN
665  w( i ) = w( j )
666  w( j ) = tmp1
667  CALL sswap( n, z( 1, i ), 1, z( 1, j ), 1 )
668  END IF
669  50 CONTINUE
670  END IF
671 *
672 * Set WORK(1) to optimal workspace size.
673 *
674  work( 1 ) = lwkopt
675  iwork( 1 ) = liwmin
676 *
677  RETURN
678 *
679 * End of SSYEVR
680 *
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
subroutine sormtr(SIDE, UPLO, TRANS, M, N, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
SORMTR
Definition: sormtr.f:174
subroutine scopy(N, SX, INCX, SY, INCY)
SCOPY
Definition: scopy.f:53
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine sstein(N, D, E, M, W, IBLOCK, ISPLIT, Z, LDZ, WORK, IWORK, IFAIL, INFO)
SSTEIN
Definition: sstein.f:176
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine ssytrd(UPLO, N, A, LDA, D, E, TAU, WORK, LWORK, INFO)
SSYTRD
Definition: ssytrd.f:194
subroutine sswap(N, SX, INCX, SY, INCY)
SSWAP
Definition: sswap.f:53
subroutine sstebz(RANGE, ORDER, N, VL, VU, IL, IU, ABSTOL, D, E, M, NSPLIT, W, IBLOCK, ISPLIT, WORK, IWORK, INFO)
SSTEBZ
Definition: sstebz.f:265
subroutine sscal(N, SA, SX, INCX)
SSCAL
Definition: sscal.f:55
real function slansy(NORM, UPLO, N, A, LDA, WORK)
SLANSY returns the value of the 1-norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a real symmetric matrix.
Definition: slansy.f:124
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
subroutine sstemr(JOBZ, RANGE, N, D, E, VL, VU, IL, IU, M, W, Z, LDZ, NZC, ISUPPZ, TRYRAC, WORK, LWORK, IWORK, LIWORK, INFO)
SSTEMR
Definition: sstemr.f:314
subroutine ssterf(N, D, E, INFO)
SSTERF
Definition: ssterf.f:88

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine ssyevx ( character  JOBZ,
character  RANGE,
character  UPLO,
integer  N,
real, dimension( lda, * )  A,
integer  LDA,
real  VL,
real  VU,
integer  IL,
integer  IU,
real  ABSTOL,
integer  M,
real, dimension( * )  W,
real, dimension( ldz, * )  Z,
integer  LDZ,
real, dimension( * )  WORK,
integer  LWORK,
integer, dimension( * )  IWORK,
integer, dimension( * )  IFAIL,
integer  INFO 
)

SSYEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for SY matrices

Download SSYEVX + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 SSYEVX computes selected eigenvalues and, optionally, eigenvectors
 of a real symmetric matrix A.  Eigenvalues and eigenvectors can be
 selected by specifying either a range of values or a range of indices
 for the desired eigenvalues.
Parameters
[in]JOBZ
          JOBZ is CHARACTER*1
          = 'N':  Compute eigenvalues only;
          = 'V':  Compute eigenvalues and eigenvectors.
[in]RANGE
          RANGE is CHARACTER*1
          = 'A': all eigenvalues will be found.
          = 'V': all eigenvalues in the half-open interval (VL,VU]
                 will be found.
          = 'I': the IL-th through IU-th eigenvalues will be found.
[in]UPLO
          UPLO is CHARACTER*1
          = 'U':  Upper triangle of A is stored;
          = 'L':  Lower triangle of A is stored.
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in,out]A
          A is REAL array, dimension (LDA, N)
          On entry, the symmetric matrix A.  If UPLO = 'U', the
          leading N-by-N upper triangular part of A contains the
          upper triangular part of the matrix A.  If UPLO = 'L',
          the leading N-by-N lower triangular part of A contains
          the lower triangular part of the matrix A.
          On exit, the lower triangle (if UPLO='L') or the upper
          triangle (if UPLO='U') of A, including the diagonal, is
          destroyed.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[in]VL
          VL is REAL
[in]VU
          VU is REAL
          If RANGE='V', the lower and upper bounds of the interval to
          be searched for eigenvalues. VL < VU.
          Not referenced if RANGE = 'A' or 'I'.
[in]IL
          IL is INTEGER
[in]IU
          IU is INTEGER
          If RANGE='I', the indices (in ascending order) of the
          smallest and largest eigenvalues to be returned.
          1 <= IL <= IU <= N, if N > 0; IL = 1 and IU = 0 if N = 0.
          Not referenced if RANGE = 'A' or 'V'.
[in]ABSTOL
          ABSTOL is REAL
          The absolute error tolerance for the eigenvalues.
          An approximate eigenvalue is accepted as converged
          when it is determined to lie in an interval [a,b]
          of width less than or equal to

                  ABSTOL + EPS *   max( |a|,|b| ) ,

          where EPS is the machine precision.  If ABSTOL is less than
          or equal to zero, then  EPS*|T|  will be used in its place,
          where |T| is the 1-norm of the tridiagonal matrix obtained
          by reducing A to tridiagonal form.

          Eigenvalues will be computed most accurately when ABSTOL is
          set to twice the underflow threshold 2*SLAMCH('S'), not zero.
          If this routine returns with INFO>0, indicating that some
          eigenvectors did not converge, try setting ABSTOL to
          2*SLAMCH('S').

          See "Computing Small Singular Values of Bidiagonal Matrices
          with Guaranteed High Relative Accuracy," by Demmel and
          Kahan, LAPACK Working Note #3.
[out]M
          M is INTEGER
          The total number of eigenvalues found.  0 <= M <= N.
          If RANGE = 'A', M = N, and if RANGE = 'I', M = IU-IL+1.
[out]W
          W is REAL array, dimension (N)
          On normal exit, the first M elements contain the selected
          eigenvalues in ascending order.
[out]Z
          Z is REAL array, dimension (LDZ, max(1,M))
          If JOBZ = 'V', then if INFO = 0, the first M columns of Z
          contain the orthonormal eigenvectors of the matrix A
          corresponding to the selected eigenvalues, with the i-th
          column of Z holding the eigenvector associated with W(i).
          If an eigenvector fails to converge, then that column of Z
          contains the latest approximation to the eigenvector, and the
          index of the eigenvector is returned in IFAIL.
          If JOBZ = 'N', then Z is not referenced.
          Note: the user must ensure that at least max(1,M) columns are
          supplied in the array Z; if RANGE = 'V', the exact value of M
          is not known in advance and an upper bound must be used.
[in]LDZ
          LDZ is INTEGER
          The leading dimension of the array Z.  LDZ >= 1, and if
          JOBZ = 'V', LDZ >= max(1,N).
[out]WORK
          WORK is REAL array, dimension (MAX(1,LWORK))
          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
[in]LWORK
          LWORK is INTEGER
          The length of the array WORK.  LWORK >= 1, when N <= 1;
          otherwise 8*N.
          For optimal efficiency, LWORK >= (NB+3)*N,
          where NB is the max of the blocksize for SSYTRD and SORMTR
          returned by ILAENV.

          If LWORK = -1, then a workspace query is assumed; the routine
          only calculates the optimal size of the WORK array, returns
          this value as the first entry of the WORK array, and no error
          message related to LWORK is issued by XERBLA.
[out]IWORK
          IWORK is INTEGER array, dimension (5*N)
[out]IFAIL
          IFAIL is INTEGER array, dimension (N)
          If JOBZ = 'V', then if INFO = 0, the first M elements of
          IFAIL are zero.  If INFO > 0, then IFAIL contains the
          indices of the eigenvectors that failed to converge.
          If JOBZ = 'N', then IFAIL is not referenced.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  if INFO = i, then i eigenvectors failed to converge.
                Their indices are stored in array IFAIL.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011

Definition at line 248 of file ssyevx.f.

248 *
249 * -- LAPACK driver routine (version 3.4.0) --
250 * -- LAPACK is a software package provided by Univ. of Tennessee, --
251 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
252 * November 2011
253 *
254 * .. Scalar Arguments ..
255  CHARACTER jobz, range, uplo
256  INTEGER il, info, iu, lda, ldz, lwork, m, n
257  REAL abstol, vl, vu
258 * ..
259 * .. Array Arguments ..
260  INTEGER ifail( * ), iwork( * )
261  REAL a( lda, * ), w( * ), work( * ), z( ldz, * )
262 * ..
263 *
264 * =====================================================================
265 *
266 * .. Parameters ..
267  REAL zero, one
268  parameter( zero = 0.0e+0, one = 1.0e+0 )
269 * ..
270 * .. Local Scalars ..
271  LOGICAL alleig, indeig, lower, lquery, test, valeig,
272  $ wantz
273  CHARACTER order
274  INTEGER i, iinfo, imax, indd, inde, indee, indibl,
275  $ indisp, indiwo, indtau, indwkn, indwrk, iscale,
276  $ itmp1, j, jj, llwork, llwrkn, lwkmin,
277  $ lwkopt, nb, nsplit
278  REAL abstll, anrm, bignum, eps, rmax, rmin, safmin,
279  $ sigma, smlnum, tmp1, vll, vuu
280 * ..
281 * .. External Functions ..
282  LOGICAL lsame
283  INTEGER ilaenv
284  REAL slamch, slansy
285  EXTERNAL lsame, ilaenv, slamch, slansy
286 * ..
287 * .. External Subroutines ..
288  EXTERNAL scopy, slacpy, sorgtr, sormtr, sscal, sstebz,
290 * ..
291 * .. Intrinsic Functions ..
292  INTRINSIC max, min, sqrt
293 * ..
294 * .. Executable Statements ..
295 *
296 * Test the input parameters.
297 *
298  lower = lsame( uplo, 'L' )
299  wantz = lsame( jobz, 'V' )
300  alleig = lsame( range, 'A' )
301  valeig = lsame( range, 'V' )
302  indeig = lsame( range, 'I' )
303  lquery = ( lwork.EQ.-1 )
304 *
305  info = 0
306  IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
307  info = -1
308  ELSE IF( .NOT.( alleig .OR. valeig .OR. indeig ) ) THEN
309  info = -2
310  ELSE IF( .NOT.( lower .OR. lsame( uplo, 'U' ) ) ) THEN
311  info = -3
312  ELSE IF( n.LT.0 ) THEN
313  info = -4
314  ELSE IF( lda.LT.max( 1, n ) ) THEN
315  info = -6
316  ELSE
317  IF( valeig ) THEN
318  IF( n.GT.0 .AND. vu.LE.vl )
319  $ info = -8
320  ELSE IF( indeig ) THEN
321  IF( il.LT.1 .OR. il.GT.max( 1, n ) ) THEN
322  info = -9
323  ELSE IF( iu.LT.min( n, il ) .OR. iu.GT.n ) THEN
324  info = -10
325  END IF
326  END IF
327  END IF
328  IF( info.EQ.0 ) THEN
329  IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) ) THEN
330  info = -15
331  END IF
332  END IF
333 *
334  IF( info.EQ.0 ) THEN
335  IF( n.LE.1 ) THEN
336  lwkmin = 1
337  work( 1 ) = lwkmin
338  ELSE
339  lwkmin = 8*n
340  nb = ilaenv( 1, 'SSYTRD', uplo, n, -1, -1, -1 )
341  nb = max( nb, ilaenv( 1, 'SORMTR', uplo, n, -1, -1, -1 ) )
342  lwkopt = max( lwkmin, ( nb + 3 )*n )
343  work( 1 ) = lwkopt
344  END IF
345 *
346  IF( lwork.LT.lwkmin .AND. .NOT.lquery )
347  $ info = -17
348  END IF
349 *
350  IF( info.NE.0 ) THEN
351  CALL xerbla( 'SSYEVX', -info )
352  RETURN
353  ELSE IF( lquery ) THEN
354  RETURN
355  END IF
356 *
357 * Quick return if possible
358 *
359  m = 0
360  IF( n.EQ.0 ) THEN
361  RETURN
362  END IF
363 *
364  IF( n.EQ.1 ) THEN
365  IF( alleig .OR. indeig ) THEN
366  m = 1
367  w( 1 ) = a( 1, 1 )
368  ELSE
369  IF( vl.LT.a( 1, 1 ) .AND. vu.GE.a( 1, 1 ) ) THEN
370  m = 1
371  w( 1 ) = a( 1, 1 )
372  END IF
373  END IF
374  IF( wantz )
375  $ z( 1, 1 ) = one
376  RETURN
377  END IF
378 *
379 * Get machine constants.
380 *
381  safmin = slamch( 'Safe minimum' )
382  eps = slamch( 'Precision' )
383  smlnum = safmin / eps
384  bignum = one / smlnum
385  rmin = sqrt( smlnum )
386  rmax = min( sqrt( bignum ), one / sqrt( sqrt( safmin ) ) )
387 *
388 * Scale matrix to allowable range, if necessary.
389 *
390  iscale = 0
391  abstll = abstol
392  IF( valeig ) THEN
393  vll = vl
394  vuu = vu
395  END IF
396  anrm = slansy( 'M', uplo, n, a, lda, work )
397  IF( anrm.GT.zero .AND. anrm.LT.rmin ) THEN
398  iscale = 1
399  sigma = rmin / anrm
400  ELSE IF( anrm.GT.rmax ) THEN
401  iscale = 1
402  sigma = rmax / anrm
403  END IF
404  IF( iscale.EQ.1 ) THEN
405  IF( lower ) THEN
406  DO 10 j = 1, n
407  CALL sscal( n-j+1, sigma, a( j, j ), 1 )
408  10 CONTINUE
409  ELSE
410  DO 20 j = 1, n
411  CALL sscal( j, sigma, a( 1, j ), 1 )
412  20 CONTINUE
413  END IF
414  IF( abstol.GT.0 )
415  $ abstll = abstol*sigma
416  IF( valeig ) THEN
417  vll = vl*sigma
418  vuu = vu*sigma
419  END IF
420  END IF
421 *
422 * Call SSYTRD to reduce symmetric matrix to tridiagonal form.
423 *
424  indtau = 1
425  inde = indtau + n
426  indd = inde + n
427  indwrk = indd + n
428  llwork = lwork - indwrk + 1
429  CALL ssytrd( uplo, n, a, lda, work( indd ), work( inde ),
430  $ work( indtau ), work( indwrk ), llwork, iinfo )
431 *
432 * If all eigenvalues are desired and ABSTOL is less than or equal to
433 * zero, then call SSTERF or SORGTR and SSTEQR. If this fails for
434 * some eigenvalue, then try SSTEBZ.
435 *
436  test = .false.
437  IF( indeig ) THEN
438  IF( il.EQ.1 .AND. iu.EQ.n ) THEN
439  test = .true.
440  END IF
441  END IF
442  IF( ( alleig .OR. test ) .AND. ( abstol.LE.zero ) ) THEN
443  CALL scopy( n, work( indd ), 1, w, 1 )
444  indee = indwrk + 2*n
445  IF( .NOT.wantz ) THEN
446  CALL scopy( n-1, work( inde ), 1, work( indee ), 1 )
447  CALL ssterf( n, w, work( indee ), info )
448  ELSE
449  CALL slacpy( 'A', n, n, a, lda, z, ldz )
450  CALL sorgtr( uplo, n, z, ldz, work( indtau ),
451  $ work( indwrk ), llwork, iinfo )
452  CALL scopy( n-1, work( inde ), 1, work( indee ), 1 )
453  CALL ssteqr( jobz, n, w, work( indee ), z, ldz,
454  $ work( indwrk ), info )
455  IF( info.EQ.0 ) THEN
456  DO 30 i = 1, n
457  ifail( i ) = 0
458  30 CONTINUE
459  END IF
460  END IF
461  IF( info.EQ.0 ) THEN
462  m = n
463  GO TO 40
464  END IF
465  info = 0
466  END IF
467 *
468 * Otherwise, call SSTEBZ and, if eigenvectors are desired, SSTEIN.
469 *
470  IF( wantz ) THEN
471  order = 'B'
472  ELSE
473  order = 'E'
474  END IF
475  indibl = 1
476  indisp = indibl + n
477  indiwo = indisp + n
478  CALL sstebz( range, order, n, vll, vuu, il, iu, abstll,
479  $ work( indd ), work( inde ), m, nsplit, w,
480  $ iwork( indibl ), iwork( indisp ), work( indwrk ),
481  $ iwork( indiwo ), info )
482 *
483  IF( wantz ) THEN
484  CALL sstein( n, work( indd ), work( inde ), m, w,
485  $ iwork( indibl ), iwork( indisp ), z, ldz,
486  $ work( indwrk ), iwork( indiwo ), ifail, info )
487 *
488 * Apply orthogonal matrix used in reduction to tridiagonal
489 * form to eigenvectors returned by SSTEIN.
490 *
491  indwkn = inde
492  llwrkn = lwork - indwkn + 1
493  CALL sormtr( 'L', uplo, 'N', n, m, a, lda, work( indtau ), z,
494  $ ldz, work( indwkn ), llwrkn, iinfo )
495  END IF
496 *
497 * If matrix was scaled, then rescale eigenvalues appropriately.
498 *
499  40 CONTINUE
500  IF( iscale.EQ.1 ) THEN
501  IF( info.EQ.0 ) THEN
502  imax = m
503  ELSE
504  imax = info - 1
505  END IF
506  CALL sscal( imax, one / sigma, w, 1 )
507  END IF
508 *
509 * If eigenvalues are not in order, then sort them, along with
510 * eigenvectors.
511 *
512  IF( wantz ) THEN
513  DO 60 j = 1, m - 1
514  i = 0
515  tmp1 = w( j )
516  DO 50 jj = j + 1, m
517  IF( w( jj ).LT.tmp1 ) THEN
518  i = jj
519  tmp1 = w( jj )
520  END IF
521  50 CONTINUE
522 *
523  IF( i.NE.0 ) THEN
524  itmp1 = iwork( indibl+i-1 )
525  w( i ) = w( j )
526  iwork( indibl+i-1 ) = iwork( indibl+j-1 )
527  w( j ) = tmp1
528  iwork( indibl+j-1 ) = itmp1
529  CALL sswap( n, z( 1, i ), 1, z( 1, j ), 1 )
530  IF( info.NE.0 ) THEN
531  itmp1 = ifail( i )
532  ifail( i ) = ifail( j )
533  ifail( j ) = itmp1
534  END IF
535  END IF
536  60 CONTINUE
537  END IF
538 *
539 * Set WORK(1) to optimal workspace size.
540 *
541  work( 1 ) = lwkopt
542 *
543  RETURN
544 *
545 * End of SSYEVX
546 *
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
subroutine sormtr(SIDE, UPLO, TRANS, M, N, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
SORMTR
Definition: sormtr.f:174
subroutine sorgtr(UPLO, N, A, LDA, TAU, WORK, LWORK, INFO)
SORGTR
Definition: sorgtr.f:125
subroutine scopy(N, SX, INCX, SY, INCY)
SCOPY
Definition: scopy.f:53
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine sstein(N, D, E, M, W, IBLOCK, ISPLIT, Z, LDZ, WORK, IWORK, IFAIL, INFO)
SSTEIN
Definition: sstein.f:176
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine ssytrd(UPLO, N, A, LDA, D, E, TAU, WORK, LWORK, INFO)
SSYTRD
Definition: ssytrd.f:194
subroutine sswap(N, SX, INCX, SY, INCY)
SSWAP
Definition: sswap.f:53
subroutine sstebz(RANGE, ORDER, N, VL, VU, IL, IU, ABSTOL, D, E, M, NSPLIT, W, IBLOCK, ISPLIT, WORK, IWORK, INFO)
SSTEBZ
Definition: sstebz.f:265
subroutine sscal(N, SA, SX, INCX)
SSCAL
Definition: sscal.f:55
subroutine slacpy(UPLO, M, N, A, LDA, B, LDB)
SLACPY copies all or part of one two-dimensional array to another.
Definition: slacpy.f:105
real function slansy(NORM, UPLO, N, A, LDA, WORK)
SLANSY returns the value of the 1-norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a real symmetric matrix.
Definition: slansy.f:124
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
subroutine ssterf(N, D, E, INFO)
SSTERF
Definition: ssterf.f:88
subroutine ssteqr(COMPZ, N, D, E, Z, LDZ, WORK, INFO)
SSTEQR
Definition: ssteqr.f:133

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine ssygv ( integer  ITYPE,
character  JOBZ,
character  UPLO,
integer  N,
real, dimension( lda, * )  A,
integer  LDA,
real, dimension( ldb, * )  B,
integer  LDB,
real, dimension( * )  W,
real, dimension( * )  WORK,
integer  LWORK,
integer  INFO 
)

SSYGV

Download SSYGV + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 SSYGV computes all the eigenvalues, and optionally, the eigenvectors
 of a real generalized symmetric-definite eigenproblem, of the form
 A*x=(lambda)*B*x,  A*Bx=(lambda)*x,  or B*A*x=(lambda)*x.
 Here A and B are assumed to be symmetric and B is also
 positive definite.
Parameters
[in]ITYPE
          ITYPE is INTEGER
          Specifies the problem type to be solved:
          = 1:  A*x = (lambda)*B*x
          = 2:  A*B*x = (lambda)*x
          = 3:  B*A*x = (lambda)*x
[in]JOBZ
          JOBZ is CHARACTER*1
          = 'N':  Compute eigenvalues only;
          = 'V':  Compute eigenvalues and eigenvectors.
[in]UPLO
          UPLO is CHARACTER*1
          = 'U':  Upper triangles of A and B are stored;
          = 'L':  Lower triangles of A and B are stored.
[in]N
          N is INTEGER
          The order of the matrices A and B.  N >= 0.
[in,out]A
          A is REAL array, dimension (LDA, N)
          On entry, the symmetric matrix A.  If UPLO = 'U', the
          leading N-by-N upper triangular part of A contains the
          upper triangular part of the matrix A.  If UPLO = 'L',
          the leading N-by-N lower triangular part of A contains
          the lower triangular part of the matrix A.

          On exit, if JOBZ = 'V', then if INFO = 0, A contains the
          matrix Z of eigenvectors.  The eigenvectors are normalized
          as follows:
          if ITYPE = 1 or 2, Z**T*B*Z = I;
          if ITYPE = 3, Z**T*inv(B)*Z = I.
          If JOBZ = 'N', then on exit the upper triangle (if UPLO='U')
          or the lower triangle (if UPLO='L') of A, including the
          diagonal, is destroyed.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[in,out]B
          B is REAL array, dimension (LDB, N)
          On entry, the symmetric positive definite matrix B.
          If UPLO = 'U', the leading N-by-N upper triangular part of B
          contains the upper triangular part of the matrix B.
          If UPLO = 'L', the leading N-by-N lower triangular part of B
          contains the lower triangular part of the matrix B.

          On exit, if INFO <= N, the part of B containing the matrix is
          overwritten by the triangular factor U or L from the Cholesky
          factorization B = U**T*U or B = L*L**T.
[in]LDB
          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).
[out]W
          W is REAL array, dimension (N)
          If INFO = 0, the eigenvalues in ascending order.
[out]WORK
          WORK is REAL array, dimension (MAX(1,LWORK))
          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
[in]LWORK
          LWORK is INTEGER
          The length of the array WORK.  LWORK >= max(1,3*N-1).
          For optimal efficiency, LWORK >= (NB+2)*N,
          where NB is the blocksize for SSYTRD returned by ILAENV.

          If LWORK = -1, then a workspace query is assumed; the routine
          only calculates the optimal size of the WORK array, returns
          this value as the first entry of the WORK array, and no error
          message related to LWORK is issued by XERBLA.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  SPOTRF or SSYEV returned an error code:
             <= N:  if INFO = i, SSYEV failed to converge;
                    i off-diagonal elements of an intermediate
                    tridiagonal form did not converge to zero;
             > N:   if INFO = N + i, for 1 <= i <= N, then the leading
                    minor of order i of B is not positive definite.
                    The factorization of B could not be completed and
                    no eigenvalues or eigenvectors were computed.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2015

Definition at line 177 of file ssygv.f.

177 *
178 * -- LAPACK driver routine (version 3.6.0) --
179 * -- LAPACK is a software package provided by Univ. of Tennessee, --
180 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
181 * November 2015
182 *
183 * .. Scalar Arguments ..
184  CHARACTER jobz, uplo
185  INTEGER info, itype, lda, ldb, lwork, n
186 * ..
187 * .. Array Arguments ..
188  REAL a( lda, * ), b( ldb, * ), w( * ), work( * )
189 * ..
190 *
191 * =====================================================================
192 *
193 * .. Parameters ..
194  REAL one
195  parameter( one = 1.0e+0 )
196 * ..
197 * .. Local Scalars ..
198  LOGICAL lquery, upper, wantz
199  CHARACTER trans
200  INTEGER lwkmin, lwkopt, nb, neig
201 * ..
202 * .. External Functions ..
203  LOGICAL lsame
204  INTEGER ilaenv
205  EXTERNAL ilaenv, lsame
206 * ..
207 * .. External Subroutines ..
208  EXTERNAL spotrf, ssyev, ssygst, strmm, strsm, xerbla
209 * ..
210 * .. Intrinsic Functions ..
211  INTRINSIC max
212 * ..
213 * .. Executable Statements ..
214 *
215 * Test the input parameters.
216 *
217  wantz = lsame( jobz, 'V' )
218  upper = lsame( uplo, 'U' )
219  lquery = ( lwork.EQ.-1 )
220 *
221  info = 0
222  IF( itype.LT.1 .OR. itype.GT.3 ) THEN
223  info = -1
224  ELSE IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
225  info = -2
226  ELSE IF( .NOT.( upper .OR. lsame( uplo, 'L' ) ) ) THEN
227  info = -3
228  ELSE IF( n.LT.0 ) THEN
229  info = -4
230  ELSE IF( lda.LT.max( 1, n ) ) THEN
231  info = -6
232  ELSE IF( ldb.LT.max( 1, n ) ) THEN
233  info = -8
234  END IF
235 *
236  IF( info.EQ.0 ) THEN
237  lwkmin = max( 1, 3*n - 1 )
238  nb = ilaenv( 1, 'SSYTRD', uplo, n, -1, -1, -1 )
239  lwkopt = max( lwkmin, ( nb + 2 )*n )
240  work( 1 ) = lwkopt
241 *
242  IF( lwork.LT.lwkmin .AND. .NOT.lquery ) THEN
243  info = -11
244  END IF
245  END IF
246 *
247  IF( info.NE.0 ) THEN
248  CALL xerbla( 'SSYGV ', -info )
249  RETURN
250  ELSE IF( lquery ) THEN
251  RETURN
252  END IF
253 *
254 * Quick return if possible
255 *
256  IF( n.EQ.0 )
257  $ RETURN
258 *
259 * Form a Cholesky factorization of B.
260 *
261  CALL spotrf( uplo, n, b, ldb, info )
262  IF( info.NE.0 ) THEN
263  info = n + info
264  RETURN
265  END IF
266 *
267 * Transform problem to standard eigenvalue problem and solve.
268 *
269  CALL ssygst( itype, uplo, n, a, lda, b, ldb, info )
270  CALL ssyev( jobz, uplo, n, a, lda, w, work, lwork, info )
271 *
272  IF( wantz ) THEN
273 *
274 * Backtransform eigenvectors to the original problem.
275 *
276  neig = n
277  IF( info.GT.0 )
278  $ neig = info - 1
279  IF( itype.EQ.1 .OR. itype.EQ.2 ) THEN
280 *
281 * For A*x=(lambda)*B*x and A*B*x=(lambda)*x;
282 * backtransform eigenvectors: x = inv(L)**T*y or inv(U)*y
283 *
284  IF( upper ) THEN
285  trans = 'N'
286  ELSE
287  trans = 'T'
288  END IF
289 *
290  CALL strsm( 'Left', uplo, trans, 'Non-unit', n, neig, one,
291  $ b, ldb, a, lda )
292 *
293  ELSE IF( itype.EQ.3 ) THEN
294 *
295 * For B*A*x=(lambda)*x;
296 * backtransform eigenvectors: x = L*y or U**T*y
297 *
298  IF( upper ) THEN
299  trans = 'T'
300  ELSE
301  trans = 'N'
302  END IF
303 *
304  CALL strmm( 'Left', uplo, trans, 'Non-unit', n, neig, one,
305  $ b, ldb, a, lda )
306  END IF
307  END IF
308 *
309  work( 1 ) = lwkopt
310  RETURN
311 *
312 * End of SSYGV
313 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine ssygst(ITYPE, UPLO, N, A, LDA, B, LDB, INFO)
SSYGST
Definition: ssygst.f:129
subroutine spotrf(UPLO, N, A, LDA, INFO)
SPOTRF
Definition: spotrf.f:109
subroutine strmm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
STRMM
Definition: strmm.f:179
subroutine strsm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
STRSM
Definition: strsm.f:183
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:134
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine ssygvd ( integer  ITYPE,
character  JOBZ,
character  UPLO,
integer  N,
real, dimension( lda, * )  A,
integer  LDA,
real, dimension( ldb, * )  B,
integer  LDB,
real, dimension( * )  W,
real, dimension( * )  WORK,
integer  LWORK,
integer, dimension( * )  IWORK,
integer  LIWORK,
integer  INFO 
)

SSYGVD

Download SSYGVD + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 SSYGVD computes all the eigenvalues, and optionally, the eigenvectors
 of a real generalized symmetric-definite eigenproblem, of the form
 A*x=(lambda)*B*x,  A*Bx=(lambda)*x,  or B*A*x=(lambda)*x.  Here A and
 B are assumed to be symmetric and B is also positive definite.
 If eigenvectors are desired, it uses a divide and conquer algorithm.

 The divide and conquer algorithm makes very mild assumptions about
 floating point arithmetic. It will work on machines with a guard
 digit in add/subtract, or on those binary machines without guard
 digits which subtract like the Cray X-MP, Cray Y-MP, Cray C-90, or
 Cray-2. It could conceivably fail on hexadecimal or decimal machines
 without guard digits, but we know of none.
Parameters
[in]ITYPE
          ITYPE is INTEGER
          Specifies the problem type to be solved:
          = 1:  A*x = (lambda)*B*x
          = 2:  A*B*x = (lambda)*x
          = 3:  B*A*x = (lambda)*x
[in]JOBZ
          JOBZ is CHARACTER*1
          = 'N':  Compute eigenvalues only;
          = 'V':  Compute eigenvalues and eigenvectors.
[in]UPLO
          UPLO is CHARACTER*1
          = 'U':  Upper triangles of A and B are stored;
          = 'L':  Lower triangles of A and B are stored.
[in]N
          N is INTEGER
          The order of the matrices A and B.  N >= 0.
[in,out]A
          A is REAL array, dimension (LDA, N)
          On entry, the symmetric matrix A.  If UPLO = 'U', the
          leading N-by-N upper triangular part of A contains the
          upper triangular part of the matrix A.  If UPLO = 'L',
          the leading N-by-N lower triangular part of A contains
          the lower triangular part of the matrix A.

          On exit, if JOBZ = 'V', then if INFO = 0, A contains the
          matrix Z of eigenvectors.  The eigenvectors are normalized
          as follows:
          if ITYPE = 1 or 2, Z**T*B*Z = I;
          if ITYPE = 3, Z**T*inv(B)*Z = I.
          If JOBZ = 'N', then on exit the upper triangle (if UPLO='U')
          or the lower triangle (if UPLO='L') of A, including the
          diagonal, is destroyed.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[in,out]B
          B is REAL array, dimension (LDB, N)
          On entry, the symmetric matrix B.  If UPLO = 'U', the
          leading N-by-N upper triangular part of B contains the
          upper triangular part of the matrix B.  If UPLO = 'L',
          the leading N-by-N lower triangular part of B contains
          the lower triangular part of the matrix B.

          On exit, if INFO <= N, the part of B containing the matrix is
          overwritten by the triangular factor U or L from the Cholesky
          factorization B = U**T*U or B = L*L**T.
[in]LDB
          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).
[out]W
          W is REAL array, dimension (N)
          If INFO = 0, the eigenvalues in ascending order.
[out]WORK
          WORK is REAL array, dimension (MAX(1,LWORK))
          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
[in]LWORK
          LWORK is INTEGER
          The dimension of the array WORK.
          If N <= 1,               LWORK >= 1.
          If JOBZ = 'N' and N > 1, LWORK >= 2*N+1.
          If JOBZ = 'V' and N > 1, LWORK >= 1 + 6*N + 2*N**2.

          If LWORK = -1, then a workspace query is assumed; the routine
          only calculates the optimal sizes of the WORK and IWORK
          arrays, returns these values as the first entries of the WORK
          and IWORK arrays, and no error message related to LWORK or
          LIWORK is issued by XERBLA.
[out]IWORK
          IWORK is INTEGER array, dimension (MAX(1,LIWORK))
          On exit, if INFO = 0, IWORK(1) returns the optimal LIWORK.
[in]LIWORK
          LIWORK is INTEGER
          The dimension of the array IWORK.
          If N <= 1,                LIWORK >= 1.
          If JOBZ  = 'N' and N > 1, LIWORK >= 1.
          If JOBZ  = 'V' and N > 1, LIWORK >= 3 + 5*N.

          If LIWORK = -1, then a workspace query is assumed; the
          routine only calculates the optimal sizes of the WORK and
          IWORK arrays, returns these values as the first entries of
          the WORK and IWORK arrays, and no error message related to
          LWORK or LIWORK is issued by XERBLA.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  SPOTRF or SSYEVD returned an error code:
             <= N:  if INFO = i and JOBZ = 'N', then the algorithm
                    failed to converge; i off-diagonal elements of an
                    intermediate tridiagonal form did not converge to
                    zero;
                    if INFO = i and JOBZ = 'V', then the algorithm
                    failed to compute an eigenvalue while working on
                    the submatrix lying in rows and columns INFO/(N+1)
                    through mod(INFO,N+1);
             > N:   if INFO = N + i, for 1 <= i <= N, then the leading
                    minor of order i of B is not positive definite.
                    The factorization of B could not be completed and
                    no eigenvalues or eigenvectors were computed.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2015
Further Details:
  Modified so that no backsubstitution is performed if SSYEVD fails to
  converge (NEIG in old code could be greater than N causing out of
  bounds reference to A - reported by Ralf Meyer).  Also corrected the
  description of INFO and the test on ITYPE. Sven, 16 Feb 05.
Contributors:
Mark Fahey, Department of Mathematics, Univ. of Kentucky, USA

Definition at line 229 of file ssygvd.f.

229 *
230 * -- LAPACK driver routine (version 3.6.0) --
231 * -- LAPACK is a software package provided by Univ. of Tennessee, --
232 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
233 * November 2015
234 *
235 * .. Scalar Arguments ..
236  CHARACTER jobz, uplo
237  INTEGER info, itype, lda, ldb, liwork, lwork, n
238 * ..
239 * .. Array Arguments ..
240  INTEGER iwork( * )
241  REAL a( lda, * ), b( ldb, * ), w( * ), work( * )
242 * ..
243 *
244 * =====================================================================
245 *
246 * .. Parameters ..
247  REAL one
248  parameter( one = 1.0e+0 )
249 * ..
250 * .. Local Scalars ..
251  LOGICAL lquery, upper, wantz
252  CHARACTER trans
253  INTEGER liopt, liwmin, lopt, lwmin
254 * ..
255 * .. External Functions ..
256  LOGICAL lsame
257  EXTERNAL lsame
258 * ..
259 * .. External Subroutines ..
260  EXTERNAL spotrf, ssyevd, ssygst, strmm, strsm, xerbla
261 * ..
262 * .. Intrinsic Functions ..
263  INTRINSIC max, real
264 * ..
265 * .. Executable Statements ..
266 *
267 * Test the input parameters.
268 *
269  wantz = lsame( jobz, 'V' )
270  upper = lsame( uplo, 'U' )
271  lquery = ( lwork.EQ.-1 .OR. liwork.EQ.-1 )
272 *
273  info = 0
274  IF( n.LE.1 ) THEN
275  liwmin = 1
276  lwmin = 1
277  ELSE IF( wantz ) THEN
278  liwmin = 3 + 5*n
279  lwmin = 1 + 6*n + 2*n**2
280  ELSE
281  liwmin = 1
282  lwmin = 2*n + 1
283  END IF
284  lopt = lwmin
285  liopt = liwmin
286  IF( itype.LT.1 .OR. itype.GT.3 ) THEN
287  info = -1
288  ELSE IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
289  info = -2
290  ELSE IF( .NOT.( upper .OR. lsame( uplo, 'L' ) ) ) THEN
291  info = -3
292  ELSE IF( n.LT.0 ) THEN
293  info = -4
294  ELSE IF( lda.LT.max( 1, n ) ) THEN
295  info = -6
296  ELSE IF( ldb.LT.max( 1, n ) ) THEN
297  info = -8
298  END IF
299 *
300  IF( info.EQ.0 ) THEN
301  work( 1 ) = lopt
302  iwork( 1 ) = liopt
303 *
304  IF( lwork.LT.lwmin .AND. .NOT.lquery ) THEN
305  info = -11
306  ELSE IF( liwork.LT.liwmin .AND. .NOT.lquery ) THEN
307  info = -13
308  END IF
309  END IF
310 *
311  IF( info.NE.0 ) THEN
312  CALL xerbla( 'SSYGVD', -info )
313  RETURN
314  ELSE IF( lquery ) THEN
315  RETURN
316  END IF
317 *
318 * Quick return if possible
319 *
320  IF( n.EQ.0 )
321  $ RETURN
322 *
323 * Form a Cholesky factorization of B.
324 *
325  CALL spotrf( uplo, n, b, ldb, info )
326  IF( info.NE.0 ) THEN
327  info = n + info
328  RETURN
329  END IF
330 *
331 * Transform problem to standard eigenvalue problem and solve.
332 *
333  CALL ssygst( itype, uplo, n, a, lda, b, ldb, info )
334  CALL ssyevd( jobz, uplo, n, a, lda, w, work, lwork, iwork, liwork,
335  $ info )
336  lopt = max( REAL( LOPT ), REAL( WORK( 1 ) ) )
337  liopt = max( REAL( LIOPT ), REAL( IWORK( 1 ) ) )
338 *
339  IF( wantz .AND. info.EQ.0 ) THEN
340 *
341 * Backtransform eigenvectors to the original problem.
342 *
343  IF( itype.EQ.1 .OR. itype.EQ.2 ) THEN
344 *
345 * For A*x=(lambda)*B*x and A*B*x=(lambda)*x;
346 * backtransform eigenvectors: x = inv(L)**T*y or inv(U)*y
347 *
348  IF( upper ) THEN
349  trans = 'N'
350  ELSE
351  trans = 'T'
352  END IF
353 *
354  CALL strsm( 'Left', uplo, trans, 'Non-unit', n, n, one,
355  $ b, ldb, a, lda )
356 *
357  ELSE IF( itype.EQ.3 ) THEN
358 *
359 * For B*A*x=(lambda)*x;
360 * backtransform eigenvectors: x = L*y or U**T*y
361 *
362  IF( upper ) THEN
363  trans = 'T'
364  ELSE
365  trans = 'N'
366  END IF
367 *
368  CALL strmm( 'Left', uplo, trans, 'Non-unit', n, n, one,
369  $ b, ldb, a, lda )
370  END IF
371  END IF
372 *
373  work( 1 ) = lopt
374  iwork( 1 ) = liopt
375 *
376  RETURN
377 *
378 * End of SSYGVD
379 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine ssygst(ITYPE, UPLO, N, A, LDA, B, LDB, INFO)
SSYGST
Definition: ssygst.f:129
subroutine spotrf(UPLO, N, A, LDA, INFO)
SPOTRF
Definition: spotrf.f:109
subroutine ssyevd(JOBZ, UPLO, N, A, LDA, W, WORK, LWORK, IWORK, LIWORK, INFO)
SSYEVD computes the eigenvalues and, optionally, the left and/or right eigenvectors for SY matrices ...
Definition: ssyevd.f:185
subroutine strmm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
STRMM
Definition: strmm.f:179
subroutine strsm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
STRSM
Definition: strsm.f:183
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine ssygvx ( integer  ITYPE,
character  JOBZ,
character  RANGE,
character  UPLO,
integer  N,
real, dimension( lda, * )  A,
integer  LDA,
real, dimension( ldb, * )  B,
integer  LDB,
real  VL,
real  VU,
integer  IL,
integer  IU,
real  ABSTOL,
integer  M,
real, dimension( * )  W,
real, dimension( ldz, * )  Z,
integer  LDZ,
real, dimension( * )  WORK,
integer  LWORK,
integer, dimension( * )  IWORK,
integer, dimension( * )  IFAIL,
integer  INFO 
)

SSYGVX

Download SSYGVX + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 SSYGVX computes selected eigenvalues, and optionally, eigenvectors
 of a real generalized symmetric-definite eigenproblem, of the form
 A*x=(lambda)*B*x,  A*Bx=(lambda)*x,  or B*A*x=(lambda)*x.  Here A
 and B are assumed to be symmetric and B is also positive definite.
 Eigenvalues and eigenvectors can be selected by specifying either a
 range of values or a range of indices for the desired eigenvalues.
Parameters
[in]ITYPE
          ITYPE is INTEGER
          Specifies the problem type to be solved:
          = 1:  A*x = (lambda)*B*x
          = 2:  A*B*x = (lambda)*x
          = 3:  B*A*x = (lambda)*x
[in]JOBZ
          JOBZ is CHARACTER*1
          = 'N':  Compute eigenvalues only;
          = 'V':  Compute eigenvalues and eigenvectors.
[in]RANGE
          RANGE is CHARACTER*1
          = 'A': all eigenvalues will be found.
          = 'V': all eigenvalues in the half-open interval (VL,VU]
                 will be found.
          = 'I': the IL-th through IU-th eigenvalues will be found.
[in]UPLO
          UPLO is CHARACTER*1
          = 'U':  Upper triangle of A and B are stored;
          = 'L':  Lower triangle of A and B are stored.
[in]N
          N is INTEGER
          The order of the matrix pencil (A,B).  N >= 0.
[in,out]A
          A is REAL array, dimension (LDA, N)
          On entry, the symmetric matrix A.  If UPLO = 'U', the
          leading N-by-N upper triangular part of A contains the
          upper triangular part of the matrix A.  If UPLO = 'L',
          the leading N-by-N lower triangular part of A contains
          the lower triangular part of the matrix A.

          On exit, the lower triangle (if UPLO='L') or the upper
          triangle (if UPLO='U') of A, including the diagonal, is
          destroyed.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[in,out]B
          B is REAL array, dimension (LDA, N)
          On entry, the symmetric matrix B.  If UPLO = 'U', the
          leading N-by-N upper triangular part of B contains the
          upper triangular part of the matrix B.  If UPLO = 'L',
          the leading N-by-N lower triangular part of B contains
          the lower triangular part of the matrix B.

          On exit, if INFO <= N, the part of B containing the matrix is
          overwritten by the triangular factor U or L from the Cholesky
          factorization B = U**T*U or B = L*L**T.
[in]LDB
          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).
[in]VL
          VL is REAL
[in]VU
          VU is REAL
          If RANGE='V', the lower and upper bounds of the interval to
          be searched for eigenvalues. VL < VU.
          Not referenced if RANGE = 'A' or 'I'.
[in]IL
          IL is INTEGER
[in]IU
          IU is INTEGER
          If RANGE='I', the indices (in ascending order) of the
          smallest and largest eigenvalues to be returned.
          1 <= IL <= IU <= N, if N > 0; IL = 1 and IU = 0 if N = 0.
          Not referenced if RANGE = 'A' or 'V'.
[in]ABSTOL
          ABSTOL is REAL
          The absolute error tolerance for the eigenvalues.
          An approximate eigenvalue is accepted as converged
          when it is determined to lie in an interval [a,b]
          of width less than or equal to

                  ABSTOL + EPS *   max( |a|,|b| ) ,

          where EPS is the machine precision.  If ABSTOL is less than
          or equal to zero, then  EPS*|T|  will be used in its place,
          where |T| is the 1-norm of the tridiagonal matrix obtained
          by reducing C to tridiagonal form, where C is the symmetric
          matrix of the standard symmetric problem to which the
          generalized problem is transformed.

          Eigenvalues will be computed most accurately when ABSTOL is
          set to twice the underflow threshold 2*DLAMCH('S'), not zero.
          If this routine returns with INFO>0, indicating that some
          eigenvectors did not converge, try setting ABSTOL to
          2*SLAMCH('S').
[out]M
          M is INTEGER
          The total number of eigenvalues found.  0 <= M <= N.
          If RANGE = 'A', M = N, and if RANGE = 'I', M = IU-IL+1.
[out]W
          W is REAL array, dimension (N)
          On normal exit, the first M elements contain the selected
          eigenvalues in ascending order.
[out]Z
          Z is REAL array, dimension (LDZ, max(1,M))
          If JOBZ = 'N', then Z is not referenced.
          If JOBZ = 'V', then if INFO = 0, the first M columns of Z
          contain the orthonormal eigenvectors of the matrix A
          corresponding to the selected eigenvalues, with the i-th
          column of Z holding the eigenvector associated with W(i).
          The eigenvectors are normalized as follows:
          if ITYPE = 1 or 2, Z**T*B*Z = I;
          if ITYPE = 3, Z**T*inv(B)*Z = I.

          If an eigenvector fails to converge, then that column of Z
          contains the latest approximation to the eigenvector, and the
          index of the eigenvector is returned in IFAIL.
          Note: the user must ensure that at least max(1,M) columns are
          supplied in the array Z; if RANGE = 'V', the exact value of M
          is not known in advance and an upper bound must be used.
[in]LDZ
          LDZ is INTEGER
          The leading dimension of the array Z.  LDZ >= 1, and if
          JOBZ = 'V', LDZ >= max(1,N).
[out]WORK
          WORK is REAL array, dimension (MAX(1,LWORK))
          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
[in]LWORK
          LWORK is INTEGER
          The length of the array WORK.  LWORK >= max(1,8*N).
          For optimal efficiency, LWORK >= (NB+3)*N,
          where NB is the blocksize for SSYTRD returned by ILAENV.

          If LWORK = -1, then a workspace query is assumed; the routine
          only calculates the optimal size of the WORK array, returns
          this value as the first entry of the WORK array, and no error
          message related to LWORK is issued by XERBLA.
[out]IWORK
          IWORK is INTEGER array, dimension (5*N)
[out]IFAIL
          IFAIL is INTEGER array, dimension (N)
          If JOBZ = 'V', then if INFO = 0, the first M elements of
          IFAIL are zero.  If INFO > 0, then IFAIL contains the
          indices of the eigenvectors that failed to converge.
          If JOBZ = 'N', then IFAIL is not referenced.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  SPOTRF or SSYEVX returned an error code:
             <= N:  if INFO = i, SSYEVX failed to converge;
                    i eigenvectors failed to converge.  Their indices
                    are stored in array IFAIL.
             > N:   if INFO = N + i, for 1 <= i <= N, then the leading
                    minor of order i of B is not positive definite.
                    The factorization of B could not be completed and
                    no eigenvalues or eigenvectors were computed.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2015
Contributors:
Mark Fahey, Department of Mathematics, Univ. of Kentucky, USA

Definition at line 292 of file ssygvx.f.

292 *
293 * -- LAPACK driver routine (version 3.6.0) --
294 * -- LAPACK is a software package provided by Univ. of Tennessee, --
295 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
296 * November 2015
297 *
298 * .. Scalar Arguments ..
299  CHARACTER jobz, range, uplo
300  INTEGER il, info, itype, iu, lda, ldb, ldz, lwork, m, n
301  REAL abstol, vl, vu
302 * ..
303 * .. Array Arguments ..
304  INTEGER ifail( * ), iwork( * )
305  REAL a( lda, * ), b( ldb, * ), w( * ), work( * ),
306  $ z( ldz, * )
307 * ..
308 *
309 * =====================================================================
310 *
311 * .. Parameters ..
312  REAL one
313  parameter( one = 1.0e+0 )
314 * ..
315 * .. Local Scalars ..
316  LOGICAL alleig, indeig, lquery, upper, valeig, wantz
317  CHARACTER trans
318  INTEGER lwkmin, lwkopt, nb
319 * ..
320 * .. External Functions ..
321  LOGICAL lsame
322  INTEGER ilaenv
323  EXTERNAL ilaenv, lsame
324 * ..
325 * .. External Subroutines ..
326  EXTERNAL spotrf, ssyevx, ssygst, strmm, strsm, xerbla
327 * ..
328 * .. Intrinsic Functions ..
329  INTRINSIC max, min
330 * ..
331 * .. Executable Statements ..
332 *
333 * Test the input parameters.
334 *
335  upper = lsame( uplo, 'U' )
336  wantz = lsame( jobz, 'V' )
337  alleig = lsame( range, 'A' )
338  valeig = lsame( range, 'V' )
339  indeig = lsame( range, 'I' )
340  lquery = ( lwork.EQ.-1 )
341 *
342  info = 0
343  IF( itype.LT.1 .OR. itype.GT.3 ) THEN
344  info = -1
345  ELSE IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
346  info = -2
347  ELSE IF( .NOT.( alleig .OR. valeig .OR. indeig ) ) THEN
348  info = -3
349  ELSE IF( .NOT.( upper .OR. lsame( uplo, 'L' ) ) ) THEN
350  info = -4
351  ELSE IF( n.LT.0 ) THEN
352  info = -5
353  ELSE IF( lda.LT.max( 1, n ) ) THEN
354  info = -7
355  ELSE IF( ldb.LT.max( 1, n ) ) THEN
356  info = -9
357  ELSE
358  IF( valeig ) THEN
359  IF( n.GT.0 .AND. vu.LE.vl )
360  $ info = -11
361  ELSE IF( indeig ) THEN
362  IF( il.LT.1 .OR. il.GT.max( 1, n ) ) THEN
363  info = -12
364  ELSE IF( iu.LT.min( n, il ) .OR. iu.GT.n ) THEN
365  info = -13
366  END IF
367  END IF
368  END IF
369  IF (info.EQ.0) THEN
370  IF (ldz.LT.1 .OR. (wantz .AND. ldz.LT.n)) THEN
371  info = -18
372  END IF
373  END IF
374 *
375  IF( info.EQ.0 ) THEN
376  lwkmin = max( 1, 8*n )
377  nb = ilaenv( 1, 'SSYTRD', uplo, n, -1, -1, -1 )
378  lwkopt = max( lwkmin, ( nb + 3 )*n )
379  work( 1 ) = lwkopt
380 *
381  IF( lwork.LT.lwkmin .AND. .NOT.lquery ) THEN
382  info = -20
383  END IF
384  END IF
385 *
386  IF( info.NE.0 ) THEN
387  CALL xerbla( 'SSYGVX', -info )
388  RETURN
389  ELSE IF( lquery ) THEN
390  RETURN
391  END IF
392 *
393 * Quick return if possible
394 *
395  m = 0
396  IF( n.EQ.0 ) THEN
397  RETURN
398  END IF
399 *
400 * Form a Cholesky factorization of B.
401 *
402  CALL spotrf( uplo, n, b, ldb, info )
403  IF( info.NE.0 ) THEN
404  info = n + info
405  RETURN
406  END IF
407 *
408 * Transform problem to standard eigenvalue problem and solve.
409 *
410  CALL ssygst( itype, uplo, n, a, lda, b, ldb, info )
411  CALL ssyevx( jobz, range, uplo, n, a, lda, vl, vu, il, iu, abstol,
412  $ m, w, z, ldz, work, lwork, iwork, ifail, info )
413 *
414  IF( wantz ) THEN
415 *
416 * Backtransform eigenvectors to the original problem.
417 *
418  IF( info.GT.0 )
419  $ m = info - 1
420  IF( itype.EQ.1 .OR. itype.EQ.2 ) THEN
421 *
422 * For A*x=(lambda)*B*x and A*B*x=(lambda)*x;
423 * backtransform eigenvectors: x = inv(L)**T*y or inv(U)*y
424 *
425  IF( upper ) THEN
426  trans = 'N'
427  ELSE
428  trans = 'T'
429  END IF
430 *
431  CALL strsm( 'Left', uplo, trans, 'Non-unit', n, m, one, b,
432  $ ldb, z, ldz )
433 *
434  ELSE IF( itype.EQ.3 ) THEN
435 *
436 * For B*A*x=(lambda)*x;
437 * backtransform eigenvectors: x = L*y or U**T*y
438 *
439  IF( upper ) THEN
440  trans = 'T'
441  ELSE
442  trans = 'N'
443  END IF
444 *
445  CALL strmm( 'Left', uplo, trans, 'Non-unit', n, m, one, b,
446  $ ldb, z, ldz )
447  END IF
448  END IF
449 *
450 * Set WORK(1) to optimal workspace size.
451 *
452  work( 1 ) = lwkopt
453 *
454  RETURN
455 *
456 * End of SSYGVX
457 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine ssygst(ITYPE, UPLO, N, A, LDA, B, LDB, INFO)
SSYGST
Definition: ssygst.f:129
subroutine spotrf(UPLO, N, A, LDA, INFO)
SPOTRF
Definition: spotrf.f:109
subroutine strmm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
STRMM
Definition: strmm.f:179
subroutine strsm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
STRSM
Definition: strsm.f:183
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine ssyevx(JOBZ, RANGE, UPLO, N, A, LDA, VL, VU, IL, IU, ABSTOL, M, W, Z, LDZ, WORK, LWORK, IWORK, IFAIL, INFO)
SSYEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for SY matrices ...
Definition: ssyevx.f:248
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83

Here is the call graph for this function:

Here is the caller graph for this function: