LAPACK 3.3.0

dsbgv.f

Go to the documentation of this file.
00001       SUBROUTINE DSBGV( JOBZ, UPLO, N, KA, KB, AB, LDAB, BB, LDBB, W, Z,
00002      $                  LDZ, WORK, 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, KA, KB, LDAB, LDBB, LDZ, N
00012 *     ..
00013 *     .. Array Arguments ..
00014       DOUBLE PRECISION   AB( LDAB, * ), BB( LDBB, * ), W( * ),
00015      $                   WORK( * ), Z( LDZ, * )
00016 *     ..
00017 *
00018 *  Purpose
00019 *  =======
00020 *
00021 *  DSBGV computes all the eigenvalues, and optionally, the eigenvectors
00022 *  of a real generalized symmetric-definite banded eigenproblem, of
00023 *  the form A*x=(lambda)*B*x. Here A and B are assumed to be symmetric
00024 *  and banded, and B is also positive definite.
00025 *
00026 *  Arguments
00027 *  =========
00028 *
00029 *  JOBZ    (input) CHARACTER*1
00030 *          = 'N':  Compute eigenvalues only;
00031 *          = 'V':  Compute eigenvalues and eigenvectors.
00032 *
00033 *  UPLO    (input) CHARACTER*1
00034 *          = 'U':  Upper triangles of A and B are stored;
00035 *          = 'L':  Lower triangles of A and B are stored.
00036 *
00037 *  N       (input) INTEGER
00038 *          The order of the matrices A and B.  N >= 0.
00039 *
00040 *  KA      (input) INTEGER
00041 *          The number of superdiagonals of the matrix A if UPLO = 'U',
00042 *          or the number of subdiagonals if UPLO = 'L'. KA >= 0.
00043 *
00044 *  KB      (input) INTEGER
00045 *          The number of superdiagonals of the matrix B if UPLO = 'U',
00046 *          or the number of subdiagonals if UPLO = 'L'. KB >= 0.
00047 *
00048 *  AB      (input/output) DOUBLE PRECISION array, dimension (LDAB, N)
00049 *          On entry, the upper or lower triangle of the symmetric band
00050 *          matrix A, stored in the first ka+1 rows of the array.  The
00051 *          j-th column of A is stored in the j-th column of the array AB
00052 *          as follows:
00053 *          if UPLO = 'U', AB(ka+1+i-j,j) = A(i,j) for max(1,j-ka)<=i<=j;
00054 *          if UPLO = 'L', AB(1+i-j,j)    = A(i,j) for j<=i<=min(n,j+ka).
00055 *
00056 *          On exit, the contents of AB are destroyed.
00057 *
00058 *  LDAB    (input) INTEGER
00059 *          The leading dimension of the array AB.  LDAB >= KA+1.
00060 *
00061 *  BB      (input/output) DOUBLE PRECISION array, dimension (LDBB, N)
00062 *          On entry, the upper or lower triangle of the symmetric band
00063 *          matrix B, stored in the first kb+1 rows of the array.  The
00064 *          j-th column of B is stored in the j-th column of the array BB
00065 *          as follows:
00066 *          if UPLO = 'U', BB(kb+1+i-j,j) = B(i,j) for max(1,j-kb)<=i<=j;
00067 *          if UPLO = 'L', BB(1+i-j,j)    = B(i,j) for j<=i<=min(n,j+kb).
00068 *
00069 *          On exit, the factor S from the split Cholesky factorization
00070 *          B = S**T*S, as returned by DPBSTF.
00071 *
00072 *  LDBB    (input) INTEGER
00073 *          The leading dimension of the array BB.  LDBB >= KB+1.
00074 *
00075 *  W       (output) DOUBLE PRECISION array, dimension (N)
00076 *          If INFO = 0, the eigenvalues in ascending order.
00077 *
00078 *  Z       (output) DOUBLE PRECISION array, dimension (LDZ, N)
00079 *          If JOBZ = 'V', then if INFO = 0, Z contains the matrix Z of
00080 *          eigenvectors, with the i-th column of Z holding the
00081 *          eigenvector associated with W(i). The eigenvectors are
00082 *          normalized so that Z**T*B*Z = I.
00083 *          If JOBZ = 'N', then Z is not referenced.
00084 *
00085 *  LDZ     (input) INTEGER
00086 *          The leading dimension of the array Z.  LDZ >= 1, and if
00087 *          JOBZ = 'V', LDZ >= N.
00088 *
00089 *  WORK    (workspace) DOUBLE PRECISION array, dimension (3*N)
00090 *
00091 *  INFO    (output) INTEGER
00092 *          = 0:  successful exit
00093 *          < 0:  if INFO = -i, the i-th argument had an illegal value
00094 *          > 0:  if INFO = i, and i is:
00095 *             <= N:  the algorithm failed to converge:
00096 *                    i off-diagonal elements of an intermediate
00097 *                    tridiagonal form did not converge to zero;
00098 *             > N:   if INFO = N + i, for 1 <= i <= N, then DPBSTF
00099 *                    returned INFO = i: B is not positive definite.
00100 *                    The factorization of B could not be completed and
00101 *                    no eigenvalues or eigenvectors were computed.
00102 *
00103 *  =====================================================================
00104 *
00105 *     .. Local Scalars ..
00106       LOGICAL            UPPER, WANTZ
00107       CHARACTER          VECT
00108       INTEGER            IINFO, INDE, INDWRK
00109 *     ..
00110 *     .. External Functions ..
00111       LOGICAL            LSAME
00112       EXTERNAL           LSAME
00113 *     ..
00114 *     .. External Subroutines ..
00115       EXTERNAL           DPBSTF, DSBGST, DSBTRD, DSTEQR, DSTERF, XERBLA
00116 *     ..
00117 *     .. Executable Statements ..
00118 *
00119 *     Test the input parameters.
00120 *
00121       WANTZ = LSAME( JOBZ, 'V' )
00122       UPPER = LSAME( UPLO, 'U' )
00123 *
00124       INFO = 0
00125       IF( .NOT.( WANTZ .OR. LSAME( JOBZ, 'N' ) ) ) THEN
00126          INFO = -1
00127       ELSE IF( .NOT.( UPPER .OR. LSAME( UPLO, 'L' ) ) ) THEN
00128          INFO = -2
00129       ELSE IF( N.LT.0 ) THEN
00130          INFO = -3
00131       ELSE IF( KA.LT.0 ) THEN
00132          INFO = -4
00133       ELSE IF( KB.LT.0 .OR. KB.GT.KA ) THEN
00134          INFO = -5
00135       ELSE IF( LDAB.LT.KA+1 ) THEN
00136          INFO = -7
00137       ELSE IF( LDBB.LT.KB+1 ) THEN
00138          INFO = -9
00139       ELSE IF( LDZ.LT.1 .OR. ( WANTZ .AND. LDZ.LT.N ) ) THEN
00140          INFO = -12
00141       END IF
00142       IF( INFO.NE.0 ) THEN
00143          CALL XERBLA( 'DSBGV ', -INFO )
00144          RETURN
00145       END IF
00146 *
00147 *     Quick return if possible
00148 *
00149       IF( N.EQ.0 )
00150      $   RETURN
00151 *
00152 *     Form a split Cholesky factorization of B.
00153 *
00154       CALL DPBSTF( UPLO, N, KB, BB, LDBB, INFO )
00155       IF( INFO.NE.0 ) THEN
00156          INFO = N + INFO
00157          RETURN
00158       END IF
00159 *
00160 *     Transform problem to standard eigenvalue problem.
00161 *
00162       INDE = 1
00163       INDWRK = INDE + N
00164       CALL DSBGST( JOBZ, UPLO, N, KA, KB, AB, LDAB, BB, LDBB, Z, LDZ,
00165      $             WORK( INDWRK ), IINFO )
00166 *
00167 *     Reduce to tridiagonal form.
00168 *
00169       IF( WANTZ ) THEN
00170          VECT = 'U'
00171       ELSE
00172          VECT = 'N'
00173       END IF
00174       CALL DSBTRD( VECT, UPLO, N, KA, AB, LDAB, W, WORK( INDE ), Z, LDZ,
00175      $             WORK( INDWRK ), IINFO )
00176 *
00177 *     For eigenvalues only, call DSTERF.  For eigenvectors, call SSTEQR.
00178 *
00179       IF( .NOT.WANTZ ) THEN
00180          CALL DSTERF( N, W, WORK( INDE ), INFO )
00181       ELSE
00182          CALL DSTEQR( JOBZ, N, W, WORK( INDE ), Z, LDZ, WORK( INDWRK ),
00183      $                INFO )
00184       END IF
00185       RETURN
00186 *
00187 *     End of DSBGV
00188 *
00189       END
 All Files Functions