zlaqr3(l) - Linux man page
Name
Synopsis
- SUBROUTINE ZLAQR3(
WANTT, WANTZ, N, KTOP, KBOT, NW, H, LDH, ILOZ, IHIZ, Z, LDZ, NS, ND, SH, 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
COMPLEX*16
H( LDH, * ), SH( * ), T( LDT, * ), V( LDV, * ), WORK( * ), 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
BETA, CDUM, S, TAU
DOUBLE
PRECISION FOO, SAFMAX, SAFMIN, SMLNUM, ULP
INTEGER
I, IFST, ILST, INFO, INFQR, J, JW, KCOL, KLN, KNT, KROW, KWTOP, LTOP, LWK1, LWK2, LWK3, LWKOPT, NMIN
DOUBLE
PRECISION DLAMCH
INTEGER
ILAENV
EXTERNAL
DLAMCH, ILAENV
EXTERNAL
DLABAD, ZCOPY, ZGEHRD, ZGEMM, ZLACPY, ZLAHQR, ZLAQR4, ZLARF, ZLARFG, ZLASET, ZTREXC, ZUNMHR
INTRINSIC
ABS, DBLE, DCMPLX, DCONJG, DIMAG, INT, MAX, MIN
DOUBLE
PRECISION CABS1
CABS1(
CDUM ) = ABS( DBLE( CDUM ) ) + ABS( DIMAG( CDUM ) )
JW
= MIN( NW, KBOT-KTOP+1 )
IF(
JW.LE.2 ) THEN
LWKOPT
= 1
ELSE
CALL
ZGEHRD( JW, 1, JW-1, T, LDT, WORK, WORK, -1, INFO )
LWK1
= INT( WORK( 1 ) )
CALL
ZUNMHR( 'R', 'N', JW, JW, 1, JW-1, T, LDT, WORK, V, LDV, WORK, -1, INFO )
LWK2
= INT( WORK( 1 ) )
CALL
ZLAQR4( .true., .true., JW, 1, JW, T, LDT, SH, 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 ) = DCMPLX( LWKOPT, 0 )
RETURN
END
IF
NS
= 0
ND
= 0
WORK(
1 ) = ONE
IF(
KTOP.GT.KBOT ) RETURN
IF(
NW.LT.1 ) RETURN
SAFMIN
= DLAMCH( 'SAFE MINIMUM' )
SAFMAX
= RONE / SAFMIN
CALL
DLABAD( SAFMIN, SAFMAX )
ULP
= DLAMCH( 'PRECISION' )
SMLNUM
= SAFMIN*( DBLE( 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
SH(
KWTOP ) = H( KWTOP, KWTOP )
NS
= 1
ND
= 0
IF(
CABS1( S ).LE.MAX( SMLNUM, ULP*CABS1( 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
ZLACPY( 'U', JW, JW, H( KWTOP, KWTOP ), LDH, T, LDT )
CALL
ZCOPY( JW-1, H( KWTOP+1, KWTOP ), LDH+1, T( 2, 1 ), LDT+1 )
CALL
ZLASET( 'A', JW, JW, ZERO, ONE, V, LDV )
NMIN
= ILAENV( 12, 'ZLAQR3', 'SV', JW, 1, JW, LWORK )
IF(
JW.GT.NMIN ) THEN
CALL
ZLAQR4( .true., .true., JW, 1, JW, T, LDT, SH( KWTOP ), 1, JW, V, LDV, WORK, LWORK, INFQR )
ELSE
CALL
ZLAHQR( .true., .true., JW, 1, JW, T, LDT, SH( KWTOP ), 1, JW, V, LDV, INFQR )
END
IF
NS
= JW
ILST
= INFQR + 1
DO
10 KNT = INFQR + 1, JW
FOO
= CABS1( T( NS, NS ) )
IF(
FOO.EQ.RZERO ) FOO = CABS1( S )
IF(
CABS1( S )*CABS1( V( 1, NS ) ).LE.MAX( SMLNUM, ULP*FOO ) ) THEN
NS
= NS - 1
ELSE
IFST
= NS
CALL
ZTREXC( 'V', JW, T, LDT, V, LDV, IFST, ILST, INFO )
ILST
= ILST + 1
END
IF
10
CONTINUE
IF(
NS.EQ.0 ) S = ZERO
IF(
NS.LT.JW ) THEN
DO
30 I = INFQR + 1, NS
IFST
= I
DO
20 J = I + 1, NS
IF(
CABS1( T( J, J ) ).GT.CABS1( T( IFST, IFST ) ) ) IFST = J
20
CONTINUE
ILST
= I
IF(
IFST.NE.ILST ) CALL ZTREXC( 'V', JW, T, LDT, V, LDV, IFST, ILST, INFO )
30
CONTINUE
END
IF
DO
40 I = INFQR + 1, JW
SH(
KWTOP+I-1 ) = T( I, I )
40
CONTINUE
IF(
NS.LT.JW .OR. S.EQ.ZERO ) THEN
IF(
NS.GT.1 .AND. S.NE.ZERO ) THEN
CALL
ZCOPY( NS, V, LDV, WORK, 1 )
DO
50 I = 1, NS
WORK(
I ) = DCONJG( WORK( I ) )
50
CONTINUE
BETA
= WORK( 1 )
CALL
ZLARFG( NS, BETA, WORK( 2 ), 1, TAU )
WORK(
1 ) = ONE
CALL
ZLASET( 'L', JW-2, JW-2, ZERO, ZERO, T( 3, 1 ), LDT )
CALL
ZLARF( 'L', NS, JW, WORK, 1, DCONJG( TAU ), T, LDT, WORK( JW+1 ) )
CALL
ZLARF( 'R', NS, NS, WORK, 1, TAU, T, LDT, WORK( JW+1 ) )
CALL
ZLARF( 'R', JW, NS, WORK, 1, TAU, V, LDV, WORK( JW+1 ) )
CALL
ZGEHRD( JW, 1, NS, T, LDT, WORK, WORK( JW+1 ), LWORK-JW, INFO )
END
IF
IF(
KWTOP.GT.1 ) H( KWTOP, KWTOP-1 ) = S*DCONJG( V( 1, 1 ) )
CALL
ZLACPY( 'U', JW, JW, T, LDT, H( KWTOP, KWTOP ), LDH )
CALL
ZCOPY( JW-1, T( 2, 1 ), LDT+1, H( KWTOP+1, KWTOP ), LDH+1 )
IF(
NS.GT.1 .AND. S.NE.ZERO ) CALL ZUNMHR( '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
60 KROW = LTOP, KWTOP - 1, NV
KLN
= MIN( NV, KWTOP-KROW )
CALL
ZGEMM( 'N', 'N', KLN, JW, JW, ONE, H( KROW, KWTOP ), LDH, V, LDV, ZERO, WV, LDWV )
CALL
ZLACPY( 'A', KLN, JW, WV, LDWV, H( KROW, KWTOP ), LDH )
60
CONTINUE
IF(
WANTT ) THEN
DO
70 KCOL = KBOT + 1, N, NH
KLN
= MIN( NH, N-KCOL+1 )
CALL
ZGEMM( 'C', 'N', JW, KLN, JW, ONE, V, LDV, H( KWTOP, KCOL ), LDH, ZERO, T, LDT )
CALL
ZLACPY( 'A', JW, KLN, T, LDT, H( KWTOP, KCOL ), LDH )
70
CONTINUE
END
IF
IF(
WANTZ ) THEN
DO
80 KROW = ILOZ, IHIZ, NV
KLN
= MIN( NV, IHIZ-KROW+1 )
CALL
ZGEMM( 'N', 'N', KLN, JW, JW, ONE, Z( KROW, KWTOP ), LDZ, V, LDV, ZERO, WV, LDWV )
CALL
ZLACPY( 'A', KLN, JW, WV, LDWV, Z( KROW, KWTOP ), LDZ )
80
CONTINUE
END
IF
END
IF
ND
= JW - NS
NS
= NS - INFQR
WORK(
1 ) = DCMPLX( LWKOPT, 0 )
END