LAPACK  3.9.0
LAPACK: Linear Algebra PACKage

◆ chpgvx()

subroutine chpgvx ( integer  ITYPE,
character  JOBZ,
character  RANGE,
character  UPLO,
integer  N,
complex, dimension( * )  AP,
complex, dimension( * )  BP,
real  VL,
real  VU,
integer  IL,
integer  IU,
real  ABSTOL,
integer  M,
real, dimension( * )  W,
complex, dimension( ldz, * )  Z,
integer  LDZ,
complex, dimension( * )  WORK,
real, dimension( * )  RWORK,
integer, dimension( * )  IWORK,
integer, dimension( * )  IFAIL,
integer  INFO 
)

CHPGVX

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

Purpose:
 CHPGVX computes selected eigenvalues and, optionally, eigenvectors
 of a complex generalized Hermitian-definite eigenproblem, of the form
 A*x=(lambda)*B*x,  A*Bx=(lambda)*x,  or B*A*x=(lambda)*x.  Here A and
 B are assumed to be Hermitian, stored in packed format, and B is also
 positive definite.  Eigenvalues and eigenvectors can be selected by
 specifying either a range of values or a range of indices for the
 desired eigenvalues.
Parameters
[in]ITYPE
          ITYPE is INTEGER
          Specifies the problem type to be solved:
          = 1:  A*x = (lambda)*B*x
          = 2:  A*B*x = (lambda)*x
          = 3:  B*A*x = (lambda)*x
[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,out]AP
          AP is COMPLEX array, dimension (N*(N+1)/2)
          On entry, the upper or lower triangle of the Hermitian matrix
          A, packed columnwise in a linear array.  The j-th column of A
          is stored in the array AP as follows:
          if UPLO = 'U', AP(i + (j-1)*j/2) = A(i,j) for 1<=i<=j;
          if UPLO = 'L', AP(i + (j-1)*(2*n-j)/2) = A(i,j) for j<=i<=n.

          On exit, the contents of AP are destroyed.
[in,out]BP
          BP is COMPLEX array, dimension (N*(N+1)/2)
          On entry, the upper or lower triangle of the Hermitian matrix
          B, packed columnwise in a linear array.  The j-th column of B
          is stored in the array BP as follows:
          if UPLO = 'U', BP(i + (j-1)*j/2) = B(i,j) for 1<=i<=j;
          if UPLO = 'L', BP(i + (j-1)*(2*n-j)/2) = B(i,j) for j<=i<=n.

          On exit, the triangular factor U or L from the Cholesky
          factorization B = U**H*U or B = L*L**H, in the same storage
          format as B.
[in]VL
          VL is REAL

          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 REAL

          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 REAL
          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 AP to tridiagonal form.

          Eigenvalues will be computed most accurately when ABSTOL is
          set to twice the underflow threshold 2*SLAMCH('S'), not zero.
          If this routine returns with INFO>0, indicating that some
          eigenvectors did not converge, try setting ABSTOL to
          2*SLAMCH('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 REAL array, dimension (N)
          On normal exit, the first M elements contain the selected
          eigenvalues in ascending order.
[out]Z
          Z is COMPLEX array, dimension (LDZ, N)
          If JOBZ = 'N', then Z is not referenced.
          If JOBZ = 'V', then if INFO = 0, the first M columns of Z
          contain the orthonormal eigenvectors of the matrix A
          corresponding to the selected eigenvalues, with the i-th
          column of Z holding the eigenvector associated with W(i).
          The eigenvectors are normalized as follows:
          if ITYPE = 1 or 2, Z**H*B*Z = I;
          if ITYPE = 3, Z**H*inv(B)*Z = I.

          If an eigenvector fails to converge, then that column of Z
          contains the latest approximation to the eigenvector, and the
          index of the eigenvector is returned in IFAIL.
          Note: the user must ensure that at least max(1,M) columns are
          supplied in the array Z; if RANGE = 'V', the exact value of M
          is not known in advance and an upper bound must be used.
[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 COMPLEX array, dimension (2*N)
[out]RWORK
          RWORK is REAL array, dimension (7*N)
[out]IWORK
          IWORK is INTEGER array, dimension (5*N)
[out]IFAIL
          IFAIL is INTEGER array, dimension (N)
          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 eigenvectors 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
          > 0:  CPPTRF or CHPEVX returned an error code:
             <= N:  if INFO = i, CHPEVX failed to converge;
                    i eigenvectors failed to converge.  Their indices
                    are stored in array IFAIL.
             > N:   if INFO = N + i, for 1 <= i <= n, then the leading
                    minor of order i of 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.
Date
June 2016
Contributors:
Mark Fahey, Department of Mathematics, Univ. of Kentucky, USA

Definition at line 279 of file chpgvx.f.

279 *
280 * -- LAPACK driver routine (version 3.7.0) --
281 * -- LAPACK is a software package provided by Univ. of Tennessee, --
282 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
283 * June 2016
284 *
285 * .. Scalar Arguments ..
286  CHARACTER JOBZ, RANGE, UPLO
287  INTEGER IL, INFO, ITYPE, IU, LDZ, M, N
288  REAL ABSTOL, VL, VU
289 * ..
290 * .. Array Arguments ..
291  INTEGER IFAIL( * ), IWORK( * )
292  REAL RWORK( * ), W( * )
293  COMPLEX AP( * ), BP( * ), WORK( * ), Z( LDZ, * )
294 * ..
295 *
296 * =====================================================================
297 *
298 * .. Local Scalars ..
299  LOGICAL ALLEIG, INDEIG, UPPER, VALEIG, WANTZ
300  CHARACTER TRANS
301  INTEGER J
302 * ..
303 * .. External Functions ..
304  LOGICAL LSAME
305  EXTERNAL lsame
306 * ..
307 * .. External Subroutines ..
308  EXTERNAL chpevx, chpgst, cpptrf, ctpmv, ctpsv, xerbla
309 * ..
310 * .. Intrinsic Functions ..
311  INTRINSIC min
312 * ..
313 * .. Executable Statements ..
314 *
315 * Test the input parameters.
316 *
317  wantz = lsame( jobz, 'V' )
318  upper = lsame( uplo, 'U' )
319  alleig = lsame( range, 'A' )
320  valeig = lsame( range, 'V' )
321  indeig = lsame( range, 'I' )
322 *
323  info = 0
324  IF( itype.LT.1 .OR. itype.GT.3 ) THEN
325  info = -1
326  ELSE IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
327  info = -2
328  ELSE IF( .NOT.( alleig .OR. valeig .OR. indeig ) ) THEN
329  info = -3
330  ELSE IF( .NOT.( upper .OR. lsame( uplo, 'L' ) ) ) THEN
331  info = -4
332  ELSE IF( n.LT.0 ) THEN
333  info = -5
334  ELSE
335  IF( valeig ) THEN
336  IF( n.GT.0 .AND. vu.LE.vl ) THEN
337  info = -9
338  END IF
339  ELSE IF( indeig ) THEN
340  IF( il.LT.1 ) THEN
341  info = -10
342  ELSE IF( iu.LT.min( n, il ) .OR. iu.GT.n ) THEN
343  info = -11
344  END IF
345  END IF
346  END IF
347  IF( info.EQ.0 ) THEN
348  IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) ) THEN
349  info = -16
350  END IF
351  END IF
352 *
353  IF( info.NE.0 ) THEN
354  CALL xerbla( 'CHPGVX', -info )
355  RETURN
356  END IF
357 *
358 * Quick return if possible
359 *
360  IF( n.EQ.0 )
361  $ RETURN
362 *
363 * Form a Cholesky factorization of B.
364 *
365  CALL cpptrf( uplo, n, bp, info )
366  IF( info.NE.0 ) THEN
367  info = n + info
368  RETURN
369  END IF
370 *
371 * Transform problem to standard eigenvalue problem and solve.
372 *
373  CALL chpgst( itype, uplo, n, ap, bp, info )
374  CALL chpevx( jobz, range, uplo, n, ap, vl, vu, il, iu, abstol, m,
375  $ w, z, ldz, work, rwork, iwork, ifail, info )
376 *
377  IF( wantz ) THEN
378 *
379 * Backtransform eigenvectors to the original problem.
380 *
381  IF( info.GT.0 )
382  $ m = info - 1
383  IF( itype.EQ.1 .OR. itype.EQ.2 ) THEN
384 *
385 * For A*x=(lambda)*B*x and A*B*x=(lambda)*x;
386 * backtransform eigenvectors: x = inv(L)**H*y or inv(U)*y
387 *
388  IF( upper ) THEN
389  trans = 'N'
390  ELSE
391  trans = 'C'
392  END IF
393 *
394  DO 10 j = 1, m
395  CALL ctpsv( uplo, trans, 'Non-unit', n, bp, z( 1, j ),
396  $ 1 )
397  10 CONTINUE
398 *
399  ELSE IF( itype.EQ.3 ) THEN
400 *
401 * For B*A*x=(lambda)*x;
402 * backtransform eigenvectors: x = L*y or U**H*y
403 *
404  IF( upper ) THEN
405  trans = 'C'
406  ELSE
407  trans = 'N'
408  END IF
409 *
410  DO 20 j = 1, m
411  CALL ctpmv( uplo, trans, 'Non-unit', n, bp, z( 1, j ),
412  $ 1 )
413  20 CONTINUE
414  END IF
415  END IF
416 *
417  RETURN
418 *
419 * End of CHPGVX
420 *
Here is the call graph for this function:
Here is the caller graph for this function:
xerbla
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
ctpsv
subroutine ctpsv(UPLO, TRANS, DIAG, N, AP, X, INCX)
CTPSV
Definition: ctpsv.f:146
lsame
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
cpptrf
subroutine cpptrf(UPLO, N, AP, INFO)
CPPTRF
Definition: cpptrf.f:121
ctpmv
subroutine ctpmv(UPLO, TRANS, DIAG, N, AP, X, INCX)
CTPMV
Definition: ctpmv.f:144
chpevx
subroutine chpevx(JOBZ, RANGE, UPLO, N, AP, VL, VU, IL, IU, ABSTOL, M, W, Z, LDZ, WORK, RWORK, IWORK, IFAIL, INFO)
CHPEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for OTHER matrice...
Definition: chpevx.f:242
chpgst
subroutine chpgst(ITYPE, UPLO, N, AP, BP, INFO)
CHPGST
Definition: chpgst.f:115