dlaqr5 (l) - Linux Manuals

NAME

SYNOPSIS

SUBROUTINE DLAQR5(

WANTT, WANTZ, KACC22, N, KTOP, KBOT, NSHFTS, SR, SI, 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

    

DOUBLE PRECISION H( LDH, * ), SI( * ), SR( * ), U( LDU, * ), V( LDV, * ), WH( LDWH, * ), WV( LDWV, * ), Z( LDZ, * )

    

DOUBLE PRECISION ZERO, ONE

    

PARAMETER ( ZERO = 0.0d0, ONE = 1.0d0 )

    

DOUBLE PRECISION ALPHA, BETA, H11, H12, H21, H22, REFSUM, SAFMAX, SAFMIN, SCL, SMLNUM, SWAP, TST1, TST2, ULP

    

INTEGER I, 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, MAX, MIN, MOD

    

DOUBLE PRECISION VT( 3 )

    

EXTERNAL DGEMM, DLABAD, DLACPY, DLAQR1, DLARFG, DLASET, DTRMM

    

IF( NSHFTS.LT.2 ) RETURN

    

IF( KTOP.GE.KBOT ) RETURN

    

DO 10 I = 1, NSHFTS - 2, 2

    

IF( SI( I ).NE.-SI( I+1 ) ) THEN

    

SWAP = SR( I )

    

SR( I ) = SR( I+1 )

    

SR( I+1 ) = SR( I+2 )

    

SR( I+2 ) = SWAP

    

SWAP = SI( I )

    

SI( I ) = SI( I+1 )

    

SI( I+1 ) = SI( I+2 )

    

SI( I+2 ) = SWAP

    

END IF

    

10 CONTINUE

    

NS = NSHFTS - MOD( NSHFTS, 2 )

    

SAFMIN = DLAMCH( aqSAFE MINIMUMaq )

    

SAFMAX = ONE / SAFMIN

    

CALL DLABAD( SAFMIN, SAFMAX )

    

ULP = DLAMCH( aqPRECISIONaq )

    

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 220 INCOL = 3*( 1-NBMPS ) + KTOP - 1, KBOT - 2, 3*NBMPS - 2

    

NDCOL = INCOL + KDU

    

IF( ACCUM ) CALL DLASET( aqALLaq, KDU, KDU, ZERO, ONE, U, LDU )

    

DO 150 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 20 M = MTOP, MBOT

    

K = KRCOL + 3*( M-1 )

    

IF( K.EQ.KTOP-1 ) THEN

    

CALL DLAQR1( 3, H( KTOP, KTOP ), LDH, SR( 2*M-1 ), SI( 2*M-1 ), SR( 2*M ), SI( 2*M ), V( 1, M ) )

    

ALPHA = V( 1, M )

    

CALL DLARFG( 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 DLARFG( 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 DLAQR1( 3, H( K+1, K+1 ), LDH, SR( 2*M-1 ), SI( 2*M-1 ), SR( 2*M ), SI( 2*M ), VT )

    

ALPHA = VT( 1 )

    

CALL DLARFG( 3, ALPHA, VT( 2 ), 1, VT( 1 ) )

    

REFSUM = VT( 1 )*( H( K+1, K )+VT( 2 )* H( K+2, K ) )

    

IF( ABS( H( K+2, K )-REFSUM*VT( 2 ) )+ ABS( REFSUM*VT( 3 ) ).GT.ULP* ( ABS( H( K, K ) )+ABS( H( K+1, K+1 ) )+ABS( 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

    

20 CONTINUE

    

K = KRCOL + 3*( M22-1 )

    

IF( BMP22 ) THEN

    

IF( K.EQ.KTOP-1 ) THEN

    

CALL DLAQR1( 2, H( K+1, K+1 ), LDH, SR( 2*M22-1 ), SI( 2*M22-1 ), SR( 2*M22 ), SI( 2*M22 ), V( 1, M22 ) )

    

BETA = V( 1, M22 )

    

CALL DLARFG( 2, BETA, V( 2, M22 ), 1, V( 1, M22 ) )

    

ELSE

    

BETA = H( K+1, K )

    

V( 2, M22 ) = H( K+2, K )

    

CALL DLARFG( 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 40 J = MAX( KTOP, KRCOL ), JBOT

    

MEND = MIN( MBOT, ( J-KRCOL+2 ) / 3 )

    

DO 30 M = MTOP, MEND

    

K = KRCOL + 3*( M-1 )

    

REFSUM = V( 1, M )*( H( K+1, J )+V( 2, M )* H( K+2, J )+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 )

    

30 CONTINUE

    

40 CONTINUE

    

IF( BMP22 ) THEN

    

K = KRCOL + 3*( M22-1 )

    

DO 50 J = MAX( K+1, KTOP ), JBOT

    

REFSUM = V( 1, M22 )*( H( K+1, J )+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 )

    

50 CONTINUE

    

END IF

    

IF( ACCUM ) THEN

    

JTOP = MAX( KTOP, INCOL )

    

ELSE IF( WANTT ) THEN

    

JTOP = 1

    

ELSE

    

JTOP = KTOP

    

END IF

    

DO 90 M = MTOP, MBOT

    

IF( V( 1, M ).NE.ZERO ) THEN

    

K = KRCOL + 3*( M-1 )

    

DO 60 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*V( 2, M )

    

H( J, K+3 ) = H( J, K+3 ) - REFSUM*V( 3, M )

    

60 CONTINUE

    

IF( ACCUM ) THEN

    

KMS = K - INCOL

    

DO 70 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*V( 2, M )

    

U( J, KMS+3 ) = U( J, KMS+3 ) - REFSUM*V( 3, M )

    

70 CONTINUE

    

ELSE IF( WANTZ ) THEN

    

DO 80 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*V( 2, M )

    

Z( J, K+3 ) = Z( J, K+3 ) - REFSUM*V( 3, M )

    

80 CONTINUE

    

END IF

    

END IF

    

90 CONTINUE

    

K = KRCOL + 3*( M22-1 )

    

IF( BMP22 .AND. ( V( 1, M22 ).NE.ZERO ) ) THEN

    

DO 100 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*V( 2, M22 )

    

100 CONTINUE

    

IF( ACCUM ) THEN

    

KMS = K - INCOL

    

DO 110 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*V( 2, M22 )

    

110 CONTINUE

    

ELSE IF( WANTZ ) THEN

    

DO 120 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*V( 2, M22 )

    

120 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 130 M = MSTART, MEND

    

K = MIN( KBOT-1, KRCOL+3*( M-1 ) )

    

IF( H( K+1, K ).NE.ZERO ) THEN

    

TST1 = ABS( H( K, K ) ) + ABS( H( K+1, K+1 ) )

    

IF( TST1.EQ.ZERO ) THEN

    

IF( K.GE.KTOP+1 ) TST1 = TST1 + ABS( H( K, K-1 ) )

    

IF( K.GE.KTOP+2 ) TST1 = TST1 + ABS( H( K, K-2 ) )

    

IF( K.GE.KTOP+3 ) TST1 = TST1 + ABS( H( K, K-3 ) )

    

IF( K.LE.KBOT-2 ) TST1 = TST1 + ABS( H( K+2, K+1 ) )

    

IF( K.LE.KBOT-3 ) TST1 = TST1 + ABS( H( K+3, K+1 ) )

    

IF( K.LE.KBOT-4 ) TST1 = TST1 + ABS( H( K+4, K+1 ) )

    

END IF

    

IF( ABS( H( K+1, K ) ).LE.MAX( SMLNUM, ULP*TST1 ) ) THEN

    

H12 = MAX( ABS( H( K+1, K ) ), ABS( H( K, K+1 ) ) )

    

H21 = MIN( ABS( H( K+1, K ) ), ABS( H( K, K+1 ) ) )

    

H11 = MAX( ABS( H( K+1, K+1 ) ), ABS( H( K, K )-H( K+1, K+1 ) ) )

    

H22 = MIN( ABS( H( K+1, K+1 ) ), ABS( H( K, K )-H( K+1, K+1 ) ) )

    

SCL = H11 + H12

    

TST2 = H22*( H11 / SCL )

    

IF( TST2.EQ.ZERO .OR. H21*( H12 / SCL ).LE. MAX( SMLNUM, ULP*TST2 ) )H( K+1, K ) = ZERO

    

END IF

    

END IF

    

130 CONTINUE

    

MEND = MIN( NBMPS, ( KBOT-KRCOL-1 ) / 3 )

    

DO 140 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*V( 2, M )

    

H( K+4, K+3 ) = H( K+4, K+3 ) - REFSUM*V( 3, M )

    

140 CONTINUE

    

150 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 160 JCOL = MIN( NDCOL, KBOT ) + 1, JBOT, NH

    

JLEN = MIN( NH, JBOT-JCOL+1 )

    

CALL DGEMM( aqCaq, aqNaq, NU, JLEN, NU, ONE, U( K1, K1 ), LDU, H( INCOL+K1, JCOL ), LDH, ZERO, WH, LDWH )

    

CALL DLACPY( aqALLaq, NU, JLEN, WH, LDWH, H( INCOL+K1, JCOL ), LDH )

    

160 CONTINUE

    

DO 170 JROW = JTOP, MAX( KTOP, INCOL ) - 1, NV

    

JLEN = MIN( NV, MAX( KTOP, INCOL )-JROW )

    

CALL DGEMM( aqNaq, aqNaq, JLEN, NU, NU, ONE, H( JROW, INCOL+K1 ), LDH, U( K1, K1 ), LDU, ZERO, WV, LDWV )

    

CALL DLACPY( aqALLaq, JLEN, NU, WV, LDWV, H( JROW, INCOL+K1 ), LDH )

    

170 CONTINUE

    

IF( WANTZ ) THEN

    

DO 180 JROW = ILOZ, IHIZ, NV

    

JLEN = MIN( NV, IHIZ-JROW+1 )

    

CALL DGEMM( aqNaq, aqNaq, JLEN, NU, NU, ONE, Z( JROW, INCOL+K1 ), LDZ, U( K1, K1 ), LDU, ZERO, WV, LDWV )

    

CALL DLACPY( aqALLaq, JLEN, NU, WV, LDWV, Z( JROW, INCOL+K1 ), LDZ )

    

180 CONTINUE

    

END IF

    

ELSE

    

I2 = ( KDU+1 ) / 2

    

I4 = KDU

    

J2 = I4 - I2

    

J4 = KDU

    

KZS = ( J4-J2 ) - ( NS+1 )

    

KNZ = NS + 1

    

DO 190 JCOL = MIN( NDCOL, KBOT ) + 1, JBOT, NH

    

JLEN = MIN( NH, JBOT-JCOL+1 )

    

CALL DLACPY( aqALLaq, KNZ, JLEN, H( INCOL+1+J2, JCOL ), LDH, WH( KZS+1, 1 ), LDWH )

    

CALL DLASET( aqALLaq, KZS, JLEN, ZERO, ZERO, WH, LDWH )

    

CALL DTRMM( aqLaq, aqUaq, aqCaq, aqNaq, KNZ, JLEN, ONE, U( J2+1, 1+KZS ), LDU, WH( KZS+1, 1 ), LDWH )

    

CALL DGEMM( aqCaq, aqNaq, I2, JLEN, J2, ONE, U, LDU, H( INCOL+1, JCOL ), LDH, ONE, WH, LDWH )

    

CALL DLACPY( aqALLaq, J2, JLEN, H( INCOL+1, JCOL ), LDH, WH( I2+1, 1 ), LDWH )

    

CALL DTRMM( aqLaq, aqLaq, aqCaq, aqNaq, J2, JLEN, ONE, U( 1, I2+1 ), LDU, WH( I2+1, 1 ), LDWH )

    

CALL DGEMM( aqCaq, aqNaq, 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 DLACPY( aqALLaq, KDU, JLEN, WH, LDWH, H( INCOL+1, JCOL ), LDH )

    

190 CONTINUE

    

DO 200 JROW = JTOP, MAX( INCOL, KTOP ) - 1, NV

    

JLEN = MIN( NV, MAX( INCOL, KTOP )-JROW )

    

CALL DLACPY( aqALLaq, JLEN, KNZ, H( JROW, INCOL+1+J2 ), LDH, WV( 1, 1+KZS ), LDWV )

    

CALL DLASET( aqALLaq, JLEN, KZS, ZERO, ZERO, WV, LDWV )

    

CALL DTRMM( aqRaq, aqUaq, aqNaq, aqNaq, JLEN, KNZ, ONE, U( J2+1, 1+KZS ), LDU, WV( 1, 1+KZS ), LDWV )

    

CALL DGEMM( aqNaq, aqNaq, JLEN, I2, J2, ONE, H( JROW, INCOL+1 ), LDH, U, LDU, ONE, WV, LDWV )

    

CALL DLACPY( aqALLaq, JLEN, J2, H( JROW, INCOL+1 ), LDH, WV( 1, 1+I2 ), LDWV )

    

CALL DTRMM( aqRaq, aqLaq, aqNaq, aqNaq, JLEN, I4-I2, ONE, U( 1, I2+1 ), LDU, WV( 1, 1+I2 ), LDWV )

    

CALL DGEMM( aqNaq, aqNaq, 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 DLACPY( aqALLaq, JLEN, KDU, WV, LDWV, H( JROW, INCOL+1 ), LDH )

    

200 CONTINUE

    

IF( WANTZ ) THEN

    

DO 210 JROW = ILOZ, IHIZ, NV

    

JLEN = MIN( NV, IHIZ-JROW+1 )

    

CALL DLACPY( aqALLaq, JLEN, KNZ, Z( JROW, INCOL+1+J2 ), LDZ, WV( 1, 1+KZS ), LDWV )

    

CALL DLASET( aqALLaq, JLEN, KZS, ZERO, ZERO, WV, LDWV )

    

CALL DTRMM( aqRaq, aqUaq, aqNaq, aqNaq, JLEN, KNZ, ONE, U( J2+1, 1+KZS ), LDU, WV( 1, 1+KZS ), LDWV )

    

CALL DGEMM( aqNaq, aqNaq, JLEN, I2, J2, ONE, Z( JROW, INCOL+1 ), LDZ, U, LDU, ONE, WV, LDWV )

    

CALL DLACPY( aqALLaq, JLEN, J2, Z( JROW, INCOL+1 ), LDZ, WV( 1, 1+I2 ), LDWV )

    

CALL DTRMM( aqRaq, aqLaq, aqNaq, aqNaq, JLEN, I4-I2, ONE, U( 1, I2+1 ), LDU, WV( 1, 1+I2 ), LDWV )

    

CALL DGEMM( aqNaq, aqNaq, 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 DLACPY( aqALLaq, JLEN, KDU, WV, LDWV, Z( JROW, INCOL+1 ), LDZ )

    

210 CONTINUE

    

END IF

    

END IF

    

END IF

    

220 CONTINUE

    

END

PURPOSE