LAPACK  3.9.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.
Date
June 2017

Definition at line 162 of file dspgv.f.

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