LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine zhpevd ( character  JOBZ,
character  UPLO,
integer  N,
complex*16, dimension( * )  AP,
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 
)

ZHPEVD computes the eigenvalues and, optionally, the left and/or right eigenvectors for OTHER matrices

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

Purpose:
 ZHPEVD computes all the eigenvalues and, optionally, eigenvectors of
 a complex Hermitian matrix A in packed storage.  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]JOBZ
          JOBZ is CHARACTER*1
          = 'N':  Compute eigenvalues only;
          = 'V':  Compute eigenvalues and eigenvectors.
[in]UPLO
          UPLO is CHARACTER*1
          = 'U':  Upper triangle of A is stored;
          = 'L':  Lower triangle of A is stored.
[in]N
          N is INTEGER
          The order of the matrix A.  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, AP is overwritten by values generated during the
          reduction to tridiagonal form.  If UPLO = 'U', the diagonal
          and first superdiagonal of the tridiagonal matrix T overwrite
          the corresponding elements of A, and if UPLO = 'L', the
          diagonal and first subdiagonal of T overwrite the
          corresponding elements of A.
[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 orthonormal
          eigenvectors of the matrix A, with the i-th column of Z
          holding the eigenvector associated with W(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 array WORK.
          If N <= 1,               LWORK must be at least 1.
          If JOBZ = 'N' and N > 1, LWORK must be at least N.
          If JOBZ = 'V' and N > 1, LWORK must be at least 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 (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 must be at least 1.
          If JOBZ = 'N' and N > 1, LRWORK must be at least N.
          If JOBZ = 'V' and N > 1, LRWORK must be at least
                    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 must be at least 1.
          If JOBZ  = 'V' and N > 1, LIWORK must be at least 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:  if INFO = i, the algorithm failed to converge; i
                off-diagonal elements of an intermediate tridiagonal
                form did not converge to zero.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011

Definition at line 203 of file zhpevd.f.

203 *
204 * -- LAPACK driver routine (version 3.4.0) --
205 * -- LAPACK is a software package provided by Univ. of Tennessee, --
206 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
207 * November 2011
208 *
209 * .. Scalar Arguments ..
210  CHARACTER jobz, uplo
211  INTEGER info, ldz, liwork, lrwork, lwork, n
212 * ..
213 * .. Array Arguments ..
214  INTEGER iwork( * )
215  DOUBLE PRECISION rwork( * ), w( * )
216  COMPLEX*16 ap( * ), work( * ), z( ldz, * )
217 * ..
218 *
219 * =====================================================================
220 *
221 * .. Parameters ..
222  DOUBLE PRECISION zero, one
223  parameter ( zero = 0.0d+0, one = 1.0d+0 )
224  COMPLEX*16 cone
225  parameter ( cone = ( 1.0d+0, 0.0d+0 ) )
226 * ..
227 * .. Local Scalars ..
228  LOGICAL lquery, wantz
229  INTEGER iinfo, imax, inde, indrwk, indtau, indwrk,
230  $ iscale, liwmin, llrwk, llwrk, lrwmin, lwmin
231  DOUBLE PRECISION anrm, bignum, eps, rmax, rmin, safmin, sigma,
232  $ smlnum
233 * ..
234 * .. External Functions ..
235  LOGICAL lsame
236  DOUBLE PRECISION dlamch, zlanhp
237  EXTERNAL lsame, dlamch, zlanhp
238 * ..
239 * .. External Subroutines ..
240  EXTERNAL dscal, dsterf, xerbla, zdscal, zhptrd, zstedc,
241  $ zupmtr
242 * ..
243 * .. Intrinsic Functions ..
244  INTRINSIC sqrt
245 * ..
246 * .. Executable Statements ..
247 *
248 * Test the input parameters.
249 *
250  wantz = lsame( jobz, 'V' )
251  lquery = ( lwork.EQ.-1 .OR. lrwork.EQ.-1 .OR. liwork.EQ.-1 )
252 *
253  info = 0
254  IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
255  info = -1
256  ELSE IF( .NOT.( lsame( uplo, 'L' ) .OR. lsame( uplo, 'U' ) ) )
257  $ THEN
258  info = -2
259  ELSE IF( n.LT.0 ) THEN
260  info = -3
261  ELSE IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) ) THEN
262  info = -7
263  END IF
264 *
265  IF( info.EQ.0 ) THEN
266  IF( n.LE.1 ) THEN
267  lwmin = 1
268  liwmin = 1
269  lrwmin = 1
270  ELSE
271  IF( wantz ) THEN
272  lwmin = 2*n
273  lrwmin = 1 + 5*n + 2*n**2
274  liwmin = 3 + 5*n
275  ELSE
276  lwmin = n
277  lrwmin = n
278  liwmin = 1
279  END IF
280  END IF
281  work( 1 ) = lwmin
282  rwork( 1 ) = lrwmin
283  iwork( 1 ) = liwmin
284 *
285  IF( lwork.LT.lwmin .AND. .NOT.lquery ) THEN
286  info = -9
287  ELSE IF( lrwork.LT.lrwmin .AND. .NOT.lquery ) THEN
288  info = -11
289  ELSE IF( liwork.LT.liwmin .AND. .NOT.lquery ) THEN
290  info = -13
291  END IF
292  END IF
293 *
294  IF( info.NE.0 ) THEN
295  CALL xerbla( 'ZHPEVD', -info )
296  RETURN
297  ELSE IF( lquery ) THEN
298  RETURN
299  END IF
300 *
301 * Quick return if possible
302 *
303  IF( n.EQ.0 )
304  $ RETURN
305 *
306  IF( n.EQ.1 ) THEN
307  w( 1 ) = ap( 1 )
308  IF( wantz )
309  $ z( 1, 1 ) = cone
310  RETURN
311  END IF
312 *
313 * Get machine constants.
314 *
315  safmin = dlamch( 'Safe minimum' )
316  eps = dlamch( 'Precision' )
317  smlnum = safmin / eps
318  bignum = one / smlnum
319  rmin = sqrt( smlnum )
320  rmax = sqrt( bignum )
321 *
322 * Scale matrix to allowable range, if necessary.
323 *
324  anrm = zlanhp( 'M', uplo, n, ap, rwork )
325  iscale = 0
326  IF( anrm.GT.zero .AND. anrm.LT.rmin ) THEN
327  iscale = 1
328  sigma = rmin / anrm
329  ELSE IF( anrm.GT.rmax ) THEN
330  iscale = 1
331  sigma = rmax / anrm
332  END IF
333  IF( iscale.EQ.1 ) THEN
334  CALL zdscal( ( n*( n+1 ) ) / 2, sigma, ap, 1 )
335  END IF
336 *
337 * Call ZHPTRD to reduce Hermitian packed matrix to tridiagonal form.
338 *
339  inde = 1
340  indtau = 1
341  indrwk = inde + n
342  indwrk = indtau + n
343  llwrk = lwork - indwrk + 1
344  llrwk = lrwork - indrwk + 1
345  CALL zhptrd( uplo, n, ap, w, rwork( inde ), work( indtau ),
346  $ iinfo )
347 *
348 * For eigenvalues only, call DSTERF. For eigenvectors, first call
349 * ZUPGTR to generate the orthogonal matrix, then call ZSTEDC.
350 *
351  IF( .NOT.wantz ) THEN
352  CALL dsterf( n, w, rwork( inde ), info )
353  ELSE
354  CALL zstedc( 'I', n, w, rwork( inde ), z, ldz, work( indwrk ),
355  $ llwrk, rwork( indrwk ), llrwk, iwork, liwork,
356  $ info )
357  CALL zupmtr( 'L', uplo, 'N', n, n, ap, work( indtau ), z, ldz,
358  $ work( indwrk ), iinfo )
359  END IF
360 *
361 * If matrix was scaled, then rescale eigenvalues appropriately.
362 *
363  IF( iscale.EQ.1 ) THEN
364  IF( info.EQ.0 ) THEN
365  imax = n
366  ELSE
367  imax = info - 1
368  END IF
369  CALL dscal( imax, one / sigma, w, 1 )
370  END IF
371 *
372  work( 1 ) = lwmin
373  rwork( 1 ) = lrwmin
374  iwork( 1 ) = liwmin
375  RETURN
376 *
377 * End of ZHPEVD
378 *
subroutine dsterf(N, D, E, INFO)
DSTERF
Definition: dsterf.f:88
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:65
double precision function zlanhp(NORM, UPLO, N, AP, WORK)
ZLANHP returns the value of the 1-norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a complex Hermitian matrix supplied in packed form.
Definition: zlanhp.f:119
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine zstedc(COMPZ, N, D, E, Z, LDZ, WORK, LWORK, RWORK, LRWORK, IWORK, LIWORK, INFO)
ZSTEDC
Definition: zstedc.f:215
subroutine zhptrd(UPLO, N, AP, D, E, TAU, INFO)
ZHPTRD
Definition: zhptrd.f:153
subroutine dscal(N, DA, DX, INCX)
DSCAL
Definition: dscal.f:55
subroutine zdscal(N, DA, ZX, INCX)
ZDSCAL
Definition: zdscal.f:54
subroutine zupmtr(SIDE, UPLO, TRANS, M, N, AP, TAU, C, LDC, WORK, INFO)
ZUPMTR
Definition: zupmtr.f:152
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: