LAPACK: SRC/ssbev.f Source File

Go to the documentation of this file.
00001       SUBROUTINE SSBEV( JOBZ, UPLO, N, KD, AB, LDAB, W, Z, LDZ, WORK,
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, KD, LDAB, LDZ, N
00012 *     ..
00013 *     .. Array Arguments ..
00014       REAL               AB( LDAB, * ), W( * ), WORK( * ), Z( LDZ, * )
00015 *     ..
00016 *
00017 *  Purpose
00018 *  =======
00019 *
00020 *  SSBEV computes all the eigenvalues and, optionally, eigenvectors of
00021 *  a real symmetric band matrix A.
00022 *
00023 *  Arguments
00024 *  =========
00025 *
00026 *  JOBZ    (input) CHARACTER*1
00027 *          = 'N':  Compute eigenvalues only;
00028 *          = 'V':  Compute eigenvalues and eigenvectors.
00029 *
00030 *  UPLO    (input) CHARACTER*1
00031 *          = 'U':  Upper triangle of A is stored;
00032 *          = 'L':  Lower triangle of A is stored.
00033 *
00034 *  N       (input) INTEGER
00035 *          The order of the matrix A.  N >= 0.
00036 *
00037 *  KD      (input) INTEGER
00038 *          The number of superdiagonals of the matrix A if UPLO = 'U',
00039 *          or the number of subdiagonals if UPLO = 'L'.  KD >= 0.
00040 *
00041 *  AB      (input/output) REAL array, dimension (LDAB, N)
00042 *          On entry, the upper or lower triangle of the symmetric band
00043 *          matrix A, stored in the first KD+1 rows of the array.  The
00044 *          j-th column of A is stored in the j-th column of the array AB
00045 *          as follows:
00046 *          if UPLO = 'U', AB(kd+1+i-j,j) = A(i,j) for max(1,j-kd)<=i<=j;
00047 *          if UPLO = 'L', AB(1+i-j,j)    = A(i,j) for j<=i<=min(n,j+kd).
00048 *
00049 *          On exit, AB is overwritten by values generated during the
00050 *          reduction to tridiagonal form.  If UPLO = 'U', the first
00051 *          superdiagonal and the diagonal of the tridiagonal matrix T
00052 *          are returned in rows KD and KD+1 of AB, and if UPLO = 'L',
00053 *          the diagonal and first subdiagonal of T are returned in the
00054 *          first two rows of AB.
00055 *
00056 *  LDAB    (input) INTEGER
00057 *          The leading dimension of the array AB.  LDAB >= KD + 1.
00058 *
00059 *  W       (output) REAL array, dimension (N)
00060 *          If INFO = 0, the eigenvalues in ascending order.
00061 *
00062 *  Z       (output) REAL array, dimension (LDZ, N)
00063 *          If JOBZ = 'V', then if INFO = 0, Z contains the orthonormal
00064 *          eigenvectors of the matrix A, with the i-th column of Z
00065 *          holding the eigenvector associated with W(i).
00066 *          If JOBZ = 'N', then Z is not referenced.
00067 *
00068 *  LDZ     (input) INTEGER
00069 *          The leading dimension of the array Z.  LDZ >= 1, and if
00070 *          JOBZ = 'V', LDZ >= max(1,N).
00071 *
00072 *  WORK    (workspace) REAL array, dimension (max(1,3*N-2))
00073 *
00074 *  INFO    (output) INTEGER
00075 *          = 0:  successful exit
00076 *          < 0:  if INFO = -i, the i-th argument had an illegal value
00077 *          > 0:  if INFO = i, the algorithm failed to converge; i
00078 *                off-diagonal elements of an intermediate tridiagonal
00079 *                form did not converge to zero.
00080 *
00081 *  =====================================================================
00082 *
00083 *     .. Parameters ..
00084       REAL               ZERO, ONE
00085       PARAMETER          ( ZERO = 0.0E0, ONE = 1.0E0 )
00086 *     ..
00087 *     .. Local Scalars ..
00088       LOGICAL            LOWER, WANTZ
00089       INTEGER            IINFO, IMAX, INDE, INDWRK, ISCALE
00090       REAL               ANRM, BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA,
00091      $                   SMLNUM
00092 *     ..
00093 *     .. External Functions ..
00094       LOGICAL            LSAME
00095       REAL               SLAMCH, SLANSB
00096       EXTERNAL           LSAME, SLAMCH, SLANSB
00097 *     ..
00098 *     .. External Subroutines ..
00099       EXTERNAL           SLASCL, SSBTRD, SSCAL, SSTEQR, SSTERF, XERBLA
00100 *     ..
00101 *     .. Intrinsic Functions ..
00102       INTRINSIC          SQRT
00103 *     ..
00104 *     .. Executable Statements ..
00105 *
00106 *     Test the input parameters.
00107 *
00108       WANTZ = LSAME( JOBZ, 'V' )
00109       LOWER = LSAME( UPLO, 'L' )
00110 *
00111       INFO = 0
00112       IF( .NOT.( WANTZ .OR. LSAME( JOBZ, 'N' ) ) ) THEN
00113          INFO = -1
00114       ELSE IF( .NOT.( LOWER .OR. LSAME( UPLO, 'U' ) ) ) THEN
00115          INFO = -2
00116       ELSE IF( N.LT.0 ) THEN
00117          INFO = -3
00118       ELSE IF( KD.LT.0 ) THEN
00119          INFO = -4
00120       ELSE IF( LDAB.LT.KD+1 ) THEN
00121          INFO = -6
00122       ELSE IF( LDZ.LT.1 .OR. ( WANTZ .AND. LDZ.LT.N ) ) THEN
00123          INFO = -9
00124       END IF
00125 *
00126       IF( INFO.NE.0 ) THEN
00127          CALL XERBLA( 'SSBEV ', -INFO )
00128          RETURN
00129       END IF
00130 *
00131 *     Quick return if possible
00132 *
00133       IF( N.EQ.0 )
00134      $   RETURN
00135 *
00136       IF( N.EQ.1 ) THEN
00137          IF( LOWER ) THEN
00138             W( 1 ) = AB( 1, 1 )
00139          ELSE
00140             W( 1 ) = AB( KD+1, 1 )
00141          END IF
00142          IF( WANTZ )
00143      $      Z( 1, 1 ) = ONE
00144          RETURN
00145       END IF
00146 *
00147 *     Get machine constants.
00148 *
00149       SAFMIN = SLAMCH( 'Safe minimum' )
00150       EPS = SLAMCH( 'Precision' )
00151       SMLNUM = SAFMIN / EPS
00152       BIGNUM = ONE / SMLNUM
00153       RMIN = SQRT( SMLNUM )
00154       RMAX = SQRT( BIGNUM )
00155 *
00156 *     Scale matrix to allowable range, if necessary.
00157 *
00158       ANRM = SLANSB( 'M', UPLO, N, KD, AB, LDAB, WORK )
00159       ISCALE = 0
00160       IF( ANRM.GT.ZERO .AND. ANRM.LT.RMIN ) THEN
00161          ISCALE = 1
00162          SIGMA = RMIN / ANRM
00163       ELSE IF( ANRM.GT.RMAX ) THEN
00164          ISCALE = 1
00165          SIGMA = RMAX / ANRM
00166       END IF
00167       IF( ISCALE.EQ.1 ) THEN
00168          IF( LOWER ) THEN
00169             CALL SLASCL( 'B', KD, KD, ONE, SIGMA, N, N, AB, LDAB, INFO )
00170          ELSE
00171             CALL SLASCL( 'Q', KD, KD, ONE, SIGMA, N, N, AB, LDAB, INFO )
00172          END IF
00173       END IF
00174 *
00175 *     Call SSBTRD to reduce symmetric band matrix to tridiagonal form.
00176 *
00177       INDE = 1
00178       INDWRK = INDE + N
00179       CALL SSBTRD( JOBZ, UPLO, N, KD, AB, LDAB, W, WORK( INDE ), Z, LDZ,
00180      $             WORK( INDWRK ), IINFO )
00181 *
00182 *     For eigenvalues only, call SSTERF.  For eigenvectors, call SSTEQR.
00183 *
00184       IF( .NOT.WANTZ ) THEN
00185          CALL SSTERF( N, W, WORK( INDE ), INFO )
00186       ELSE
00187          CALL SSTEQR( JOBZ, N, W, WORK( INDE ), Z, LDZ, WORK( INDWRK ),
00188      $                INFO )
00189       END IF
00190 *
00191 *     If matrix was scaled, then rescale eigenvalues appropriately.
00192 *
00193       IF( ISCALE.EQ.1 ) THEN
00194          IF( INFO.EQ.0 ) THEN
00195             IMAX = N
00196          ELSE
00197             IMAX = INFO - 1
00198          END IF
00199          CALL SSCAL( IMAX, ONE / SIGMA, W, 1 )
00200       END IF
00201 *
00202       RETURN
00203 *
00204 *     End of SSBEV
00205 *
00206       END