LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ dspgv()

subroutine dspgv ( integer  ITYPE,
character  JOBZ,
character  UPLO,
integer  N,
double precision, dimension( * )  AP,
double precision, dimension( * )  BP,
double precision, dimension( * )  W,
double precision, dimension( ldz, * )  Z,
integer  LDZ,
double precision, dimension( * )  WORK,
integer  INFO 
)

DSPGV

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

Purpose:
 DSPGV computes all the eigenvalues and, optionally, the eigenvectors
 of a real generalized symmetric-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 symmetric, stored in packed format,
 and B is also positive definite.
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]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 DOUBLE PRECISION array, dimension (N*(N+1)/2)
          On entry, the upper or lower triangle of the symmetric 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 DOUBLE PRECISION array, dimension (N*(N+1)/2)
          On entry, the upper or lower triangle of the symmetric 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.
[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.  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.
[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 (3*N)
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  DPPTRF or DSPEV returned an error code:
             <= N:  if INFO = i, DSPEV 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 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.

Definition at line 158 of file dspgv.f.

160 *
161 * -- LAPACK driver routine --
162 * -- LAPACK is a software package provided by Univ. of Tennessee, --
163 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
164 *
165 * .. Scalar Arguments ..
166  CHARACTER JOBZ, UPLO
167  INTEGER INFO, ITYPE, LDZ, N
168 * ..
169 * .. Array Arguments ..
170  DOUBLE PRECISION AP( * ), BP( * ), W( * ), WORK( * ),
171  $ Z( LDZ, * )
172 * ..
173 *
174 * =====================================================================
175 *
176 * .. Local Scalars ..
177  LOGICAL UPPER, WANTZ
178  CHARACTER TRANS
179  INTEGER J, NEIG
180 * ..
181 * .. External Functions ..
182  LOGICAL LSAME
183  EXTERNAL lsame
184 * ..
185 * .. External Subroutines ..
186  EXTERNAL dpptrf, dspev, dspgst, dtpmv, dtpsv, xerbla
187 * ..
188 * .. Executable Statements ..
189 *
190 * Test the input parameters.
191 *
192  wantz = lsame( jobz, 'V' )
193  upper = lsame( uplo, 'U' )
194 *
195  info = 0
196  IF( itype.LT.1 .OR. itype.GT.3 ) THEN
197  info = -1
198  ELSE IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
199  info = -2
200  ELSE IF( .NOT.( upper .OR. lsame( uplo, 'L' ) ) ) THEN
201  info = -3
202  ELSE IF( n.LT.0 ) THEN
203  info = -4
204  ELSE IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) ) THEN
205  info = -9
206  END IF
207  IF( info.NE.0 ) THEN
208  CALL xerbla( 'DSPGV ', -info )
209  RETURN
210  END IF
211 *
212 * Quick return if possible
213 *
214  IF( n.EQ.0 )
215  $ RETURN
216 *
217 * Form a Cholesky factorization of B.
218 *
219  CALL dpptrf( uplo, n, bp, info )
220  IF( info.NE.0 ) THEN
221  info = n + info
222  RETURN
223  END IF
224 *
225 * Transform problem to standard eigenvalue problem and solve.
226 *
227  CALL dspgst( itype, uplo, n, ap, bp, info )
228  CALL dspev( jobz, uplo, n, ap, w, z, ldz, work, info )
229 *
230  IF( wantz ) THEN
231 *
232 * Backtransform eigenvectors to the original problem.
233 *
234  neig = n
235  IF( info.GT.0 )
236  $ neig = info - 1
237  IF( itype.EQ.1 .OR. itype.EQ.2 ) THEN
238 *
239 * For A*x=(lambda)*B*x and A*B*x=(lambda)*x;
240 * backtransform eigenvectors: x = inv(L)**T*y or inv(U)*y
241 *
242  IF( upper ) THEN
243  trans = 'N'
244  ELSE
245  trans = 'T'
246  END IF
247 *
248  DO 10 j = 1, neig
249  CALL dtpsv( uplo, trans, 'Non-unit', n, bp, z( 1, j ),
250  $ 1 )
251  10 CONTINUE
252 *
253  ELSE IF( itype.EQ.3 ) THEN
254 *
255 * For B*A*x=(lambda)*x;
256 * backtransform eigenvectors: x = L*y or U**T*y
257 *
258  IF( upper ) THEN
259  trans = 'T'
260  ELSE
261  trans = 'N'
262  END IF
263 *
264  DO 20 j = 1, neig
265  CALL dtpmv( uplo, trans, 'Non-unit', n, bp, z( 1, j ),
266  $ 1 )
267  20 CONTINUE
268  END IF
269  END IF
270  RETURN
271 *
272 * End of DSPGV
273 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
subroutine dtpsv(UPLO, TRANS, DIAG, N, AP, X, INCX)
DTPSV
Definition: dtpsv.f:144
subroutine dtpmv(UPLO, TRANS, DIAG, N, AP, X, INCX)
DTPMV
Definition: dtpmv.f:142
subroutine dpptrf(UPLO, N, AP, INFO)
DPPTRF
Definition: dpptrf.f:119
subroutine dspgst(ITYPE, UPLO, N, AP, BP, INFO)
DSPGST
Definition: dspgst.f:113
subroutine dspev(JOBZ, UPLO, N, AP, W, Z, LDZ, WORK, INFO)
DSPEV computes the eigenvalues and, optionally, the left and/or right eigenvectors for OTHER matrices
Definition: dspev.f:130
Here is the call graph for this function:
Here is the caller graph for this function: