 LAPACK  3.10.1 LAPACK: Linear Algebra PACKage

## ◆ dstemr()

 subroutine dstemr ( character JOBZ, character RANGE, integer N, double precision, dimension( * ) D, double precision, dimension( * ) E, double precision VL, double precision VU, integer IL, integer IU, integer M, double precision, dimension( * ) W, double precision, dimension( ldz, * ) Z, integer LDZ, integer NZC, integer, dimension( * ) ISUPPZ, logical TRYRAC, double precision, dimension( * ) WORK, integer LWORK, integer, dimension( * ) IWORK, integer LIWORK, integer INFO )

DSTEMR

Purpose:
``` DSTEMR computes selected eigenvalues and, optionally, eigenvectors
of a real symmetric tridiagonal matrix T. Any such unreduced matrix has
a well defined set of pairwise different real eigenvalues, the corresponding
real eigenvectors are pairwise orthogonal.

The spectrum may be computed either completely or partially by specifying
either an interval (VL,VU] or a range of indices IL:IU for the desired
eigenvalues.

Depending on the number of desired eigenvalues, these are computed either
by bisection or the dqds algorithm. Numerically orthogonal eigenvectors are
computed by the use of various suitable L D L^T factorizations near clusters
of close eigenvalues (referred to as RRRs, Relatively Robust
Representations). An informal sketch of the algorithm 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.

For more details, see:
- 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.

Further Details
1.DSTEMR works only on machines which follow IEEE-754
floating-point standard in their handling of infinities and NaNs.
This permits the use of efficient inner loops avoiding a check for
zero divisors.```
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 T. On exit, D is overwritten.``` [in,out] E ``` E is DOUBLE PRECISION array, dimension (N) On entry, the (N-1) subdiagonal elements of the tridiagonal matrix T in elements 1 to N-1 of E. E(N) need not be set on input, but is used internally as workspace. On exit, E is overwritten.``` [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. 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. Not referenced if RANGE = 'A' or 'V'.``` [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', and if INFO = 0, then the first M columns of Z contain the orthonormal eigenvectors of the matrix T 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 can be computed with a workspace query by setting NZC = -1, see below.``` [in] LDZ ``` LDZ is INTEGER The leading dimension of the array Z. LDZ >= 1, and if JOBZ = 'V', then LDZ >= max(1,N).``` [in] NZC ``` NZC is INTEGER The number of eigenvectors to be held in the array Z. If RANGE = 'A', then NZC >= max(1,N). If RANGE = 'V', then NZC >= the number of eigenvalues in (VL,VU]. If RANGE = 'I', then NZC >= IU-IL+1. If NZC = -1, then a workspace query is assumed; the routine calculates the number of columns of the array Z that are needed to hold the eigenvectors. This value is returned as the first entry of the Z array, and no error message related to NZC is issued by XERBLA.``` [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 computed eigenvector is nonzero only in elements ISUPPZ( 2*i-1 ) through ISUPPZ( 2*i ). This is relevant in the case when the matrix is split. ISUPPZ is only accessed when JOBZ is 'V' and N > 0.``` [in,out] TRYRAC ``` TRYRAC is LOGICAL If TRYRAC = .TRUE., indicates that the code should check whether the tridiagonal matrix defines its eigenvalues to high relative accuracy. If so, the code uses relative-accuracy preserving algorithms that might be (a bit) slower depending on the matrix. If the matrix does not define its eigenvalues to high relative accuracy, the code can uses possibly faster algorithms. If TRYRAC = .FALSE., the code is not required to guarantee relatively accurate eigenvalues and can use the fastest possible techniques. On exit, a .TRUE. TRYRAC will be set to .FALSE. if the matrix does not define its eigenvalues to high relative accuracy.``` [out] WORK ``` WORK is DOUBLE PRECISION array, dimension (LWORK) On exit, if INFO = 0, WORK(1) returns the optimal (and minimal) LWORK.``` [in] LWORK ``` LWORK is INTEGER The dimension of the array WORK. LWORK >= max(1,18*N) if JOBZ = 'V', and LWORK >= max(1,12*N) if JOBZ = 'N'. 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 (LIWORK) On exit, if INFO = 0, IWORK(1) returns the optimal LIWORK.``` [in] LIWORK ``` LIWORK is INTEGER The dimension of the array IWORK. LIWORK >= max(1,10*N) if the eigenvectors are desired, and LIWORK >= max(1,8*N) if only the eigenvalues are to be computed. 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 On exit, INFO = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value > 0: if INFO = 1X, internal error in DLARRE, if INFO = 2X, internal error in DLARRV. Here, the digit X = ABS( IINFO ) < 10, where IINFO is the nonzero error code returned by DLARRE or DLARRV, respectively.```
Contributors:
Beresford Parlett, University of California, Berkeley, USA
Jim Demmel, University of California, Berkeley, USA
Inderjit Dhillon, University of Texas, Austin, USA
Osni Marques, LBNL/NERSC, USA
Christof Voemel, University of California, Berkeley, USA

Definition at line 318 of file dstemr.f.

321 *
322 * -- LAPACK computational routine --
323 * -- LAPACK is a software package provided by Univ. of Tennessee, --
324 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
325 *
326 * .. Scalar Arguments ..
327  CHARACTER JOBZ, RANGE
328  LOGICAL TRYRAC
329  INTEGER IL, INFO, IU, LDZ, NZC, LIWORK, LWORK, M, N
330  DOUBLE PRECISION VL, VU
331 * ..
332 * .. Array Arguments ..
333  INTEGER ISUPPZ( * ), IWORK( * )
334  DOUBLE PRECISION D( * ), E( * ), W( * ), WORK( * )
335  DOUBLE PRECISION Z( LDZ, * )
336 * ..
337 *
338 * =====================================================================
339 *
340 * .. Parameters ..
341  DOUBLE PRECISION ZERO, ONE, FOUR, MINRGP
342  parameter( zero = 0.0d0, one = 1.0d0,
343  \$ four = 4.0d0,
344  \$ minrgp = 1.0d-3 )
345 * ..
346 * .. Local Scalars ..
347  LOGICAL ALLEIG, INDEIG, LQUERY, VALEIG, WANTZ, ZQUERY
348  INTEGER I, IBEGIN, IEND, IFIRST, IIL, IINDBL, IINDW,
349  \$ IINDWK, IINFO, IINSPL, IIU, ILAST, IN, INDD,
350  \$ INDE2, INDERR, INDGP, INDGRS, INDWRK, ITMP,
351  \$ ITMP2, J, JBLK, JJ, LIWMIN, LWMIN, NSPLIT,
352  \$ NZCMIN, OFFSET, WBEGIN, WEND
353  DOUBLE PRECISION BIGNUM, CS, EPS, PIVMIN, R1, R2, RMAX, RMIN,
354  \$ RTOL1, RTOL2, SAFMIN, SCALE, SMLNUM, SN,
355  \$ THRESH, TMP, TNRM, WL, WU
356 * ..
357 * ..
358 * .. External Functions ..
359  LOGICAL LSAME
360  DOUBLE PRECISION DLAMCH, DLANST
361  EXTERNAL lsame, dlamch, dlanst
362 * ..
363 * .. External Subroutines ..
364  EXTERNAL dcopy, dlae2, dlaev2, dlarrc, dlarre, dlarrj,
366 * ..
367 * .. Intrinsic Functions ..
368  INTRINSIC max, min, sqrt
369
370
371 * ..
372 * .. Executable Statements ..
373 *
374 * Test the input parameters.
375 *
376  wantz = lsame( jobz, 'V' )
377  alleig = lsame( range, 'A' )
378  valeig = lsame( range, 'V' )
379  indeig = lsame( range, 'I' )
380 *
381  lquery = ( ( lwork.EQ.-1 ).OR.( liwork.EQ.-1 ) )
382  zquery = ( nzc.EQ.-1 )
383
384 * DSTEMR needs WORK of size 6*N, IWORK of size 3*N.
385 * In addition, DLARRE needs WORK of size 6*N, IWORK of size 5*N.
386 * Furthermore, DLARRV needs WORK of size 12*N, IWORK of size 7*N.
387  IF( wantz ) THEN
388  lwmin = 18*n
389  liwmin = 10*n
390  ELSE
391 * need less workspace if only the eigenvalues are wanted
392  lwmin = 12*n
393  liwmin = 8*n
394  ENDIF
395
396  wl = zero
397  wu = zero
398  iil = 0
399  iiu = 0
400  nsplit = 0
401
402  IF( valeig ) THEN
403 * We do not reference VL, VU in the cases RANGE = 'I','A'
404 * The interval (WL, WU] contains all the wanted eigenvalues.
405 * It is either given by the user or computed in DLARRE.
406  wl = vl
407  wu = vu
408  ELSEIF( indeig ) THEN
409 * We do not reference IL, IU in the cases RANGE = 'V','A'
410  iil = il
411  iiu = iu
412  ENDIF
413 *
414  info = 0
415  IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
416  info = -1
417  ELSE IF( .NOT.( alleig .OR. valeig .OR. indeig ) ) THEN
418  info = -2
419  ELSE IF( n.LT.0 ) THEN
420  info = -3
421  ELSE IF( valeig .AND. n.GT.0 .AND. wu.LE.wl ) THEN
422  info = -7
423  ELSE IF( indeig .AND. ( iil.LT.1 .OR. iil.GT.n ) ) THEN
424  info = -8
425  ELSE IF( indeig .AND. ( iiu.LT.iil .OR. iiu.GT.n ) ) THEN
426  info = -9
427  ELSE IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) ) THEN
428  info = -13
429  ELSE IF( lwork.LT.lwmin .AND. .NOT.lquery ) THEN
430  info = -17
431  ELSE IF( liwork.LT.liwmin .AND. .NOT.lquery ) THEN
432  info = -19
433  END IF
434 *
435 * Get machine constants.
436 *
437  safmin = dlamch( 'Safe minimum' )
438  eps = dlamch( 'Precision' )
439  smlnum = safmin / eps
440  bignum = one / smlnum
441  rmin = sqrt( smlnum )
442  rmax = min( sqrt( bignum ), one / sqrt( sqrt( safmin ) ) )
443 *
444  IF( info.EQ.0 ) THEN
445  work( 1 ) = lwmin
446  iwork( 1 ) = liwmin
447 *
448  IF( wantz .AND. alleig ) THEN
449  nzcmin = n
450  ELSE IF( wantz .AND. valeig ) THEN
451  CALL dlarrc( 'T', n, vl, vu, d, e, safmin,
452  \$ nzcmin, itmp, itmp2, info )
453  ELSE IF( wantz .AND. indeig ) THEN
454  nzcmin = iiu-iil+1
455  ELSE
456 * WANTZ .EQ. FALSE.
457  nzcmin = 0
458  ENDIF
459  IF( zquery .AND. info.EQ.0 ) THEN
460  z( 1,1 ) = nzcmin
461  ELSE IF( nzc.LT.nzcmin .AND. .NOT.zquery ) THEN
462  info = -14
463  END IF
464  END IF
465
466  IF( info.NE.0 ) THEN
467 *
468  CALL xerbla( 'DSTEMR', -info )
469 *
470  RETURN
471  ELSE IF( lquery .OR. zquery ) THEN
472  RETURN
473  END IF
474 *
475 * Handle N = 0, 1, and 2 cases immediately
476 *
477  m = 0
478  IF( n.EQ.0 )
479  \$ RETURN
480 *
481  IF( n.EQ.1 ) THEN
482  IF( alleig .OR. indeig ) THEN
483  m = 1
484  w( 1 ) = d( 1 )
485  ELSE
486  IF( wl.LT.d( 1 ) .AND. wu.GE.d( 1 ) ) THEN
487  m = 1
488  w( 1 ) = d( 1 )
489  END IF
490  END IF
491  IF( wantz.AND.(.NOT.zquery) ) THEN
492  z( 1, 1 ) = one
493  isuppz(1) = 1
494  isuppz(2) = 1
495  END IF
496  RETURN
497  END IF
498 *
499  IF( n.EQ.2 ) THEN
500  IF( .NOT.wantz ) THEN
501  CALL dlae2( d(1), e(1), d(2), r1, r2 )
502  ELSE IF( wantz.AND.(.NOT.zquery) ) THEN
503  CALL dlaev2( d(1), e(1), d(2), r1, r2, cs, sn )
504  END IF
505  IF( alleig.OR.
506  \$ (valeig.AND.(r2.GT.wl).AND.
507  \$ (r2.LE.wu)).OR.
508  \$ (indeig.AND.(iil.EQ.1)) ) THEN
509  m = m+1
510  w( m ) = r2
511  IF( wantz.AND.(.NOT.zquery) ) THEN
512  z( 1, m ) = -sn
513  z( 2, m ) = cs
514 * Note: At most one of SN and CS can be zero.
515  IF (sn.NE.zero) THEN
516  IF (cs.NE.zero) THEN
517  isuppz(2*m-1) = 1
518  isuppz(2*m) = 2
519  ELSE
520  isuppz(2*m-1) = 1
521  isuppz(2*m) = 1
522  END IF
523  ELSE
524  isuppz(2*m-1) = 2
525  isuppz(2*m) = 2
526  END IF
527  ENDIF
528  ENDIF
529  IF( alleig.OR.
530  \$ (valeig.AND.(r1.GT.wl).AND.
531  \$ (r1.LE.wu)).OR.
532  \$ (indeig.AND.(iiu.EQ.2)) ) THEN
533  m = m+1
534  w( m ) = r1
535  IF( wantz.AND.(.NOT.zquery) ) THEN
536  z( 1, m ) = cs
537  z( 2, m ) = sn
538 * Note: At most one of SN and CS can be zero.
539  IF (sn.NE.zero) THEN
540  IF (cs.NE.zero) THEN
541  isuppz(2*m-1) = 1
542  isuppz(2*m) = 2
543  ELSE
544  isuppz(2*m-1) = 1
545  isuppz(2*m) = 1
546  END IF
547  ELSE
548  isuppz(2*m-1) = 2
549  isuppz(2*m) = 2
550  END IF
551  ENDIF
552  ENDIF
553
554  ELSE
555
556 * Continue with general N
557
558  indgrs = 1
559  inderr = 2*n + 1
560  indgp = 3*n + 1
561  indd = 4*n + 1
562  inde2 = 5*n + 1
563  indwrk = 6*n + 1
564 *
565  iinspl = 1
566  iindbl = n + 1
567  iindw = 2*n + 1
568  iindwk = 3*n + 1
569 *
570 * Scale matrix to allowable range, if necessary.
571 * The allowable range is related to the PIVMIN parameter; see the
572 * comments in DLARRD. The preference for scaling small values
573 * up is heuristic; we expect users' matrices not to be close to the
574 * RMAX threshold.
575 *
576  scale = one
577  tnrm = dlanst( 'M', n, d, e )
578  IF( tnrm.GT.zero .AND. tnrm.LT.rmin ) THEN
579  scale = rmin / tnrm
580  ELSE IF( tnrm.GT.rmax ) THEN
581  scale = rmax / tnrm
582  END IF
583  IF( scale.NE.one ) THEN
584  CALL dscal( n, scale, d, 1 )
585  CALL dscal( n-1, scale, e, 1 )
586  tnrm = tnrm*scale
587  IF( valeig ) THEN
588 * If eigenvalues in interval have to be found,
589 * scale (WL, WU] accordingly
590  wl = wl*scale
591  wu = wu*scale
592  ENDIF
593  END IF
594 *
595 * Compute the desired eigenvalues of the tridiagonal after splitting
596 * into smaller subblocks if the corresponding off-diagonal elements
597 * are small
598 * THRESH is the splitting parameter for DLARRE
599 * A negative THRESH forces the old splitting criterion based on the
600 * size of the off-diagonal. A positive THRESH switches to splitting
601 * which preserves relative accuracy.
602 *
603  IF( tryrac ) THEN
604 * Test whether the matrix warrants the more expensive relative approach.
605  CALL dlarrr( n, d, e, iinfo )
606  ELSE
607 * The user does not care about relative accurately eigenvalues
608  iinfo = -1
609  ENDIF
610 * Set the splitting criterion
611  IF (iinfo.EQ.0) THEN
612  thresh = eps
613  ELSE
614  thresh = -eps
615 * relative accuracy is desired but T does not guarantee it
616  tryrac = .false.
617  ENDIF
618 *
619  IF( tryrac ) THEN
620 * Copy original diagonal, needed to guarantee relative accuracy
621  CALL dcopy(n,d,1,work(indd),1)
622  ENDIF
623 * Store the squares of the offdiagonal values of T
624  DO 5 j = 1, n-1
625  work( inde2+j-1 ) = e(j)**2
626  5 CONTINUE
627
628 * Set the tolerance parameters for bisection
629  IF( .NOT.wantz ) THEN
630 * DLARRE computes the eigenvalues to full precision.
631  rtol1 = four * eps
632  rtol2 = four * eps
633  ELSE
634 * DLARRE computes the eigenvalues to less than full precision.
635 * DLARRV will refine the eigenvalue approximations, and we can
636 * need less accurate initial bisection in DLARRE.
637 * Note: these settings do only affect the subset case and DLARRE
638  rtol1 = sqrt(eps)
639  rtol2 = max( sqrt(eps)*5.0d-3, four * eps )
640  ENDIF
641  CALL dlarre( range, n, wl, wu, iil, iiu, d, e,
642  \$ work(inde2), rtol1, rtol2, thresh, nsplit,
643  \$ iwork( iinspl ), m, w, work( inderr ),
644  \$ work( indgp ), iwork( iindbl ),
645  \$ iwork( iindw ), work( indgrs ), pivmin,
646  \$ work( indwrk ), iwork( iindwk ), iinfo )
647  IF( iinfo.NE.0 ) THEN
648  info = 10 + abs( iinfo )
649  RETURN
650  END IF
651 * Note that if RANGE .NE. 'V', DLARRE computes bounds on the desired
652 * part of the spectrum. All desired eigenvalues are contained in
653 * (WL,WU]
654
655
656  IF( wantz ) THEN
657 *
658 * Compute the desired eigenvectors corresponding to the computed
659 * eigenvalues
660 *
661  CALL dlarrv( n, wl, wu, d, e,
662  \$ pivmin, iwork( iinspl ), m,
663  \$ 1, m, minrgp, rtol1, rtol2,
664  \$ w, work( inderr ), work( indgp ), iwork( iindbl ),
665  \$ iwork( iindw ), work( indgrs ), z, ldz,
666  \$ isuppz, work( indwrk ), iwork( iindwk ), iinfo )
667  IF( iinfo.NE.0 ) THEN
668  info = 20 + abs( iinfo )
669  RETURN
670  END IF
671  ELSE
672 * DLARRE computes eigenvalues of the (shifted) root representation
673 * DLARRV returns the eigenvalues of the unshifted matrix.
674 * However, if the eigenvectors are not desired by the user, we need
675 * to apply the corresponding shifts from DLARRE to obtain the
676 * eigenvalues of the original matrix.
677  DO 20 j = 1, m
678  itmp = iwork( iindbl+j-1 )
679  w( j ) = w( j ) + e( iwork( iinspl+itmp-1 ) )
680  20 CONTINUE
681  END IF
682 *
683
684  IF ( tryrac ) THEN
685 * Refine computed eigenvalues so that they are relatively accurate
686 * with respect to the original matrix T.
687  ibegin = 1
688  wbegin = 1
689  DO 39 jblk = 1, iwork( iindbl+m-1 )
690  iend = iwork( iinspl+jblk-1 )
691  in = iend - ibegin + 1
692  wend = wbegin - 1
693 * check if any eigenvalues have to be refined in this block
694  36 CONTINUE
695  IF( wend.LT.m ) THEN
696  IF( iwork( iindbl+wend ).EQ.jblk ) THEN
697  wend = wend + 1
698  GO TO 36
699  END IF
700  END IF
701  IF( wend.LT.wbegin ) THEN
702  ibegin = iend + 1
703  GO TO 39
704  END IF
705
706  offset = iwork(iindw+wbegin-1)-1
707  ifirst = iwork(iindw+wbegin-1)
708  ilast = iwork(iindw+wend-1)
709  rtol2 = four * eps
710  CALL dlarrj( in,
711  \$ work(indd+ibegin-1), work(inde2+ibegin-1),
712  \$ ifirst, ilast, rtol2, offset, w(wbegin),
713  \$ work( inderr+wbegin-1 ),
714  \$ work( indwrk ), iwork( iindwk ), pivmin,
715  \$ tnrm, iinfo )
716  ibegin = iend + 1
717  wbegin = wend + 1
718  39 CONTINUE
719  ENDIF
720 *
721 * If matrix was scaled, then rescale eigenvalues appropriately.
722 *
723  IF( scale.NE.one ) THEN
724  CALL dscal( m, one / scale, w, 1 )
725  END IF
726
727  END IF
728
729 *
730 * If eigenvalues are not in increasing order, then sort them,
731 * possibly along with eigenvectors.
732 *
733  IF( nsplit.GT.1 .OR. n.EQ.2 ) THEN
734  IF( .NOT. wantz ) THEN
735  CALL dlasrt( 'I', m, w, iinfo )
736  IF( iinfo.NE.0 ) THEN
737  info = 3
738  RETURN
739  END IF
740  ELSE
741  DO 60 j = 1, m - 1
742  i = 0
743  tmp = w( j )
744  DO 50 jj = j + 1, m
745  IF( w( jj ).LT.tmp ) THEN
746  i = jj
747  tmp = w( jj )
748  END IF
749  50 CONTINUE
750  IF( i.NE.0 ) THEN
751  w( i ) = w( j )
752  w( j ) = tmp
753  IF( wantz ) THEN
754  CALL dswap( n, z( 1, i ), 1, z( 1, j ), 1 )
755  itmp = isuppz( 2*i-1 )
756  isuppz( 2*i-1 ) = isuppz( 2*j-1 )
757  isuppz( 2*j-1 ) = itmp
758  itmp = isuppz( 2*i )
759  isuppz( 2*i ) = isuppz( 2*j )
760  isuppz( 2*j ) = itmp
761  END IF
762  END IF
763  60 CONTINUE
764  END IF
765  ENDIF
766 *
767 *
768  work( 1 ) = lwmin
769  iwork( 1 ) = liwmin
770  RETURN
771 *
772 * End of DSTEMR
773 *
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 dlaev2(A, B, C, RT1, RT2, CS1, SN1)
DLAEV2 computes the eigenvalues and eigenvectors of a 2-by-2 symmetric/Hermitian matrix.
Definition: dlaev2.f:120
subroutine dlarrj(N, D, E2, IFIRST, ILAST, RTOL, OFFSET, W, WERR, WORK, IWORK, PIVMIN, SPDIAM, INFO)
DLARRJ performs refinement of the initial estimates of the eigenvalues of the matrix T.
Definition: dlarrj.f:168
subroutine dlae2(A, B, C, RT1, RT2)
DLAE2 computes the eigenvalues of a 2-by-2 symmetric matrix.
Definition: dlae2.f:102
subroutine dlarrc(JOBT, N, VL, VU, D, E, PIVMIN, EIGCNT, LCNT, RCNT, INFO)
DLARRC computes the number of eigenvalues of the symmetric tridiagonal matrix.
Definition: dlarrc.f:137
subroutine dlarre(RANGE, N, VL, VU, IL, IU, D, E, E2, RTOL1, RTOL2, SPLTOL, NSPLIT, ISPLIT, M, W, WERR, WGAP, IBLOCK, INDEXW, GERS, PIVMIN, WORK, IWORK, INFO)
DLARRE given the tridiagonal matrix T, sets small off-diagonal elements to zero and for each unreduce...
Definition: dlarre.f:305
subroutine dlarrr(N, D, E, INFO)
DLARRR performs tests to decide whether the symmetric tridiagonal matrix T warrants expensive computa...
Definition: dlarrr.f:94
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
subroutine dlasrt(ID, N, D, INFO)
DLASRT sorts numbers in increasing or decreasing order.
Definition: dlasrt.f:88
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 dlarrv(N, VL, VU, D, L, PIVMIN, ISPLIT, M, DOL, DOU, MINRGP, RTOL1, RTOL2, W, WERR, WGAP, IBLOCK, INDEXW, GERS, Z, LDZ, ISUPPZ, WORK, IWORK, INFO)
DLARRV computes the eigenvectors of the tridiagonal matrix T = L D LT given L, D and the eigenvalues ...
Definition: dlarrv.f:292
Here is the call graph for this function:
Here is the caller graph for this function: