LAPACK: SRC/cheev.f Source File

Go to the documentation of this file.
00001       SUBROUTINE CHEEV( JOBZ, UPLO, N, A, LDA, W, WORK, LWORK, RWORK,
00002      $                  INFO )
00003 *
00004 *  -- LAPACK driver routine (version 3.2) --
00005 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
00006 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
00007 *     November 2006
00008 *
00009 *     .. Scalar Arguments ..
00010       CHARACTER          JOBZ, UPLO
00011       INTEGER            INFO, LDA, LWORK, N
00012 *     ..
00013 *     .. Array Arguments ..
00014       REAL               RWORK( * ), W( * )
00015       COMPLEX            A( LDA, * ), WORK( * )
00016 *     ..
00017 *
00018 *  Purpose
00019 *  =======
00020 *
00021 *  CHEEV computes all eigenvalues and, optionally, eigenvectors of a
00022 *  complex Hermitian matrix A.
00023 *
00024 *  Arguments
00025 *  =========
00026 *
00027 *  JOBZ    (input) CHARACTER*1
00028 *          = 'N':  Compute eigenvalues only;
00029 *          = 'V':  Compute eigenvalues and eigenvectors.
00030 *
00031 *  UPLO    (input) CHARACTER*1
00032 *          = 'U':  Upper triangle of A is stored;
00033 *          = 'L':  Lower triangle of A is stored.
00034 *
00035 *  N       (input) INTEGER
00036 *          The order of the matrix A.  N >= 0.
00037 *
00038 *  A       (input/output) COMPLEX array, dimension (LDA, N)
00039 *          On entry, the Hermitian matrix A.  If UPLO = 'U', the
00040 *          leading N-by-N upper triangular part of A contains the
00041 *          upper triangular part of the matrix A.  If UPLO = 'L',
00042 *          the leading N-by-N lower triangular part of A contains
00043 *          the lower triangular part of the matrix A.
00044 *          On exit, if JOBZ = 'V', then if INFO = 0, A contains the
00045 *          orthonormal eigenvectors of the matrix A.
00046 *          If JOBZ = 'N', then on exit the lower triangle (if UPLO='L')
00047 *          or the upper triangle (if UPLO='U') of A, including the
00048 *          diagonal, is destroyed.
00049 *
00050 *  LDA     (input) INTEGER
00051 *          The leading dimension of the array A.  LDA >= max(1,N).
00052 *
00053 *  W       (output) REAL array, dimension (N)
00054 *          If INFO = 0, the eigenvalues in ascending order.
00055 *
00056 *  WORK    (workspace/output) COMPLEX array, dimension (MAX(1,LWORK))
00057 *          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
00058 *
00059 *  LWORK   (input) INTEGER
00060 *          The length of the array WORK.  LWORK >= max(1,2*N-1).
00061 *          For optimal efficiency, LWORK >= (NB+1)*N,
00062 *          where NB is the blocksize for CHETRD returned by ILAENV.
00063 *
00064 *          If LWORK = -1, then a workspace query is assumed; the routine
00065 *          only calculates the optimal size of the WORK array, returns
00066 *          this value as the first entry of the WORK array, and no error
00067 *          message related to LWORK is issued by XERBLA.
00068 *
00069 *  RWORK   (workspace) REAL array, dimension (max(1, 3*N-2))
00070 *
00071 *  INFO    (output) INTEGER
00072 *          = 0:  successful exit
00073 *          < 0:  if INFO = -i, the i-th argument had an illegal value
00074 *          > 0:  if INFO = i, the algorithm failed to converge; i
00075 *                off-diagonal elements of an intermediate tridiagonal
00076 *                form did not converge to zero.
00077 *
00078 *  =====================================================================
00079 *
00080 *     .. Parameters ..
00081       REAL               ZERO, ONE
00082       PARAMETER          ( ZERO = 0.0E0, ONE = 1.0E0 )
00083       COMPLEX            CONE
00084       PARAMETER          ( CONE = ( 1.0E0, 0.0E0 ) )
00085 *     ..
00086 *     .. Local Scalars ..
00087       LOGICAL            LOWER, LQUERY, WANTZ
00088       INTEGER            IINFO, IMAX, INDE, INDTAU, INDWRK, ISCALE,
00089      $                   LLWORK, LWKOPT, NB
00090       REAL               ANRM, BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA,
00091      $                   SMLNUM
00092 *     ..
00093 *     .. External Functions ..
00094       LOGICAL            LSAME
00095       INTEGER            ILAENV
00096       REAL               CLANHE, SLAMCH
00097       EXTERNAL           ILAENV, LSAME, CLANHE, SLAMCH
00098 *     ..
00099 *     .. External Subroutines ..
00100       EXTERNAL           CHETRD, CLASCL, CSTEQR, CUNGTR, SSCAL, SSTERF,
00101      $                   XERBLA
00102 *     ..
00103 *     .. Intrinsic Functions ..
00104       INTRINSIC          MAX, SQRT
00105 *     ..
00106 *     .. Executable Statements ..
00107 *
00108 *     Test the input parameters.
00109 *
00110       WANTZ = LSAME( JOBZ, 'V' )
00111       LOWER = LSAME( UPLO, 'L' )
00112       LQUERY = ( LWORK.EQ.-1 )
00113 *
00114       INFO = 0
00115       IF( .NOT.( WANTZ .OR. LSAME( JOBZ, 'N' ) ) ) THEN
00116          INFO = -1
00117       ELSE IF( .NOT.( LOWER .OR. LSAME( UPLO, 'U' ) ) ) THEN
00118          INFO = -2
00119       ELSE IF( N.LT.0 ) THEN
00120          INFO = -3
00121       ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
00122          INFO = -5
00123       END IF
00124 *
00125       IF( INFO.EQ.0 ) THEN
00126          NB = ILAENV( 1, 'CHETRD', UPLO, N, -1, -1, -1 )
00127          LWKOPT = MAX( 1, ( NB+1 )*N )
00128          WORK( 1 ) = LWKOPT
00129 *
00130          IF( LWORK.LT.MAX( 1, 2*N-1 ) .AND. .NOT.LQUERY )
00131      $      INFO = -8
00132       END IF
00133 *
00134       IF( INFO.NE.0 ) THEN
00135          CALL XERBLA( 'CHEEV ', -INFO )
00136          RETURN
00137       ELSE IF( LQUERY ) THEN
00138          RETURN
00139       END IF
00140 *
00141 *     Quick return if possible
00142 *
00143       IF( N.EQ.0 ) THEN
00144          RETURN
00145       END IF
00146 *
00147       IF( N.EQ.1 ) THEN
00148          W( 1 ) = A( 1, 1 )
00149          WORK( 1 ) = 1
00150          IF( WANTZ )
00151      $      A( 1, 1 ) = CONE
00152          RETURN
00153       END IF
00154 *
00155 *     Get machine constants.
00156 *
00157       SAFMIN = SLAMCH( 'Safe minimum' )
00158       EPS = SLAMCH( 'Precision' )
00159       SMLNUM = SAFMIN / EPS
00160       BIGNUM = ONE / SMLNUM
00161       RMIN = SQRT( SMLNUM )
00162       RMAX = SQRT( BIGNUM )
00163 *
00164 *     Scale matrix to allowable range, if necessary.
00165 *
00166       ANRM = CLANHE( 'M', UPLO, N, A, LDA, RWORK )
00167       ISCALE = 0
00168       IF( ANRM.GT.ZERO .AND. ANRM.LT.RMIN ) THEN
00169          ISCALE = 1
00170          SIGMA = RMIN / ANRM
00171       ELSE IF( ANRM.GT.RMAX ) THEN
00172          ISCALE = 1
00173          SIGMA = RMAX / ANRM
00174       END IF
00175       IF( ISCALE.EQ.1 )
00176      $   CALL CLASCL( UPLO, 0, 0, ONE, SIGMA, N, N, A, LDA, INFO )
00177 *
00178 *     Call CHETRD to reduce Hermitian matrix to tridiagonal form.
00179 *
00180       INDE = 1
00181       INDTAU = 1
00182       INDWRK = INDTAU + N
00183       LLWORK = LWORK - INDWRK + 1
00184       CALL CHETRD( UPLO, N, A, LDA, W, RWORK( INDE ), WORK( INDTAU ),
00185      $             WORK( INDWRK ), LLWORK, IINFO )
00186 *
00187 *     For eigenvalues only, call SSTERF.  For eigenvectors, first call
00188 *     CUNGTR to generate the unitary matrix, then call CSTEQR.
00189 *
00190       IF( .NOT.WANTZ ) THEN
00191          CALL SSTERF( N, W, RWORK( INDE ), INFO )
00192       ELSE
00193          CALL CUNGTR( UPLO, N, A, LDA, WORK( INDTAU ), WORK( INDWRK ),
00194      $                LLWORK, IINFO )
00195          INDWRK = INDE + N
00196          CALL CSTEQR( JOBZ, N, W, RWORK( INDE ), A, LDA,
00197      $                RWORK( INDWRK ), INFO )
00198       END IF
00199 *
00200 *     If matrix was scaled, then rescale eigenvalues appropriately.
00201 *
00202       IF( ISCALE.EQ.1 ) THEN
00203          IF( INFO.EQ.0 ) THEN
00204             IMAX = N
00205          ELSE
00206             IMAX = INFO - 1
00207          END IF
00208          CALL SSCAL( IMAX, ONE / SIGMA, W, 1 )
00209       END IF
00210 *
00211 *     Set WORK(1) to optimal complex workspace size.
00212 *
00213       WORK( 1 ) = LWKOPT
00214 *
00215       RETURN
00216 *
00217 *     End of CHEEV
00218 *
00219       END