slaqr2 man page on Scientific
[printable version]
SLAQR2(1) LAPACK auxiliary routine (version 3.2) SLAQR2(1)
NAME
SYNOPSIS
SUBROUTINE SLAQR2( 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, LWKOPT
LOGICAL BULGE, SORTED
REAL SLAMCH
EXTERNAL SLAMCH
EXTERNAL SCOPY, SGEHRD, SGEMM, SLABAD, SLACPY, SLAHQR,
SLANV2, 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 ) )
LWKOPT = JW + MAX( LWK1, LWK2 )
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 )
CALL SLAHQR( .true., .true., JW, 1, JW, T, LDT, SR( KWTOP
), SI( KWTOP ), 1, JW, V, LDV, INFQR )
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
PURPOSE
LAPACK auxiliary routine (versioNovember 2008 SLAQR2(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
|
