SRC/cheev.f

   1 *> \brief <b> CHEEV computes the eigenvalues and, optionally, the left and/or right eigenvectors for HE 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 CHEEV + dependencies
  10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/cheev.f">
  11 *> [TGZ]</a>
  12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/cheev.f">
  13 *> [ZIP]</a>
  14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/cheev.f">
  15 *> [TXT]</a>
  16 *> \endhtmlonly
  17 *
  18 *  Definition:
  19 *  ===========
  20 *
  21 *       SUBROUTINE CHEEV( JOBZ, UPLO, N, A, LDA, W, WORK, LWORK, RWORK,
  22 *                         INFO )
  23 *
  24 *       .. Scalar Arguments ..
  25 *       CHARACTER          JOBZ, UPLO
  26 *       INTEGER            INFO, LDA, LWORK, N
  27 *       ..
  28 *       .. Array Arguments ..
  29 *       REAL               RWORK( * ), W( * )
  30 *       COMPLEX            A( LDA, * ), WORK( * )
  31 *       ..
  32 *
  33 *
  34 *> \par Purpose:
  35 *  =============
  36 *>
  37 *> \verbatim
  38 *>
  39 *> CHEEV computes all eigenvalues and, optionally, eigenvectors of a
  40 *> complex Hermitian matrix A.
  41 *> \endverbatim
  42 *
  43 *  Arguments:
  44 *  ==========
  45 *
  46 *> \param[in] JOBZ
  47 *> \verbatim
  48 *>          JOBZ is CHARACTER*1
  49 *>          = 'N':  Compute eigenvalues only;
  50 *>          = 'V':  Compute eigenvalues and eigenvectors.
  51 *> \endverbatim
  52 *>
  53 *> \param[in] UPLO
  54 *> \verbatim
  55 *>          UPLO is CHARACTER*1
  56 *>          = 'U':  Upper triangle of A is stored;
  57 *>          = 'L':  Lower triangle of A is stored.
  58 *> \endverbatim
  59 *>
  60 *> \param[in] N
  61 *> \verbatim
  62 *>          N is INTEGER
  63 *>          The order of the matrix A.  N >= 0.
  64 *> \endverbatim
  65 *>
  66 *> \param[in,out] A
  67 *> \verbatim
  68 *>          A is COMPLEX array, dimension (LDA, N)
  69 *>          On entry, the Hermitian matrix A.  If UPLO = 'U', the
  70 *>          leading N-by-N upper triangular part of A contains the
  71 *>          upper triangular part of the matrix A.  If UPLO = 'L',
  72 *>          the leading N-by-N lower triangular part of A contains
  73 *>          the lower triangular part of the matrix A.
  74 *>          On exit, if JOBZ = 'V', then if INFO = 0, A contains the
  75 *>          orthonormal eigenvectors of the matrix A.
  76 *>          If JOBZ = 'N', then on exit the lower triangle (if UPLO='L')
  77 *>          or the upper triangle (if UPLO='U') of A, including the
  78 *>          diagonal, is destroyed.
  79 *> \endverbatim
  80 *>
  81 *> \param[in] LDA
  82 *> \verbatim
  83 *>          LDA is INTEGER
  84 *>          The leading dimension of the array A.  LDA >= max(1,N).
  85 *> \endverbatim
  86 *>
  87 *> \param[out] W
  88 *> \verbatim
  89 *>          W is REAL array, dimension (N)
  90 *>          If INFO = 0, the eigenvalues in ascending order.
  91 *> \endverbatim
  92 *>
  93 *> \param[out] WORK
  94 *> \verbatim
  95 *>          WORK is COMPLEX array, dimension (MAX(1,LWORK))
  96 *>          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
  97 *> \endverbatim
  98 *>
  99 *> \param[in] LWORK
 100 *> \verbatim
 101 *>          LWORK is INTEGER
 102 *>          The length of the array WORK.  LWORK >= max(1,2*N-1).
 103 *>          For optimal efficiency, LWORK >= (NB+1)*N,
 104 *>          where NB is the blocksize for CHETRD returned by ILAENV.
 105 *>
 106 *>          If LWORK = -1, then a workspace query is assumed; the routine
 107 *>          only calculates the optimal size of the WORK array, returns
 108 *>          this value as the first entry of the WORK array, and no error
 109 *>          message related to LWORK is issued by XERBLA.
 110 *> \endverbatim
 111 *>
 112 *> \param[out] RWORK
 113 *> \verbatim
 114 *>          RWORK is REAL array, dimension (max(1, 3*N-2))
 115 *> \endverbatim
 116 *>
 117 *> \param[out] INFO
 118 *> \verbatim
 119 *>          INFO is INTEGER
 120 *>          = 0:  successful exit
 121 *>          < 0:  if INFO = -i, the i-th argument had an illegal value
 122 *>          > 0:  if INFO = i, the algorithm failed to converge; i
 123 *>                off-diagonal elements of an intermediate tridiagonal
 124 *>                form did not converge to zero.
 125 *> \endverbatim
 126 *
 127 *  Authors:
 128 *  ========
 129 *
 130 *> \author Univ. of Tennessee
 131 *> \author Univ. of California Berkeley
 132 *> \author Univ. of Colorado Denver
 133 *> \author NAG Ltd.
 134 *
 135 *> \date November 2011
 136 *
 137 *> \ingroup complexHEeigen
 138 *
 139 *  =====================================================================
 140       SUBROUTINE CHEEV( JOBZ, UPLO, N, A, LDA, W, WORK, LWORK, RWORK,
 141      $                  INFO )
 142 *
 143 *  -- LAPACK driver routine (version 3.4.0) --
 144 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 145 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 146 *     November 2011
 147 *
 148 *     .. Scalar Arguments ..
 149       CHARACTER          JOBZ, UPLO
 150       INTEGER            INFO, LDA, LWORK, N
 151 *     ..
 152 *     .. Array Arguments ..
 153       REAL               RWORK( * ), W( * )
 154       COMPLEX            A( LDA, * ), WORK( * )
 155 *     ..
 156 *
 157 *  =====================================================================
 158 *
 159 *     .. Parameters ..
 160       REAL               ZERO, ONE
 161       PARAMETER          ( ZERO = 0.0E0, ONE = 1.0E0 )
 162       COMPLEX            CONE
 163       PARAMETER          ( CONE = ( 1.0E0, 0.0E0 ) )
 164 *     ..
 165 *     .. Local Scalars ..
 166       LOGICAL            LOWER, LQUERY, WANTZ
 167       INTEGER            IINFO, IMAX, INDE, INDTAU, INDWRK, ISCALE,
 168      $                   LLWORK, LWKOPT, NB
 169       REAL               ANRM, BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA,
 170      $                   SMLNUM
 171 *     ..
 172 *     .. External Functions ..
 173       LOGICAL            LSAME
 174       INTEGER            ILAENV
 175       REAL               CLANHE, SLAMCH
 176       EXTERNAL           ILAENV, LSAME, CLANHE, SLAMCH
 177 *     ..
 178 *     .. External Subroutines ..
 179       EXTERNAL           CHETRD, CLASCL, CSTEQR, CUNGTR, SSCAL, SSTERF,
 180      $                   XERBLA
 181 *     ..
 182 *     .. Intrinsic Functions ..
 183       INTRINSIC          MAX, SQRT
 184 *     ..
 185 *     .. Executable Statements ..
 186 *
 187 *     Test the input parameters.
 188 *
 189       WANTZ = LSAME( JOBZ, 'V' )
 190       LOWER = LSAME( UPLO, 'L' )
 191       LQUERY = ( LWORK.EQ.-1 )
 192 *
 193       INFO = 0
 194       IF( .NOT.( WANTZ .OR. LSAME( JOBZ, 'N' ) ) ) THEN
 195          INFO = -1
 196       ELSE IF( .NOT.( LOWER .OR. LSAME( UPLO, 'U' ) ) ) THEN
 197          INFO = -2
 198       ELSE IF( N.LT.0 ) THEN
 199          INFO = -3
 200       ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
 201          INFO = -5
 202       END IF
 203 *
 204       IF( INFO.EQ.0 ) THEN
 205          NB = ILAENV( 1, 'CHETRD', UPLO, N, -1, -1, -1 )
 206          LWKOPT = MAX( 1, ( NB+1 )*N )
 207          WORK( 1 ) = LWKOPT
 208 *
 209          IF( LWORK.LT.MAX( 1, 2*N-1 ) .AND. .NOT.LQUERY )
 210      $      INFO = -8
 211       END IF
 212 *
 213       IF( INFO.NE.0 ) THEN
 214          CALL XERBLA( 'CHEEV ', -INFO )
 215          RETURN
 216       ELSE IF( LQUERY ) THEN
 217          RETURN
 218       END IF
 219 *
 220 *     Quick return if possible
 221 *
 222       IF( N.EQ.0 ) THEN
 223          RETURN
 224       END IF
 225 *
 226       IF( N.EQ.1 ) THEN
 227          W( 1 ) = A( 1, 1 )
 228          WORK( 1 ) = 1
 229          IF( WANTZ )
 230      $      A( 1, 1 ) = CONE
 231          RETURN
 232       END IF
 233 *
 234 *     Get machine constants.
 235 *
 236       SAFMIN = SLAMCH( 'Safe minimum' )
 237       EPS = SLAMCH( 'Precision' )
 238       SMLNUM = SAFMIN / EPS
 239       BIGNUM = ONE / SMLNUM
 240       RMIN = SQRT( SMLNUM )
 241       RMAX = SQRT( BIGNUM )
 242 *
 243 *     Scale matrix to allowable range, if necessary.
 244 *
 245       ANRM = CLANHE( 'M', UPLO, N, A, LDA, RWORK )
 246       ISCALE = 0
 247       IF( ANRM.GT.ZERO .AND. ANRM.LT.RMIN ) THEN
 248          ISCALE = 1
 249          SIGMA = RMIN / ANRM
 250       ELSE IF( ANRM.GT.RMAX ) THEN
 251          ISCALE = 1
 252          SIGMA = RMAX / ANRM
 253       END IF
 254       IF( ISCALE.EQ.1 )
 255      $   CALL CLASCL( UPLO, 0, 0, ONE, SIGMA, N, N, A, LDA, INFO )
 256 *
 257 *     Call CHETRD to reduce Hermitian matrix to tridiagonal form.
 258 *
 259       INDE = 1
 260       INDTAU = 1
 261       INDWRK = INDTAU + N
 262       LLWORK = LWORK - INDWRK + 1
 263       CALL CHETRD( UPLO, N, A, LDA, W, RWORK( INDE ), WORK( INDTAU ),
 264      $             WORK( INDWRK ), LLWORK, IINFO )
 265 *
 266 *     For eigenvalues only, call SSTERF.  For eigenvectors, first call
 267 *     CUNGTR to generate the unitary matrix, then call CSTEQR.
 268 *
 269       IF( .NOT.WANTZ ) THEN
 270          CALL SSTERF( N, W, RWORK( INDE ), INFO )
 271       ELSE
 272          CALL CUNGTR( UPLO, N, A, LDA, WORK( INDTAU ), WORK( INDWRK ),
 273      $                LLWORK, IINFO )
 274          INDWRK = INDE + N
 275          CALL CSTEQR( JOBZ, N, W, RWORK( INDE ), A, LDA,
 276      $                RWORK( INDWRK ), INFO )
 277       END IF
 278 *
 279 *     If matrix was scaled, then rescale eigenvalues appropriately.
 280 *
 281       IF( ISCALE.EQ.1 ) THEN
 282          IF( INFO.EQ.0 ) THEN
 283             IMAX = N
 284          ELSE
 285             IMAX = INFO - 1
 286          END IF
 287          CALL SSCAL( IMAX, ONE / SIGMA, W, 1 )
 288       END IF
 289 *
 290 *     Set WORK(1) to optimal complex workspace size.
 291 *
 292       WORK( 1 ) = LWKOPT
 293 *
 294       RETURN
 295 *
 296 *     End of CHEEV
 297 *
 298       END