Annotation of rpl/lapack/lapack/dstevd.f, revision 1.8

1.8     ! bertrand    1: *> \brief <b> DSTEVD computes the eigenvalues and, optionally, the left and/or right eigenvectors for OTHER matrices</b>
        !             2: *
        !             3: *  =========== DOCUMENTATION ===========
        !             4: *
        !             5: * Online html documentation available at 
        !             6: *            http://www.netlib.org/lapack/explore-html/ 
        !             7: *
        !             8: *> \htmlonly
        !             9: *> Download DSTEVD + dependencies 
        !            10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dstevd.f"> 
        !            11: *> [TGZ]</a> 
        !            12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dstevd.f"> 
        !            13: *> [ZIP]</a> 
        !            14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dstevd.f"> 
        !            15: *> [TXT]</a>
        !            16: *> \endhtmlonly 
        !            17: *
        !            18: *  Definition:
        !            19: *  ===========
        !            20: *
        !            21: *       SUBROUTINE DSTEVD( JOBZ, N, D, E, Z, LDZ, WORK, LWORK, IWORK,
        !            22: *                          LIWORK, INFO )
        !            23: * 
        !            24: *       .. Scalar Arguments ..
        !            25: *       CHARACTER          JOBZ
        !            26: *       INTEGER            INFO, LDZ, LIWORK, LWORK, N
        !            27: *       ..
        !            28: *       .. Array Arguments ..
        !            29: *       INTEGER            IWORK( * )
        !            30: *       DOUBLE PRECISION   D( * ), E( * ), WORK( * ), Z( LDZ, * )
        !            31: *       ..
        !            32: *  
        !            33: *
        !            34: *> \par Purpose:
        !            35: *  =============
        !            36: *>
        !            37: *> \verbatim
        !            38: *>
        !            39: *> DSTEVD computes all eigenvalues and, optionally, eigenvectors of a
        !            40: *> real symmetric tridiagonal matrix. If eigenvectors are desired, it
        !            41: *> uses a divide and conquer algorithm.
        !            42: *>
        !            43: *> The divide and conquer algorithm makes very mild assumptions about
        !            44: *> floating point arithmetic. It will work on machines with a guard
        !            45: *> digit in add/subtract, or on those binary machines without guard
        !            46: *> digits which subtract like the Cray X-MP, Cray Y-MP, Cray C-90, or
        !            47: *> Cray-2. It could conceivably fail on hexadecimal or decimal machines
        !            48: *> without guard digits, but we know of none.
        !            49: *> \endverbatim
        !            50: *
        !            51: *  Arguments:
        !            52: *  ==========
        !            53: *
        !            54: *> \param[in] JOBZ
        !            55: *> \verbatim
        !            56: *>          JOBZ is CHARACTER*1
        !            57: *>          = 'N':  Compute eigenvalues only;
        !            58: *>          = 'V':  Compute eigenvalues and eigenvectors.
        !            59: *> \endverbatim
        !            60: *>
        !            61: *> \param[in] N
        !            62: *> \verbatim
        !            63: *>          N is INTEGER
        !            64: *>          The order of the matrix.  N >= 0.
        !            65: *> \endverbatim
        !            66: *>
        !            67: *> \param[in,out] D
        !            68: *> \verbatim
        !            69: *>          D is DOUBLE PRECISION array, dimension (N)
        !            70: *>          On entry, the n diagonal elements of the tridiagonal matrix
        !            71: *>          A.
        !            72: *>          On exit, if INFO = 0, the eigenvalues in ascending order.
        !            73: *> \endverbatim
        !            74: *>
        !            75: *> \param[in,out] E
        !            76: *> \verbatim
        !            77: *>          E is DOUBLE PRECISION array, dimension (N-1)
        !            78: *>          On entry, the (n-1) subdiagonal elements of the tridiagonal
        !            79: *>          matrix A, stored in elements 1 to N-1 of E.
        !            80: *>          On exit, the contents of E are destroyed.
        !            81: *> \endverbatim
        !            82: *>
        !            83: *> \param[out] Z
        !            84: *> \verbatim
        !            85: *>          Z is DOUBLE PRECISION array, dimension (LDZ, N)
        !            86: *>          If JOBZ = 'V', then if INFO = 0, Z contains the orthonormal
        !            87: *>          eigenvectors of the matrix A, with the i-th column of Z
        !            88: *>          holding the eigenvector associated with D(i).
        !            89: *>          If JOBZ = 'N', then Z is not referenced.
        !            90: *> \endverbatim
        !            91: *>
        !            92: *> \param[in] LDZ
        !            93: *> \verbatim
        !            94: *>          LDZ is INTEGER
        !            95: *>          The leading dimension of the array Z.  LDZ >= 1, and if
        !            96: *>          JOBZ = 'V', LDZ >= max(1,N).
        !            97: *> \endverbatim
        !            98: *>
        !            99: *> \param[out] WORK
        !           100: *> \verbatim
        !           101: *>          WORK is DOUBLE PRECISION array,
        !           102: *>                                         dimension (LWORK)
        !           103: *>          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
        !           104: *> \endverbatim
        !           105: *>
        !           106: *> \param[in] LWORK
        !           107: *> \verbatim
        !           108: *>          LWORK is INTEGER
        !           109: *>          The dimension of the array WORK.
        !           110: *>          If JOBZ  = 'N' or N <= 1 then LWORK must be at least 1.
        !           111: *>          If JOBZ  = 'V' and N > 1 then LWORK must be at least
        !           112: *>                         ( 1 + 4*N + N**2 ).
        !           113: *>
        !           114: *>          If LWORK = -1, then a workspace query is assumed; the routine
        !           115: *>          only calculates the optimal sizes of the WORK and IWORK
        !           116: *>          arrays, returns these values as the first entries of the WORK
        !           117: *>          and IWORK arrays, and no error message related to LWORK or
        !           118: *>          LIWORK is issued by XERBLA.
        !           119: *> \endverbatim
        !           120: *>
        !           121: *> \param[out] IWORK
        !           122: *> \verbatim
        !           123: *>          IWORK is INTEGER array, dimension (MAX(1,LIWORK))
        !           124: *>          On exit, if INFO = 0, IWORK(1) returns the optimal LIWORK.
        !           125: *> \endverbatim
        !           126: *>
        !           127: *> \param[in] LIWORK
        !           128: *> \verbatim
        !           129: *>          LIWORK is INTEGER
        !           130: *>          The dimension of the array IWORK.
        !           131: *>          If JOBZ  = 'N' or N <= 1 then LIWORK must be at least 1.
        !           132: *>          If JOBZ  = 'V' and N > 1 then LIWORK must be at least 3+5*N.
        !           133: *>
        !           134: *>          If LIWORK = -1, then a workspace query is assumed; the
        !           135: *>          routine only calculates the optimal sizes of the WORK and
        !           136: *>          IWORK arrays, returns these values as the first entries of
        !           137: *>          the WORK and IWORK arrays, and no error message related to
        !           138: *>          LWORK or LIWORK is issued by XERBLA.
        !           139: *> \endverbatim
        !           140: *>
        !           141: *> \param[out] INFO
        !           142: *> \verbatim
        !           143: *>          INFO is INTEGER
        !           144: *>          = 0:  successful exit
        !           145: *>          < 0:  if INFO = -i, the i-th argument had an illegal value
        !           146: *>          > 0:  if INFO = i, the algorithm failed to converge; i
        !           147: *>                off-diagonal elements of E did not converge to zero.
        !           148: *> \endverbatim
        !           149: *
        !           150: *  Authors:
        !           151: *  ========
        !           152: *
        !           153: *> \author Univ. of Tennessee 
        !           154: *> \author Univ. of California Berkeley 
        !           155: *> \author Univ. of Colorado Denver 
        !           156: *> \author NAG Ltd. 
        !           157: *
        !           158: *> \date November 2011
        !           159: *
        !           160: *> \ingroup doubleOTHEReigen
        !           161: *
        !           162: *  =====================================================================
