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

Functions

subroutine dsyev (JOBZ, UPLO, N, A, LDA, W, WORK, LWORK, INFO)
  DSYEV computes the eigenvalues and, optionally, the left and/or right eigenvectors for SY matrices More...
 
subroutine dsyevd (JOBZ, UPLO, N, A, LDA, W, WORK, LWORK, IWORK, LIWORK, INFO)
  DSYEVD computes the eigenvalues and, optionally, the left and/or right eigenvectors for SY matrices More...
 
subroutine dsyevr (JOBZ, RANGE, UPLO, N, A, LDA, VL, VU, IL, IU, ABSTOL, M, W, Z, LDZ, ISUPPZ, WORK, LWORK, IWORK, LIWORK, INFO)
  DSYEVR computes the eigenvalues and, optionally, the left and/or right eigenvectors for SY matrices More...
 
subroutine dsyevx (JOBZ, RANGE, UPLO, N, A, LDA, VL, VU, IL, IU, ABSTOL, M, W, Z, LDZ, WORK, LWORK, IWORK, IFAIL, INFO)
  DSYEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for SY matrices More...
 
subroutine dsygv (ITYPE, JOBZ, UPLO, N, A, LDA, B, LDB, W, WORK, LWORK, INFO)
 DSYGV More...
 
subroutine dsygvd (ITYPE, JOBZ, UPLO, N, A, LDA, B, LDB, W, WORK, LWORK, IWORK, LIWORK, INFO)
 DSYGVD More...
 
subroutine dsygvx (ITYPE, JOBZ, RANGE, UPLO, N, A, LDA, B, LDB, VL, VU, IL, IU, ABSTOL, M, W, Z, LDZ, WORK, LWORK, IWORK, IFAIL, INFO)
 DSYGVX More...
 

Detailed Description

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

Function Documentation

subroutine dsyev ( character  JOBZ,
character  UPLO,
integer  N,
double precision, dimension( lda, * )  A,
integer  LDA,
double precision, dimension( * )  W,
double precision, dimension( * )  WORK,
integer  LWORK,
integer  INFO 
)

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

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

Purpose:
 DSYEV 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 DOUBLE PRECISION 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 DOUBLE PRECISION array, dimension (N)
          If INFO = 0, the eigenvalues in ascending order.
[out]WORK
          WORK is DOUBLE PRECISION 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 DSYTRD 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 dsyev.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  DOUBLE PRECISION a( lda, * ), w( * ), work( * )
