claqr3 (l) - Linux Manuals

NAME

SYNOPSIS

SUBROUTINE CLAQR3(

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 H( LDH, * ), SH( * ), T( LDT, * ), V( LDV, * ), WORK( * ), WV( LDWV, * ), Z( LDZ, * )

    

COMPLEX ZERO, ONE

    

PARAMETER ( ZERO = ( 0.0e0, 0.0e0 ), ONE = ( 1.0e0, 0.0e0 ) )

    

REAL RZERO, RONE

    

PARAMETER ( RZERO = 0.0e0, RONE = 1.0e0 )

    

COMPLEX BETA, CDUM, S, TAU

    

REAL FOO, SAFMAX, SAFMIN, SMLNUM, ULP

    

INTEGER I, IFST, ILST, INFO, INFQR, J, JW, KCOL, KLN, KNT, KROW, KWTOP, LTOP, LWK1, LWK2, LWK3, LWKOPT, NMIN

    

REAL SLAMCH

    

INTEGER ILAENV

    

EXTERNAL SLAMCH, ILAENV

    

EXTERNAL CCOPY, CGEHRD, CGEMM, CLACPY, CLAHQR, CLAQR4, CLARF, CLARFG, CLASET, CTREXC, CUNMHR, SLABAD

    

INTRINSIC ABS, AIMAG, CMPLX, CONJG, INT, MAX, MIN, REAL

    

REAL CABS1

    

CABS1( CDUM ) = ABS( REAL( CDUM ) ) + ABS( AIMAG( CDUM ) )

    

JW = MIN( NW, KBOT-KTOP+1 )

    

IF( JW.LE.2 ) THEN

    

LWKOPT = 1

    

ELSE

    

CALL CGEHRD( JW, 1, JW-1, T, LDT, WORK, WORK, -1, INFO )

    

LWK1 = INT( WORK( 1 ) )

    

CALL CUNMHR( aqRaq, aqNaq, JW, JW, 1, JW-1, T, LDT, WORK, V, LDV, WORK, -1, INFO )

    

LWK2 = INT( WORK( 1 ) )

    

CALL CLAQR4( .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 ) = CMPLX( LWKOPT, 0 )

    

RETURN

    

END IF

    

NS = 0

    

ND = 0

    

WORK( 1 ) = ONE

    

IF( KTOP.GT.KBOT ) RETURN

    

IF( NW.LT.1 ) RETURN

    

SAFMIN = SLAMCH( aqSAFE MINIMUMaq )

    

SAFMAX = RONE / SAFMIN

    

CALL SLABAD( SAFMIN, SAFMAX )

    

ULP = SLAMCH( aqPRECISIONaq )

    

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

    

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 CLACPY( aqUaq, JW, JW, H( KWTOP, KWTOP ), LDH, T, LDT )

    

CALL CCOPY( JW-1, H( KWTOP+1, KWTOP ), LDH+1, T( 2, 1 ), LDT+1 )

    

CALL CLASET( aqAaq, JW, JW, ZERO, ONE, V, LDV )

    

NMIN = ILAENV( 12, aqCLAQR3aq, aqSVaq, JW, 1, JW, LWORK )

    

IF( JW.GT.NMIN ) THEN

    

CALL CLAQR4( .true., .true., JW, 1, JW, T, LDT, SH( KWTOP ), 1, JW, V, LDV, WORK, LWORK, INFQR )

    

ELSE

    

CALL CLAHQR( .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 CTREXC( 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 CTREXC( 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 CCOPY( NS, V, LDV, WORK, 1 )

    

DO 50 I = 1, NS

    

WORK( I ) = CONJG( WORK( I ) )

    

50 CONTINUE

    

BETA = WORK( 1 )

    

CALL CLARFG( NS, BETA, WORK( 2 ), 1, TAU )

    

WORK( 1 ) = ONE

    

CALL CLASET( aqLaq, JW-2, JW-2, ZERO, ZERO, T( 3, 1 ), LDT )

    

CALL CLARF( aqLaq, NS, JW, WORK, 1, CONJG( TAU ), T, LDT, WORK( JW+1 ) )

    

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

    

CALL CLARF( aqRaq, JW, NS, WORK, 1, TAU, V, LDV, WORK( JW+1 ) )

    

CALL CGEHRD( JW, 1, NS, T, LDT, WORK, WORK( JW+1 ), LWORK-JW, INFO )

    

END IF

    

IF( KWTOP.GT.1 ) H( KWTOP, KWTOP-1 ) = S*CONJG( V( 1, 1 ) )

    

CALL CLACPY( aqUaq, JW, JW, T, LDT, H( KWTOP, KWTOP ), LDH )

    

CALL CCOPY( JW-1, T( 2, 1 ), LDT+1, H( KWTOP+1, KWTOP ), LDH+1 )

    

IF( NS.GT.1 .AND. S.NE.ZERO ) CALL CUNMHR( 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 CGEMM( aqNaq, aqNaq, KLN, JW, JW, ONE, H( KROW, KWTOP ), LDH, V, LDV, ZERO, WV, LDWV )

    

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

    

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

    

CALL CLACPY( aqAaq, KLN, JW, WV, LDWV, Z( KROW, KWTOP ), LDZ )

    

80 CONTINUE

    

END IF

    

END IF

    

ND = JW - NS

    

NS = NS - INFQR

    

WORK( 1 ) = CMPLX( LWKOPT, 0 )

    

END

PURPOSE