1.1       bertrand  163:       SUBROUTINE DSTEVD( JOBZ, N, D, E, Z, LDZ, WORK, LWORK, IWORK,
                    164:      $                   LIWORK, INFO )
                    165: *
1.8     ! bertrand  166: *  -- LAPACK driver routine (version 3.4.0) --
1.1       bertrand  167: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
                    168: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
1.8     ! bertrand  169: *     November 2011
1.1       bertrand  170: *
                    171: *     .. Scalar Arguments ..
                    172:       CHARACTER          JOBZ
                    173:       INTEGER            INFO, LDZ, LIWORK, LWORK, N
                    174: *     ..
                    175: *     .. Array Arguments ..
                    176:       INTEGER            IWORK( * )
                    177:       DOUBLE PRECISION   D( * ), E( * ), WORK( * ), Z( LDZ, * )
                    178: *     ..
                    179: *
                    180: *  =====================================================================
                    181: *
                    182: *     .. Parameters ..
                    183:       DOUBLE PRECISION   ZERO, ONE
                    184:       PARAMETER          ( ZERO = 0.0D0, ONE = 1.0D0 )
                    185: *     ..
                    186: *     .. Local Scalars ..
                    187:       LOGICAL            LQUERY, WANTZ
                    188:       INTEGER            ISCALE, LIWMIN, LWMIN
                    189:       DOUBLE PRECISION   BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA, SMLNUM,
                    190:      $                   TNRM
                    191: *     ..
                    192: *     .. External Functions ..
                    193:       LOGICAL            LSAME
                    194:       DOUBLE PRECISION   DLAMCH, DLANST
                    195:       EXTERNAL           LSAME, DLAMCH, DLANST
                    196: *     ..
                    197: *     .. External Subroutines ..
                    198:       EXTERNAL           DSCAL, DSTEDC, DSTERF, XERBLA
                    199: *     ..
                    200: *     .. Intrinsic Functions ..
                    201:       INTRINSIC          SQRT
                    202: *     ..
                    203: *     .. Executable Statements ..
                    204: *
                    205: *     Test the input parameters.
                    206: *
                    207:       WANTZ = LSAME( JOBZ, 'V' )
                    208:       LQUERY = ( LWORK.EQ.-1 .OR. LIWORK.EQ.-1 )
                    209: *
                    210:       INFO = 0
                    211:       LIWMIN = 1
                    212:       LWMIN = 1
                    213:       IF( N.GT.1 .AND. WANTZ ) THEN
                    214:          LWMIN = 1 + 4*N + N**2
                    215:          LIWMIN = 3 + 5*N
                    216:       END IF
                    217: *
                    218:       IF( .NOT.( WANTZ .OR. LSAME( JOBZ, 'N' ) ) ) THEN
                    219:          INFO = -1
                    220:       ELSE IF( N.LT.0 ) THEN
                    221:          INFO = -2
                    222:       ELSE IF( LDZ.LT.1 .OR. ( WANTZ .AND. LDZ.LT.N ) ) THEN
                    223:          INFO = -6
                    224:       END IF
                    225: *
                    226:       IF( INFO.EQ.0 ) THEN
                    227:          WORK( 1 ) = LWMIN
                    228:          IWORK( 1 ) = LIWMIN
                    229: *
                    230:          IF( LWORK.LT.LWMIN .AND. .NOT.LQUERY ) THEN
                    231:             INFO = -8
                    232:          ELSE IF( LIWORK.LT.LIWMIN .AND. .NOT.LQUERY ) THEN
                    233:             INFO = -10
                    234:          END IF
                    235:       END IF
                    236: *
                    237:       IF( INFO.NE.0 ) THEN
                    238:          CALL XERBLA( 'DSTEVD', -INFO )
                    239:          RETURN
                    240:       ELSE IF( LQUERY ) THEN
                    241:          RETURN
                    242:       END IF
                    243: *
                    244: *     Quick return if possible
                    245: *
                    246:       IF( N.EQ.0 )
                    247:      $   RETURN
                    248: *
                    249:       IF( N.EQ.1 ) THEN
                    250:          IF( WANTZ )
                    251:      $      Z( 1, 1 ) = ONE
                    252:          RETURN
                    253:       END IF
                    254: *
                    255: *     Get machine constants.
                    256: *
                    257:       SAFMIN = DLAMCH( 'Safe minimum' )
                    258:       EPS = DLAMCH( 'Precision' )
                    259:       SMLNUM = SAFMIN / EPS
                    260:       BIGNUM = ONE / SMLNUM
                    261:       RMIN = SQRT( SMLNUM )
                    262:       RMAX = SQRT( BIGNUM )
                    263: *
                    264: *     Scale matrix to allowable range, if necessary.
                    265: *
                    266:       ISCALE = 0
                    267:       TNRM = DLANST( 'M', N, D, E )
                    268:       IF( TNRM.GT.ZERO .AND. TNRM.LT.RMIN ) THEN
                    269:          ISCALE = 1
                    270:          SIGMA = RMIN / TNRM
                    271:       ELSE IF( TNRM.GT.RMAX ) THEN
                    272:          ISCALE = 1
                    273:          SIGMA = RMAX / TNRM
                    274:       END IF
                    275:       IF( ISCALE.EQ.1 ) THEN
                    276:          CALL DSCAL( N, SIGMA, D, 1 )
                    277:          CALL DSCAL( N-1, SIGMA, E( 1 ), 1 )
                    278:       END IF
                    279: *
                    280: *     For eigenvalues only, call DSTERF.  For eigenvalues and
                    281: *     eigenvectors, call DSTEDC.
                    282: *
                    283:       IF( .NOT.WANTZ ) THEN
                    284:          CALL DSTERF( N, D, E, INFO )
                    285:       ELSE
                    286:          CALL DSTEDC( 'I', N, D, E, Z, LDZ, WORK, LWORK, IWORK, LIWORK,
                    287:      $                INFO )
                    288:       END IF
                    289: *
                    290: *     If matrix was scaled, then rescale eigenvalues appropriately.
                    291: *
                    292:       IF( ISCALE.EQ.1 )
                    293:      $   CALL DSCAL( N, ONE / SIGMA, D, 1 )
                    294: *
                    295:       WORK( 1 ) = LWMIN
                    296:       IWORK( 1 ) = LIWMIN
                    297: *
                    298:       RETURN
                    299: *
                    300: *     End of DSTEVD
                    301: *
                    302:       END

CVSweb interface <joel.bertrand@systella.fr>