Annotation of rpl/lapack/lapack/zlargv.f, revision 1.1
1.1 ! bertrand 1: SUBROUTINE ZLARGV( N, X, INCX, Y, INCY, C, INCC )
! 2: *
! 3: * -- LAPACK auxiliary routine (version 3.2) --
! 4: * -- LAPACK is a software package provided by Univ. of Tennessee, --
! 5: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
! 6: * November 2006
! 7: *
! 8: * .. Scalar Arguments ..
! 9: INTEGER INCC, INCX, INCY, N
! 10: * ..
! 11: * .. Array Arguments ..
! 12: DOUBLE PRECISION C( * )
! 13: COMPLEX*16 X( * ), Y( * )
! 14: * ..
! 15: *
! 16: * Purpose
! 17: * =======
! 18: *
! 19: * ZLARGV generates a vector of complex plane rotations with real
! 20: * cosines, determined by elements of the complex vectors x and y.
! 21: * For i = 1,2,...,n
! 22: *
! 23: * ( c(i) s(i) ) ( x(i) ) = ( r(i) )
! 24: * ( -conjg(s(i)) c(i) ) ( y(i) ) = ( 0 )
! 25: *
! 26: * where c(i)**2 + ABS(s(i))**2 = 1
! 27: *
! 28: * The following conventions are used (these are the same as in ZLARTG,
! 29: * but differ from the BLAS1 routine ZROTG):
! 30: * If y(i)=0, then c(i)=1 and s(i)=0.
! 31: * If x(i)=0, then c(i)=0 and s(i) is chosen so that r(i) is real.
! 32: *
! 33: * Arguments
! 34: * =========
! 35: *
! 36: * N (input) INTEGER
! 37: * The number of plane rotations to be generated.
! 38: *
! 39: * X (input/output) COMPLEX*16 array, dimension (1+(N-1)*INCX)
! 40: * On entry, the vector x.
! 41: * On exit, x(i) is overwritten by r(i), for i = 1,...,n.
! 42: *
! 43: * INCX (input) INTEGER
! 44: * The increment between elements of X. INCX > 0.
! 45: *
! 46: * Y (input/output) COMPLEX*16 array, dimension (1+(N-1)*INCY)
! 47: * On entry, the vector y.
! 48: * On exit, the sines of the plane rotations.
! 49: *
! 50: * INCY (input) INTEGER
! 51: * The increment between elements of Y. INCY > 0.
! 52: *
! 53: * C (output) DOUBLE PRECISION array, dimension (1+(N-1)*INCC)
! 54: * The cosines of the plane rotations.
! 55: *
! 56: * INCC (input) INTEGER
! 57: * The increment between elements of C. INCC > 0.
! 58: *
! 59: * Further Details
! 60: * ======= =======
! 61: *
! 62: * 6-6-96 - Modified with a new algorithm by W. Kahan and J. Demmel
! 63: *
! 64: * This version has a few statements commented out for thread safety
! 65: * (machine parameters are computed on each entry). 10 feb 03, SJH.
! 66: *
! 67: * =====================================================================
! 68: *
! 69: * .. Parameters ..
! 70: DOUBLE PRECISION TWO, ONE, ZERO
! 71: PARAMETER ( TWO = 2.0D+0, ONE = 1.0D+0, ZERO = 0.0D+0 )
! 72: COMPLEX*16 CZERO
! 73: PARAMETER ( CZERO = ( 0.0D+0, 0.0D+0 ) )
! 74: * ..
! 75: * .. Local Scalars ..
! 76: * LOGICAL FIRST
! 77:
! 78: INTEGER COUNT, I, IC, IX, IY, J
! 79: DOUBLE PRECISION CS, D, DI, DR, EPS, F2, F2S, G2, G2S, SAFMIN,
! 80: $ SAFMN2, SAFMX2, SCALE
! 81: COMPLEX*16 F, FF, FS, G, GS, R, SN
! 82: * ..
! 83: * .. External Functions ..
! 84: DOUBLE PRECISION DLAMCH, DLAPY2
! 85: EXTERNAL DLAMCH, DLAPY2
! 86: * ..
! 87: * .. Intrinsic Functions ..
! 88: INTRINSIC ABS, DBLE, DCMPLX, DCONJG, DIMAG, INT, LOG,
! 89: $ MAX, SQRT
! 90: * ..
! 91: * .. Statement Functions ..
! 92: DOUBLE PRECISION ABS1, ABSSQ
! 93: * ..
! 94: * .. Save statement ..
! 95: * SAVE FIRST, SAFMX2, SAFMIN, SAFMN2
! 96: * ..
! 97: * .. Data statements ..
! 98: * DATA FIRST / .TRUE. /
! 99: * ..
! 100: * .. Statement Function definitions ..
! 101: ABS1( FF ) = MAX( ABS( DBLE( FF ) ), ABS( DIMAG( FF ) ) )
! 102: ABSSQ( FF ) = DBLE( FF )**2 + DIMAG( FF )**2
! 103: * ..
! 104: * .. Executable Statements ..
! 105: *
! 106: * IF( FIRST ) THEN
! 107: * FIRST = .FALSE.
! 108: SAFMIN = DLAMCH( 'S' )
! 109: EPS = DLAMCH( 'E' )
! 110: SAFMN2 = DLAMCH( 'B' )**INT( LOG( SAFMIN / EPS ) /
! 111: $ LOG( DLAMCH( 'B' ) ) / TWO )
! 112: SAFMX2 = ONE / SAFMN2
! 113: * END IF
! 114: IX = 1
! 115: IY = 1
! 116: IC = 1
! 117: DO 60 I = 1, N
! 118: F = X( IX )
! 119: G = Y( IY )
! 120: *
! 121: * Use identical algorithm as in ZLARTG
! 122: *
! 123: SCALE = MAX( ABS1( F ), ABS1( G ) )
! 124: FS = F
! 125: GS = G
! 126: COUNT = 0
! 127: IF( SCALE.GE.SAFMX2 ) THEN
! 128: 10 CONTINUE
! 129: COUNT = COUNT + 1
! 130: FS = FS*SAFMN2
! 131: GS = GS*SAFMN2
! 132: SCALE = SCALE*SAFMN2
! 133: IF( SCALE.GE.SAFMX2 )
! 134: $ GO TO 10
! 135: ELSE IF( SCALE.LE.SAFMN2 ) THEN
! 136: IF( G.EQ.CZERO ) THEN
! 137: CS = ONE
! 138: SN = CZERO
! 139: R = F
! 140: GO TO 50
! 141: END IF
! 142: 20 CONTINUE
! 143: COUNT = COUNT - 1
! 144: FS = FS*SAFMX2
! 145: GS = GS*SAFMX2
! 146: SCALE = SCALE*SAFMX2
! 147: IF( SCALE.LE.SAFMN2 )
! 148: $ GO TO 20
! 149: END IF
! 150: F2 = ABSSQ( FS )
! 151: G2 = ABSSQ( GS )
! 152: IF( F2.LE.MAX( G2, ONE )*SAFMIN ) THEN
! 153: *
! 154: * This is a rare case: F is very small.
! 155: *
! 156: IF( F.EQ.CZERO ) THEN
! 157: CS = ZERO
! 158: R = DLAPY2( DBLE( G ), DIMAG( G ) )
! 159: * Do complex/real division explicitly with two real
! 160: * divisions
! 161: D = DLAPY2( DBLE( GS ), DIMAG( GS ) )
! 162: SN = DCMPLX( DBLE( GS ) / D, -DIMAG( GS ) / D )
! 163: GO TO 50
! 164: END IF
! 165: F2S = DLAPY2( DBLE( FS ), DIMAG( FS ) )
! 166: * G2 and G2S are accurate
! 167: * G2 is at least SAFMIN, and G2S is at least SAFMN2
! 168: G2S = SQRT( G2 )
! 169: * Error in CS from underflow in F2S is at most
! 170: * UNFL / SAFMN2 .lt. sqrt(UNFL*EPS) .lt. EPS
! 171: * If MAX(G2,ONE)=G2, then F2 .lt. G2*SAFMIN,
! 172: * and so CS .lt. sqrt(SAFMIN)
! 173: * If MAX(G2,ONE)=ONE, then F2 .lt. SAFMIN
! 174: * and so CS .lt. sqrt(SAFMIN)/SAFMN2 = sqrt(EPS)
! 175: * Therefore, CS = F2S/G2S / sqrt( 1 + (F2S/G2S)**2 ) = F2S/G2S
! 176: CS = F2S / G2S
! 177: * Make sure abs(FF) = 1
! 178: * Do complex/real division explicitly with 2 real divisions
! 179: IF( ABS1( F ).GT.ONE ) THEN
! 180: D = DLAPY2( DBLE( F ), DIMAG( F ) )
! 181: FF = DCMPLX( DBLE( F ) / D, DIMAG( F ) / D )
! 182: ELSE
! 183: DR = SAFMX2*DBLE( F )
! 184: DI = SAFMX2*DIMAG( F )
! 185: D = DLAPY2( DR, DI )
! 186: FF = DCMPLX( DR / D, DI / D )
! 187: END IF
! 188: SN = FF*DCMPLX( DBLE( GS ) / G2S, -DIMAG( GS ) / G2S )
! 189: R = CS*F + SN*G
! 190: ELSE
! 191: *
! 192: * This is the most common case.
! 193: * Neither F2 nor F2/G2 are less than SAFMIN
! 194: * F2S cannot overflow, and it is accurate
! 195: *
! 196: F2S = SQRT( ONE+G2 / F2 )
! 197: * Do the F2S(real)*FS(complex) multiply with two real
! 198: * multiplies
! 199: R = DCMPLX( F2S*DBLE( FS ), F2S*DIMAG( FS ) )
! 200: CS = ONE / F2S
! 201: D = F2 + G2
! 202: * Do complex/real division explicitly with two real divisions
! 203: SN = DCMPLX( DBLE( R ) / D, DIMAG( R ) / D )
! 204: SN = SN*DCONJG( GS )
! 205: IF( COUNT.NE.0 ) THEN
! 206: IF( COUNT.GT.0 ) THEN
! 207: DO 30 J = 1, COUNT
! 208: R = R*SAFMX2
! 209: 30 CONTINUE
! 210: ELSE
! 211: DO 40 J = 1, -COUNT
! 212: R = R*SAFMN2
! 213: 40 CONTINUE
! 214: END IF
! 215: END IF
! 216: END IF
! 217: 50 CONTINUE
! 218: C( IC ) = CS
! 219: Y( IY ) = SN
! 220: X( IX ) = R
! 221: IC = IC + INCC
! 222: IY = IY + INCY
! 223: IX = IX + INCX
! 224: 60 CONTINUE
! 225: RETURN
! 226: *
! 227: * End of ZLARGV
! 228: *
! 229: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to zlargv.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack