SYNOPSIS¶
- SUBROUTINE
SLAQR3(
- WANTT, WANTZ, N, KTOP, KBOT, NW, H, LDH, ILOZ, IHIZ, Z, LDZ, NS, ND, SR,
SI, V, LDV, NH, T, LDT, NV, WV, LDWV, WORK, LWORK )
INTEGER IHIZ, ILOZ, KBOT, KTOP, LDH, LDT, LDV, LDWV, LDZ, LWORK,
N, ND, NH, NS, NV, NW LOGICAL WANTT, WANTZ REAL H( LDH, * ), SI( * ), SR( *
), T( LDT, * ), V( LDV, * ), WORK( * ), WV( LDWV, * ), Z( LDZ, * ) REAL
ZERO, ONE PARAMETER ( ZERO = 0.0e0, ONE = 1.0e0 ) REAL AA, BB, BETA, CC, CS,
DD, EVI, EVK, FOO, S, SAFMAX, SAFMIN, SMLNUM, SN, TAU, ULP INTEGER I, IFST,
ILST, INFO, INFQR, J, JW, K, KCOL, KEND, KLN, KROW, KWTOP, LTOP, LWK1, LWK2,
LWK3, LWKOPT, NMIN LOGICAL BULGE, SORTED REAL SLAMCH INTEGER ILAENV EXTERNAL
SLAMCH, ILAENV EXTERNAL SCOPY, SGEHRD, SGEMM, SLABAD, SLACPY, SLAHQR,
SLANV2, SLAQR4, SLARF, SLARFG, SLASET, SORMHR, STREXC INTRINSIC ABS, INT,
MAX, MIN, REAL, SQRT JW = MIN( NW, KBOT-KTOP+1 ) IF( JW.LE.2 ) THEN LWKOPT =
1 ELSE CALL SGEHRD( JW, 1, JW-1, T, LDT, WORK, WORK, -1, INFO ) LWK1 = INT(
WORK( 1 ) ) CALL SORMHR( 'R', 'N', JW, JW, 1, JW-1, T, LDT, WORK, V, LDV,
WORK, -1, INFO ) LWK2 = INT( WORK( 1 ) ) CALL SLAQR4( .true., .true., JW, 1,
JW, T, LDT, SR, SI, 1, JW, V, LDV, WORK, -1, INFQR ) LWK3 = INT( WORK( 1 ) )
LWKOPT = MAX( JW+MAX( LWK1, LWK2 ), LWK3 ) END IF IF( LWORK.EQ.-1 ) THEN
WORK( 1 ) = REAL( LWKOPT ) RETURN END IF NS = 0 ND = 0 WORK( 1 ) = ONE IF(
KTOP.GT.KBOT ) RETURN IF( NW.LT.1 ) RETURN SAFMIN = SLAMCH( 'SAFE MINIMUM' )
SAFMAX = ONE / SAFMIN CALL SLABAD( SAFMIN, SAFMAX ) ULP = SLAMCH(
'PRECISION' ) SMLNUM = SAFMIN*( REAL( N ) / ULP ) JW = MIN( NW, KBOT-KTOP+1
) KWTOP = KBOT - JW + 1 IF( KWTOP.EQ.KTOP ) THEN S = ZERO ELSE S = H( KWTOP,
KWTOP-1 ) END IF IF( KBOT.EQ.KWTOP ) THEN SR( KWTOP ) = H( KWTOP, KWTOP )
SI( KWTOP ) = ZERO NS = 1 ND = 0 IF( ABS( S ).LE.MAX( SMLNUM, ULP*ABS( H(
KWTOP, KWTOP ) ) ) ) THEN NS = 0 ND = 1 IF( KWTOP.GT.KTOP ) H( KWTOP,
KWTOP-1 ) = ZERO END IF WORK( 1 ) = ONE RETURN END IF CALL SLACPY( 'U', JW,
JW, H( KWTOP, KWTOP ), LDH, T, LDT ) CALL SCOPY( JW-1, H( KWTOP+1, KWTOP ),
LDH+1, T( 2, 1 ), LDT+1 ) CALL SLASET( 'A', JW, JW, ZERO, ONE, V, LDV ) NMIN
= ILAENV( 12, 'SLAQR3', 'SV', JW, 1, JW, LWORK ) IF( JW.GT.NMIN ) THEN CALL
SLAQR4( .true., .true., JW, 1, JW, T, LDT, SR( KWTOP ), SI( KWTOP ), 1, JW,
V, LDV, WORK, LWORK, INFQR ) ELSE CALL SLAHQR( .true., .true., JW, 1, JW, T,
LDT, SR( KWTOP ), SI( KWTOP ), 1, JW, V, LDV, INFQR ) END IF DO 10 J = 1, JW
- 3 T( J+2, J ) = ZERO T( J+3, J ) = ZERO 10 CONTINUE IF( JW.GT.2 ) T( JW,
JW-2 ) = ZERO NS = JW ILST = INFQR + 1 20 CONTINUE IF( ILST.LE.NS ) THEN IF(
NS.EQ.1 ) THEN BULGE = .FALSE. ELSE BULGE = T( NS, NS-1 ).NE.ZERO END IF IF(
FOO = ABS( T( NS, NS ) ) IF( FOO.EQ.ZERO ) FOO = ABS( S ) IF( ABS( S*V( 1,
NS ) ).LE.MAX( SMLNUM, ULP*FOO ) ) THEN NS = NS - 1 ELSE IFST = NS CALL
STREXC( 'V', JW, T, LDT, V, LDV, IFST, ILST, WORK, INFO ) ILST = ILST + 1
END IF ELSE FOO = ABS( T( NS, NS ) ) + SQRT( ABS( T( NS, NS-1 ) ) )* SQRT(
ABS( T( NS-1, NS ) ) ) IF( FOO.EQ.ZERO ) FOO = ABS( S ) IF( MAX( ABS( S*V(
1, NS ) ), ABS( S*V( 1, NS-1 ) ) ).LE. MAX( SMLNUM, ULP*FOO ) ) THEN NS = NS
- 2 ELSE IFST = NS CALL STREXC( 'V', JW, T, LDT, V, LDV, IFST, ILST, WORK,
INFO ) ILST = ILST + 2 END IF END IF GO TO 20 END IF IF( NS.EQ.0 ) S = ZERO
IF( NS.LT.JW ) THEN SORTED = .false. I = NS + 1 30 CONTINUE IF( SORTED ) GO
TO 50 SORTED = .true. KEND = I - 1 I = INFQR + 1 IF( I.EQ.NS ) THEN K = I +
1 ELSE IF( T( I+1, I ).EQ.ZERO ) THEN K = I + 1 ELSE K = I + 2 END IF 40
CONTINUE IF( K.LE.KEND ) THEN IF( K.EQ.I+1 ) THEN EVI = ABS( T( I, I ) )
ELSE EVI = ABS( T( I, I ) ) + SQRT( ABS( T( I+1, I ) ) )* SQRT( ABS( T( I,
I+1 ) ) ) END IF IF( K.EQ.KEND ) THEN EVK = ABS( T( K, K ) ) ELSE IF( T(
K+1, K ).EQ.ZERO ) THEN EVK = ABS( T( K, K ) ) ELSE EVK = ABS( T( K, K ) ) +
SQRT( ABS( T( K+1, K ) ) )* SQRT( ABS( T( K, K+1 ) ) ) END IF IF( EVI.GE.EVK
) THEN I = K ELSE SORTED = .false. IFST = I ILST = K CALL STREXC( 'V', JW,
T, LDT, V, LDV, IFST, ILST, WORK, INFO ) IF( INFO.EQ.0 ) THEN I = ILST ELSE
I = K END IF END IF IF( I.EQ.KEND ) THEN K = I + 1 ELSE IF( T( I+1, I
).EQ.ZERO ) THEN K = I + 1 ELSE K = I + 2 END IF GO TO 40 END IF GO TO 30 50
CONTINUE END IF I = JW 60 CONTINUE IF( I.GE.INFQR+1 ) THEN IF( I.EQ.INFQR+1
) THEN SR( KWTOP+I-1 ) = T( I, I ) SI( KWTOP+I-1 ) = ZERO I = I - 1 ELSE IF(
T( I, I-1 ).EQ.ZERO ) THEN SR( KWTOP+I-1 ) = T( I, I ) SI( KWTOP+I-1 ) =
ZERO I = I - 1 ELSE AA = T( I-1, I-1 ) CC = T( I, I-1 ) BB = T( I-1, I ) DD
= T( I, I ) CALL SLANV2( AA, BB, CC, DD, SR( KWTOP+I-2 ), SI( KWTOP+I-2 ),
SR( KWTOP+I-1 ), SI( KWTOP+I-1 ), CS, SN ) I = I - 2 END IF GO TO 60 END IF
IF( NS.LT.JW .OR. S.EQ.ZERO ) THEN IF( NS.GT.1 .AND. S.NE.ZERO ) THEN CALL
SCOPY( NS, V, LDV, WORK, 1 ) BETA = WORK( 1 ) CALL SLARFG( NS, BETA, WORK( 2
), 1, TAU ) WORK( 1 ) = ONE CALL SLASET( 'L', JW-2, JW-2, ZERO, ZERO, T( 3,
1 ), LDT ) CALL SLARF( 'L', NS, JW, WORK, 1, TAU, T, LDT, WORK( JW+1 ) )
CALL SLARF( 'R', NS, NS, WORK, 1, TAU, T, LDT, WORK( JW+1 ) ) CALL SLARF(
'R', JW, NS, WORK, 1, TAU, V, LDV, WORK( JW+1 ) ) CALL SGEHRD( JW, 1, NS, T,
LDT, WORK, WORK( JW+1 ), LWORK-JW, INFO ) END IF IF( KWTOP.GT.1 ) H( KWTOP,
KWTOP-1 ) = S*V( 1, 1 ) CALL SLACPY( 'U', JW, JW, T, LDT, H( KWTOP, KWTOP ),
LDH ) CALL SCOPY( JW-1, T( 2, 1 ), LDT+1, H( KWTOP+1, KWTOP ), LDH+1 ) IF(
NS.GT.1 .AND. S.NE.ZERO ) CALL SORMHR( 'R', 'N', JW, NS, 1, NS, T, LDT,
WORK, V, LDV, WORK( JW+1 ), LWORK-JW, INFO ) IF( WANTT ) THEN LTOP = 1 ELSE
LTOP = KTOP END IF DO 70 KROW = LTOP, KWTOP - 1, NV KLN = MIN( NV,
KWTOP-KROW ) CALL SGEMM( 'N', 'N', KLN, JW, JW, ONE, H( KROW, KWTOP ), LDH,
V, LDV, ZERO, WV, LDWV ) CALL SLACPY( 'A', KLN, JW, WV, LDWV, H( KROW, KWTOP
), LDH ) 70 CONTINUE IF( WANTT ) THEN DO 80 KCOL = KBOT + 1, N, NH KLN =
MIN( NH, N-KCOL+1 ) CALL SGEMM( 'C', 'N', JW, KLN, JW, ONE, V, LDV, H(
KWTOP, KCOL ), LDH, ZERO, T, LDT ) CALL SLACPY( 'A', JW, KLN, T, LDT, H(
KWTOP, KCOL ), LDH ) 80 CONTINUE END IF IF( WANTZ ) THEN DO 90 KROW = ILOZ,
IHIZ, NV KLN = MIN( NV, IHIZ-KROW+1 ) CALL SGEMM( 'N', 'N', KLN, JW, JW,
ONE, Z( KROW, KWTOP ), LDZ, V, LDV, ZERO, WV, LDWV ) CALL SLACPY( 'A', KLN,
JW, WV, LDWV, Z( KROW, KWTOP ), LDZ ) 90 CONTINUE END IF END IF ND = JW - NS
NS = NS - INFQR WORK( 1 ) = REAL( LWKOPT ) END