146 * ..
147 *
148 * =====================================================================
149 *
150 * .. Parameters ..
151  DOUBLE PRECISION zero, one
152  parameter( zero = 0.0d0, one = 1.0d0 )
153 * ..
154 * .. Local Scalars ..
155  LOGICAL lower, lquery, wantz
156  INTEGER iinfo, imax, inde, indtau, indwrk, iscale,
157  $ llwork, lwkopt, nb
158  DOUBLE PRECISION anrm, bignum, eps, rmax, rmin, safmin, sigma,
159  $ smlnum
160 * ..
161 * .. External Functions ..
162  LOGICAL lsame
163  INTEGER ilaenv
164  DOUBLE PRECISION dlamch, dlansy
165  EXTERNAL lsame, ilaenv, dlamch, dlansy
166 * ..
167 * .. External Subroutines ..
168  EXTERNAL dlascl, dorgtr, dscal, dsteqr, dsterf, dsytrd,
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, 'DSYTRD', 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( 'DSYEV ', -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 = dlamch( 'Safe minimum' )
226  eps = dlamch( '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 = dlansy( '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 dlascl( uplo, 0, 0, one, sigma, n, n, a, lda, info )
245 *
246 * Call DSYTRD 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 dsytrd( uplo, n, a, lda, w, work( inde ), work( indtau ),
253  $ work( indwrk ), llwork, iinfo )
254 *
255 * For eigenvalues only, call DSTERF. For eigenvectors, first call
256 * DORGTR to generate the orthogonal matrix, then call DSTEQR.
257 *
258  IF( .NOT.wantz ) THEN
259  CALL dsterf( n, w, work( inde ), info )
260  ELSE
261  CALL dorgtr( uplo, n, a, lda, work( indtau ), work( indwrk ),
262  $ llwork, iinfo )
263  CALL dsteqr( 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 dscal( 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 DSYEV
285 *
subroutine dsteqr(COMPZ, N, D, E, Z, LDZ, WORK, INFO)
DSTEQR
Definition: dsteqr.f:133
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine dsterf(N, D, E, INFO)
DSTERF
Definition: dsterf.f:88
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:65
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine dsytrd(UPLO, N, A, LDA, D, E, TAU, WORK, LWORK, INFO)
DSYTRD
Definition: dsytrd.f:194
subroutine dlascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
DLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: dlascl.f:141
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
subroutine dorgtr(UPLO, N, A, LDA, TAU, WORK, LWORK, INFO)
DORGTR
Definition: dorgtr.f:125
subroutine dscal(N, DA, DX, INCX)
DSCAL
Definition: dscal.f:55
double precision function dlansy(NORM, UPLO, N, A, LDA, WORK)
DLANSY 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: dlansy.f:124

Here is the call graph for this function:

Here is the caller graph for this function:

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

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

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

Purpose:
 DSYEVD 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, DSYEVD needs N**2 more
 workspace than DSYEVX.
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 DOUBLE PRECISION 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 DOUBLE PRECISION array, dimension (N)
          If INFO = 0, the eigenvalues in ascending order.
[out]WORK
          WORK is DOUBLE PRECISION 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
September 2012
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 187 of file dsyevd.f.

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

Here is the call graph for this function:

Here is the caller graph for this function:

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

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

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

Purpose:
 DSYEVR 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.

 DSYEVR first reduces the matrix A to tridiagonal form T with a call
 to DSYTRD.  Then, whenever possible, DSYEVR calls DSTEMR to compute
 the eigenspectrum using Relatively Robust Representations.  DSTEMR
 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 DSTEMR'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 : DSYEVR calls DSTEMR when the full spectrum is requested
 on machines which conform to the ieee-754 floating point standard.
 DSYEVR calls DSTEBZ and SSTEIN on non-ieee machines and
 when partial spectrum requests are made.

 Normal execution of DSTEMR 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, DSTEBZ and
          DSTEIN 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 DOUBLE PRECISION 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 DOUBLE PRECISION
[in]VU
          VU is DOUBLE PRECISION
          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 DOUBLE PRECISION
          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
          DLAMCH( '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 DOUBLE PRECISION array, dimension (N)
          The first M elements contain the selected eigenvalues in
          ascending order.
[out]Z
          Z is DOUBLE PRECISION 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 DOUBLE PRECISION 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 DSYTRD and DORMTR
          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 (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 size of the IWORK array,
          returns this value as the first entry of the IWORK array, and
          no error message related to 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 327 of file dsyevr.f.

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

Here is the call graph for this function:

Here is the caller graph for this function:

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

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

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

Purpose:
 DSYEVX 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 DOUBLE PRECISION 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 DOUBLE PRECISION
[in]VU
          VU is DOUBLE PRECISION
          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 DOUBLE PRECISION
          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*DLAMCH('S'), not zero.
          If this routine returns with INFO>0, indicating that some
          eigenvectors did not converge, try setting ABSTOL to
          2*DLAMCH('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 DOUBLE PRECISION array, dimension (N)
          On normal exit, the first M elements contain the selected
          eigenvalues in ascending order.
[out]Z
          Z is DOUBLE PRECISION 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 DOUBLE PRECISION 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 DSYTRD and DORMTR
          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 dsyevx.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  DOUBLE PRECISION abstol, vl, vu
258 * ..
259 * .. Array Arguments ..
260  INTEGER ifail( * ), iwork( * )
261  DOUBLE PRECISION a( lda, * ), w( * ), work( * ), z( ldz, * )
262 * ..
263 *
264 * =====================================================================
265 *
266 * .. Parameters ..
267  DOUBLE PRECISION zero, one
268  parameter( zero = 0.0d+0, one = 1.0d+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  DOUBLE PRECISION abstll, anrm, bignum, eps, rmax, rmin, safmin,
279  $ sigma, smlnum, tmp1, vll, vuu
280 * ..
281 * .. External Functions ..
282  LOGICAL lsame
283  INTEGER ilaenv
284  DOUBLE PRECISION dlamch, dlansy
285  EXTERNAL lsame, ilaenv, dlamch, dlansy
286 * ..
287 * .. External Subroutines ..
288  EXTERNAL dcopy, dlacpy, dorgtr, dormtr, dscal, dstebz,
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, 'DSYTRD', uplo, n, -1, -1, -1 )
341  nb = max( nb, ilaenv( 1, 'DORMTR', 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( 'DSYEVX', -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 = dlamch( 'Safe minimum' )
382  eps = dlamch( '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 = dlansy( '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 dscal( n-j+1, sigma, a( j, j ), 1 )
408  10 CONTINUE
409  ELSE
410  DO 20 j = 1, n
411  CALL dscal( 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 DSYTRD 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 dsytrd( 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 DSTERF or DORGTR and SSTEQR. If this fails for
434 * some eigenvalue, then try DSTEBZ.
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 dcopy( n, work( indd ), 1, w, 1 )
444  indee = indwrk + 2*n
445  IF( .NOT.wantz ) THEN
446  CALL dcopy( n-1, work( inde ), 1, work( indee ), 1 )
447  CALL dsterf( n, w, work( indee ), info )
448  ELSE
449  CALL dlacpy( 'A', n, n, a, lda, z, ldz )
450  CALL dorgtr( uplo, n, z, ldz, work( indtau ),
451  $ work( indwrk ), llwork, iinfo )
452  CALL dcopy( n-1, work( inde ), 1, work( indee ), 1 )
453  CALL dsteqr( 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 DSTEBZ 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 dstebz( 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 dstein( 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 DSTEIN.
490 *
491  indwkn = inde
492  llwrkn = lwork - indwkn + 1
493  CALL dormtr( '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 dscal( 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 dswap( 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 DSYEVX
546 *
subroutine dsteqr(COMPZ, N, D, E, Z, LDZ, WORK, INFO)
DSTEQR
Definition: dsteqr.f:133
subroutine dstein(N, D, E, M, W, IBLOCK, ISPLIT, Z, LDZ, WORK, IWORK, IFAIL, INFO)
DSTEIN
Definition: dstein.f:176
subroutine dormtr(SIDE, UPLO, TRANS, M, N, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
DORMTR
Definition: dormtr.f:173
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
Definition: dlacpy.f:105
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine dsterf(N, D, E, INFO)
DSTERF
Definition: dsterf.f:88
subroutine dswap(N, DX, INCX, DY, INCY)
DSWAP
Definition: dswap.f:53
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:65
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine dsytrd(UPLO, N, A, LDA, D, E, TAU, WORK, LWORK, INFO)
DSYTRD
Definition: dsytrd.f:194
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
subroutine dorgtr(UPLO, N, A, LDA, TAU, WORK, LWORK, INFO)
DORGTR
Definition: dorgtr.f:125
subroutine dstebz(RANGE, ORDER, N, VL, VU, IL, IU, ABSTOL, D, E, M, NSPLIT, W, IBLOCK, ISPLIT, WORK, IWORK, INFO)
DSTEBZ
Definition: dstebz.f:265
subroutine dscal(N, DA, DX, INCX)
DSCAL
Definition: dscal.f:55
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
Definition: dcopy.f:53
double precision function dlansy(NORM, UPLO, N, A, LDA, WORK)
DLANSY 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: dlansy.f:124

Here is the call graph for this function:

Here is the caller graph for this function:

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

DSYGV

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

Purpose:
 DSYGV 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 DOUBLE PRECISION 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 DOUBLE PRECISION 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 DOUBLE PRECISION array, dimension (N)
          If INFO = 0, the eigenvalues in ascending order.
[out]WORK
          WORK is DOUBLE PRECISION 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 DSYTRD 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:  DPOTRF or DSYEV returned an error code:
             <= N:  if INFO = i, DSYEV 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 dsygv.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  DOUBLE PRECISION a( lda, * ), b( ldb, * ), w( * ), work( * )
189 * ..
190 *
191 * =====================================================================
192 *
193 * .. Parameters ..
194  DOUBLE PRECISION one
195  parameter( one = 1.0d+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 lsame, ilaenv
206 * ..
207 * .. External Subroutines ..
208  EXTERNAL dpotrf, dsyev, dsygst, dtrmm, dtrsm, 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, 'DSYTRD', 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( 'DSYGV ', -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 dpotrf( 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 dsygst( itype, uplo, n, a, lda, b, ldb, info )
270  CALL dsyev( 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 dtrsm( '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 dtrmm( '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 DSYGV
313 *
subroutine dsygst(ITYPE, UPLO, N, A, LDA, B, LDB, INFO)
DSYGST
Definition: dsygst.f:129
subroutine dtrmm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
DTRMM
Definition: dtrmm.f:179
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine dtrsm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
DTRSM
Definition: dtrsm.f:183
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine dpotrf(UPLO, N, A, LDA, INFO)
DPOTRF
Definition: dpotrf.f:109
subroutine dsyev(JOBZ, UPLO, N, A, LDA, W, WORK, LWORK, INFO)
DSYEV computes the eigenvalues and, optionally, the left and/or right eigenvectors for SY matrices ...
Definition: dsyev.f:134
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 dsygvd ( integer  ITYPE,
character  JOBZ,
character  UPLO,
integer  N,
double precision, dimension( lda, * )  A,
integer  LDA,
double precision, dimension( ldb, * )  B,
integer  LDB,
double precision, dimension( * )  W,
double precision, dimension( * )  WORK,
integer  LWORK,
integer, dimension( * )  IWORK,
integer  LIWORK,
integer  INFO 
)

DSYGVD

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

Purpose:
 DSYGVD 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 DOUBLE PRECISION 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 DOUBLE PRECISION 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 DOUBLE PRECISION array, dimension (N)
          If INFO = 0, the eigenvalues in ascending order.
[out]WORK
          WORK is DOUBLE PRECISION 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:  DPOTRF or DSYEVD 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 DSYEVD 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 dsygvd.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  DOUBLE PRECISION a( lda, * ), b( ldb, * ), w( * ), work( * )
242 * ..
243 *
244 * =====================================================================
245 *
246 * .. Parameters ..
247  DOUBLE PRECISION one
248  parameter( one = 1.0d+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 dpotrf, dsyevd, dsygst, dtrmm, dtrsm, xerbla
261 * ..
262 * .. Intrinsic Functions ..
263  INTRINSIC dble, max
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( 'DSYGVD', -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 dpotrf( 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 dsygst( itype, uplo, n, a, lda, b, ldb, info )
334  CALL dsyevd( jobz, uplo, n, a, lda, w, work, lwork, iwork, liwork,
335  $ info )
336  lopt = max( dble( lopt ), dble( work( 1 ) ) )
337  liopt = max( dble( liopt ), dble( 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 dtrsm( '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 dtrmm( '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 DSYGVD
379 *
subroutine dsygst(ITYPE, UPLO, N, A, LDA, B, LDB, INFO)
DSYGST
Definition: dsygst.f:129
subroutine dsyevd(JOBZ, UPLO, N, A, LDA, W, WORK, LWORK, IWORK, LIWORK, INFO)
DSYEVD computes the eigenvalues and, optionally, the left and/or right eigenvectors for SY matrices ...
Definition: dsyevd.f:187
subroutine dtrmm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
DTRMM
Definition: dtrmm.f:179
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine dtrsm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
DTRSM
Definition: dtrsm.f:183
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine dpotrf(UPLO, N, A, LDA, INFO)
DPOTRF
Definition: dpotrf.f:109

Here is the call graph for this function:

Here is the caller graph for this function:

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

DSYGVX

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

Purpose:
 DSYGVX 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 DOUBLE PRECISION 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 DOUBLE PRECISION 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).
[in]VL
          VL is DOUBLE PRECISION
[in]VU
          VU is DOUBLE PRECISION
          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 DOUBLE PRECISION
          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*DLAMCH('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 DOUBLE PRECISION array, dimension (N)
          On normal exit, the first M elements contain the selected
          eigenvalues in ascending order.
[out]Z
          Z is DOUBLE PRECISION 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 DOUBLE PRECISION 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 DSYTRD 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:  DPOTRF or DSYEVX returned an error code:
             <= N:  if INFO = i, DSYEVX 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 dsygvx.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  DOUBLE PRECISION abstol, vl, vu
302 * ..
303 * .. Array Arguments ..
304  INTEGER ifail( * ), iwork( * )
305  DOUBLE PRECISION a( lda, * ), b( ldb, * ), w( * ), work( * ),
306  $ z( ldz, * )
307 * ..
308 *
309 * =====================================================================
310 *
311 * .. Parameters ..
312  DOUBLE PRECISION one
313  parameter( one = 1.0d+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 lsame, ilaenv
324 * ..
325 * .. External Subroutines ..
326  EXTERNAL dpotrf, dsyevx, dsygst, dtrmm, dtrsm, 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, 'DSYTRD', 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( 'DSYGVX', -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 dpotrf( 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 dsygst( itype, uplo, n, a, lda, b, ldb, info )
411  CALL dsyevx( 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 dtrsm( '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 dtrmm( '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 DSYGVX
457 *
subroutine dsyevx(JOBZ, RANGE, UPLO, N, A, LDA, VL, VU, IL, IU, ABSTOL, M, W, Z, LDZ, WORK, LWORK, IWORK, IFAIL, INFO)
DSYEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for SY matrices ...
Definition: dsyevx.f:248
subroutine dsygst(ITYPE, UPLO, N, A, LDA, B, LDB, INFO)
DSYGST
Definition: dsygst.f:129
subroutine dtrmm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
DTRMM
Definition: dtrmm.f:179
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine dtrsm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
DTRSM
Definition: dtrsm.f:183
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine dpotrf(UPLO, N, A, LDA, INFO)
DPOTRF
Definition: dpotrf.f:109
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: