LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ zhbevd()

subroutine zhbevd ( character  JOBZ,
character  UPLO,
integer  N,
integer  KD,
complex*16, dimension( ldab, * )  AB,
integer  LDAB,
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 
)

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

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

Purpose:
 ZHBEVD computes all the eigenvalues and, optionally, eigenvectors of
 a complex Hermitian band matrix A.  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]KD
          KD is INTEGER
          The number of superdiagonals of the matrix A if UPLO = 'U',
          or the number of subdiagonals if UPLO = 'L'.  KD >= 0.
[in,out]AB
          AB is COMPLEX*16 array, dimension (LDAB, N)
          On entry, the upper or lower triangle of the Hermitian band
          matrix A, stored in the first KD+1 rows of the array.  The
          j-th column of A is stored in the j-th column of the array AB
          as follows:
          if UPLO = 'U', AB(kd+1+i-j,j) = A(i,j) for max(1,j-kd)<=i<=j;
          if UPLO = 'L', AB(1+i-j,j)    = A(i,j) for j<=i<=min(n,j+kd).

          On exit, AB is overwritten by values generated during the
          reduction to tridiagonal form.  If UPLO = 'U', the first
          superdiagonal and the diagonal of the tridiagonal matrix T
          are returned in rows KD and KD+1 of AB, and if UPLO = 'L',
          the diagonal and first subdiagonal of T are returned in the
          first two rows of AB.
[in]LDAB
          LDAB is INTEGER
          The leading dimension of the array AB.  LDAB >= KD + 1.
