LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ dsbgv()

subroutine dsbgv ( character  jobz,
character  uplo,
integer  n,
integer  ka,
integer  kb,
double precision, dimension( ldab, * )  ab,
integer  ldab,
double precision, dimension( ldbb, * )  bb,
integer  ldbb,
double precision, dimension( * )  w,
double precision, dimension( ldz, * )  z,
integer  ldz,
double precision, dimension( * )  work,
integer  info 
)

DSBGV

Download DSBGV + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 DSBGV computes all the eigenvalues, and optionally, the eigenvectors
 of a real generalized symmetric-definite banded eigenproblem, of
 the form A*x=(lambda)*B*x. Here A and B are assumed to be symmetric
 and banded, and B is also positive definite.
Parameters
[in]JOBZ
          JOBZ is CHARACTER*1
          = 'N':  Compute eigenvalues only;
          = 'V':  Compute eigenvalues and eigenvectors.
[in]UPLO
          UPLO is CHARACTER*1
          = 'U':  Upper triangles of A and B are stored;
          = 'L':  Lower triangles of A and B are stored.
[in]N
          N is INTEGER
          The order of the matrices A and B.  N >= 0.
[in]KA
          KA is INTEGER
          The number of superdiagonals of the matrix A if UPLO = 'U',
          or the number of subdiagonals if UPLO = 'L'. KA >= 0.
[in]KB
          KB is INTEGER
          The number of superdiagonals of the matrix B if UPLO = 'U',
          or the number of subdiagonals if UPLO = 'L'. KB >= 0.
[in,out]AB
          AB is DOUBLE PRECISION array, dimension (LDAB, N)
          On entry, the upper or lower triangle of the symmetric band
          matrix A, stored in the first ka+1 rows of the array.  The
          j-th column of A is stored in the j-th column of the array AB
          as follows:
          if UPLO = 'U', AB(ka+1+i-j,j) = A(i,j) for max(1,j-ka)<=i<=j;
          if UPLO = 'L', AB(1+i-j,j)    = A(i,j) for j<=i<=min(n,j+ka).

          On exit, the contents of AB are destroyed.
[in]LDAB
          LDAB is INTEGER
          The leading dimension of the array AB.  LDAB >= KA+1.
[in,out]BB
          BB is DOUBLE PRECISION array, dimension (LDBB, N)
          On entry, the upper or lower triangle of the symmetric band
          matrix B, stored in the first kb+1 rows of the array.  The
          j-th column of B is stored in the j-th column of the array BB
          as follows:
          if UPLO = 'U', BB(kb+1+i-j,j) = B(i,j) for max(1,j-kb)<=i<=j;
          if UPLO = 'L', BB(1+i-j,j)    = B(i,j) for j<=i<=min(n,j+kb).

          On exit, the factor S from the split Cholesky factorization
          B = S**T*S, as returned by DPBSTF.
[in]LDBB
          LDBB is INTEGER
          The leading dimension of the array BB.  LDBB >= KB+1.
[out]W
          W is DOUBLE PRECISION array, dimension (N)
          If INFO = 0, the eigenvalues in ascending order.
[out]Z
          Z is DOUBLE PRECISION array, dimension (LDZ, N)
          If JOBZ = 'V', then if INFO = 0, Z contains the matrix Z of
          eigenvectors, with the i-th column of Z holding the
          eigenvector associated with W(i). The eigenvectors are
          normalized so that Z**T*B*Z = I.
          If JOBZ = 'N', then Z is not referenced.
[in]LDZ
          LDZ is INTEGER
          The leading dimension of the array Z.  LDZ >= 1, and if
          JOBZ = 'V', LDZ >= N.
[out]WORK
          WORK is DOUBLE PRECISION array, dimension (3*N)
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  if INFO = i, and i is:
             <= N:  the algorithm failed to converge:
                    i off-diagonal elements of an intermediate
                    tridiagonal form did not converge to zero;
             > N:   if INFO = N + i, for 1 <= i <= N, then DPBSTF
                    returned INFO = i: B is not positive definite.
                    The factorization of B could not be completed and
                    no eigenvalues or eigenvectors were computed.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 175 of file dsbgv.f.

177*
178* -- LAPACK driver routine --
179* -- LAPACK is a software package provided by Univ. of Tennessee, --
180* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
181*
182* .. Scalar Arguments ..
183 CHARACTER JOBZ, UPLO
184 INTEGER INFO, KA, KB, LDAB, LDBB, LDZ, N
185* ..
186* .. Array Arguments ..
187 DOUBLE PRECISION AB( LDAB, * ), BB( LDBB, * ), W( * ),
188 $ WORK( * ), Z( LDZ, * )
189* ..
190*
191* =====================================================================
192*
193* .. Local Scalars ..
194 LOGICAL UPPER, WANTZ
195 CHARACTER VECT
196 INTEGER IINFO, INDE, INDWRK
197* ..
198* .. External Functions ..
199 LOGICAL LSAME
200 EXTERNAL lsame
201* ..
202* .. External Subroutines ..
203 EXTERNAL dpbstf, dsbgst, dsbtrd, dsteqr, dsterf, xerbla
204* ..
205* .. Executable Statements ..
206*
207* Test the input parameters.
208*
209 wantz = lsame( jobz, 'V' )
210 upper = lsame( uplo, 'U' )
211*
212 info = 0
213 IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
214 info = -1
215 ELSE IF( .NOT.( upper .OR. lsame( uplo, 'L' ) ) ) THEN
216 info = -2
217 ELSE IF( n.LT.0 ) THEN
218 info = -3
219 ELSE IF( ka.LT.0 ) THEN
220 info = -4
221 ELSE IF( kb.LT.0 .OR. kb.GT.ka ) THEN
222 info = -5
223 ELSE IF( ldab.LT.ka+1 ) THEN
224 info = -7
225 ELSE IF( ldbb.LT.kb+1 ) THEN
226 info = -9
227 ELSE IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) ) THEN
228 info = -12
229 END IF
230 IF( info.NE.0 ) THEN
231 CALL xerbla( 'DSBGV', -info )
232 RETURN
233 END IF
234*
235* Quick return if possible
236*
237 IF( n.EQ.0 )
238 $ RETURN
239*
240* Form a split Cholesky factorization of B.
241*
242 CALL dpbstf( uplo, n, kb, bb, ldbb, info )
243 IF( info.NE.0 ) THEN
244 info = n + info
245 RETURN
246 END IF
247*
248* Transform problem to standard eigenvalue problem.
249*
250 inde = 1
251 indwrk = inde + n
252 CALL dsbgst( jobz, uplo, n, ka, kb, ab, ldab, bb, ldbb, z, ldz,
253 $ work( indwrk ), iinfo )
254*
255* Reduce to tridiagonal form.
256*
257 IF( wantz ) THEN
258 vect = 'U'
259 ELSE
260 vect = 'N'
261 END IF
262 CALL dsbtrd( vect, uplo, n, ka, ab, ldab, w, work( inde ), z, ldz,
263 $ work( indwrk ), iinfo )
264*
265* For eigenvalues only, call DSTERF. For eigenvectors, call SSTEQR.
266*
267 IF( .NOT.wantz ) THEN
268 CALL dsterf( n, w, work( inde ), info )
269 ELSE
270 CALL dsteqr( jobz, n, w, work( inde ), z, ldz, work( indwrk ),
271 $ info )
272 END IF
273 RETURN
274*
275* End of DSBGV
276*
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine dsbgst(vect, uplo, n, ka, kb, ab, ldab, bb, ldbb, x, ldx, work, info)
DSBGST
Definition dsbgst.f:159
subroutine dsbtrd(vect, uplo, n, kd, ab, ldab, d, e, q, ldq, work, info)
DSBTRD
Definition dsbtrd.f:163
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
subroutine dpbstf(uplo, n, kd, ab, ldab, info)
DPBSTF
Definition dpbstf.f:152
subroutine dsteqr(compz, n, d, e, z, ldz, work, info)
DSTEQR
Definition dsteqr.f:131
subroutine dsterf(n, d, e, info)
DSTERF
Definition dsterf.f:86
Here is the call graph for this function:
Here is the caller graph for this function: