Contents


NAME

     dspgvd - compute all the eigenvalues,  and  optionally,  the
     eigenvectors of a real generalized symmetric-definite eigen-
     problem, of the form A*x=(lambda)*B*x,  A*Bx=(lambda)*x,  or
     B*A*x=(lambda)*x

SYNOPSIS

     SUBROUTINE DSPGVD(ITYPE, JOBZ, UPLO, N, AP, BP, W, Z, LDZ, WORK,
           LWORK, IWORK, LIWORK, INFO)

     CHARACTER * 1 JOBZ, UPLO
     INTEGER ITYPE, N, LDZ, LWORK, LIWORK, INFO
     INTEGER IWORK(*)
     DOUBLE PRECISION AP(*), BP(*), W(*), Z(LDZ,*), WORK(*)

     SUBROUTINE DSPGVD_64(ITYPE, JOBZ, UPLO, N, AP, BP, W, Z, LDZ, WORK,
           LWORK, IWORK, LIWORK, INFO)

     CHARACTER * 1 JOBZ, UPLO
     INTEGER*8 ITYPE, N, LDZ, LWORK, LIWORK, INFO
     INTEGER*8 IWORK(*)
     DOUBLE PRECISION AP(*), BP(*), W(*), Z(LDZ,*), WORK(*)

  F95 INTERFACE
     SUBROUTINE SPGVD(ITYPE, JOBZ, UPLO, [N], AP, BP, W, Z, [LDZ], [WORK],
            [LWORK], [IWORK], [LIWORK], [INFO])

     CHARACTER(LEN=1) :: JOBZ, UPLO
     INTEGER :: ITYPE, N, LDZ, LWORK, LIWORK, INFO
     INTEGER, DIMENSION(:) :: IWORK
     REAL(8), DIMENSION(:) :: AP, BP, W, WORK
     REAL(8), DIMENSION(:,:) :: Z

     SUBROUTINE SPGVD_64(ITYPE, JOBZ, UPLO, [N], AP, BP, W, Z, [LDZ],
            [WORK], [LWORK], [IWORK], [LIWORK], [INFO])

     CHARACTER(LEN=1) :: JOBZ, UPLO
     INTEGER(8) :: ITYPE, N, LDZ, LWORK, LIWORK, INFO
     INTEGER(8), DIMENSION(:) :: IWORK
     REAL(8), DIMENSION(:) :: AP, BP, W, WORK
     REAL(8), DIMENSION(:,:) :: Z

  C INTERFACE
     #include <sunperf.h>

     void dspgvd(int itype, char jobz, char uplo, int  n,  double
               *ap,  double  *bp,  double *w, double *z, int ldz,
               int *info);

     void dspgvd_64(long itype, char jobz,  char  uplo,  long  n,
               double *ap, double *bp, double *w, double *z, long
               ldz, long *info);

PURPOSE

     dspgvd computes all the  eigenvalues,  and  optionally,  the
     eigenvectors of a real generalized symmetric-definite eigen-
     problem, 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 symmetric,
     stored in packed format, and B is also positive definite.
     If eigenvectors are desired, it uses a  divide  and  conquer
     algorithm.

     The divide and conquer algorithm makes very mild assumptions
     about  floating  point  arithmetic. It will work on machines
     with a guard digit  in  add/subtract,  or  on  those  binary
     machines  without  guard digits which subtract like the Cray
     X-MP, Cray Y-MP, Cray C-90, or Cray-2. It could  conceivably
     fail  on  hexadecimal  or  decimal  machines  without  guard
     digits, but we know of none.

ARGUMENTS

     ITYPE (input)
               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

     JOBZ (input)
               = 'N':  Compute eigenvalues only;
               = 'V':  Compute eigenvalues and eigenvectors.

     UPLO (input)
               = 'U':  Upper triangles of A and B are stored;
               = 'L':  Lower triangles of A and B are stored.

     N (input) The order of the matrices A and B.  N >= 0.

     AP (input/output)
               Double precision array, dimension  (N*(N+1)/2)  On
               entry,  the  upper  or  lower triangle of the sym-
               metric 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.

     BP (input/output)
               Double precision array, dimension  (N*(N+1)/2)  On
               entry,  the  upper  or  lower triangle of the sym-
               metric 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**T*U or B = L*L**T,
               in the same storage format as B.

     W (output)
               Double precision array, dimension (N) If INFO = 0,
               the eigenvalues in ascending order.

     Z (output)
               Double precision array, dimension (LDZ, N) If JOBZ
               =  'V',  then if INFO = 0, Z contains the matrix Z
               of eigenvectors.  The eigenvectors are  normalized
               as  follows:   if ITYPE = 1 or 2, Z**T*B*Z = I; if
               ITYPE = 3, Z**T*inv(B)*Z = I.  If JOBZ = 'N', then
               Z is not referenced.

     LDZ (input)
               The leading dimension of the array Z.  LDZ  >=  1,
               and if JOBZ = 'V', LDZ >= max(1,N).

     WORK (workspace/output)
               Double precision array, dimension (LWORK) On exit,
               if INFO = 0, WORK(1) returns the optimal LWORK.

     LWORK (input)
               The dimension of the  array  WORK.   If  N  <=  1,
               LWORK  >=  1.   If  JOBZ = 'N' and N > 1, LWORK >=
               2*N.  If JOBZ = 'V' and N > 1, LWORK >= 1 + 6*N  +
               2*N**2.

               If LWORK = -1, then a workspace query is  assumed;
               the  routine  only  calculates the optimal size of
               the WORK array, returns this value  as  the  first
               entry  of  the  WORK  array,  and no error message
               related to LWORK is issued by XERBLA.

     IWORK (workspace/output)
               Integer array, dimension (LIWORK) On exit, if INFO
               = 0, IWORK(1) returns the optimal LIWORK.

     LIWORK (input)
               The dimension of the array IWORK.  If JOBZ  =  'N'
               or N <= 1, LIWORK >= 1.  If JOBZ  = 'V' and N > 1,
               LIWORK >= 3 + 5*N.

               If LIWORK = -1, then a workspace query is assumed;
               the  routine  only  calculates the optimal size of
               the IWORK array, returns this value as  the  first
               entry  of  the  IWORK  array, and no error message
               related to LIWORK is issued by XERBLA.

     INFO (output)
               = 0:  successful exit
               < 0:  if INFO = -i, the i-th argument had an ille-
               gal value
               > 0:  DPPTRF or DSPEVD returned an error code:
               <= N:  if INFO = i, DSPEVD failed to  converge;  i
               off-diagonal elements of an intermediate tridiago-
               nal form did not converge to zero; > 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.

FURTHER DETAILS

     Based on contributions by
        Mark Fahey, Department of Mathematics, Univ. of Kentucky,
     USA