LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
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.
Date
June 2016

Definition at line 229 of file dstevx.f.

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

Here is the call graph for this function:

Here is the caller graph for this function: