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

1.8     ! bertrand    1: *> \brief <b> DSPEV 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 DSPEV + dependencies 
        !            10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dspev.f"> 
        !            11: *> [TGZ]</a> 
        !            12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dspev.f"> 
        !            13: *> [ZIP]</a> 
        !            14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dspev.f"> 
        !            15: *> [TXT]</a>
        !            16: *> \endhtmlonly 
        !            17: *
        !            18: *  Definition:
        !            19: *  ===========
        !            20: *
        !            21: *       SUBROUTINE DSPEV( JOBZ, UPLO, N, AP, W, Z, LDZ, WORK, INFO )
        !            22: * 
        !            23: *       .. Scalar Arguments ..
        !            24: *       CHARACTER          JOBZ, UPLO
        !            25: *       INTEGER            INFO, LDZ, N
        !            26: *       ..
        !            27: *       .. Array Arguments ..
        !            28: *       DOUBLE PRECISION   AP( * ), W( * ), WORK( * ), Z( LDZ, * )
        !            29: *       ..
        !            30: *  
        !            31: *
        !            32: *> \par Purpose:
        !            33: *  =============
        !            34: *>
        !            35: *> \verbatim
        !            36: *>
        !            37: *> DSPEV computes all the eigenvalues and, optionally, eigenvectors of a
        !            38: *> real symmetric matrix A in packed storage.
        !            39: *> \endverbatim
        !            40: *
        !            41: *  Arguments:
        !            42: *  ==========
        !            43: *
        !            44: *> \param[in] JOBZ
        !            45: *> \verbatim
        !            46: *>          JOBZ is CHARACTER*1
        !            47: *>          = 'N':  Compute eigenvalues only;
        !            48: *>          = 'V':  Compute eigenvalues and eigenvectors.
        !            49: *> \endverbatim
        !            50: *>
        !            51: *> \param[in] UPLO
        !            52: *> \verbatim
        !            53: *>          UPLO is CHARACTER*1
        !            54: *>          = 'U':  Upper triangle of A is stored;
        !            55: *>          = 'L':  Lower triangle of A is stored.
        !            56: *> \endverbatim
        !            57: *>
        !            58: *> \param[in] N
        !            59: *> \verbatim
        !            60: *>          N is INTEGER
        !            61: *>          The order of the matrix A.  N >= 0.
        !            62: *> \endverbatim
        !            63: *>
        !            64: *> \param[in,out] AP
        !            65: *> \verbatim
        !            66: *>          AP is DOUBLE PRECISION array, dimension (N*(N+1)/2)
        !            67: *>          On entry, the upper or lower triangle of the symmetric matrix
        !            68: *>          A, packed columnwise in a linear array.  The j-th column of A
        !            69: *>          is stored in the array AP as follows:
        !            70: *>          if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j;
        !            71: *>          if UPLO = 'L', AP(i + (j-1)*(2*n-j)/2) = A(i,j) for j<=i<=n.
        !            72: *>
        !            73: *>          On exit, AP is overwritten by values generated during the
        !            74: *>          reduction to tridiagonal form.  If UPLO = 'U', the diagonal
        !            75: *>          and first superdiagonal of the tridiagonal matrix T overwrite
        !            76: *>          the corresponding elements of A, and if UPLO = 'L', the
        !            77: *>          diagonal and first subdiagonal of T overwrite the
        !            78: *>          corresponding elements of A.
        !            79: *> \endverbatim
        !            80: *>
        !            81: *> \param[out] W
        !            82: *> \verbatim
        !            83: *>          W is DOUBLE PRECISION array, dimension (N)
        !            84: *>          If INFO = 0, the eigenvalues in ascending order.
        !            85: *> \endverbatim
        !            86: *>
        !            87: *> \param[out] Z
        !            88: *> \verbatim
        !            89: *>          Z is DOUBLE PRECISION array, dimension (LDZ, N)
        !            90: *>          If JOBZ = 'V', then if INFO = 0, Z contains the orthonormal
        !            91: *>          eigenvectors of the matrix A, with the i-th column of Z
        !            92: *>          holding the eigenvector associated with W(i).
        !            93: *>          If JOBZ = 'N', then Z is not referenced.
        !            94: *> \endverbatim
        !            95: *>
        !            96: *> \param[in] LDZ
        !            97: *> \verbatim
        !            98: *>          LDZ is INTEGER
        !            99: *>          The leading dimension of the array Z.  LDZ >= 1, and if
        !           100: *>          JOBZ = 'V', LDZ >= max(1,N).
        !           101: *> \endverbatim
        !           102: *>
        !           103: *> \param[out] WORK
        !           104: *> \verbatim
        !           105: *>          WORK is DOUBLE PRECISION array, dimension (3*N)
        !           106: *> \endverbatim
        !           107: *>
        !           108: *> \param[out] INFO
        !           109: *> \verbatim
        !           110: *>          INFO is INTEGER
        !           111: *>          = 0:  successful exit.
        !           112: *>          < 0:  if INFO = -i, the i-th argument had an illegal value.
        !           113: *>          > 0:  if INFO = i, the algorithm failed to converge; i
        !           114: *>                off-diagonal elements of an intermediate tridiagonal
        !           115: *>                form did not converge to zero.
        !           116: *> \endverbatim
        !           117: *
        !           118: *  Authors:
        !           119: *  ========
        !           120: *
        !           121: *> \author Univ. of Tennessee 
        !           122: *> \author Univ. of California Berkeley 
        !           123: *> \author Univ. of Colorado Denver 
        !           124: *> \author NAG Ltd. 
        !           125: *
        !           126: *> \date November 2011
        !           127: *
        !           128: *> \ingroup doubleOTHEReigen
        !           129: *
        !           130: *  =====================================================================
1.1       bertrand  131:       SUBROUTINE DSPEV( JOBZ, UPLO, N, AP, W, Z, LDZ, WORK, INFO )
                    132: *
1.8     ! bertrand  133: *  -- LAPACK driver routine (version 3.4.0) --
1.1       bertrand  134: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
                    135: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
1.8     ! bertrand  136: *     November 2011
1.1       bertrand  137: *
                    138: *     .. Scalar Arguments ..
                    139:       CHARACTER          JOBZ, UPLO
                    140:       INTEGER            INFO, LDZ, N
                    141: *     ..
                    142: *     .. Array Arguments ..
                    143:       DOUBLE PRECISION   AP( * ), W( * ), WORK( * ), Z( LDZ, * )
                    144: *     ..
                    145: *
                    146: *  =====================================================================
                    147: *
                    148: *     .. Parameters ..
                    149:       DOUBLE PRECISION   ZERO, ONE
                    150:       PARAMETER          ( ZERO = 0.0D0, ONE = 1.0D0 )
                    151: *     ..
                    152: *     .. Local Scalars ..
                    153:       LOGICAL            WANTZ
                    154:       INTEGER            IINFO, IMAX, INDE, INDTAU, INDWRK, ISCALE
                    155:       DOUBLE PRECISION   ANRM, BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA,
                    156:      $                   SMLNUM
                    157: *     ..
                    158: *     .. External Functions ..
                    159:       LOGICAL            LSAME
                    160:       DOUBLE PRECISION   DLAMCH, DLANSP
                    161:       EXTERNAL           LSAME, DLAMCH, DLANSP
                    162: *     ..
                    163: *     .. External Subroutines ..
                    164:       EXTERNAL           DOPGTR, DSCAL, DSPTRD, DSTEQR, DSTERF, XERBLA
                    165: *     ..
                    166: *     .. Intrinsic Functions ..
                    167:       INTRINSIC          SQRT
                    168: *     ..
                    169: *     .. Executable Statements ..
                    170: *
                    171: *     Test the input parameters.
                    172: *
                    173:       WANTZ = LSAME( JOBZ, 'V' )
                    174: *
                    175:       INFO = 0
                    176:       IF( .NOT.( WANTZ .OR. LSAME( JOBZ, 'N' ) ) ) THEN
                    177:          INFO = -1
                    178:       ELSE IF( .NOT.( LSAME( UPLO, 'U' ) .OR. LSAME( UPLO, 'L' ) ) )
                    179:      $          THEN
                    180:          INFO = -2
                    181:       ELSE IF( N.LT.0 ) THEN
                    182:          INFO = -3
                    183:       ELSE IF( LDZ.LT.1 .OR. ( WANTZ .AND. LDZ.LT.N ) ) THEN
                    184:          INFO = -7
                    185:       END IF
                    186: *
                    187:       IF( INFO.NE.0 ) THEN
                    188:          CALL XERBLA( 'DSPEV ', -INFO )
                    189:          RETURN
                    190:       END IF
                    191: *
                    192: *     Quick return if possible
                    193: *
                    194:       IF( N.EQ.0 )
                    195:      $   RETURN
                    196: *
                    197:       IF( N.EQ.1 ) THEN
                    198:          W( 1 ) = AP( 1 )
                    199:          IF( WANTZ )
                    200:      $      Z( 1, 1 ) = ONE
                    201:          RETURN
                    202:       END IF
                    203: *
                    204: *     Get machine constants.
                    205: *
                    206:       SAFMIN = DLAMCH( 'Safe minimum' )
                    207:       EPS = DLAMCH( 'Precision' )
                    208:       SMLNUM = SAFMIN / EPS
                    209:       BIGNUM = ONE / SMLNUM
                    210:       RMIN = SQRT( SMLNUM )
                    211:       RMAX = SQRT( BIGNUM )
                    212: *
                    213: *     Scale matrix to allowable range, if necessary.
                    214: *
                    215:       ANRM = DLANSP( 'M', UPLO, N, AP, WORK )
                    216:       ISCALE = 0
                    217:       IF( ANRM.GT.ZERO .AND. ANRM.LT.RMIN ) THEN
                    218:          ISCALE = 1
                    219:          SIGMA = RMIN / ANRM
                    220:       ELSE IF( ANRM.GT.RMAX ) THEN
                    221:          ISCALE = 1
                    222:          SIGMA = RMAX / ANRM
                    223:       END IF
                    224:       IF( ISCALE.EQ.1 ) THEN
                    225:          CALL DSCAL( ( N*( N+1 ) ) / 2, SIGMA, AP, 1 )
                    226:       END IF
                    227: *
                    228: *     Call DSPTRD to reduce symmetric packed matrix to tridiagonal form.
                    229: *
                    230:       INDE = 1
                    231:       INDTAU = INDE + N
                    232:       CALL DSPTRD( UPLO, N, AP, W, WORK( INDE ), WORK( INDTAU ), IINFO )
                    233: *
                    234: *     For eigenvalues only, call DSTERF.  For eigenvectors, first call
                    235: *     DOPGTR to generate the orthogonal matrix, then call DSTEQR.
                    236: *
                    237:       IF( .NOT.WANTZ ) THEN
                    238:          CALL DSTERF( N, W, WORK( INDE ), INFO )
                    239:       ELSE
                    240:          INDWRK = INDTAU + N
                    241:          CALL DOPGTR( UPLO, N, AP, WORK( INDTAU ), Z, LDZ,
                    242:      $                WORK( INDWRK ), IINFO )
                    243:          CALL DSTEQR( JOBZ, N, W, WORK( INDE ), Z, LDZ, WORK( INDTAU ),
                    244:      $                INFO )
                    245:       END IF
                    246: *
                    247: *     If matrix was scaled, then rescale eigenvalues appropriately.
                    248: *
                    249:       IF( ISCALE.EQ.1 ) THEN
                    250:          IF( INFO.EQ.0 ) THEN
                    251:             IMAX = N
                    252:          ELSE
                    253:             IMAX = INFO - 1
                    254:          END IF
                    255:          CALL DSCAL( IMAX, ONE / SIGMA, W, 1 )
                    256:       END IF
                    257: *
                    258:       RETURN
                    259: *
                    260: *     End of DSPEV
                    261: *
                    262:       END

CVSweb interface <joel.bertrand@systella.fr>