LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ dstevx()

subroutine dstevx ( character  JOBZ,
character  RANGE,
integer  N,
double precision, dimension( * )  D,
double precision, dimension( * )  E,
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, dimension( * )  IWORK,
integer, dimension( * )  IFAIL,
integer  INFO 
)

DSTEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for OTHER matrices

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

Purpose:
 DSTEVX computes selected eigenvalues and, optionally, eigenvectors
 of a real symmetric tridiagonal 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]N
          N is INTEGER
          The order of the matrix.  N >= 0.
[in,out]D
          D is DOUBLE PRECISION array, dimension (N)
          On entry, the n diagonal elements of the tridiagonal matrix
          A.
          On exit, D may be multiplied by a constant factor chosen
          to avoid over/underflow in computing the eigenvalues.
[in,out]E
          E is DOUBLE PRECISION array, dimension (max(1,N-1))
          On entry, the (n-1) subdiagonal elements of the tridiagonal
          matrix A in elements 1 to N-1 of E.
          On exit, E may be multiplied by a constant factor chosen
          to avoid over/underflow in computing the eigenvalues.
[in]VL
          VL is DOUBLE PRECISION
          If RANGE='V', the lower bound of the interval to
          be searched for eigenvalues. VL < VU.
          Not referenced if RANGE = 'A' or 'I'.
[in]VU
          VU is DOUBLE PRECISION
          If RANGE='V', the upper bound of the interval to
          be searched for eigenvalues. VL < VU.
          Not referenced if RANGE = 'A' or 'I'.
[in]IL
          IL is INTEGER
          If RANGE='I', the index of the
          smallest eigenvalue 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]IU
          IU is INTEGER
          If RANGE='I', the index of the
          largest eigenvalue 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.

          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)
          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 (INFO > 0), 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 (5*N)
[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.

Definition at line 225 of file dstevx.f.

227 *
228 * -- LAPACK driver routine --
229 * -- LAPACK is a software package provided by Univ. of Tennessee, --
230 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
231 *
232 * .. Scalar Arguments ..
233  CHARACTER JOBZ, RANGE
234  INTEGER IL, INFO, IU, LDZ, M, N
235  DOUBLE PRECISION ABSTOL, VL, VU
236 * ..
237 * .. Array Arguments ..
238  INTEGER IFAIL( * ), IWORK( * )
239  DOUBLE PRECISION D( * ), E( * ), W( * ), WORK( * ), Z( LDZ, * )
240 * ..
241 *
242 * =====================================================================
243 *
244 * .. Parameters ..
245  DOUBLE PRECISION ZERO, ONE
246  parameter( zero = 0.0d0, one = 1.0d0 )
247 * ..
248 * .. Local Scalars ..
249  LOGICAL ALLEIG, INDEIG, TEST, VALEIG, WANTZ
250  CHARACTER ORDER
251  INTEGER I, IMAX, INDIBL, INDISP, INDIWO, INDWRK,
252  $ ISCALE, ITMP1, J, JJ, NSPLIT
253  DOUBLE PRECISION BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA, SMLNUM,
254  $ TMP1, TNRM, VLL, VUU
255 * ..
256 * .. External Functions ..
257  LOGICAL LSAME
258  DOUBLE PRECISION DLAMCH, DLANST
259  EXTERNAL lsame, dlamch, dlanst
260 * ..
261 * .. External Subroutines ..
262  EXTERNAL dcopy, dscal, dstebz, dstein, dsteqr, dsterf,
263  $ dswap, xerbla
264 * ..
265 * .. Intrinsic Functions ..
266  INTRINSIC max, min, sqrt
267 * ..
268 * .. Executable Statements ..
269 *
270 * Test the input parameters.
271 *
272  wantz = lsame( jobz, 'V' )
273  alleig = lsame( range, 'A' )
274  valeig = lsame( range, 'V' )
275  indeig = lsame( range, 'I' )
276 *
277  info = 0
278  IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
279  info = -1
280  ELSE IF( .NOT.( alleig .OR. valeig .OR. indeig ) ) THEN
281  info = -2
282  ELSE IF( n.LT.0 ) THEN
283  info = -3
284  ELSE
285  IF( valeig ) THEN
286  IF( n.GT.0 .AND. vu.LE.vl )
287  $ info = -7
288  ELSE IF( indeig ) THEN
289  IF( il.LT.1 .OR. il.GT.max( 1, n ) ) THEN
290  info = -8
291  ELSE IF( iu.LT.min( n, il ) .OR. iu.GT.n ) THEN
292  info = -9
293  END IF
294  END IF
295  END IF
296  IF( info.EQ.0 ) THEN
297  IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) )
298  $ info = -14
299  END IF
300 *
301  IF( info.NE.0 ) THEN
302  CALL xerbla( 'DSTEVX', -info )
303  RETURN
304  END IF
305 *
306 * Quick return if possible
307 *
308  m = 0
309  IF( n.EQ.0 )
310  $ RETURN
311 *
312  IF( n.EQ.1 ) THEN
313  IF( alleig .OR. indeig ) THEN
314  m = 1
315  w( 1 ) = d( 1 )
316  ELSE
317  IF( vl.LT.d( 1 ) .AND. vu.GE.d( 1 ) ) THEN
318  m = 1
319  w( 1 ) = d( 1 )
320  END IF
321  END IF
322  IF( wantz )
323  $ z( 1, 1 ) = one
324  RETURN
325  END IF
326 *
327 * Get machine constants.
328 *
329  safmin = dlamch( 'Safe minimum' )
330  eps = dlamch( 'Precision' )
331  smlnum = safmin / eps
332  bignum = one / smlnum
333  rmin = sqrt( smlnum )
334  rmax = min( sqrt( bignum ), one / sqrt( sqrt( safmin ) ) )
335 *
336 * Scale matrix to allowable range, if necessary.
337 *
338  iscale = 0
339  IF( valeig ) THEN
340  vll = vl
341  vuu = vu
342  ELSE
343  vll = zero
344  vuu = zero
345  END IF
346  tnrm = dlanst( 'M', n, d, e )
347  IF( tnrm.GT.zero .AND. tnrm.LT.rmin ) THEN
348  iscale = 1
349  sigma = rmin / tnrm
350  ELSE IF( tnrm.GT.rmax ) THEN
351  iscale = 1
352  sigma = rmax / tnrm
353  END IF
354  IF( iscale.EQ.1 ) THEN
355  CALL dscal( n, sigma, d, 1 )
356  CALL dscal( n-1, sigma, e( 1 ), 1 )
357  IF( valeig ) THEN
358  vll = vl*sigma
359  vuu = vu*sigma
360  END IF
361  END IF
362 *
363 * If all eigenvalues are desired and ABSTOL is less than zero, then
364 * call DSTERF or SSTEQR. If this fails for some eigenvalue, then
365 * try DSTEBZ.
366 *
367  test = .false.
368  IF( indeig ) THEN
369  IF( il.EQ.1 .AND. iu.EQ.n ) THEN
370  test = .true.
371  END IF
372  END IF
373  IF( ( alleig .OR. test ) .AND. ( abstol.LE.zero ) ) THEN
374  CALL dcopy( n, d, 1, w, 1 )
375  CALL dcopy( n-1, e( 1 ), 1, work( 1 ), 1 )
376  indwrk = n + 1
377  IF( .NOT.wantz ) THEN
378  CALL dsterf( n, w, work, info )
379  ELSE
380  CALL dsteqr( 'I', n, w, work, z, ldz, work( indwrk ), info )
381  IF( info.EQ.0 ) THEN
382  DO 10 i = 1, n
383  ifail( i ) = 0
384  10 CONTINUE
385  END IF
386  END IF
387  IF( info.EQ.0 ) THEN
388  m = n
389  GO TO 20
390  END IF
391  info = 0
392  END IF
393 *
394 * Otherwise, call DSTEBZ and, if eigenvectors are desired, SSTEIN.
395 *
396  IF( wantz ) THEN
397  order = 'B'
398  ELSE
399  order = 'E'
400  END IF
401  indwrk = 1
402  indibl = 1
403  indisp = indibl + n
404  indiwo = indisp + n
405  CALL dstebz( range, order, n, vll, vuu, il, iu, abstol, d, e, m,
406  $ nsplit, w, iwork( indibl ), iwork( indisp ),
407  $ work( indwrk ), iwork( indiwo ), info )
408 *
409  IF( wantz ) THEN
410  CALL dstein( n, d, e, m, w, iwork( indibl ), iwork( indisp ),
411  $ z, ldz, work( indwrk ), iwork( indiwo ), ifail,
412  $ info )
413  END IF
414 *
415 * If matrix was scaled, then rescale eigenvalues appropriately.
416 *
417  20 CONTINUE
418  IF( iscale.EQ.1 ) THEN
419  IF( info.EQ.0 ) THEN
420  imax = m
421  ELSE
422  imax = info - 1
423  END IF
424  CALL dscal( imax, one / sigma, w, 1 )
425  END IF
426 *
427 * If eigenvalues are not in order, then sort them, along with
428 * eigenvectors.
429 *
430  IF( wantz ) THEN
431  DO 40 j = 1, m - 1
432  i = 0
433  tmp1 = w( j )
434  DO 30 jj = j + 1, m
435  IF( w( jj ).LT.tmp1 ) THEN
436  i = jj
437  tmp1 = w( jj )
438  END IF
439  30 CONTINUE
440 *
441  IF( i.NE.0 ) THEN
442  itmp1 = iwork( indibl+i-1 )
443  w( i ) = w( j )
444  iwork( indibl+i-1 ) = iwork( indibl+j-1 )
445  w( j ) = tmp1
446  iwork( indibl+j-1 ) = itmp1
447  CALL dswap( n, z( 1, i ), 1, z( 1, j ), 1 )
448  IF( info.NE.0 ) THEN
449  itmp1 = ifail( i )
450  ifail( i ) = ifail( j )
451  ifail( j ) = itmp1
452  END IF
453  END IF
454  40 CONTINUE
455  END IF
456 *
457  RETURN
458 *
459 * End of DSTEVX
460 *
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:69
double precision function dlanst(NORM, N, D, E)
DLANST returns the value of the 1-norm, or the Frobenius norm, or the infinity norm,...
Definition: dlanst.f:100
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
subroutine dstebz(RANGE, ORDER, N, VL, VU, IL, IU, ABSTOL, D, E, M, NSPLIT, W, IBLOCK, ISPLIT, WORK, IWORK, INFO)
DSTEBZ
Definition: dstebz.f:273
subroutine dsteqr(COMPZ, N, D, E, Z, LDZ, WORK, INFO)
DSTEQR
Definition: dsteqr.f:131
subroutine dsterf(N, D, E, INFO)
DSTERF
Definition: dsterf.f:86
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
Definition: dcopy.f:82
subroutine dscal(N, DA, DX, INCX)
DSCAL
Definition: dscal.f:79
subroutine dswap(N, DX, INCX, DY, INCY)
DSWAP
Definition: dswap.f:82
subroutine dstein(N, D, E, M, W, IBLOCK, ISPLIT, Z, LDZ, WORK, IWORK, IFAIL, INFO)
DSTEIN
Definition: dstein.f:174
Here is the call graph for this function:
Here is the caller graph for this function: