zlaqr3 (l) - Linux Manuals

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( aqRaq, aqNaq, 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( aqSAFE MINIMUMaq )

    
SAFMAX = RONE / SAFMIN

    
CALL DLABAD( SAFMIN, SAFMAX )

    
ULP = DLAMCH( aqPRECISIONaq )

    
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( aqUaq, 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( aqAaq, JW, JW, ZERO, ONE, V, LDV )

    
NMIN = ILAENV( 12, aqZLAQR3aq, aqSVaq, 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( aqVaq, 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( aqVaq, 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( aqLaq, JW-2, JW-2, ZERO, ZERO, T( 3, 1 ), LDT )

    
CALL ZLARF( aqLaq, NS, JW, WORK, 1, DCONJG( TAU ), T, LDT, WORK( JW+1 ) )

    
CALL ZLARF( aqRaq, NS, NS, WORK, 1, TAU, T, LDT, WORK( JW+1 ) )

    
CALL ZLARF( aqRaq, 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( aqUaq, 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( aqRaq, aqNaq, 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( aqNaq, aqNaq, KLN, JW, JW, ONE, H( KROW, KWTOP ), LDH, V, LDV, ZERO, WV, LDWV )

    
CALL ZLACPY( aqAaq, 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( aqCaq, aqNaq, JW, KLN, JW, ONE, V, LDV, H( KWTOP, KCOL ), LDH, ZERO, T, LDT )

    
CALL ZLACPY( aqAaq, 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( aqNaq, aqNaq, KLN, JW, JW, ONE, Z( KROW, KWTOP ), LDZ, V, LDV, ZERO, WV, LDWV )

    
CALL ZLACPY( aqAaq, 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

PURPOSE