240 SUBROUTINE claqr0( WANTT, WANTZ, N, ILO, IHI, H, LDH, W, ILOZ,
241 $ ihiz, z, ldz, work, lwork, info )
249 INTEGER IHI, IHIZ, ILO, ILOZ, INFO, LDH, LDZ, LWORK, N
253 COMPLEX H( ldh, * ), W( * ), WORK( * ), Z( ldz, * )
263 parameter ( ntiny = 11 )
269 parameter ( kexnw = 5 )
275 parameter ( kexsh = 6 )
280 parameter ( wilk1 = 0.75e0 )
282 parameter ( zero = ( 0.0e0, 0.0e0 ),
283 $ one = ( 1.0e0, 0.0e0 ) )
285 parameter ( two = 2.0e0 )
288 COMPLEX AA, BB, CC, CDUM, DD, DET, RTDISC, SWAP, TR2
290 INTEGER I, INF, IT, ITMAX, K, KACC22, KBOT, KDU, KS,
291 $ kt, ktop, ku, kv, kwh, kwtop, kwv, ld, ls,
292 $ lwkopt, ndec, ndfl, nh, nho, nibble, nmin, ns,
293 $ nsmax, nsr, nve, nw, nwmax, nwr, nwupbd
308 INTRINSIC abs, aimag, cmplx, int, max, min, mod,
REAL,
315 cabs1( cdum ) = abs(
REAL( CDUM ) ) + abs( AIMAG( cdum ) )
327 IF( n.LE.ntiny )
THEN
333 $
CALL clahqr( wantt, wantz, n, ilo, ihi, h, ldh, w, iloz,
334 $ ihiz, z, ldz, info )
363 nwr = ilaenv( 13,
'CLAQR0', jbcmpz, n, ilo, ihi, lwork )
365 nwr = min( ihi-ilo+1, ( n-1 ) / 3, nwr )
372 nsr = ilaenv( 15,
'CLAQR0', jbcmpz, n, ilo, ihi, lwork )
373 nsr = min( nsr, ( n+6 ) / 9, ihi-ilo )
374 nsr = max( 2, nsr-mod( nsr, 2 ) )
380 CALL claqr3( wantt, wantz, n, ilo, ihi, nwr+1, h, ldh, iloz,
381 $ ihiz, z, ldz, ls, ld, w, h, ldh, n, h, ldh, n, h,
386 lwkopt = max( 3*nsr / 2, int( work( 1 ) ) )
390 IF( lwork.EQ.-1 )
THEN
391 work( 1 ) = cmplx( lwkopt, 0 )
397 nmin = ilaenv( 12,
'CLAQR0', jbcmpz, n, ilo, ihi, lwork )
398 nmin = max( ntiny, nmin )
402 nibble = ilaenv( 14,
'CLAQR0', jbcmpz, n, ilo, ihi, lwork )
403 nibble = max( 0, nibble )
408 kacc22 = ilaenv( 16,
'CLAQR0', jbcmpz, n, ilo, ihi, lwork )
409 kacc22 = max( 0, kacc22 )
410 kacc22 = min( 2, kacc22 )
415 nwmax = min( ( n-1 ) / 3, lwork / 2 )
421 nsmax = min( ( n+6 ) / 9, 2*lwork / 3 )
422 nsmax = nsmax - mod( nsmax, 2 )
430 itmax = max( 30, 2*kexsh )*max( 10, ( ihi-ilo+1 ) )
447 DO 10 k = kbot, ilo + 1, -1
448 IF( h( k, k-1 ).EQ.zero )
472 nwupbd = min( nh, nwmax )
473 IF( ndfl.LT.kexnw )
THEN
474 nw = min( nwupbd, nwr )
476 nw = min( nwupbd, 2*nw )
478 IF( nw.LT.nwmax )
THEN
479 IF( nw.GE.nh-1 )
THEN
482 kwtop = kbot - nw + 1
483 IF( cabs1( h( kwtop, kwtop-1 ) ).GT.
484 $ cabs1( h( kwtop-1, kwtop-2 ) ) )nw = nw + 1
487 IF( ndfl.LT.kexnw )
THEN
489 ELSE IF( ndec.GE.0 .OR. nw.GE.nwupbd )
THEN
509 nho = ( n-nw-1 ) - kt + 1
511 nve = ( n-nw ) - kwv + 1
515 CALL claqr3( wantt, wantz, n, ktop, kbot, nw, h, ldh, iloz,
516 $ ihiz, z, ldz, ls, ld, w, h( kv, 1 ), ldh, nho,
517 $ h( kv, kt ), ldh, nve, h( kwv, 1 ), ldh, work,
534 IF( ( ld.EQ.0 ) .OR. ( ( 100*ld.LE.nw*nibble ) .AND. ( kbot-
535 $ ktop+1.GT.min( nmin, nwmax ) ) ) )
THEN
541 ns = min( nsmax, nsr, max( 2, kbot-ktop ) )
542 ns = ns - mod( ns, 2 )
551 IF( mod( ndfl, kexsh ).EQ.0 )
THEN
553 DO 30 i = kbot, ks + 1, -2
554 w( i ) = h( i, i ) + wilk1*cabs1( h( i, i-1 ) )
565 IF( kbot-ks+1.LE.ns / 2 )
THEN
568 CALL clacpy(
'A', ns, ns, h( ks, ks ), ldh,
570 IF( ns.GT.nmin )
THEN
571 CALL claqr4( .false., .false., ns, 1, ns,
572 $ h( kt, 1 ), ldh, w( ks ), 1, 1,
573 $ zdum, 1, work, lwork, inf )
575 CALL clahqr( .false., .false., ns, 1, ns,
576 $ h( kt, 1 ), ldh, w( ks ), 1, 1,
588 IF( ks.GE.kbot )
THEN
589 s = cabs1( h( kbot-1, kbot-1 ) ) +
590 $ cabs1( h( kbot, kbot-1 ) ) +
591 $ cabs1( h( kbot-1, kbot ) ) +
592 $ cabs1( h( kbot, kbot ) )
593 aa = h( kbot-1, kbot-1 ) / s
594 cc = h( kbot, kbot-1 ) / s
595 bb = h( kbot-1, kbot ) / s
596 dd = h( kbot, kbot ) / s
597 tr2 = ( aa+dd ) / two
598 det = ( aa-tr2 )*( dd-tr2 ) - bb*cc
599 rtdisc = sqrt( -det )
600 w( kbot-1 ) = ( tr2+rtdisc )*s
601 w( kbot ) = ( tr2-rtdisc )*s
607 IF( kbot-ks+1.GT.ns )
THEN
612 DO 50 k = kbot, ks + 1, -1
617 IF( cabs1( w( i ) ).LT.cabs1( w( i+1 ) ) )
633 IF( kbot-ks+1.EQ.2 )
THEN
634 IF( cabs1( w( kbot )-h( kbot, kbot ) ).LT.
635 $ cabs1( w( kbot-1 )-h( kbot, kbot ) ) )
THEN
636 w( kbot-1 ) = w( kbot )
638 w( kbot ) = w( kbot-1 )
647 ns = min( ns, kbot-ks+1 )
648 ns = ns - mod( ns, 2 )
665 nho = ( n-kdu+1-4 ) - ( kdu+1 ) + 1
667 nve = n - kdu - kwv + 1
671 CALL claqr5( wantt, wantz, kacc22, n, ktop, kbot, ns,
672 $ w( ks ), h, ldh, iloz, ihiz, z, ldz, work,
673 $ 3, h( ku, 1 ), ldh, nve, h( kwv, 1 ), ldh,
674 $ nho, h( ku, kwh ), ldh )
697 work( 1 ) = cmplx( lwkopt, 0 )
subroutine clahqr(WANTT, WANTZ, N, ILO, IHI, H, LDH, W, ILOZ, IHIZ, Z, LDZ, INFO)
CLAHQR computes the eigenvalues and Schur factorization of an upper Hessenberg matrix, using the double-shift/single-shift QR algorithm.
subroutine claqr5(WANTT, WANTZ, KACC22, N, KTOP, KBOT, NSHFTS, S, H, LDH, ILOZ, IHIZ, Z, LDZ, V, LDV, U, LDU, NV, WV, LDWV, NH, WH, LDWH)
CLAQR5 performs a single small-bulge multi-shift QR sweep.
subroutine claqr4(WANTT, WANTZ, N, ILO, IHI, H, LDH, W, ILOZ, IHIZ, Z, LDZ, WORK, LWORK, INFO)
CLAQR4 computes the eigenvalues of a Hessenberg matrix, and optionally the matrices from the Schur de...
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
subroutine claqr3(WANTT, WANTZ, N, KTOP, KBOT, NW, H, LDH, ILOZ, IHIZ, Z, LDZ, NS, ND, SH, V, LDV, NH, T, LDT, NV, WV, LDWV, WORK, LWORK)
CLAQR3 performs the unitary similarity transformation of a Hessenberg matrix to detect and deflate fu...
subroutine claqr0(WANTT, WANTZ, N, ILO, IHI, H, LDH, W, ILOZ, IHIZ, Z, LDZ, WORK, LWORK, INFO)
CLAQR0 computes the eigenvalues of a Hessenberg matrix, and optionally the matrices from the Schur de...