[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 optimal LWORK.
[in]LWORK
          LWORK is INTEGER
          The dimension of the 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**2.

          If LWORK = -1, then a workspace query is assumed; the routine
          only calculates the optimal 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 optimal 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 optimal 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 optimal 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 optimal 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.

Definition at line 213 of file zhbevd.f.

215 *
216 * -- LAPACK driver routine --
217 * -- LAPACK is a software package provided by Univ. of Tennessee, --
218 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
219 *
220 * .. Scalar Arguments ..
221  CHARACTER JOBZ, UPLO
222  INTEGER INFO, KD, LDAB, LDZ, LIWORK, LRWORK, LWORK, N
223 * ..
224 * .. Array Arguments ..
225  INTEGER IWORK( * )
226  DOUBLE PRECISION RWORK( * ), W( * )
227  COMPLEX*16 AB( LDAB, * ), WORK( * ), Z( LDZ, * )
228 * ..
229 *
230 * =====================================================================
231 *
232 * .. Parameters ..
233  DOUBLE PRECISION ZERO, ONE
234  parameter( zero = 0.0d0, one = 1.0d0 )
235  COMPLEX*16 CZERO, CONE
236  parameter( czero = ( 0.0d0, 0.0d0 ),
237  $ cone = ( 1.0d0, 0.0d0 ) )
238 * ..
239 * .. Local Scalars ..
240  LOGICAL LOWER, LQUERY, WANTZ
241  INTEGER IINFO, IMAX, INDE, INDWK2, INDWRK, ISCALE,
242  $ LIWMIN, LLRWK, LLWK2, LRWMIN, LWMIN
243  DOUBLE PRECISION ANRM, BIGNUM, EPS, RMAX, RMIN, SAFMIN, SIGMA,
244  $ SMLNUM
245 * ..
246 * .. External Functions ..
247  LOGICAL LSAME
248  DOUBLE PRECISION DLAMCH, ZLANHB
249  EXTERNAL lsame, dlamch, zlanhb
250 * ..
251 * .. External Subroutines ..
252  EXTERNAL dscal, dsterf, xerbla, zgemm, zhbtrd, zlacpy,
253  $ zlascl, zstedc
254 * ..
255 * .. Intrinsic Functions ..
256  INTRINSIC sqrt
257 * ..
258 * .. Executable Statements ..
259 *
260 * Test the input parameters.
261 *
262  wantz = lsame( jobz, 'V' )
263  lower = lsame( uplo, 'L' )
264  lquery = ( lwork.EQ.-1 .OR. liwork.EQ.-1 .OR. lrwork.EQ.-1 )
265 *
266  info = 0
267  IF( n.LE.1 ) THEN
268  lwmin = 1
269  lrwmin = 1
270  liwmin = 1
271  ELSE
272  IF( wantz ) THEN
273  lwmin = 2*n**2
274  lrwmin = 1 + 5*n + 2*n**2
275  liwmin = 3 + 5*n
276  ELSE
277  lwmin = n
278  lrwmin = n
279  liwmin = 1
280  END IF
281  END IF
282  IF( .NOT.( wantz .OR. lsame( jobz, 'N' ) ) ) THEN
283  info = -1
284  ELSE IF( .NOT.( lower .OR. lsame( uplo, 'U' ) ) ) THEN
285  info = -2
286  ELSE IF( n.LT.0 ) THEN
287  info = -3
288  ELSE IF( kd.LT.0 ) THEN
289  info = -4
290  ELSE IF( ldab.LT.kd+1 ) THEN
291  info = -6
292  ELSE IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.n ) ) THEN
293  info = -9
294  END IF
295 *
296  IF( info.EQ.0 ) THEN
297  work( 1 ) = lwmin
298  rwork( 1 ) = lrwmin
299  iwork( 1 ) = liwmin
300 *
301  IF( lwork.LT.lwmin .AND. .NOT.lquery ) THEN
302  info = -11
303  ELSE IF( lrwork.LT.lrwmin .AND. .NOT.lquery ) THEN
304  info = -13
305  ELSE IF( liwork.LT.liwmin .AND. .NOT.lquery ) THEN
306  info = -15
307  END IF
308  END IF
309 *
310  IF( info.NE.0 ) THEN
311  CALL xerbla( 'ZHBEVD', -info )
312  RETURN
313  ELSE IF( lquery ) THEN
314  RETURN
315  END IF
316 *
317 * Quick return if possible
318 *
319  IF( n.EQ.0 )
320  $ RETURN
321 *
322  IF( n.EQ.1 ) THEN
323  w( 1 ) = dble( ab( 1, 1 ) )
324  IF( wantz )
325  $ z( 1, 1 ) = cone
326  RETURN
327  END IF
328 *
329 * Get machine constants.
330 *
331  safmin = dlamch( 'Safe minimum' )
332  eps = dlamch( 'Precision' )
333  smlnum = safmin / eps
334  bignum = one / smlnum
335  rmin = sqrt( smlnum )
336  rmax = sqrt( bignum )
337 *
338 * Scale matrix to allowable range, if necessary.
339 *
340  anrm = zlanhb( 'M', uplo, n, kd, ab, ldab, rwork )
341  iscale = 0
342  IF( anrm.GT.zero .AND. anrm.LT.rmin ) THEN
343  iscale = 1
344  sigma = rmin / anrm
345  ELSE IF( anrm.GT.rmax ) THEN
346  iscale = 1
347  sigma = rmax / anrm
348  END IF
349  IF( iscale.EQ.1 ) THEN
350  IF( lower ) THEN
351  CALL zlascl( 'B', kd, kd, one, sigma, n, n, ab, ldab, info )
352  ELSE
353  CALL zlascl( 'Q', kd, kd, one, sigma, n, n, ab, ldab, info )
354  END IF
355  END IF
356 *
357 * Call ZHBTRD to reduce Hermitian band matrix to tridiagonal form.
358 *
359  inde = 1
360  indwrk = inde + n
361  indwk2 = 1 + n*n
362  llwk2 = lwork - indwk2 + 1
363  llrwk = lrwork - indwrk + 1
364  CALL zhbtrd( jobz, uplo, n, kd, ab, ldab, w, rwork( inde ), z,
365  $ ldz, work, iinfo )
366 *
367 * For eigenvalues only, call DSTERF. For eigenvectors, call ZSTEDC.
368 *
369  IF( .NOT.wantz ) THEN
370  CALL dsterf( n, w, rwork( inde ), info )
371  ELSE
372  CALL zstedc( 'I', n, w, rwork( inde ), work, n, work( indwk2 ),
373  $ llwk2, rwork( indwrk ), llrwk, iwork, liwork,
374  $ info )
375  CALL zgemm( 'N', 'N', n, n, n, cone, z, ldz, work, n, czero,
376  $ work( indwk2 ), n )
377  CALL zlacpy( 'A', n, n, work( indwk2 ), n, z, ldz )
378  END IF
379 *
380 * If matrix was scaled, then rescale eigenvalues appropriately.
381 *
382  IF( iscale.EQ.1 ) THEN
383  IF( info.EQ.0 ) THEN
384  imax = n
385  ELSE
386  imax = info - 1
387  END IF
388  CALL dscal( imax, one / sigma, w, 1 )
389  END IF
390 *
391  work( 1 ) = lwmin
392  rwork( 1 ) = lrwmin
393  iwork( 1 ) = liwmin
394  RETURN
395 *
396 * End of ZHBEVD
397 *
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:69
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
subroutine dsterf(N, D, E, INFO)
DSTERF
Definition: dsterf.f:86
subroutine zgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
ZGEMM
Definition: zgemm.f:187
subroutine zlascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
ZLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: zlascl.f:143
subroutine zlacpy(UPLO, M, N, A, LDA, B, LDB)
ZLACPY copies all or part of one two-dimensional array to another.
Definition: zlacpy.f:103
double precision function zlanhb(NORM, UPLO, N, K, AB, LDAB, WORK)
ZLANHB returns the value of the 1-norm, or the Frobenius norm, or the infinity norm,...
Definition: zlanhb.f:132
subroutine zstedc(COMPZ, N, D, E, Z, LDZ, WORK, LWORK, RWORK, LRWORK, IWORK, LIWORK, INFO)
ZSTEDC
Definition: zstedc.f:212
subroutine zhbtrd(VECT, UPLO, N, KD, AB, LDAB, D, E, Q, LDQ, WORK, INFO)
ZHBTRD
Definition: zhbtrd.f:163
subroutine dscal(N, DA, DX, INCX)
DSCAL
Definition: dscal.f:79
Here is the call graph for this function:
Here is the caller graph for this function: