dlaqr5 man page on Scientific
[printable version]
DLAQR5(1) LAPACK auxiliary routine (version 3.2) DLAQR5(1)
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( 'SAFE MINIMUM' )
SAFMAX = ONE / 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 220 INCOL = 3*( 1-NBMPS ) + KTOP - 1, KBOT - 2,
3*NBMPS - 2
NDCOL = INCOL + KDU
IF( ACCUM ) CALL DLASET( 'ALL', 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( SML‐
NUM, 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( 'C', 'N', NU, JLEN, NU, ONE, U( K1, K1 ),
LDU, H( INCOL+K1, JCOL ), LDH, ZERO, WH, LDWH )
CALL DLACPY( 'ALL', 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( 'N', 'N', JLEN, NU, NU, ONE, H( JROW,
INCOL+K1 ), LDH, U( K1, K1 ), LDU, ZERO, WV, LDWV )
CALL DLACPY( 'ALL', 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( 'N', 'N', JLEN, NU, NU, ONE, Z( JROW,
INCOL+K1 ), LDZ, U( K1, K1 ), LDU, ZERO, WV, LDWV )
CALL DLACPY( 'ALL', 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( 'ALL', KNZ, JLEN, H( INCOL+1+J2, JCOL ),
LDH, WH( KZS+1, 1 ), LDWH )
CALL DLASET( 'ALL', KZS, JLEN, ZERO, ZERO, WH, LDWH )
CALL DTRMM( 'L', 'U', 'C', 'N', KNZ, JLEN, ONE, U( J2+1,
1+KZS ), LDU, WH( KZS+1, 1 ), LDWH )
CALL DGEMM( 'C', 'N', I2, JLEN, J2, ONE, U, LDU, H(
INCOL+1, JCOL ), LDH, ONE, WH, LDWH )
CALL DLACPY( 'ALL', J2, JLEN, H( INCOL+1, JCOL ), LDH,
WH( I2+1, 1 ), LDWH )
CALL DTRMM( 'L', 'L', 'C', 'N', J2, JLEN, ONE, U( 1, I2+1
), LDU, WH( I2+1, 1 ), LDWH )
CALL DGEMM( '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 DLACPY( 'ALL', 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( 'ALL', JLEN, KNZ, H( JROW, INCOL+1+J2 ),
LDH, WV( 1, 1+KZS ), LDWV )
CALL DLASET( 'ALL', JLEN, KZS, ZERO, ZERO, WV, LDWV )
CALL DTRMM( 'R', 'U', 'N', 'N', JLEN, KNZ, ONE, U( J2+1,
1+KZS ), LDU, WV( 1, 1+KZS ), LDWV )
CALL DGEMM( 'N', 'N', JLEN, I2, J2, ONE, H( JROW, INCOL+1
), LDH, U, LDU, ONE, WV, LDWV )
CALL DLACPY( 'ALL', JLEN, J2, H( JROW, INCOL+1 ), LDH,
WV( 1, 1+I2 ), LDWV )
CALL DTRMM( 'R', 'L', 'N', 'N', JLEN, I4-I2, ONE, U( 1,
I2+1 ), LDU, WV( 1, 1+I2 ), LDWV )
CALL DGEMM( '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 DLACPY( 'ALL', 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( 'ALL', JLEN, KNZ, Z( JROW, INCOL+1+J2 ),
LDZ, WV( 1, 1+KZS ), LDWV )
CALL DLASET( 'ALL', JLEN, KZS, ZERO, ZERO, WV, LDWV )
CALL DTRMM( 'R', 'U', 'N', 'N', JLEN, KNZ, ONE, U( J2+1,
1+KZS ), LDU, WV( 1, 1+KZS ), LDWV )
CALL DGEMM( 'N', 'N', JLEN, I2, J2, ONE, Z( JROW, INCOL+1
), LDZ, U, LDU, ONE, WV, LDWV )
CALL DLACPY( 'ALL', JLEN, J2, Z( JROW, INCOL+1 ), LDZ,
WV( 1, 1+I2 ), LDWV )
CALL DTRMM( 'R', 'L', 'N', 'N', JLEN, I4-I2, ONE, U( 1,
I2+1 ), LDU, WV( 1, 1+I2 ), LDWV )
CALL DGEMM( '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 DLACPY( 'ALL', JLEN, KDU, WV, LDWV, Z( JROW, INCOL+1
), LDZ )
210 CONTINUE
END IF
END IF
END IF
220 CONTINUE
END
PURPOSE
LAPACK auxiliary routine (versioNovember 2008 DLAQR5(1)
[top]
List of man pages available for Scientific
Copyright (c) for man pages and the logo by the respective OS vendor.
For those who want to learn more, the polarhome community provides shell access and support.
[legal]
[privacy]
[GNU]
[policy]
[cookies]
[netiquette]
[sponsors]
[FAQ]
Polarhome, production since 1999.
Member of Polarhome portal.
Based on Fawad Halim's script.
....................................................................
|
|
Vote for polarhome
|
