SRC/spteqr.f

   1 *> \brief \b SPTEQR
   2 *
   3 *  =========== DOCUMENTATION ===========
   4 *
   5 * Online html documentation available at
   6 *            http://www.netlib.org/lapack/explore-html/
   7 *
   8 *> \htmlonly
   9 *> Download SPTEQR + dependencies
  10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/spteqr.f">
  11 *> [TGZ]</a>
  12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/spteqr.f">
  13 *> [ZIP]</a>
  14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/spteqr.f">
  15 *> [TXT]</a>
  16 *> \endhtmlonly
  17 *
  18 *  Definition:
  19 *  ===========
  20 *
  21 *       SUBROUTINE SPTEQR( COMPZ, N, D, E, Z, LDZ, WORK, INFO )
  22 *
  23 *       .. Scalar Arguments ..
  24 *       CHARACTER          COMPZ
  25 *       INTEGER            INFO, LDZ, N
  26 *       ..
  27 *       .. Array Arguments ..
  28 *       REAL               D( * ), E( * ), WORK( * ), Z( LDZ, * )
  29 *       ..
  30 *
  31 *
  32 *> \par Purpose:
  33 *  =============
  34 *>
  35 *> \verbatim
  36 *>
  37 *> SPTEQR computes all eigenvalues and, optionally, eigenvectors of a
  38 *> symmetric positive definite tridiagonal matrix by first factoring the
  39 *> matrix using SPTTRF, and then calling SBDSQR to compute the singular
  40 *> values of the bidiagonal factor.
  41 *>
  42 *> This routine computes the eigenvalues of the positive definite
  43 *> tridiagonal matrix to high relative accuracy.  This means that if the
  44 *> eigenvalues range over many orders of magnitude in size, then the
  45 *> small eigenvalues and corresponding eigenvectors will be computed
  46 *> more accurately than, for example, with the standard QR method.
  47 *>
  48 *> The eigenvectors of a full or band symmetric positive definite matrix
  49 *> can also be found if SSYTRD, SSPTRD, or SSBTRD has been used to
  50 *> reduce this matrix to tridiagonal form. (The reduction to tridiagonal
  51 *> form, however, may preclude the possibility of obtaining high
  52 *> relative accuracy in the small eigenvalues of the original matrix, if
  53 *> these eigenvalues range over many orders of magnitude.)
  54 *> \endverbatim
  55 *
  56 *  Arguments:
  57 *  ==========
  58 *
  59 *> \param[in] COMPZ
  60 *> \verbatim
  61 *>          COMPZ is CHARACTER*1
  62 *>          = 'N':  Compute eigenvalues only.
  63 *>          = 'V':  Compute eigenvectors of original symmetric
  64 *>                  matrix also.  Array Z contains the orthogonal
  65 *>                  matrix used to reduce the original matrix to
  66 *>                  tridiagonal form.
  67 *>          = 'I':  Compute eigenvectors of tridiagonal matrix also.
  68 *> \endverbatim
  69 *>
  70 *> \param[in] N
  71 *> \verbatim
  72 *>          N is INTEGER
  73 *>          The order of the matrix.  N >= 0.
  74 *> \endverbatim
  75 *>
  76 *> \param[in,out] D
  77 *> \verbatim
  78 *>          D is REAL array, dimension (N)
  79 *>          On entry, the n diagonal elements of the tridiagonal
  80 *>          matrix.
  81 *>          On normal exit, D contains the eigenvalues, in descending
  82 *>          order.
  83 *> \endverbatim
  84 *>
  85 *> \param[in,out] E
  86 *> \verbatim
  87 *>          E is REAL array, dimension (N-1)
  88 *>          On entry, the (n-1) subdiagonal elements of the tridiagonal
  89 *>          matrix.
  90 *>          On exit, E has been destroyed.
  91 *> \endverbatim
  92 *>
  93 *> \param[in,out] Z
  94 *> \verbatim
  95 *>          Z is REAL array, dimension (LDZ, N)
  96 *>          On entry, if COMPZ = 'V', the orthogonal matrix used in the
  97 *>          reduction to tridiagonal form.
  98 *>          On exit, if COMPZ = 'V', the orthonormal eigenvectors of the
  99 *>          original symmetric matrix;
 100 *>          if COMPZ = 'I', the orthonormal eigenvectors of the
 101 *>          tridiagonal matrix.
 102 *>          If INFO > 0 on exit, Z contains the eigenvectors associated
 103 *>          with only the stored eigenvalues.
 104 *>          If  COMPZ = 'N', then Z is not referenced.
 105 *> \endverbatim
 106 *>
 107 *> \param[in] LDZ
 108 *> \verbatim
 109 *>          LDZ is INTEGER
 110 *>          The leading dimension of the array Z.  LDZ >= 1, and if
 111 *>          COMPZ = 'V' or 'I', LDZ >= max(1,N).
 112 *> \endverbatim
 113 *>
 114 *> \param[out] WORK
 115 *> \verbatim
 116 *>          WORK is REAL array, dimension (4*N)
 117 *> \endverbatim
 118 *>
 119 *> \param[out] INFO
 120 *> \verbatim
 121 *>          INFO is INTEGER
 122 *>          = 0:  successful exit.
 123 *>          < 0:  if INFO = -i, the i-th argument had an illegal value.
 124 *>          > 0:  if INFO = i, and i is:
 125 *>                <= N  the Cholesky factorization of the matrix could
 126 *>                      not be performed because the i-th principal minor
 127 *>                      was not positive definite.
 128 *>                > N   the SVD algorithm failed to converge;
 129 *>                      if INFO = N+i, i off-diagonal elements of the
 130 *>                      bidiagonal factor did not converge to zero.
 131 *> \endverbatim
 132 *
 133 *  Authors:
 134 *  ========
 135 *
 136 *> \author Univ. of Tennessee
 137 *> \author Univ. of California Berkeley
 138 *> \author Univ. of Colorado Denver
 139 *> \author NAG Ltd.
 140 *
 141 *> \date September 2012
 142 *
 143 *> \ingroup realPTcomputational
 144 *
 145 *  =====================================================================
 146       SUBROUTINE SPTEQR( COMPZ, N, D, E, Z, LDZ, WORK, INFO )
 147 *
 148 *  -- LAPACK computational routine (version 3.4.2) --
 149 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 150 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 151 *     September 2012
 152 *
 153 *     .. Scalar Arguments ..
 154       CHARACTER          COMPZ
 155       INTEGER            INFO, LDZ, N
 156 *     ..
 157 *     .. Array Arguments ..
 158       REAL               D( * ), E( * ), WORK( * ), Z( LDZ, * )
 159 *     ..
 160 *
 161 *  =====================================================================
 162 *
 163 *     .. Parameters ..
 164       REAL               ZERO, ONE
 165       PARAMETER          ( ZERO = 0.0E0, ONE = 1.0E0 )
 166 *     ..
 167 *     .. External Functions ..
 168       LOGICAL            LSAME
 169       EXTERNAL           LSAME
 170 *     ..
 171 *     .. External Subroutines ..
 172       EXTERNAL           SBDSQR, SLASET, SPTTRF, XERBLA
 173 *     ..
 174 *     .. Local Arrays ..
 175       REAL               C( 1, 1 ), VT( 1, 1 )
 176 *     ..
 177 *     .. Local Scalars ..
 178       INTEGER            I, ICOMPZ, NRU
 179 *     ..
 180 *     .. Intrinsic Functions ..
 181       INTRINSIC          MAX, SQRT
 182 *     ..
 183 *     .. Executable Statements ..
 184 *
 185 *     Test the input parameters.
 186 *
 187       INFO = 0
 188 *
 189       IF( LSAME( COMPZ, 'N' ) ) THEN
 190          ICOMPZ = 0
 191       ELSE IF( LSAME( COMPZ, 'V' ) ) THEN
 192          ICOMPZ = 1
 193       ELSE IF( LSAME( COMPZ, 'I' ) ) THEN
 194          ICOMPZ = 2
 195       ELSE
 196          ICOMPZ = -1
 197       END IF
 198       IF( ICOMPZ.LT.0 ) THEN
 199          INFO = -1
 200       ELSE IF( N.LT.0 ) THEN
 201          INFO = -2
 202       ELSE IF( ( LDZ.LT.1 ) .OR. ( ICOMPZ.GT.0 .AND. LDZ.LT.MAX( 1,
 203      $         N ) ) ) THEN
 204          INFO = -6
 205       END IF
 206       IF( INFO.NE.0 ) THEN
 207          CALL XERBLA( 'SPTEQR', -INFO )
 208          RETURN
 209       END IF
 210 *
 211 *     Quick return if possible
 212 *
 213       IF( N.EQ.0 )
 214      $   RETURN
 215 *
 216       IF( N.EQ.1 ) THEN
 217          IF( ICOMPZ.GT.0 )
 218      $      Z( 1, 1 ) = ONE
 219          RETURN
 220       END IF
 221       IF( ICOMPZ.EQ.2 )
 222      $   CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDZ )
 223 *
 224 *     Call SPTTRF to factor the matrix.
 225 *
 226       CALL SPTTRF( N, D, E, INFO )
 227       IF( INFO.NE.0 )
 228      $   RETURN
 229       DO 10 I = 1, N
 230          D( I ) = SQRT( D( I ) )
 231    10 CONTINUE
 232       DO 20 I = 1, N - 1
 233          E( I ) = E( I )*D( I )
 234    20 CONTINUE
 235 *
 236 *     Call SBDSQR to compute the singular values/vectors of the
 237 *     bidiagonal factor.
 238 *
 239       IF( ICOMPZ.GT.0 ) THEN
 240          NRU = N
 241       ELSE
 242          NRU = 0
 243       END IF
 244       CALL SBDSQR( 'Lower', N, 0, NRU, 0, D, E, VT, 1, Z, LDZ, C, 1,
 245      $             WORK, INFO )
 246 *
 247 *     Square the singular values.
 248 *
 249       IF( INFO.EQ.0 ) THEN
 250          DO 30 I = 1, N
 251             D( I ) = D( I )*D( I )
 252    30    CONTINUE
 253       ELSE
 254          INFO = N + INFO
 255       END IF
 256 *
 257       RETURN
 258 *
 259 *     End of SPTEQR
 260 *
 261       END