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

◆ dsbgvx()

subroutine dsbgvx ( character jobz,
character range,
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( ldq, * ) q,
integer ldq,
double precision vl,
double precision vu,
integer il,
integer iu,
double precision abstol,
integer m,
double precision, dimension( * ) w,
double precision, dimension( ldz, * ) z,
integer ldz,
double precision, dimension( * ) work,
integer, dimension( * ) iwork,
integer, dimension( * ) ifail,
integer info )

DSBGVX

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

Purpose:
!>
!> DSBGVX computes selected eigenvalues, and optionally, 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.  Eigenvalues and
!> eigenvectors can be selected by specifying either all eigenvalues,
!> a range of values or a range of indices for the desired eigenvalues.
!>
Parameters
[in]JOBZ
!>          JOBZ is CHARACTER*1
!>          = 'N':  Compute eigenvalues only;
!>          = 'V':  Compute eigenvalues and eigenvectors.
!>
[in]RANGE
!>          RANGE is CHARACTER*1
!>          = 'A': all eigenvalues will be found.
!>          = 'V': all eigenvalues in the half-open interval (VL,VU]
!>                 will be found.
!>          = 'I': the IL-th through IU-th eigenvalues will be found.
!>
[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(ka+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]Q
!>          Q is DOUBLE PRECISION array, dimension (LDQ, N)
!>          If JOBZ = 'V', the n-by-n matrix used in the reduction of
!>          A*x = (lambda)*B*x to standard form, i.e. C*x = (lambda)*x,
!>          and consequently C to tridiagonal form.
!>          If JOBZ = 'N', the array Q is not referenced.
!>
[in]LDQ
!>          LDQ is INTEGER
!>          The leading dimension of the array Q.  If JOBZ = 'N',
!>          LDQ >= 1. If JOBZ = 'V', LDQ >= max(1,N).
!>
[in]VL
!>          VL is DOUBLE PRECISION
!>
!>          If RANGE='V', the lower bound of the interval to
!>          be searched for eigenvalues. VL < VU.
!>          Not referenced if RANGE = 'A' or 'I'.
!>
[in]VU
!>          VU is DOUBLE PRECISION
!>
!>          If RANGE='V', the upper bound of the interval to
!>          be searched for eigenvalues. VL < VU.
!>          Not referenced if RANGE = 'A' or 'I'.
!>
[in]IL
!>          IL is INTEGER
!>
!>          If RANGE='I', the index of the
!>          smallest eigenvalue to be returned.
!>          1 <= IL <= IU <= N, if N > 0; IL = 1 and IU = 0 if N = 0.
!>          Not referenced if RANGE = 'A' or 'V'.
!>
[in]IU
!>          IU is INTEGER
!>
!>          If RANGE='I', the index of the
!>          largest eigenvalue to be returned.
!>          1 <= IL <= IU <= N, if N > 0; IL = 1 and IU = 0 if N = 0.
!>          Not referenced if RANGE = 'A' or 'V'.
!>
[in]ABSTOL
!>          ABSTOL is DOUBLE PRECISION
!>          The absolute error tolerance for the eigenvalues.
!>          An approximate eigenvalue is accepted as converged
!>          when it is determined to lie in an interval [a,b]
!>          of width less than or equal to
!>
!>                  ABSTOL + EPS *   max( |a|,|b| ) ,
!>
!>          where EPS is the machine precision.  If ABSTOL is less than
!>          or equal to zero, then  EPS*|T|  will be used in its place,
!>          where |T| is the 1-norm of the tridiagonal matrix obtained
!>          by reducing A to tridiagonal form.
!>
!>          Eigenvalues will be computed most accurately when ABSTOL is
!>          set to twice the underflow threshold 2*DLAMCH('S'), not zero.
!>          If this routine returns with INFO>0, indicating that some
!>          eigenvectors did not converge, try setting ABSTOL to
!>          2*DLAMCH('S').
!>
[out]M
!>          M is INTEGER
!>          The total number of eigenvalues found.  0 <= M <= N.
!>          If RANGE = 'A', M = N, and if RANGE = 'I', M = IU-IL+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 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 >= max(1,N).
!>
[out]WORK
!>          WORK is DOUBLE PRECISION array, dimension (7*N)
!>
[out]IWORK
!>          IWORK is INTEGER array, dimension (5*N)
!>
[out]IFAIL
!>          IFAIL is INTEGER array, dimension (M)
!>          If JOBZ = 'V', then if INFO = 0, the first M elements of
!>          IFAIL are zero.  If INFO > 0, then IFAIL contains the
!>          indices of the eigenvalues that failed to converge.
!>          If JOBZ = 'N', then IFAIL is not referenced.
!>
[out]INFO
!>          INFO is INTEGER
!>          = 0:  successful exit
!>          < 0:  if INFO = -i, the i-th argument had an illegal value
!>          <= N: if INFO = i, then i eigenvectors failed to converge.
!>                  Their indices are stored in IFAIL.
!>          > N:  DPBSTF returned an error code; i.e.,
!>                if INFO = N + i, for 1 <= i <= N, then the leading
!>                principal minor of order i of B is not positive.
!>                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.
Contributors:
Mark Fahey, Department of Mathematics, Univ. of Kentucky, USA

Definition at line 289 of file dsbgvx.f.

292*
293* -- LAPACK driver routine --
294* -- LAPACK is a software package provided by Univ. of Tennessee, --
295* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
296*
297* .. Scalar Arguments ..
298 CHARACTER JOBZ, RANGE, UPLO
299 INTEGER IL, INFO, IU, KA, KB, LDAB, LDBB, LDQ, LDZ, M,
300 $ N
301 DOUBLE PRECISION ABSTOL, VL, VU
302* ..
303* .. Array Arguments ..
304 INTEGER IFAIL( * ), IWORK( * )
305 DOUBLE PRECISION AB( LDAB, * ), BB( LDBB, * ), Q( LDQ, * ),
306 $ W( * ), WORK( * ), Z( LDZ, * )
307* ..
308*
309* =====================================================================
310*
311* .. Parameters ..
312 DOUBLE PRECISION ZERO, ONE
313 parameter( zero = 0.0d+0, one = 1.0d+0 )
314* ..
315* .. Local Scalars ..
316 LOGICAL ALLEIG, INDEIG, TEST, UPPER, VALEIG, WANTZ
317 CHARACTER ORDER, VECT
318 INTEGER I, IINFO, INDD, INDE, INDEE, INDISP,
319 $ INDIWO, INDWRK, ITMP1, J, JJ, NSPLIT
320 DOUBLE PRECISION TMP1
321* ..
322* .. External Functions ..
323 LOGICAL LSAME
324 EXTERNAL lsame
325* ..
326* .. External Subroutines ..
327 EXTERNAL dcopy, dgemv, dlacpy, dpbstf, dsbgst,
328 $ dsbtrd,
329 $ dstebz, dstein, dsteqr, dsterf, dswap, xerbla
330* ..
331* .. Intrinsic Functions ..
332 INTRINSIC min
333* ..
334* .. Executable Statements ..
335*
336* Test the input parameters.
337*
338 wantz = lsame( jobz, 'V' )
339 upper = lsame( uplo, 'U' )
340 alleig = lsame( range, 'A' )
341 valeig = lsame( range, 'V' )
342 indeig = lsame( range, 'I' )
343*
344 info = 0
345 IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
346 info = -1
347 ELSE IF( .NOT.( alleig .OR. valeig .OR. indeig ) ) THEN
348 info = -2
349 ELSE IF( .NOT.( upper .OR. lsame( uplo, 'L' ) ) ) THEN
350 info = -3
351 ELSE IF( n.LT.0 ) THEN
352 info = -4
353 ELSE IF( ka.LT.0 ) THEN
354 info = -5
355 ELSE IF( kb.LT.0 .OR. kb.GT.ka ) THEN
356 info = -6
357 ELSE IF( ldab.LT.ka+1 ) THEN
358 info = -8
359 ELSE IF( ldbb.LT.kb+1 ) THEN
360 info = -10
361 ELSE IF( ldq.LT.1 .OR. ( wantz .AND. ldq.LT.n ) ) THEN
362 info = -12
363 ELSE
364 IF( valeig ) THEN
365 IF( n.GT.0 .AND. vu.LE.vl )
366 $ info = -14
367 ELSE IF( indeig ) THEN
368 IF( il.LT.1 .OR. il.GT.max( 1, n ) ) THEN
369 info = -15
370 ELSE IF ( iu.LT.min( n, il ) .OR. iu.GT.n ) THEN
371 info = -16
372 END IF
373 END IF
374 END IF
375 IF( info.EQ.0) THEN
376 IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) ) THEN
377 info = -21
378 END IF
379 END IF
380*
381 IF( info.NE.0 ) THEN
382 CALL xerbla( 'DSBGVX', -info )
383 RETURN
384 END IF
385*
386* Quick return if possible
387*
388 m = 0
389 IF( n.EQ.0 )
390 $ RETURN
391*
392* Form a split Cholesky factorization of B.
393*
394 CALL dpbstf( uplo, n, kb, bb, ldbb, info )
395 IF( info.NE.0 ) THEN
396 info = n + info
397 RETURN
398 END IF
399*
400* Transform problem to standard eigenvalue problem.
401*
402 CALL dsbgst( jobz, uplo, n, ka, kb, ab, ldab, bb, ldbb, q, ldq,
403 $ work, iinfo )
404*
405* Reduce symmetric band matrix to tridiagonal form.
406*
407 indd = 1
408 inde = indd + n
409 indwrk = inde + n
410 IF( wantz ) THEN
411 vect = 'U'
412 ELSE
413 vect = 'N'
414 END IF
415 CALL dsbtrd( vect, uplo, n, ka, ab, ldab, work( indd ),
416 $ work( inde ), q, ldq, work( indwrk ), iinfo )
417*
418* If all eigenvalues are desired and ABSTOL is less than or equal
419* to zero, then call DSTERF or SSTEQR. If this fails for some
420* eigenvalue, then try DSTEBZ.
421*
422 test = .false.
423 IF( indeig ) THEN
424 IF( il.EQ.1 .AND. iu.EQ.n ) THEN
425 test = .true.
426 END IF
427 END IF
428 IF( ( alleig .OR. test ) .AND. ( abstol.LE.zero ) ) THEN
429 CALL dcopy( n, work( indd ), 1, w, 1 )
430 indee = indwrk + 2*n
431 CALL dcopy( n-1, work( inde ), 1, work( indee ), 1 )
432 IF( .NOT.wantz ) THEN
433 CALL dsterf( n, w, work( indee ), info )
434 ELSE
435 CALL dlacpy( 'A', n, n, q, ldq, z, ldz )
436 CALL dsteqr( jobz, n, w, work( indee ), z, ldz,
437 $ work( indwrk ), info )
438 IF( info.EQ.0 ) THEN
439 DO 10 i = 1, n
440 ifail( i ) = 0
441 10 CONTINUE
442 END IF
443 END IF
444 IF( info.EQ.0 ) THEN
445 m = n
446 GO TO 30
447 END IF
448 info = 0
449 END IF
450*
451* Otherwise, call DSTEBZ and, if eigenvectors are desired,
452* call DSTEIN.
453*
454 IF( wantz ) THEN
455 order = 'B'
456 ELSE
457 order = 'E'
458 END IF
459 indisp = 1 + n
460 indiwo = indisp + n
461 CALL dstebz( range, order, n, vl, vu, il, iu, abstol,
462 $ work( indd ), work( inde ), m, nsplit, w,
463 $ iwork( 1 ), iwork( indisp ), work( indwrk ),
464 $ iwork( indiwo ), info )
465*
466 IF( wantz ) THEN
467 CALL dstein( n, work( indd ), work( inde ), m, w,
468 $ iwork( 1 ), iwork( indisp ), z, ldz,
469 $ work( indwrk ), iwork( indiwo ), ifail, info )
470*
471* Apply transformation matrix used in reduction to tridiagonal
472* form to eigenvectors returned by DSTEIN.
473*
474 DO 20 j = 1, m
475 CALL dcopy( n, z( 1, j ), 1, work( 1 ), 1 )
476 CALL dgemv( 'N', n, n, one, q, ldq, work, 1, zero,
477 $ z( 1, j ), 1 )
478 20 CONTINUE
479 END IF
480*
481 30 CONTINUE
482*
483* If eigenvalues are not in order, then sort them, along with
484* eigenvectors.
485*
486 IF( wantz ) THEN
487 DO 50 j = 1, m - 1
488 i = 0
489 tmp1 = w( j )
490 DO 40 jj = j + 1, m
491 IF( w( jj ).LT.tmp1 ) THEN
492 i = jj
493 tmp1 = w( jj )
494 END IF
495 40 CONTINUE
496*
497 IF( i.NE.0 ) THEN
498 itmp1 = iwork( 1 + i-1 )
499 w( i ) = w( j )
500 iwork( 1 + i-1 ) = iwork( 1 + j-1 )
501 w( j ) = tmp1
502 iwork( 1 + j-1 ) = itmp1
503 CALL dswap( n, z( 1, i ), 1, z( 1, j ), 1 )
504 IF( info.NE.0 ) THEN
505 itmp1 = ifail( i )
506 ifail( i ) = ifail( j )
507 ifail( j ) = itmp1
508 END IF
509 END IF
510 50 CONTINUE
511 END IF
512*
513 RETURN
514*
515* End of DSBGVX
516*
xerbla
subroutine xerbla(srname, info)
Definition cblat2.f:3285
dcopy
subroutine dcopy(n, dx, incx, dy, incy)
DCOPY
Definition dcopy.f:82
dgemv
subroutine dgemv(trans, m, n, alpha, a, lda, x, incx, beta, y, incy)
DGEMV
Definition dgemv.f:158
dsbgst
subroutine dsbgst(vect, uplo, n, ka, kb, ab, ldab, bb, ldbb, x, ldx, work, info)
DSBGST
Definition dsbgst.f:158
dsbtrd
subroutine dsbtrd(vect, uplo, n, kd, ab, ldab, d, e, q, ldq, work, info)
DSBTRD
Definition dsbtrd.f:161
dlacpy
subroutine dlacpy(uplo, m, n, a, lda, b, ldb)
DLACPY copies all or part of one two-dimensional array to another.
Definition dlacpy.f:101
lsame
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
dpbstf
subroutine dpbstf(uplo, n, kd, ab, ldab, info)
DPBSTF
Definition dpbstf.f:150
dstebz
subroutine dstebz(range, order, n, vl, vu, il, iu, abstol, d, e, m, nsplit, w, iblock, isplit, work, iwork, info)
DSTEBZ
Definition dstebz.f:272
dstein
subroutine dstein(n, d, e, m, w, iblock, isplit, z, ldz, work, iwork, ifail, info)
DSTEIN
Definition dstein.f:172
dsteqr
subroutine dsteqr(compz, n, d, e, z, ldz, work, info)
DSTEQR
Definition dsteqr.f:129
dsterf
subroutine dsterf(n, d, e, info)
DSTERF
Definition dsterf.f:84
dswap
subroutine dswap(n, dx, incx, dy, incy)
DSWAP
Definition dswap.f:82
Here is the call graph for this function:
Here is the caller graph for this function: