LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine sspgvd ( integer  ITYPE,
character  JOBZ,
character  UPLO,
integer  N,
real, dimension( * )  AP,
real, dimension( * )  BP,
real, dimension( * )  W,
real, dimension( ldz, * )  Z,
integer  LDZ,
real, dimension( * )  WORK,
integer  LWORK,
integer, dimension( * )  IWORK,
integer  LIWORK,
integer  INFO 
)

SSPGVD

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

Purpose:
 SSPGVD 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.
 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.
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 REAL 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 REAL 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 REAL array, dimension (N)
          If INFO = 0, the eigenvalues in ascending order.
[out]Z
          Z is REAL 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 REAL array, dimension (MAX(1,LWORK))
          On exit, if INFO = 0, WORK(1) returns the required LWORK.
[in]LWORK
          LWORK is INTEGER
          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 required sizes of the WORK and IWORK
          arrays, returns these values as the first entries of the WORK
          and IWORK arrays, and no error message related to LWORK or
          LIWORK is issued by XERBLA.
[out]IWORK
          IWORK is INTEGER array, dimension (MAX(1,LIWORK))
          On exit, if INFO = 0, IWORK(1) returns the required LIWORK.
[in]LIWORK
          LIWORK is INTEGER
          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 required sizes of the WORK and
          IWORK arrays, returns these values as the first entries of
          the WORK and IWORK arrays, and no error message related to
          LWORK or LIWORK is issued by XERBLA.
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  SPPTRF or SSPEVD returned an error code:
             <= N:  if INFO = i, SSPEVD 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
November 2015
Contributors:
Mark Fahey, Department of Mathematics, Univ. of Kentucky, USA

Definition at line 212 of file sspgvd.f.

212 *
213 * -- LAPACK driver routine (version 3.6.0) --
214 * -- LAPACK is a software package provided by Univ. of Tennessee, --
215 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
216 * November 2015
217 *
218 * .. Scalar Arguments ..
219  CHARACTER jobz, uplo
220  INTEGER info, itype, ldz, liwork, lwork, n
221 * ..
222 * .. Array Arguments ..
223  INTEGER iwork( * )
224  REAL ap( * ), bp( * ), w( * ), work( * ),
225  $ z( ldz, * )
226 * ..
227 *
228 * =====================================================================
229 *
230 * .. Local Scalars ..
231  LOGICAL lquery, upper, wantz
232  CHARACTER trans
233  INTEGER j, liwmin, lwmin, neig
234 * ..
235 * .. External Functions ..
236  LOGICAL lsame
237  EXTERNAL lsame
238 * ..
239 * .. External Subroutines ..
240  EXTERNAL spptrf, sspevd, sspgst, stpmv, stpsv, xerbla
241 * ..
242 * .. Intrinsic Functions ..
243  INTRINSIC max, real
244 * ..
245 * .. Executable Statements ..
246 *
247 * Test the input parameters.
248 *
249  wantz = lsame( jobz, 'V' )
250  upper = lsame( uplo, 'U' )
251  lquery = ( lwork.EQ.-1 .OR. liwork.EQ.-1 )
252 *
253  info = 0
254  IF( itype.LT.1 .OR. itype.GT.3 ) THEN
255  info = -1
256  ELSE IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
257  info = -2
258  ELSE IF( .NOT.( upper .OR. lsame( uplo, 'L' ) ) ) THEN
259  info = -3
260  ELSE IF( n.LT.0 ) THEN
261  info = -4
262  ELSE IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) ) THEN
263  info = -9
264  END IF
265 *
266  IF( info.EQ.0 ) THEN
267  IF( n.LE.1 ) THEN
268  liwmin = 1
269  lwmin = 1
270  ELSE
271  IF( wantz ) THEN
272  liwmin = 3 + 5*n
273  lwmin = 1 + 6*n + 2*n**2
274  ELSE
275  liwmin = 1
276  lwmin = 2*n
277  END IF
278  END IF
279  work( 1 ) = lwmin
280  iwork( 1 ) = liwmin
281  IF( lwork.LT.lwmin .AND. .NOT.lquery ) THEN
282  info = -11
283  ELSE IF( liwork.LT.liwmin .AND. .NOT.lquery ) THEN
284  info = -13
285  END IF
286  END IF
287 *
288  IF( info.NE.0 ) THEN
289  CALL xerbla( 'SSPGVD', -info )
290  RETURN
291  ELSE IF( lquery ) THEN
292  RETURN
293  END IF
294 *
295 * Quick return if possible
296 *
297  IF( n.EQ.0 )
298  $ RETURN
299 *
300 * Form a Cholesky factorization of BP.
301 *
302  CALL spptrf( uplo, n, bp, info )
303  IF( info.NE.0 ) THEN
304  info = n + info
305  RETURN
306  END IF
307 *
308 * Transform problem to standard eigenvalue problem and solve.
309 *
310  CALL sspgst( itype, uplo, n, ap, bp, info )
311  CALL sspevd( jobz, uplo, n, ap, w, z, ldz, work, lwork, iwork,
312  $ liwork, info )
313  lwmin = max( REAL( LWMIN ), REAL( WORK( 1 ) ) )
314  liwmin = max( REAL( LIWMIN ), REAL( IWORK( 1 ) ) )
315 *
316  IF( wantz ) THEN
317 *
318 * Backtransform eigenvectors to the original problem.
319 *
320  neig = n
321  IF( info.GT.0 )
322  $ neig = info - 1
323  IF( itype.EQ.1 .OR. itype.EQ.2 ) THEN
324 *
325 * For A*x=(lambda)*B*x and A*B*x=(lambda)*x;
326 * backtransform eigenvectors: x = inv(L)**T *y or inv(U)*y
327 *
328  IF( upper ) THEN
329  trans = 'N'
330  ELSE
331  trans = 'T'
332  END IF
333 *
334  DO 10 j = 1, neig
335  CALL stpsv( uplo, trans, 'Non-unit', n, bp, z( 1, j ),
336  $ 1 )
337  10 CONTINUE
338 *
339  ELSE IF( itype.EQ.3 ) THEN
340 *
341 * For B*A*x=(lambda)*x;
342 * backtransform eigenvectors: x = L*y or U**T *y
343 *
344  IF( upper ) THEN
345  trans = 'T'
346  ELSE
347  trans = 'N'
348  END IF
349 *
350  DO 20 j = 1, neig
351  CALL stpmv( uplo, trans, 'Non-unit', n, bp, z( 1, j ),
352  $ 1 )
353  20 CONTINUE
354  END IF
355  END IF
356 *
357  work( 1 ) = lwmin
358  iwork( 1 ) = liwmin
359 *
360  RETURN
361 *
362 * End of SSPGVD
363 *
subroutine sspevd(JOBZ, UPLO, N, AP, W, Z, LDZ, WORK, LWORK, IWORK, LIWORK, INFO)
SSPEVD computes the eigenvalues and, optionally, the left and/or right eigenvectors for OTHER matric...
Definition: sspevd.f:180
subroutine stpmv(UPLO, TRANS, DIAG, N, AP, X, INCX)
STPMV
Definition: stpmv.f:144
subroutine sspgst(ITYPE, UPLO, N, AP, BP, INFO)
SSPGST
Definition: sspgst.f:115
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine stpsv(UPLO, TRANS, DIAG, N, AP, X, INCX)
STPSV
Definition: stpsv.f:146
subroutine spptrf(UPLO, N, AP, INFO)
SPPTRF
Definition: spptrf.f:121
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55

Here is the call graph for this function:

Here is the caller graph for this function: