 LAPACK  3.10.1 LAPACK: Linear Algebra PACKage

## ◆ dhseqr()

 subroutine dhseqr ( character JOB, character COMPZ, integer N, integer ILO, integer IHI, double precision, dimension( ldh, * ) H, integer LDH, double precision, dimension( * ) WR, double precision, dimension( * ) WI, double precision, dimension( ldz, * ) Z, integer LDZ, double precision, dimension( * ) WORK, integer LWORK, integer INFO )

DHSEQR

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

Purpose:
```    DHSEQR computes the eigenvalues of a Hessenberg matrix H
and, optionally, the matrices T and Z from the Schur decomposition
H = Z T Z**T, where T is an upper quasi-triangular matrix (the
Schur form), and Z is the orthogonal matrix of Schur vectors.

Optionally Z may be postmultiplied into an input orthogonal
matrix Q so that this routine can give the Schur factorization
of a matrix A which has been reduced to the Hessenberg form H
by the orthogonal matrix Q:  A = Q*H*Q**T = (QZ)*T*(QZ)**T.```
Parameters
 [in] JOB ``` JOB is CHARACTER*1 = 'E': compute eigenvalues only; = 'S': compute eigenvalues and the Schur form T.``` [in] COMPZ ``` COMPZ is CHARACTER*1 = 'N': no Schur vectors are computed; = 'I': Z is initialized to the unit matrix and the matrix Z of Schur vectors of H is returned; = 'V': Z must contain an orthogonal matrix Q on entry, and the product Q*Z is returned.``` [in] N ``` N is INTEGER The order of the matrix H. N >= 0.``` [in] ILO ` ILO is INTEGER` [in] IHI ``` IHI is INTEGER It is assumed that H is already upper triangular in rows and columns 1:ILO-1 and IHI+1:N. ILO and IHI are normally set by a previous call to DGEBAL, and then passed to ZGEHRD when the matrix output by DGEBAL is reduced to Hessenberg form. Otherwise ILO and IHI should be set to 1 and N respectively. If N > 0, then 1 <= ILO <= IHI <= N. If N = 0, then ILO = 1 and IHI = 0.``` [in,out] H ``` H is DOUBLE PRECISION array, dimension (LDH,N) On entry, the upper Hessenberg matrix H. On exit, if INFO = 0 and JOB = 'S', then H contains the upper quasi-triangular matrix T from the Schur decomposition (the Schur form); 2-by-2 diagonal blocks (corresponding to complex conjugate pairs of eigenvalues) are returned in standard form, with H(i,i) = H(i+1,i+1) and H(i+1,i)*H(i,i+1) < 0. If INFO = 0 and JOB = 'E', the contents of H are unspecified on exit. (The output value of H when INFO > 0 is given under the description of INFO below.) Unlike earlier versions of DHSEQR, this subroutine may explicitly H(i,j) = 0 for i > j and j = 1, 2, ... ILO-1 or j = IHI+1, IHI+2, ... N.``` [in] LDH ``` LDH is INTEGER The leading dimension of the array H. LDH >= max(1,N).``` [out] WR ` WR is DOUBLE PRECISION array, dimension (N)` [out] WI ``` WI is DOUBLE PRECISION array, dimension (N) The real and imaginary parts, respectively, of the computed eigenvalues. If two eigenvalues are computed as a complex conjugate pair, they are stored in consecutive elements of WR and WI, say the i-th and (i+1)th, with WI(i) > 0 and WI(i+1) < 0. If JOB = 'S', the eigenvalues are stored in the same order as on the diagonal of the Schur form returned in H, with WR(i) = H(i,i) and, if H(i:i+1,i:i+1) is a 2-by-2 diagonal block, WI(i) = sqrt(-H(i+1,i)*H(i,i+1)) and WI(i+1) = -WI(i).``` [in,out] Z ``` Z is DOUBLE PRECISION array, dimension (LDZ,N) If COMPZ = 'N', Z is not referenced. If COMPZ = 'I', on entry Z need not be set and on exit, if INFO = 0, Z contains the orthogonal matrix Z of the Schur vectors of H. If COMPZ = 'V', on entry Z must contain an N-by-N matrix Q, which is assumed to be equal to the unit matrix except for the submatrix Z(ILO:IHI,ILO:IHI). On exit, if INFO = 0, Z contains Q*Z. Normally Q is the orthogonal matrix generated by DORGHR after the call to DGEHRD which formed the Hessenberg matrix H. (The output value of Z when INFO > 0 is given under the description of INFO below.)``` [in] LDZ ``` LDZ is INTEGER The leading dimension of the array Z. if COMPZ = 'I' or COMPZ = 'V', then LDZ >= MAX(1,N). Otherwise, LDZ >= 1.``` [out] WORK ``` WORK is DOUBLE PRECISION array, dimension (LWORK) On exit, if INFO = 0, WORK(1) returns an estimate of the optimal value for LWORK.``` [in] LWORK ``` LWORK is INTEGER The dimension of the array WORK. LWORK >= max(1,N) is sufficient and delivers very good and sometimes optimal performance. However, LWORK as large as 11*N may be required for optimal performance. A workspace query is recommended to determine the optimal workspace size. If LWORK = -1, then DHSEQR does a workspace query. In this case, DHSEQR checks the input parameters and estimates the optimal workspace size for the given values of N, ILO and IHI. The estimate is returned in WORK(1). No error message related to LWORK is issued by XERBLA. Neither H nor Z are accessed.``` [out] INFO ``` INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value > 0: if INFO = i, DHSEQR failed to compute all of the eigenvalues. Elements 1:ilo-1 and i+1:n of WR and WI contain those eigenvalues which have been successfully computed. (Failures are rare.) If INFO > 0 and JOB = 'E', then on exit, the remaining unconverged eigenvalues are the eigen- values of the upper Hessenberg matrix rows and columns ILO through INFO of the final, output value of H. If INFO > 0 and JOB = 'S', then on exit (*) (initial value of H)*U = U*(final value of H) where U is an orthogonal matrix. The final value of H is upper Hessenberg and quasi-triangular in rows and columns INFO+1 through IHI. If INFO > 0 and COMPZ = 'V', then on exit (final value of Z) = (initial value of Z)*U where U is the orthogonal matrix in (*) (regard- less of the value of JOB.) If INFO > 0 and COMPZ = 'I', then on exit (final value of Z) = U where U is the orthogonal matrix in (*) (regard- less of the value of JOB.) If INFO > 0 and COMPZ = 'N', then Z is not accessed.```
Contributors:
Karen Braman and Ralph Byers, Department of Mathematics, University of Kansas, USA
Further Details:
```             Default values supplied by
ILAENV(ISPEC,'DHSEQR',JOB(:1)//COMPZ(:1),N,ILO,IHI,LWORK).
It is suggested that these defaults be adjusted in order
to attain best performance in each particular
computational environment.

ISPEC=12: The DLAHQR vs DLAQR0 crossover point.
Default: 75. (Must be at least 11.)

ISPEC=13: Recommended deflation window size.
This depends on ILO, IHI and NS.  NS is the
number of simultaneous shifts returned
by ILAENV(ISPEC=15).  (See ISPEC=15 below.)
The default for (IHI-ILO+1) <= 500 is NS.
The default for (IHI-ILO+1) >  500 is 3*NS/2.

ISPEC=14: Nibble crossover point. (See IPARMQ for
details.)  Default: 14% of deflation window
size.

ISPEC=15: Number of simultaneous shifts in a multishift
QR iteration.

If IHI-ILO+1 is ...

greater than      ...but less    ... the
or equal to ...      than        default is

1               30          NS =   2(+)
30               60          NS =   4(+)
60              150          NS =  10(+)
150              590          NS =  **
590             3000          NS =  64
3000             6000          NS = 128
6000             infinity      NS = 256

(+)  By default some or all matrices of this order
are passed to the implicit double shift routine
DLAHQR and this parameter is ignored.  See
ISPEC=12 above and comments in IPARMQ for
details.

(**)  The asterisks (**) indicate an ad-hoc
function of N increasing from 10 to 64.

ISPEC=16: Select structured matrix multiply.
If the number of simultaneous shifts (specified
by ISPEC=15) is less than 14, then the default
for ISPEC=16 is 0.  Otherwise the default for
ISPEC=16 is 2.```
References:
```  K. Braman, R. Byers and R. Mathias, The Multi-Shift QR
Algorithm Part I: Maintaining Well Focused Shifts, and Level 3
Performance, SIAM Journal of Matrix Analysis, volume 23, pages
929--947, 2002.
```

K. Braman, R. Byers and R. Mathias, The Multi-Shift QR Algorithm Part II: Aggressive Early Deflation, SIAM Journal of Matrix Analysis, volume 23, pages 948–973, 2002.

Definition at line 314 of file dhseqr.f.

316 *
317 * -- LAPACK computational routine --
318 * -- LAPACK is a software package provided by Univ. of Tennessee, --
319 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
320 *
321 * .. Scalar Arguments ..
322  INTEGER IHI, ILO, INFO, LDH, LDZ, LWORK, N
323  CHARACTER COMPZ, JOB
324 * ..
325 * .. Array Arguments ..
326  DOUBLE PRECISION H( LDH, * ), WI( * ), WORK( * ), WR( * ),
327  \$ Z( LDZ, * )
328 * ..
329 *
330 * =====================================================================
331 *
332 * .. Parameters ..
333 *
334 * ==== Matrices of order NTINY or smaller must be processed by
335 * . DLAHQR because of insufficient subdiagonal scratch space.
336 * . (This is a hard limit.) ====
337  INTEGER NTINY
338  parameter( ntiny = 15 )
339 *
340 * ==== NL allocates some local workspace to help small matrices
341 * . through a rare DLAHQR failure. NL > NTINY = 15 is
342 * . required and NL <= NMIN = ILAENV(ISPEC=12,...) is recom-
343 * . mended. (The default value of NMIN is 75.) Using NL = 49
344 * . allows up to six simultaneous shifts and a 16-by-16
345 * . deflation window. ====
346  INTEGER NL
347  parameter( nl = 49 )
348  DOUBLE PRECISION ZERO, ONE
349  parameter( zero = 0.0d0, one = 1.0d0 )
350 * ..
351 * .. Local Arrays ..
352  DOUBLE PRECISION HL( NL, NL ), WORKL( NL )
353 * ..
354 * .. Local Scalars ..
355  INTEGER I, KBOT, NMIN
356  LOGICAL INITZ, LQUERY, WANTT, WANTZ
357 * ..
358 * .. External Functions ..
359  INTEGER ILAENV
360  LOGICAL LSAME
361  EXTERNAL ilaenv, lsame
362 * ..
363 * .. External Subroutines ..
364  EXTERNAL dlacpy, dlahqr, dlaqr0, dlaset, xerbla
365 * ..
366 * .. Intrinsic Functions ..
367  INTRINSIC dble, max, min
368 * ..
369 * .. Executable Statements ..
370 *
371 * ==== Decode and check the input parameters. ====
372 *
373  wantt = lsame( job, 'S' )
374  initz = lsame( compz, 'I' )
375  wantz = initz .OR. lsame( compz, 'V' )
376  work( 1 ) = dble( max( 1, n ) )
377  lquery = lwork.EQ.-1
378 *
379  info = 0
380  IF( .NOT.lsame( job, 'E' ) .AND. .NOT.wantt ) THEN
381  info = -1
382  ELSE IF( .NOT.lsame( compz, 'N' ) .AND. .NOT.wantz ) THEN
383  info = -2
384  ELSE IF( n.LT.0 ) THEN
385  info = -3
386  ELSE IF( ilo.LT.1 .OR. ilo.GT.max( 1, n ) ) THEN
387  info = -4
388  ELSE IF( ihi.LT.min( ilo, n ) .OR. ihi.GT.n ) THEN
389  info = -5
390  ELSE IF( ldh.LT.max( 1, n ) ) THEN
391  info = -7
392  ELSE IF( ldz.LT.1 .OR. ( wantz .AND. ldz.LT.max( 1, n ) ) ) THEN
393  info = -11
394  ELSE IF( lwork.LT.max( 1, n ) .AND. .NOT.lquery ) THEN
395  info = -13
396  END IF
397 *
398  IF( info.NE.0 ) THEN
399 *
400 * ==== Quick return in case of invalid argument. ====
401 *
402  CALL xerbla( 'DHSEQR', -info )
403  RETURN
404 *
405  ELSE IF( n.EQ.0 ) THEN
406 *
407 * ==== Quick return in case N = 0; nothing to do. ====
408 *
409  RETURN
410 *
411  ELSE IF( lquery ) THEN
412 *
413 * ==== Quick return in case of a workspace query ====
414 *
415  CALL dlaqr0( wantt, wantz, n, ilo, ihi, h, ldh, wr, wi, ilo,
416  \$ ihi, z, ldz, work, lwork, info )
417 * ==== Ensure reported workspace size is backward-compatible with
418 * . previous LAPACK versions. ====
419  work( 1 ) = max( dble( max( 1, n ) ), work( 1 ) )
420  RETURN
421 *
422  ELSE
423 *
424 * ==== copy eigenvalues isolated by DGEBAL ====
425 *
426  DO 10 i = 1, ilo - 1
427  wr( i ) = h( i, i )
428  wi( i ) = zero
429  10 CONTINUE
430  DO 20 i = ihi + 1, n
431  wr( i ) = h( i, i )
432  wi( i ) = zero
433  20 CONTINUE
434 *
435 * ==== Initialize Z, if requested ====
436 *
437  IF( initz )
438  \$ CALL dlaset( 'A', n, n, zero, one, z, ldz )
439 *
440 * ==== Quick return if possible ====
441 *
442  IF( ilo.EQ.ihi ) THEN
443  wr( ilo ) = h( ilo, ilo )
444  wi( ilo ) = zero
445  RETURN
446  END IF
447 *
448 * ==== DLAHQR/DLAQR0 crossover point ====
449 *
450  nmin = ilaenv( 12, 'DHSEQR', job( : 1 ) // compz( : 1 ), n,
451  \$ ilo, ihi, lwork )
452  nmin = max( ntiny, nmin )
453 *
454 * ==== DLAQR0 for big matrices; DLAHQR for small ones ====
455 *
456  IF( n.GT.nmin ) THEN
457  CALL dlaqr0( wantt, wantz, n, ilo, ihi, h, ldh, wr, wi, ilo,
458  \$ ihi, z, ldz, work, lwork, info )
459  ELSE
460 *
461 * ==== Small matrix ====
462 *
463  CALL dlahqr( wantt, wantz, n, ilo, ihi, h, ldh, wr, wi, ilo,
464  \$ ihi, z, ldz, info )
465 *
466  IF( info.GT.0 ) THEN
467 *
468 * ==== A rare DLAHQR failure! DLAQR0 sometimes succeeds
469 * . when DLAHQR fails. ====
470 *
471  kbot = info
472 *
473  IF( n.GE.nl ) THEN
474 *
475 * ==== Larger matrices have enough subdiagonal scratch
476 * . space to call DLAQR0 directly. ====
477 *
478  CALL dlaqr0( wantt, wantz, n, ilo, kbot, h, ldh, wr,
479  \$ wi, ilo, ihi, z, ldz, work, lwork, info )
480 *
481  ELSE
482 *
483 * ==== Tiny matrices don't have enough subdiagonal
484 * . scratch space to benefit from DLAQR0. Hence,
485 * . tiny matrices must be copied into a larger
486 * . array before calling DLAQR0. ====
487 *
488  CALL dlacpy( 'A', n, n, h, ldh, hl, nl )
489  hl( n+1, n ) = zero
490  CALL dlaset( 'A', nl, nl-n, zero, zero, hl( 1, n+1 ),
491  \$ nl )
492  CALL dlaqr0( wantt, wantz, nl, ilo, kbot, hl, nl, wr,
493  \$ wi, ilo, ihi, z, ldz, workl, nl, info )
494  IF( wantt .OR. info.NE.0 )
495  \$ CALL dlacpy( 'A', n, n, hl, nl, h, ldh )
496  END IF
497  END IF
498  END IF
499 *
500 * ==== Clear out the trash, if necessary. ====
501 *
502  IF( ( wantt .OR. info.NE.0 ) .AND. n.GT.2 )
503  \$ CALL dlaset( 'L', n-2, n-2, zero, zero, h( 3, 1 ), ldh )
504 *
505 * ==== Ensure reported workspace size is backward-compatible with
506 * . previous LAPACK versions. ====
507 *
508  work( 1 ) = max( dble( max( 1, n ) ), work( 1 ) )
509  END IF
510 *
511 * ==== End of DHSEQR ====
512 *
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
Definition: dlacpy.f:103
subroutine dlaset(UPLO, M, N, ALPHA, BETA, A, LDA)
DLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition: dlaset.f:110
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
ILAENV
Definition: ilaenv.f:162
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
subroutine dlaqr0(WANTT, WANTZ, N, ILO, IHI, H, LDH, WR, WI, ILOZ, IHIZ, Z, LDZ, WORK, LWORK, INFO)
DLAQR0 computes the eigenvalues of a Hessenberg matrix, and optionally the matrices from the Schur de...
Definition: dlaqr0.f:256
subroutine dlahqr(WANTT, WANTZ, N, ILO, IHI, H, LDH, WR, WI, ILOZ, IHIZ, Z, LDZ, INFO)
DLAHQR computes the eigenvalues and Schur factorization of an upper Hessenberg matrix,...
Definition: dlahqr.f:207
Here is the call graph for this function:
Here is the caller graph for this function: