SYNOPSIS¶
- SUBROUTINE
ZLAQR5(
- WANTT, WANTZ, KACC22, N, KTOP, KBOT, NSHFTS, S, H, LDH, ILOZ, IHIZ, Z,
LDZ, V, LDV, U, LDU, NV, WV, LDWV, NH, WH, LDWH )
INTEGER IHIZ, ILOZ, KACC22, KBOT, KTOP, LDH, LDU, LDV, LDWH, LDWV,
LDZ, N, NH, NSHFTS, NV LOGICAL WANTT, WANTZ COMPLEX*16 H( LDH, * ), S( * ),
U( LDU, * ), V( LDV, * ), WH( LDWH, * ), WV( LDWV, * ), Z( LDZ, * )
COMPLEX*16 ZERO, ONE PARAMETER ( ZERO = ( 0.0d0, 0.0d0 ), ONE = ( 1.0d0,
0.0d0 ) ) DOUBLE PRECISION RZERO, RONE PARAMETER ( RZERO = 0.0d0, RONE =
1.0d0 ) COMPLEX*16 ALPHA, BETA, CDUM, REFSUM DOUBLE PRECISION H11, H12, H21,
H22, SAFMAX, SAFMIN, SCL, SMLNUM, TST1, TST2, ULP INTEGER I2, I4, INCOL, J,
J2, J4, JBOT, JCOL, JLEN, JROW, JTOP, K, K1, KDU, KMS, KNZ, KRCOL, KZS, M,
M22, MBOT, MEND, MSTART, MTOP, NBMPS, NDCOL, NS, NU LOGICAL ACCUM, BLK22,
BMP22 DOUBLE PRECISION DLAMCH EXTERNAL DLAMCH INTRINSIC ABS, DBLE, DCONJG,
DIMAG, MAX, MIN, MOD COMPLEX*16 VT( 3 ) EXTERNAL DLABAD, ZGEMM, ZLACPY,
ZLAQR1, ZLARFG, ZLASET, ZTRMM DOUBLE PRECISION CABS1 CABS1( CDUM ) = ABS(
DBLE( CDUM ) ) + ABS( DIMAG( CDUM ) ) IF( NSHFTS.LT.2 ) RETURN IF(
KTOP.GE.KBOT ) RETURN NS = NSHFTS - MOD( NSHFTS, 2 ) SAFMIN = DLAMCH( 'SAFE
MINIMUM' ) SAFMAX = RONE / SAFMIN CALL DLABAD( SAFMIN, SAFMAX ) ULP =
DLAMCH( 'PRECISION' ) SMLNUM = SAFMIN*( DBLE( N ) / ULP ) ACCUM = (
KACC22.EQ.1 ) .OR. ( KACC22.EQ.2 ) BLK22 = ( NS.GT.2 ) .AND. ( KACC22.EQ.2 )
IF( KTOP+2.LE.KBOT ) H( KTOP+2, KTOP ) = ZERO NBMPS = NS / 2 KDU = 6*NBMPS -
3 DO 210 INCOL = 3*( 1-NBMPS ) + KTOP - 1, KBOT - 2, 3*NBMPS - 2 NDCOL =
INCOL + KDU IF( ACCUM ) CALL ZLASET( 'ALL', KDU, KDU, ZERO, ONE, U, LDU ) DO
140 KRCOL = INCOL, MIN( INCOL+3*NBMPS-3, KBOT-2 ) MTOP = MAX( 1, ( ( KTOP-1
)-KRCOL+2 ) / 3+1 ) MBOT = MIN( NBMPS, ( KBOT-KRCOL ) / 3 ) M22 = MBOT + 1
BMP22 = ( MBOT.LT.NBMPS ) .AND. ( KRCOL+3*( M22-1 ) ).EQ. ( KBOT-2 ) DO 10 M
= MTOP, MBOT K = KRCOL + 3*( M-1 ) IF( K.EQ.KTOP-1 ) THEN CALL ZLAQR1( 3, H(
KTOP, KTOP ), LDH, S( 2*M-1 ), S( 2*M ), V( 1, M ) ) ALPHA = V( 1, M ) CALL
ZLARFG( 3, ALPHA, V( 2, M ), 1, V( 1, M ) ) ELSE BETA = H( K+1, K ) V( 2, M
) = H( K+2, K ) V( 3, M ) = H( K+3, K ) CALL ZLARFG( 3, BETA, V( 2, M ), 1,
V( 1, M ) ) IF( H( K+3, K ).NE.ZERO .OR. H( K+3, K+1 ).NE. ZERO .OR. H( K+3,
K+2 ).EQ.ZERO ) THEN H( K+1, K ) = BETA H( K+2, K ) = ZERO H( K+3, K ) =
ZERO ELSE CALL ZLAQR1( 3, H( K+1, K+1 ), LDH, S( 2*M-1 ), S( 2*M ), VT )
ALPHA = VT( 1 ) CALL ZLARFG( 3, ALPHA, VT( 2 ), 1, VT( 1 ) ) REFSUM =
DCONJG( VT( 1 ) )* ( H( K+1, K )+DCONJG( VT( 2 ) )* H( K+2, K ) ) IF( CABS1(
H( K+2, K )-REFSUM*VT( 2 ) )+ CABS1( REFSUM*VT( 3 ) ).GT.ULP* ( CABS1( H( K,
K ) )+CABS1( H( K+1, K+1 ) )+CABS1( H( K+2, K+2 ) ) ) ) THEN H( K+1, K ) =
BETA H( K+2, K ) = ZERO H( K+3, K ) = ZERO ELSE H( K+1, K ) = H( K+1, K ) -
REFSUM H( K+2, K ) = ZERO H( K+3, K ) = ZERO V( 1, M ) = VT( 1 ) V( 2, M ) =
VT( 2 ) V( 3, M ) = VT( 3 ) END IF END IF END IF 10 CONTINUE K = KRCOL + 3*(
M22-1 ) IF( BMP22 ) THEN IF( K.EQ.KTOP-1 ) THEN CALL ZLAQR1( 2, H( K+1, K+1
), LDH, S( 2*M22-1 ), S( 2*M22 ), V( 1, M22 ) ) BETA = V( 1, M22 ) CALL
ZLARFG( 2, BETA, V( 2, M22 ), 1, V( 1, M22 ) ) ELSE BETA = H( K+1, K ) V( 2,
M22 ) = H( K+2, K ) CALL ZLARFG( 2, BETA, V( 2, M22 ), 1, V( 1, M22 ) ) H(
K+1, K ) = BETA H( K+2, K ) = ZERO END IF END IF IF( ACCUM ) THEN JBOT =
MIN( NDCOL, KBOT ) ELSE IF( WANTT ) THEN JBOT = N ELSE JBOT = KBOT END IF DO
30 J = MAX( KTOP, KRCOL ), JBOT MEND = MIN( MBOT, ( J-KRCOL+2 ) / 3 ) DO 20
M = MTOP, MEND K = KRCOL + 3*( M-1 ) REFSUM = DCONJG( V( 1, M ) )* ( H( K+1,
J )+DCONJG( V( 2, M ) )* H( K+2, J )+DCONJG( V( 3, M ) )*H( K+3, J ) ) H(
K+1, J ) = H( K+1, J ) - REFSUM H( K+2, J ) = H( K+2, J ) - REFSUM*V( 2, M )
H( K+3, J ) = H( K+3, J ) - REFSUM*V( 3, M ) 20 CONTINUE 30 CONTINUE IF(
BMP22 ) THEN K = KRCOL + 3*( M22-1 ) DO 40 J = MAX( K+1, KTOP ), JBOT REFSUM
= DCONJG( V( 1, M22 ) )* ( H( K+1, J )+DCONJG( V( 2, M22 ) )* H( K+2, J ) )
H( K+1, J ) = H( K+1, J ) - REFSUM H( K+2, J ) = H( K+2, J ) - REFSUM*V( 2,
M22 ) 40 CONTINUE END IF IF( ACCUM ) THEN JTOP = MAX( KTOP, INCOL ) ELSE IF(
WANTT ) THEN JTOP = 1 ELSE JTOP = KTOP END IF DO 80 M = MTOP, MBOT IF( V( 1,
M ).NE.ZERO ) THEN K = KRCOL + 3*( M-1 ) DO 50 J = JTOP, MIN( KBOT, K+3 )
REFSUM = V( 1, M )*( H( J, K+1 )+V( 2, M )* H( J, K+2 )+V( 3, M )*H( J, K+3
) ) H( J, K+1 ) = H( J, K+1 ) - REFSUM H( J, K+2 ) = H( J, K+2 ) -
REFSUM*DCONJG( V( 2, M ) ) H( J, K+3 ) = H( J, K+3 ) - REFSUM*DCONJG( V( 3,
M ) ) 50 CONTINUE IF( ACCUM ) THEN KMS = K - INCOL DO 60 J = MAX( 1,
KTOP-INCOL ), KDU REFSUM = V( 1, M )*( U( J, KMS+1 )+V( 2, M )* U( J, KMS+2
)+V( 3, M )*U( J, KMS+3 ) ) U( J, KMS+1 ) = U( J, KMS+1 ) - REFSUM U( J,
KMS+2 ) = U( J, KMS+2 ) - REFSUM*DCONJG( V( 2, M ) ) U( J, KMS+3 ) = U( J,
KMS+3 ) - REFSUM*DCONJG( V( 3, M ) ) 60 CONTINUE ELSE IF( WANTZ ) THEN DO 70
J = ILOZ, IHIZ REFSUM = V( 1, M )*( Z( J, K+1 )+V( 2, M )* Z( J, K+2 )+V( 3,
M )*Z( J, K+3 ) ) Z( J, K+1 ) = Z( J, K+1 ) - REFSUM Z( J, K+2 ) = Z( J, K+2
) - REFSUM*DCONJG( V( 2, M ) ) Z( J, K+3 ) = Z( J, K+3 ) - REFSUM*DCONJG( V(
3, M ) ) 70 CONTINUE END IF END IF 80 CONTINUE K = KRCOL + 3*( M22-1 ) IF(
BMP22 .AND. ( V( 1, M22 ).NE.ZERO ) ) THEN DO 90 J = JTOP, MIN( KBOT, K+3 )
REFSUM = V( 1, M22 )*( H( J, K+1 )+V( 2, M22 )* H( J, K+2 ) ) H( J, K+1 ) =
H( J, K+1 ) - REFSUM H( J, K+2 ) = H( J, K+2 ) - REFSUM*DCONJG( V( 2, M22 )
) 90 CONTINUE IF( ACCUM ) THEN KMS = K - INCOL DO 100 J = MAX( 1, KTOP-INCOL
), KDU REFSUM = V( 1, M22 )*( U( J, KMS+1 )+V( 2, M22 )* U( J, KMS+2 ) ) U(
J, KMS+1 ) = U( J, KMS+1 ) - REFSUM U( J, KMS+2 ) = U( J, KMS+2 ) -
REFSUM*DCONJG( V( 2, M22 ) ) 100 CONTINUE ELSE IF( WANTZ ) THEN DO 110 J =
ILOZ, IHIZ REFSUM = V( 1, M22 )*( Z( J, K+1 )+V( 2, M22 )* Z( J, K+2 ) ) Z(
J, K+1 ) = Z( J, K+1 ) - REFSUM Z( J, K+2 ) = Z( J, K+2 ) - REFSUM*DCONJG(
V( 2, M22 ) ) 110 CONTINUE END IF END IF MSTART = MTOP IF( KRCOL+3*(
MSTART-1 ).LT.KTOP ) MSTART = MSTART + 1 MEND = MBOT IF( BMP22 ) MEND = MEND
+ 1 IF( KRCOL.EQ.KBOT-2 ) MEND = MEND + 1 DO 120 M = MSTART, MEND K = MIN(
KBOT-1, KRCOL+3*( M-1 ) ) IF( H( K+1, K ).NE.ZERO ) THEN TST1 = CABS1( H( K,
K ) ) + CABS1( H( K+1, K+1 ) ) IF( TST1.EQ.RZERO ) THEN IF( K.GE.KTOP+1 )
TST1 = TST1 + CABS1( H( K, K-1 ) ) IF( K.GE.KTOP+2 ) TST1 = TST1 + CABS1( H(
K, K-2 ) ) IF( K.GE.KTOP+3 ) TST1 = TST1 + CABS1( H( K, K-3 ) ) IF(
K.LE.KBOT-2 ) TST1 = TST1 + CABS1( H( K+2, K+1 ) ) IF( K.LE.KBOT-3 ) TST1 =
TST1 + CABS1( H( K+3, K+1 ) ) IF( K.LE.KBOT-4 ) TST1 = TST1 + CABS1( H( K+4,
K+1 ) ) END IF IF( CABS1( H( K+1, K ) ).LE.MAX( SMLNUM, ULP*TST1 ) ) THEN
H12 = MAX( CABS1( H( K+1, K ) ), CABS1( H( K, K+1 ) ) ) H21 = MIN( CABS1( H(
K+1, K ) ), CABS1( H( K, K+1 ) ) ) H11 = MAX( CABS1( H( K+1, K+1 ) ), CABS1(
H( K, K )-H( K+1, K+1 ) ) ) H22 = MIN( CABS1( H( K+1, K+1 ) ), CABS1( H( K,
K )-H( K+1, K+1 ) ) ) SCL = H11 + H12 TST2 = H22*( H11 / SCL ) IF(
TST2.EQ.RZERO .OR. H21*( H12 / SCL ).LE. MAX( SMLNUM, ULP*TST2 ) )H( K+1, K
) = ZERO END IF END IF 120 CONTINUE MEND = MIN( NBMPS, ( KBOT-KRCOL-1 ) / 3
) DO 130 M = MTOP, MEND K = KRCOL + 3*( M-1 ) REFSUM = V( 1, M )*V( 3, M
)*H( K+4, K+3 ) H( K+4, K+1 ) = -REFSUM H( K+4, K+2 ) = -REFSUM*DCONJG( V(
2, M ) ) H( K+4, K+3 ) = H( K+4, K+3 ) - REFSUM*DCONJG( V( 3, M ) ) 130
CONTINUE 140 CONTINUE IF( ACCUM ) THEN IF( WANTT ) THEN JTOP = 1 JBOT = N
ELSE JTOP = KTOP JBOT = KBOT END IF IF( ( .NOT.BLK22 ) .OR. ( INCOL.LT.KTOP
) .OR. ( NDCOL.GT.KBOT ) .OR. ( NS.LE.2 ) ) THEN K1 = MAX( 1, KTOP-INCOL )
NU = ( KDU-MAX( 0, NDCOL-KBOT ) ) - K1 + 1 DO 150 JCOL = MIN( NDCOL, KBOT )
+ 1, JBOT, NH JLEN = MIN( NH, JBOT-JCOL+1 ) CALL ZGEMM( 'C', 'N', NU, JLEN,
NU, ONE, U( K1, K1 ), LDU, H( INCOL+K1, JCOL ), LDH, ZERO, WH, LDWH ) CALL
ZLACPY( 'ALL', NU, JLEN, WH, LDWH, H( INCOL+K1, JCOL ), LDH ) 150 CONTINUE
DO 160 JROW = JTOP, MAX( KTOP, INCOL ) - 1, NV JLEN = MIN( NV, MAX( KTOP,
INCOL )-JROW ) CALL ZGEMM( 'N', 'N', JLEN, NU, NU, ONE, H( JROW, INCOL+K1 ),
LDH, U( K1, K1 ), LDU, ZERO, WV, LDWV ) CALL ZLACPY( 'ALL', JLEN, NU, WV,
LDWV, H( JROW, INCOL+K1 ), LDH ) 160 CONTINUE IF( WANTZ ) THEN DO 170 JROW =
ILOZ, IHIZ, NV JLEN = MIN( NV, IHIZ-JROW+1 ) CALL ZGEMM( 'N', 'N', JLEN, NU,
NU, ONE, Z( JROW, INCOL+K1 ), LDZ, U( K1, K1 ), LDU, ZERO, WV, LDWV ) CALL
ZLACPY( 'ALL', JLEN, NU, WV, LDWV, Z( JROW, INCOL+K1 ), LDZ ) 170 CONTINUE
END IF ELSE I2 = ( KDU+1 ) / 2 I4 = KDU J2 = I4 - I2 J4 = KDU KZS = ( J4-J2
) - ( NS+1 ) KNZ = NS + 1 DO 180 JCOL = MIN( NDCOL, KBOT ) + 1, JBOT, NH
JLEN = MIN( NH, JBOT-JCOL+1 ) CALL ZLACPY( 'ALL', KNZ, JLEN, H( INCOL+1+J2,
JCOL ), LDH, WH( KZS+1, 1 ), LDWH ) CALL ZLASET( 'ALL', KZS, JLEN, ZERO,
ZERO, WH, LDWH ) CALL ZTRMM( 'L', 'U', 'C', 'N', KNZ, JLEN, ONE, U( J2+1,
1+KZS ), LDU, WH( KZS+1, 1 ), LDWH ) CALL ZGEMM( 'C', 'N', I2, JLEN, J2,
ONE, U, LDU, H( INCOL+1, JCOL ), LDH, ONE, WH, LDWH ) CALL ZLACPY( 'ALL',
J2, JLEN, H( INCOL+1, JCOL ), LDH, WH( I2+1, 1 ), LDWH ) CALL ZTRMM( 'L',
'L', 'C', 'N', J2, JLEN, ONE, U( 1, I2+1 ), LDU, WH( I2+1, 1 ), LDWH ) CALL
ZGEMM( 'C', 'N', I4-I2, JLEN, J4-J2, ONE, U( J2+1, I2+1 ), LDU, H(
INCOL+1+J2, JCOL ), LDH, ONE, WH( I2+1, 1 ), LDWH ) CALL ZLACPY( 'ALL', KDU,
JLEN, WH, LDWH, H( INCOL+1, JCOL ), LDH ) 180 CONTINUE DO 190 JROW = JTOP,
MAX( INCOL, KTOP ) - 1, NV JLEN = MIN( NV, MAX( INCOL, KTOP )-JROW ) CALL
ZLACPY( 'ALL', JLEN, KNZ, H( JROW, INCOL+1+J2 ), LDH, WV( 1, 1+KZS ), LDWV )
CALL ZLASET( 'ALL', JLEN, KZS, ZERO, ZERO, WV, LDWV ) CALL ZTRMM( 'R', 'U',
'N', 'N', JLEN, KNZ, ONE, U( J2+1, 1+KZS ), LDU, WV( 1, 1+KZS ), LDWV ) CALL
ZGEMM( 'N', 'N', JLEN, I2, J2, ONE, H( JROW, INCOL+1 ), LDH, U, LDU, ONE,
WV, LDWV ) CALL ZLACPY( 'ALL', JLEN, J2, H( JROW, INCOL+1 ), LDH, WV( 1,
1+I2 ), LDWV ) CALL ZTRMM( 'R', 'L', 'N', 'N', JLEN, I4-I2, ONE, U( 1, I2+1
), LDU, WV( 1, 1+I2 ), LDWV ) CALL ZGEMM( 'N', 'N', JLEN, I4-I2, J4-J2, ONE,
H( JROW, INCOL+1+J2 ), LDH, U( J2+1, I2+1 ), LDU, ONE, WV( 1, 1+I2 ), LDWV )
CALL ZLACPY( 'ALL', JLEN, KDU, WV, LDWV, H( JROW, INCOL+1 ), LDH ) 190
CONTINUE IF( WANTZ ) THEN DO 200 JROW = ILOZ, IHIZ, NV JLEN = MIN( NV,
IHIZ-JROW+1 ) CALL ZLACPY( 'ALL', JLEN, KNZ, Z( JROW, INCOL+1+J2 ), LDZ, WV(
1, 1+KZS ), LDWV ) CALL ZLASET( 'ALL', JLEN, KZS, ZERO, ZERO, WV, LDWV )
CALL ZTRMM( 'R', 'U', 'N', 'N', JLEN, KNZ, ONE, U( J2+1, 1+KZS ), LDU, WV(
1, 1+KZS ), LDWV ) CALL ZGEMM( 'N', 'N', JLEN, I2, J2, ONE, Z( JROW, INCOL+1
), LDZ, U, LDU, ONE, WV, LDWV ) CALL ZLACPY( 'ALL', JLEN, J2, Z( JROW,
INCOL+1 ), LDZ, WV( 1, 1+I2 ), LDWV ) CALL ZTRMM( 'R', 'L', 'N', 'N', JLEN,
I4-I2, ONE, U( 1, I2+1 ), LDU, WV( 1, 1+I2 ), LDWV ) CALL ZGEMM( 'N', 'N',
JLEN, I4-I2, J4-J2, ONE, Z( JROW, INCOL+1+J2 ), LDZ, U( J2+1, I2+1 ), LDU,
ONE, WV( 1, 1+I2 ), LDWV ) CALL ZLACPY( 'ALL', JLEN, KDU, WV, LDWV, Z( JROW,
INCOL+1 ), LDZ ) 200 CONTINUE END IF END IF END IF 210 CONTINUE END