LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ zlaqz0()

recursive subroutine zlaqz0 ( character, intent(in)  WANTS,
character, intent(in)  WANTQ,
character, intent(in)  WANTZ,
integer, intent(in)  N,
integer, intent(in)  ILO,
integer, intent(in)  IHI,
complex*16, dimension( lda, * ), intent(inout)  A,
integer, intent(in)  LDA,
complex*16, dimension( ldb, * ), intent(inout)  B,
integer, intent(in)  LDB,
complex*16, dimension( * ), intent(inout)  ALPHA,
complex*16, dimension( * ), intent(inout)  BETA,
complex*16, dimension( ldq, * ), intent(inout)  Q,
integer, intent(in)  LDQ,
complex*16, dimension( ldz, * ), intent(inout)  Z,
integer, intent(in)  LDZ,
complex*16, dimension( * ), intent(inout)  WORK,
integer, intent(in)  LWORK,
double precision, dimension( * ), intent(out)  RWORK,
integer, intent(in)  REC,
integer, intent(out)  INFO 
)

ZLAQZ0

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

Purpose:
 ZLAQZ0 computes the eigenvalues of a real matrix pair (H,T),
 where H is an upper Hessenberg matrix and T is upper triangular,
 using the double-shift QZ method.
 Matrix pairs of this type are produced by the reduction to
 generalized upper Hessenberg form of a real matrix pair (A,B):

    A = Q1*H*Z1**H,  B = Q1*T*Z1**H,

 as computed by ZGGHRD.

 If JOB='S', then the Hessenberg-triangular pair (H,T) is
 also reduced to generalized Schur form,

    H = Q*S*Z**H,  T = Q*P*Z**H,

 where Q and Z are unitary matrices, P and S are an upper triangular
 matrices.

 Optionally, the unitary matrix Q from the generalized Schur
 factorization may be postmultiplied into an input matrix Q1, and the
 unitary matrix Z may be postmultiplied into an input matrix Z1.
 If Q1 and Z1 are the unitary matrices from ZGGHRD that reduced
 the matrix pair (A,B) to generalized upper Hessenberg form, then the
 output matrices Q1*Q and Z1*Z are the unitary factors from the
 generalized Schur factorization of (A,B):

    A = (Q1*Q)*S*(Z1*Z)**H,  B = (Q1*Q)*P*(Z1*Z)**H.

 To avoid overflow, eigenvalues of the matrix pair (H,T) (equivalently,
 of (A,B)) are computed as a pair of values (alpha,beta), where alpha is
 complex and beta real.
 If beta is nonzero, lambda = alpha / beta is an eigenvalue of the
 generalized nonsymmetric eigenvalue problem (GNEP)
    A*x = lambda*B*x
 and if alpha is nonzero, mu = beta / alpha is an eigenvalue of the
 alternate form of the GNEP
    mu*A*y = B*y.
 Eigenvalues can be read directly from the generalized Schur
 form:
   alpha = S(i,i), beta = P(i,i).

 Ref: C.B. Moler & G.W. Stewart, "An Algorithm for Generalized Matrix
      Eigenvalue Problems", SIAM J. Numer. Anal., 10(1973),
      pp. 241--256.

 Ref: B. Kagstrom, D. Kressner, "Multishift Variants of the QZ
      Algorithm with Aggressive Early Deflation", SIAM J. Numer.
      Anal., 29(2006), pp. 199--227.

 Ref: T. Steel, D. Camps, K. Meerbergen, R. Vandebril "A multishift,
      multipole rational QZ method with agressive early deflation"
Parameters
[in]WANTS
          WANTS is CHARACTER*1
          = 'E': Compute eigenvalues only;
          = 'S': Compute eigenvalues and the Schur form.
[in]WANTQ
          WANTQ is CHARACTER*1
          = 'N': Left Schur vectors (Q) are not computed;
          = 'I': Q is initialized to the unit matrix and the matrix Q
                 of left Schur vectors of (A,B) is returned;
          = 'V': Q must contain an unitary matrix Q1 on entry and
                 the product Q1*Q is returned.
[in]WANTZ
          WANTZ is CHARACTER*1
          = 'N': Right Schur vectors (Z) are not computed;
          = 'I': Z is initialized to the unit matrix and the matrix Z
                 of right Schur vectors of (A,B) is returned;
          = 'V': Z must contain an unitary matrix Z1 on entry and
                 the product Z1*Z is returned.
[in]N
          N is INTEGER
          The order of the matrices A, B, Q, and Z.  N >= 0.
[in]ILO
          ILO is INTEGER
[in]IHI
          IHI is INTEGER
          ILO and IHI mark the rows and columns of A which are in
          Hessenberg form.  It is assumed that A is already upper
          triangular in rows and columns 1:ILO-1 and IHI+1:N.
          If N > 0, 1 <= ILO <= IHI <= N; if N = 0, ILO=1 and IHI=0.
[in,out]A
          A is COMPLEX*16 array, dimension (LDA, N)
          On entry, the N-by-N upper Hessenberg matrix A.
          On exit, if JOB = 'S', A contains the upper triangular
          matrix S from the generalized Schur factorization.
          If JOB = 'E', the diagonal blocks of A match those of S, but
          the rest of A is unspecified.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max( 1, N ).
[in,out]B
          B is COMPLEX*16 array, dimension (LDB, N)
          On entry, the N-by-N upper triangular matrix B.
          On exit, if JOB = 'S', B contains the upper triangular
          matrix P from the generalized Schur factorization;
          If JOB = 'E', the diagonal blocks of B match those of P, but
          the rest of B is unspecified.
[in]LDB
          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max( 1, N ).
[out]ALPHA
          ALPHA is COMPLEX*16 array, dimension (N)
          Each scalar alpha defining an eigenvalue
          of GNEP.
[out]BETA
          BETA is COMPLEX*16 array, dimension (N)
          The scalars beta that define the eigenvalues of GNEP.
          Together, the quantities alpha = ALPHA(j) and
          beta = BETA(j) represent the j-th eigenvalue of the matrix
          pair (A,B), in one of the forms lambda = alpha/beta or
          mu = beta/alpha.  Since either lambda or mu may overflow,
          they should not, in general, be computed.
[in,out]Q
          Q is COMPLEX*16 array, dimension (LDQ, N)
          On entry, if COMPQ = 'V', the unitary matrix Q1 used in
          the reduction of (A,B) to generalized Hessenberg form.
          On exit, if COMPQ = 'I', the unitary matrix of left Schur
          vectors of (A,B), and if COMPQ = 'V', the unitary matrix
          of left Schur vectors of (A,B).
          Not referenced if COMPQ = 'N'.
[in]LDQ
          LDQ is INTEGER
          The leading dimension of the array Q.  LDQ >= 1.
          If COMPQ='V' or 'I', then LDQ >= N.
[in,out]Z
          Z is COMPLEX*16 array, dimension (LDZ, N)
          On entry, if COMPZ = 'V', the unitary matrix Z1 used in
          the reduction of (A,B) to generalized Hessenberg form.
          On exit, if COMPZ = 'I', the unitary matrix of
          right Schur vectors of (H,T), and if COMPZ = 'V', the
          unitary matrix of right Schur vectors of (A,B).
          Not referenced if COMPZ = 'N'.
[in]LDZ
          LDZ is INTEGER
          The leading dimension of the array Z.  LDZ >= 1.
          If COMPZ='V' or 'I', then LDZ >= N.
[out]WORK
          WORK is COMPLEX*16 array, dimension (MAX(1,LWORK))
          On exit, if INFO >= 0, WORK(1) returns the optimal LWORK.
[out]RWORK
          RWORK is DOUBLE PRECISION array, dimension (N)
[in]LWORK
          LWORK is INTEGER
          The dimension of the array WORK.  LWORK >= max(1,N).

          If LWORK = -1, then a workspace query is assumed; the routine
          only calculates the optimal size of the WORK array, returns
          this value as the first entry of the WORK array, and no error
          message related to LWORK is issued by XERBLA.
[in]REC
          REC is INTEGER
             REC indicates the current recursion level. Should be set
             to 0 on first call.
[out]INFO
          INFO is INTEGER
          = 0: successful exit
          < 0: if INFO = -i, the i-th argument had an illegal value
          = 1,...,N: the QZ iteration did not converge.  (A,B) is not
                     in Schur form, but ALPHA(i) and
                     BETA(i), i=INFO+1,...,N should be correct.
Author
Thijs Steel, KU Leuven
Date
May 2020

Definition at line 280 of file zlaqz0.f.

284  IMPLICIT NONE
285 
286 * Arguments
287  CHARACTER, INTENT( IN ) :: WANTS, WANTQ, WANTZ
288  INTEGER, INTENT( IN ) :: N, ILO, IHI, LDA, LDB, LDQ, LDZ, LWORK,
289  $ REC
290  INTEGER, INTENT( OUT ) :: INFO
291  COMPLEX*16, INTENT( INOUT ) :: A( LDA, * ), B( LDB, * ), Q( LDQ,
292  $ * ), Z( LDZ, * ), ALPHA( * ), BETA( * ), WORK( * )
293  DOUBLE PRECISION, INTENT( OUT ) :: RWORK( * )
294 
295 * Parameters
296  COMPLEX*16 CZERO, CONE
297  parameter( czero = ( 0.0d+0, 0.0d+0 ), cone = ( 1.0d+0,
298  $ 0.0d+0 ) )
299  DOUBLE PRECISION :: ZERO, ONE, HALF
300  parameter( zero = 0.0d0, one = 1.0d0, half = 0.5d0 )
301 
302 * Local scalars
303  DOUBLE PRECISION :: SMLNUM, ULP, SAFMIN, SAFMAX, C1, TEMPR
304  COMPLEX*16 :: ESHIFT, S1, TEMP
305  INTEGER :: ISTART, ISTOP, IITER, MAXIT, ISTART2, K, LD, NSHIFTS,
306  $ NBLOCK, NW, NMIN, NIBBLE, N_UNDEFLATED, N_DEFLATED,
307  $ NS, SWEEP_INFO, SHIFTPOS, LWORKREQ, K2, ISTARTM,
308  $ ISTOPM, IWANTS, IWANTQ, IWANTZ, NORM_INFO, AED_INFO,
309  $ NWR, NBR, NSR, ITEMP1, ITEMP2, RCOST
310  LOGICAL :: ILSCHUR, ILQ, ILZ
311  CHARACTER :: JBCMPZ*3
312 
313 * External Functions
314  EXTERNAL :: xerbla, zhgeqz, zlaqz2, zlaqz3, zlaset, dlabad,
315  $ zlartg, zrot
316  DOUBLE PRECISION, EXTERNAL :: DLAMCH
317  LOGICAL, EXTERNAL :: LSAME
318  INTEGER, EXTERNAL :: ILAENV
319 
320 *
321 * Decode wantS,wantQ,wantZ
322 *
323  IF( lsame( wants, 'E' ) ) THEN
324  ilschur = .false.
325  iwants = 1
326  ELSE IF( lsame( wants, 'S' ) ) THEN
327  ilschur = .true.
328  iwants = 2
329  ELSE
330  iwants = 0
331  END IF
332 
333  IF( lsame( wantq, 'N' ) ) THEN
334  ilq = .false.
335  iwantq = 1
336  ELSE IF( lsame( wantq, 'V' ) ) THEN
337  ilq = .true.
338  iwantq = 2
339  ELSE IF( lsame( wantq, 'I' ) ) THEN
340  ilq = .true.
341  iwantq = 3
342  ELSE
343  iwantq = 0
344  END IF
345 
346  IF( lsame( wantz, 'N' ) ) THEN
347  ilz = .false.
348  iwantz = 1
349  ELSE IF( lsame( wantz, 'V' ) ) THEN
350  ilz = .true.
351  iwantz = 2
352  ELSE IF( lsame( wantz, 'I' ) ) THEN
353  ilz = .true.
354  iwantz = 3
355  ELSE
356  iwantz = 0
357  END IF
358 *
359 * Check Argument Values
360 *
361  info = 0
362  IF( iwants.EQ.0 ) THEN
363  info = -1
364  ELSE IF( iwantq.EQ.0 ) THEN
365  info = -2
366  ELSE IF( iwantz.EQ.0 ) THEN
367  info = -3
368  ELSE IF( n.LT.0 ) THEN
369  info = -4
370  ELSE IF( ilo.LT.1 ) THEN
371  info = -5
372  ELSE IF( ihi.GT.n .OR. ihi.LT.ilo-1 ) THEN
373  info = -6
374  ELSE IF( lda.LT.n ) THEN
375  info = -8
376  ELSE IF( ldb.LT.n ) THEN
377  info = -10
378  ELSE IF( ldq.LT.1 .OR. ( ilq .AND. ldq.LT.n ) ) THEN
379  info = -15
380  ELSE IF( ldz.LT.1 .OR. ( ilz .AND. ldz.LT.n ) ) THEN
381  info = -17
382  END IF
383  IF( info.NE.0 ) THEN
384  CALL xerbla( 'ZLAQZ0', -info )
385  RETURN
386  END IF
387 
388 *
389 * Quick return if possible
390 *
391  IF( n.LE.0 ) THEN
392  work( 1 ) = dble( 1 )
393  RETURN
394  END IF
395 
396 *
397 * Get the parameters
398 *
399  jbcmpz( 1:1 ) = wants
400  jbcmpz( 2:2 ) = wantq
401  jbcmpz( 3:3 ) = wantz
402 
403  nmin = ilaenv( 12, 'ZLAQZ0', jbcmpz, n, ilo, ihi, lwork )
404 
405  nwr = ilaenv( 13, 'ZLAQZ0', jbcmpz, n, ilo, ihi, lwork )
406  nwr = max( 2, nwr )
407  nwr = min( ihi-ilo+1, ( n-1 ) / 3, nwr )
408 
409  nibble = ilaenv( 14, 'ZLAQZ0', jbcmpz, n, ilo, ihi, lwork )
410 
411  nsr = ilaenv( 15, 'ZLAQZ0', jbcmpz, n, ilo, ihi, lwork )
412  nsr = min( nsr, ( n+6 ) / 9, ihi-ilo )
413  nsr = max( 2, nsr-mod( nsr, 2 ) )
414 
415  rcost = ilaenv( 17, 'ZLAQZ0', jbcmpz, n, ilo, ihi, lwork )
416  itemp1 = int( nsr/sqrt( 1+2*nsr/( dble( rcost )/100*n ) ) )
417  itemp1 = ( ( itemp1-1 )/4 )*4+4
418  nbr = nsr+itemp1
419 
420  IF( n .LT. nmin .OR. rec .GE. 2 ) THEN
421  CALL zhgeqz( wants, wantq, wantz, n, ilo, ihi, a, lda, b, ldb,
422  $ alpha, beta, q, ldq, z, ldz, work, lwork, rwork,
423  $ info )
424  RETURN
425  END IF
426 
427 *
428 * Find out required workspace
429 *
430 
431 * Workspace query to ZLAQZ2
432  nw = max( nwr, nmin )
433  CALL zlaqz2( ilschur, ilq, ilz, n, ilo, ihi, nw, a, lda, b, ldb,
434  $ q, ldq, z, ldz, n_undeflated, n_deflated, alpha,
435  $ beta, work, nw, work, nw, work, -1, rwork, rec,
436  $ aed_info )
437  itemp1 = int( work( 1 ) )
438 * Workspace query to ZLAQZ3
439  CALL zlaqz3( ilschur, ilq, ilz, n, ilo, ihi, nsr, nbr, alpha,
440  $ beta, a, lda, b, ldb, q, ldq, z, ldz, work, nbr,
441  $ work, nbr, work, -1, sweep_info )
442  itemp2 = int( work( 1 ) )
443 
444  lworkreq = max( itemp1+2*nw**2, itemp2+2*nbr**2 )
445  IF ( lwork .EQ.-1 ) THEN
446  work( 1 ) = dble( lworkreq )
447  RETURN
448  ELSE IF ( lwork .LT. lworkreq ) THEN
449  info = -19
450  END IF
451  IF( info.NE.0 ) THEN
452  CALL xerbla( 'ZLAQZ0', info )
453  RETURN
454  END IF
455 *
456 * Initialize Q and Z
457 *
458  IF( iwantq.EQ.3 ) CALL zlaset( 'FULL', n, n, czero, cone, q,
459  $ ldq )
460  IF( iwantz.EQ.3 ) CALL zlaset( 'FULL', n, n, czero, cone, z,
461  $ ldz )
462 
463 * Get machine constants
464  safmin = dlamch( 'SAFE MINIMUM' )
465  safmax = one/safmin
466  CALL dlabad( safmin, safmax )
467  ulp = dlamch( 'PRECISION' )
468  smlnum = safmin*( dble( n )/ulp )
469 
470  istart = ilo
471  istop = ihi
472  maxit = 30*( ihi-ilo+1 )
473  ld = 0
474 
475  DO iiter = 1, maxit
476  IF( iiter .GE. maxit ) THEN
477  info = istop+1
478  GOTO 80
479  END IF
480  IF ( istart+1 .GE. istop ) THEN
481  istop = istart
482  EXIT
483  END IF
484 
485 * Check deflations at the end
486  IF ( abs( a( istop, istop-1 ) ) .LE. max( smlnum,
487  $ ulp*( abs( a( istop, istop ) )+abs( a( istop-1,
488  $ istop-1 ) ) ) ) ) THEN
489  a( istop, istop-1 ) = czero
490  istop = istop-1
491  ld = 0
492  eshift = czero
493  END IF
494 * Check deflations at the start
495  IF ( abs( a( istart+1, istart ) ) .LE. max( smlnum,
496  $ ulp*( abs( a( istart, istart ) )+abs( a( istart+1,
497  $ istart+1 ) ) ) ) ) THEN
498  a( istart+1, istart ) = czero
499  istart = istart+1
500  ld = 0
501  eshift = czero
502  END IF
503 
504  IF ( istart+1 .GE. istop ) THEN
505  EXIT
506  END IF
507 
508 * Check interior deflations
509  istart2 = istart
510  DO k = istop, istart+1, -1
511  IF ( abs( a( k, k-1 ) ) .LE. max( smlnum, ulp*( abs( a( k,
512  $ k ) )+abs( a( k-1, k-1 ) ) ) ) ) THEN
513  a( k, k-1 ) = czero
514  istart2 = k
515  EXIT
516  END IF
517  END DO
518 
519 * Get range to apply rotations to
520  IF ( ilschur ) THEN
521  istartm = 1
522  istopm = n
523  ELSE
524  istartm = istart2
525  istopm = istop
526  END IF
527 
528 * Check infinite eigenvalues, this is done without blocking so might
529 * slow down the method when many infinite eigenvalues are present
530  k = istop
531  DO WHILE ( k.GE.istart2 )
532  tempr = zero
533  IF( k .LT. istop ) THEN
534  tempr = tempr+abs( b( k, k+1 ) )
535  END IF
536  IF( k .GT. istart2 ) THEN
537  tempr = tempr+abs( b( k-1, k ) )
538  END IF
539 
540  IF( abs( b( k, k ) ) .LT. max( smlnum, ulp*tempr ) ) THEN
541 * A diagonal element of B is negligable, move it
542 * to the top and deflate it
543 
544  DO k2 = k, istart2+1, -1
545  CALL zlartg( b( k2-1, k2 ), b( k2-1, k2-1 ), c1, s1,
546  $ temp )
547  b( k2-1, k2 ) = temp
548  b( k2-1, k2-1 ) = czero
549 
550  CALL zrot( k2-2-istartm+1, b( istartm, k2 ), 1,
551  $ b( istartm, k2-1 ), 1, c1, s1 )
552  CALL zrot( min( k2+1, istop )-istartm+1, a( istartm,
553  $ k2 ), 1, a( istartm, k2-1 ), 1, c1, s1 )
554  IF ( ilz ) THEN
555  CALL zrot( n, z( 1, k2 ), 1, z( 1, k2-1 ), 1, c1,
556  $ s1 )
557  END IF
558 
559  IF( k2.LT.istop ) THEN
560  CALL zlartg( a( k2, k2-1 ), a( k2+1, k2-1 ), c1,
561  $ s1, temp )
562  a( k2, k2-1 ) = temp
563  a( k2+1, k2-1 ) = czero
564 
565  CALL zrot( istopm-k2+1, a( k2, k2 ), lda, a( k2+1,
566  $ k2 ), lda, c1, s1 )
567  CALL zrot( istopm-k2+1, b( k2, k2 ), ldb, b( k2+1,
568  $ k2 ), ldb, c1, s1 )
569  IF( ilq ) THEN
570  CALL zrot( n, q( 1, k2 ), 1, q( 1, k2+1 ), 1,
571  $ c1, dconjg( s1 ) )
572  END IF
573  END IF
574 
575  END DO
576 
577  IF( istart2.LT.istop )THEN
578  CALL zlartg( a( istart2, istart2 ), a( istart2+1,
579  $ istart2 ), c1, s1, temp )
580  a( istart2, istart2 ) = temp
581  a( istart2+1, istart2 ) = czero
582 
583  CALL zrot( istopm-( istart2+1 )+1, a( istart2,
584  $ istart2+1 ), lda, a( istart2+1,
585  $ istart2+1 ), lda, c1, s1 )
586  CALL zrot( istopm-( istart2+1 )+1, b( istart2,
587  $ istart2+1 ), ldb, b( istart2+1,
588  $ istart2+1 ), ldb, c1, s1 )
589  IF( ilq ) THEN
590  CALL zrot( n, q( 1, istart2 ), 1, q( 1,
591  $ istart2+1 ), 1, c1, dconjg( s1 ) )
592  END IF
593  END IF
594 
595  istart2 = istart2+1
596 
597  END IF
598  k = k-1
599  END DO
600 
601 * istart2 now points to the top of the bottom right
602 * unreduced Hessenberg block
603  IF ( istart2 .GE. istop ) THEN
604  istop = istart2-1
605  ld = 0
606  eshift = czero
607  cycle
608  END IF
609 
610  nw = nwr
611  nshifts = nsr
612  nblock = nbr
613 
614  IF ( istop-istart2+1 .LT. nmin ) THEN
615 * Setting nw to the size of the subblock will make AED deflate
616 * all the eigenvalues. This is slightly more efficient than just
617 * using qz_small because the off diagonal part gets updated via BLAS.
618  IF ( istop-istart+1 .LT. nmin ) THEN
619  nw = istop-istart+1
620  istart2 = istart
621  ELSE
622  nw = istop-istart2+1
623  END IF
624  END IF
625 
626 *
627 * Time for AED
628 *
629  CALL zlaqz2( ilschur, ilq, ilz, n, istart2, istop, nw, a, lda,
630  $ b, ldb, q, ldq, z, ldz, n_undeflated, n_deflated,
631  $ alpha, beta, work, nw, work( nw**2+1 ), nw,
632  $ work( 2*nw**2+1 ), lwork-2*nw**2, rwork, rec,
633  $ aed_info )
634 
635  IF ( n_deflated > 0 ) THEN
636  istop = istop-n_deflated
637  ld = 0
638  eshift = czero
639  END IF
640 
641  IF ( 100*n_deflated > nibble*( n_deflated+n_undeflated ) .OR.
642  $ istop-istart2+1 .LT. nmin ) THEN
643 * AED has uncovered many eigenvalues. Skip a QZ sweep and run
644 * AED again.
645  cycle
646  END IF
647 
648  ld = ld+1
649 
650  ns = min( nshifts, istop-istart2 )
651  ns = min( ns, n_undeflated )
652  shiftpos = istop-n_deflated-n_undeflated+1
653 
654  IF ( mod( ld, 6 ) .EQ. 0 ) THEN
655 *
656 * Exceptional shift. Chosen for no particularly good reason.
657 *
658  IF( ( dble( maxit )*safmin )*abs( a( istop,
659  $ istop-1 ) ).LT.abs( a( istop-1, istop-1 ) ) ) THEN
660  eshift = a( istop, istop-1 )/b( istop-1, istop-1 )
661  ELSE
662  eshift = eshift+cone/( safmin*dble( maxit ) )
663  END IF
664  alpha( shiftpos ) = cone
665  beta( shiftpos ) = eshift
666  ns = 1
667  END IF
668 
669 *
670 * Time for a QZ sweep
671 *
672  CALL zlaqz3( ilschur, ilq, ilz, n, istart2, istop, ns, nblock,
673  $ alpha( shiftpos ), beta( shiftpos ), a, lda, b,
674  $ ldb, q, ldq, z, ldz, work, nblock, work( nblock**
675  $ 2+1 ), nblock, work( 2*nblock**2+1 ),
676  $ lwork-2*nblock**2, sweep_info )
677 
678  END DO
679 
680 *
681 * Call ZHGEQZ to normalize the eigenvalue blocks and set the eigenvalues
682 * If all the eigenvalues have been found, ZHGEQZ will not do any iterations
683 * and only normalize the blocks. In case of a rare convergence failure,
684 * the single shift might perform better.
685 *
686  80 CALL zhgeqz( wants, wantq, wantz, n, ilo, ihi, a, lda, b, ldb,
687  $ alpha, beta, q, ldq, z, ldz, work, lwork, rwork,
688  $ norm_info )
689 
690  info = norm_info
691 
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:69
subroutine zlartg(f, g, c, s, r)
ZLARTG generates a plane rotation with real cosine and complex sine.
Definition: zlartg.f90:118
subroutine dlabad(SMALL, LARGE)
DLABAD
Definition: dlabad.f:74
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 zlaqz3(ILSCHUR, ILQ, ILZ, N, ILO, IHI, NSHIFTS, NBLOCK_DESIRED, ALPHA, BETA, A, LDA, B, LDB, Q, LDQ, Z, LDZ, QC, LDQC, ZC, LDZC, WORK, LWORK, INFO)
ZLAQZ3
Definition: zlaqz3.f:208
subroutine zhgeqz(JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT, ALPHA, BETA, Q, LDQ, Z, LDZ, WORK, LWORK, RWORK, INFO)
ZHGEQZ
Definition: zhgeqz.f:284
recursive subroutine zlaqz2(ILSCHUR, ILQ, ILZ, N, ILO, IHI, NW, A, LDA, B, LDB, Q, LDQ, Z, LDZ, NS, ND, ALPHA, BETA, QC, LDQC, ZC, LDZC, WORK, LWORK, RWORK, REC, INFO)
ZLAQZ2
Definition: zlaqz2.f:234
subroutine zrot(N, CX, INCX, CY, INCY, C, S)
ZROT applies a plane rotation with real cosine and complex sine to a pair of complex vectors.
Definition: zrot.f:103
subroutine zlaset(UPLO, M, N, ALPHA, BETA, A, LDA)
ZLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition: zlaset.f:106
Here is the call graph for this function:
Here is the caller graph for this function: