LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ zhpgvd()

subroutine zhpgvd ( integer  itype,
character  jobz,
character  uplo,
integer  n,
complex*16, dimension( * )  ap,
complex*16, dimension( * )  bp,
double precision, dimension( * )  w,
complex*16, dimension( ldz, * )  z,
integer  ldz,
complex*16, dimension( * )  work,
integer  lwork,
double precision, dimension( * )  rwork,
integer  lrwork,
integer, dimension( * )  iwork,
integer  liwork,
integer  info 
)

ZHPGVD

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

Purpose:
 ZHPGVD computes all the eigenvalues and, optionally, the 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.
 If eigenvectors are desired, it uses a divide and conquer algorithm.
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 COMPLEX*16 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*16 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.
[out]W
          W is DOUBLE PRECISION array, dimension (N)
          If INFO = 0, the eigenvalues in ascending order.
[out]Z
          Z is COMPLEX*16 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**H*B*Z = I;
          if ITYPE = 3, Z**H*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 COMPLEX*16 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 >= N.
          If JOBZ = 'V' and N > 1, LWORK >= 2*N.

          If LWORK = -1, then a workspace query is assumed; the routine
          only calculates the required sizes of the WORK, RWORK and
          IWORK arrays, returns these values as the first entries of
          the WORK, RWORK and IWORK arrays, and no error message
          related to LWORK or LRWORK or LIWORK is issued by XERBLA.
[out]RWORK
          RWORK is DOUBLE PRECISION array, dimension (MAX(1,LRWORK))
          On exit, if INFO = 0, RWORK(1) returns the required LRWORK.
[in]LRWORK
          LRWORK is INTEGER
          The dimension of array RWORK.
          If N <= 1,               LRWORK >= 1.
          If JOBZ = 'N' and N > 1, LRWORK >= N.
          If JOBZ = 'V' and N > 1, LRWORK >= 1 + 5*N + 2*N**2.

          If LRWORK = -1, then a workspace query is assumed; the
          routine only calculates the required sizes of the WORK, RWORK
          and IWORK arrays, returns these values as the first entries
          of the WORK, RWORK and IWORK arrays, and no error message
          related to LWORK or LRWORK 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 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, RWORK
          and IWORK arrays, returns these values as the first entries
          of the WORK, RWORK and IWORK arrays, and no error message
          related to LWORK or LRWORK 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:  ZPPTRF or ZHPEVD returned an error code:
             <= N:  if INFO = i, ZHPEVD failed to converge;
                    i off-diagonal elements of an intermediate
                    tridiagonal form did not convergeto zero;
             > N:   if INFO = N + i, for 1 <= i <= n, then the leading
                    principal minor of order i of B is not positive.
                    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.
Contributors:
Mark Fahey, Department of Mathematics, Univ. of Kentucky, USA

Definition at line 223 of file zhpgvd.f.

225*
226* -- LAPACK driver routine --
227* -- LAPACK is a software package provided by Univ. of Tennessee, --
228* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
229*
230* .. Scalar Arguments ..
231 CHARACTER JOBZ, UPLO
232 INTEGER INFO, ITYPE, LDZ, LIWORK, LRWORK, LWORK, N
233* ..
234* .. Array Arguments ..
235 INTEGER IWORK( * )
236 DOUBLE PRECISION RWORK( * ), W( * )
237 COMPLEX*16 AP( * ), BP( * ), WORK( * ), Z( LDZ, * )
238* ..
239*
240* =====================================================================
241*
242* .. Local Scalars ..
243 LOGICAL LQUERY, UPPER, WANTZ
244 CHARACTER TRANS
245 INTEGER J, LIWMIN, LRWMIN, LWMIN, NEIG
246* ..
247* .. External Functions ..
248 LOGICAL LSAME
249 EXTERNAL lsame
250* ..
251* .. External Subroutines ..
252 EXTERNAL xerbla, zhpevd, zhpgst, zpptrf, ztpmv, ztpsv
253* ..
254* .. Intrinsic Functions ..
255 INTRINSIC dble, max
256* ..
257* .. Executable Statements ..
258*
259* Test the input parameters.
260*
261 wantz = lsame( jobz, 'V' )
262 upper = lsame( uplo, 'U' )
263 lquery = ( lwork.EQ.-1 .OR. lrwork.EQ.-1 .OR. liwork.EQ.-1 )
264*
265 info = 0
266 IF( itype.LT.1 .OR. itype.GT.3 ) THEN
267 info = -1
268 ELSE IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
269 info = -2
270 ELSE IF( .NOT.( upper .OR. lsame( uplo, 'L' ) ) ) THEN
271 info = -3
272 ELSE IF( n.LT.0 ) THEN
273 info = -4
274 ELSE IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) ) THEN
275 info = -9
276 END IF
277*
278 IF( info.EQ.0 ) THEN
279 IF( n.LE.1 ) THEN
280 lwmin = 1
281 liwmin = 1
282 lrwmin = 1
283 ELSE
284 IF( wantz ) THEN
285 lwmin = 2*n
286 lrwmin = 1 + 5*n + 2*n**2
287 liwmin = 3 + 5*n
288 ELSE
289 lwmin = n
290 lrwmin = n
291 liwmin = 1
292 END IF
293 END IF
294*
295 work( 1 ) = lwmin
296 rwork( 1 ) = lrwmin
297 iwork( 1 ) = liwmin
298 IF( lwork.LT.lwmin .AND. .NOT.lquery ) THEN
299 info = -11
300 ELSE IF( lrwork.LT.lrwmin .AND. .NOT.lquery ) THEN
301 info = -13
302 ELSE IF( liwork.LT.liwmin .AND. .NOT.lquery ) THEN
303 info = -15
304 END IF
305 END IF
306*
307 IF( info.NE.0 ) THEN
308 CALL xerbla( 'ZHPGVD', -info )
309 RETURN
310 ELSE IF( lquery ) THEN
311 RETURN
312 END IF
313*
314* Quick return if possible
315*
316 IF( n.EQ.0 )
317 $ RETURN
318*
319* Form a Cholesky factorization of B.
320*
321 CALL zpptrf( uplo, n, bp, info )
322 IF( info.NE.0 ) THEN
323 info = n + info
324 RETURN
325 END IF
326*
327* Transform problem to standard eigenvalue problem and solve.
328*
329 CALL zhpgst( itype, uplo, n, ap, bp, info )
330 CALL zhpevd( jobz, uplo, n, ap, w, z, ldz, work, lwork, rwork,
331 $ lrwork, iwork, liwork, info )
332 lwmin = int( max( dble( lwmin ), dble( work( 1 ) ) ) )
333 lrwmin = int( max( dble( lrwmin ), dble( rwork( 1 ) ) ) )
334 liwmin = int( max( dble( liwmin ), dble( iwork( 1 ) ) ) )
335*
336 IF( wantz ) THEN
337*
338* Backtransform eigenvectors to the original problem.
339*
340 neig = n
341 IF( info.GT.0 )
342 $ neig = info - 1
343 IF( itype.EQ.1 .OR. itype.EQ.2 ) THEN
344*
345* For A*x=(lambda)*B*x and A*B*x=(lambda)*x;
346* backtransform eigenvectors: x = inv(L)**H *y or inv(U)*y
347*
348 IF( upper ) THEN
349 trans = 'N'
350 ELSE
351 trans = 'C'
352 END IF
353*
354 DO 10 j = 1, neig
355 CALL ztpsv( uplo, trans, 'Non-unit', n, bp, z( 1, j ),
356 $ 1 )
357 10 CONTINUE
358*
359 ELSE IF( itype.EQ.3 ) THEN
360*
361* For B*A*x=(lambda)*x;
362* backtransform eigenvectors: x = L*y or U**H *y
363*
364 IF( upper ) THEN
365 trans = 'C'
366 ELSE
367 trans = 'N'
368 END IF
369*
370 DO 20 j = 1, neig
371 CALL ztpmv( uplo, trans, 'Non-unit', n, bp, z( 1, j ),
372 $ 1 )
373 20 CONTINUE
374 END IF
375 END IF
376*
377 work( 1 ) = lwmin
378 rwork( 1 ) = lrwmin
379 iwork( 1 ) = liwmin
380 RETURN
381*
382* End of ZHPGVD
383*
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine zhpevd(jobz, uplo, n, ap, w, z, ldz, work, lwork, rwork, lrwork, iwork, liwork, info)
ZHPEVD computes the eigenvalues and, optionally, the left and/or right eigenvectors for OTHER matrice...
Definition zhpevd.f:194
subroutine zhpgst(itype, uplo, n, ap, bp, info)
ZHPGST
Definition zhpgst.f:113
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
subroutine zpptrf(uplo, n, ap, info)
ZPPTRF
Definition zpptrf.f:119
subroutine ztpmv(uplo, trans, diag, n, ap, x, incx)
ZTPMV
Definition ztpmv.f:142
subroutine ztpsv(uplo, trans, diag, n, ap, x, incx)
ZTPSV
Definition ztpsv.f:144
Here is the call graph for this function:
Here is the caller graph for this function: