Annotation of rpl/lapack/lapack/zlartg.f, revision 1.8
1.8 ! bertrand 1: *> \brief \b ZLARTG
! 2: *
! 3: * =========== DOCUMENTATION ===========
! 4: *
! 5: * Online html documentation available at
! 6: * http://www.netlib.org/lapack/explore-html/
! 7: *
! 8: *> \htmlonly
! 9: *> Download ZLARTG + dependencies
! 10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zlartg.f">
! 11: *> [TGZ]</a>
! 12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zlartg.f">
! 13: *> [ZIP]</a>
! 14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zlartg.f">
! 15: *> [TXT]</a>
! 16: *> \endhtmlonly
! 17: *
! 18: * Definition:
! 19: * ===========
! 20: *
! 21: * SUBROUTINE ZLARTG( F, G, CS, SN, R )
! 22: *
! 23: * .. Scalar Arguments ..
! 24: * DOUBLE PRECISION CS
! 25: * COMPLEX*16 F, G, R, SN
! 26: * ..
! 27: *
! 28: *
! 29: *> \par Purpose:
! 30: * =============
! 31: *>
! 32: *> \verbatim
! 33: *>
! 34: *> ZLARTG generates a plane rotation so that
! 35: *>
! 36: *> [ CS SN ] [ F ] [ R ]
! 37: *> [ __ ] . [ ] = [ ] where CS**2 + |SN|**2 = 1.
! 38: *> [ -SN CS ] [ G ] [ 0 ]
! 39: *>
! 40: *> This is a faster version of the BLAS1 routine ZROTG, except for
! 41: *> the following differences:
! 42: *> F and G are unchanged on return.
! 43: *> If G=0, then CS=1 and SN=0.
! 44: *> If F=0, then CS=0 and SN is chosen so that R is real.
! 45: *> \endverbatim
! 46: *
! 47: * Arguments:
! 48: * ==========
! 49: *
! 50: *> \param[in] F
! 51: *> \verbatim
! 52: *> F is COMPLEX*16
! 53: *> The first component of vector to be rotated.
! 54: *> \endverbatim
! 55: *>
! 56: *> \param[in] G
! 57: *> \verbatim
! 58: *> G is COMPLEX*16
! 59: *> The second component of vector to be rotated.
! 60: *> \endverbatim
! 61: *>
! 62: *> \param[out] CS
! 63: *> \verbatim
! 64: *> CS is DOUBLE PRECISION
! 65: *> The cosine of the rotation.
! 66: *> \endverbatim
! 67: *>
! 68: *> \param[out] SN
! 69: *> \verbatim
! 70: *> SN is COMPLEX*16
! 71: *> The sine of the rotation.
! 72: *> \endverbatim
! 73: *>
! 74: *> \param[out] R
! 75: *> \verbatim
! 76: *> R is COMPLEX*16
! 77: *> The nonzero component of the rotated vector.
! 78: *> \endverbatim
! 79: *
! 80: * Authors:
! 81: * ========
! 82: *
! 83: *> \author Univ. of Tennessee
! 84: *> \author Univ. of California Berkeley
! 85: *> \author Univ. of Colorado Denver
! 86: *> \author NAG Ltd.
! 87: *
! 88: *> \date November 2011
! 89: *
! 90: *> \ingroup complex16OTHERauxiliary
! 91: *
! 92: *> \par Further Details:
! 93: * =====================
! 94: *>
! 95: *> \verbatim
! 96: *>
! 97: *> 3-5-96 - Modified with a new algorithm by W. Kahan and J. Demmel
! 98: *>
! 99: *> This version has a few statements commented out for thread safety
! 100: *> (machine parameters are computed on each entry). 10 feb 03, SJH.
! 101: *> \endverbatim
! 102: *>
! 103: * =====================================================================
1.1 bertrand 104: SUBROUTINE ZLARTG( F, G, CS, SN, R )
105: *
1.8 ! bertrand 106: * -- LAPACK auxiliary routine (version 3.4.0) --
1.1 bertrand 107: * -- LAPACK is a software package provided by Univ. of Tennessee, --
108: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
1.8 ! bertrand 109: * November 2011
1.1 bertrand 110: *
111: * .. Scalar Arguments ..
112: DOUBLE PRECISION CS
113: COMPLEX*16 F, G, R, SN
114: * ..
115: *
116: * =====================================================================
117: *
118: * .. Parameters ..
119: DOUBLE PRECISION TWO, ONE, ZERO
120: PARAMETER ( TWO = 2.0D+0, ONE = 1.0D+0, ZERO = 0.0D+0 )
121: COMPLEX*16 CZERO
122: PARAMETER ( CZERO = ( 0.0D+0, 0.0D+0 ) )
123: * ..
124: * .. Local Scalars ..
125: * LOGICAL FIRST
126: INTEGER COUNT, I
127: DOUBLE PRECISION D, DI, DR, EPS, F2, F2S, G2, G2S, SAFMIN,
128: $ SAFMN2, SAFMX2, SCALE
129: COMPLEX*16 FF, FS, GS
130: * ..
131: * .. External Functions ..
132: DOUBLE PRECISION DLAMCH, DLAPY2
133: EXTERNAL DLAMCH, DLAPY2
134: * ..
135: * .. Intrinsic Functions ..
136: INTRINSIC ABS, DBLE, DCMPLX, DCONJG, DIMAG, INT, LOG,
137: $ MAX, SQRT
138: * ..
139: * .. Statement Functions ..
140: DOUBLE PRECISION ABS1, ABSSQ
141: * ..
142: * .. Save statement ..
143: * SAVE FIRST, SAFMX2, SAFMIN, SAFMN2
144: * ..
145: * .. Data statements ..
146: * DATA FIRST / .TRUE. /
147: * ..
148: * .. Statement Function definitions ..
149: ABS1( FF ) = MAX( ABS( DBLE( FF ) ), ABS( DIMAG( FF ) ) )
150: ABSSQ( FF ) = DBLE( FF )**2 + DIMAG( FF )**2
151: * ..
152: * .. Executable Statements ..
153: *
154: * IF( FIRST ) THEN
155: SAFMIN = DLAMCH( 'S' )
156: EPS = DLAMCH( 'E' )
157: SAFMN2 = DLAMCH( 'B' )**INT( LOG( SAFMIN / EPS ) /
158: $ LOG( DLAMCH( 'B' ) ) / TWO )
159: SAFMX2 = ONE / SAFMN2
160: * FIRST = .FALSE.
161: * END IF
162: SCALE = MAX( ABS1( F ), ABS1( G ) )
163: FS = F
164: GS = G
165: COUNT = 0
166: IF( SCALE.GE.SAFMX2 ) THEN
167: 10 CONTINUE
168: COUNT = COUNT + 1
169: FS = FS*SAFMN2
170: GS = GS*SAFMN2
171: SCALE = SCALE*SAFMN2
172: IF( SCALE.GE.SAFMX2 )
173: $ GO TO 10
174: ELSE IF( SCALE.LE.SAFMN2 ) THEN
175: IF( G.EQ.CZERO ) THEN
176: CS = ONE
177: SN = CZERO
178: R = F
179: RETURN
180: END IF
181: 20 CONTINUE
182: COUNT = COUNT - 1
183: FS = FS*SAFMX2
184: GS = GS*SAFMX2
185: SCALE = SCALE*SAFMX2
186: IF( SCALE.LE.SAFMN2 )
187: $ GO TO 20
188: END IF
189: F2 = ABSSQ( FS )
190: G2 = ABSSQ( GS )
191: IF( F2.LE.MAX( G2, ONE )*SAFMIN ) THEN
192: *
193: * This is a rare case: F is very small.
194: *
195: IF( F.EQ.CZERO ) THEN
196: CS = ZERO
197: R = DLAPY2( DBLE( G ), DIMAG( G ) )
198: * Do complex/real division explicitly with two real divisions
199: D = DLAPY2( DBLE( GS ), DIMAG( GS ) )
200: SN = DCMPLX( DBLE( GS ) / D, -DIMAG( GS ) / D )
201: RETURN
202: END IF
203: F2S = DLAPY2( DBLE( FS ), DIMAG( FS ) )
204: * G2 and G2S are accurate
205: * G2 is at least SAFMIN, and G2S is at least SAFMN2
206: G2S = SQRT( G2 )
207: * Error in CS from underflow in F2S is at most
208: * UNFL / SAFMN2 .lt. sqrt(UNFL*EPS) .lt. EPS
209: * If MAX(G2,ONE)=G2, then F2 .lt. G2*SAFMIN,
210: * and so CS .lt. sqrt(SAFMIN)
211: * If MAX(G2,ONE)=ONE, then F2 .lt. SAFMIN
212: * and so CS .lt. sqrt(SAFMIN)/SAFMN2 = sqrt(EPS)
213: * Therefore, CS = F2S/G2S / sqrt( 1 + (F2S/G2S)**2 ) = F2S/G2S
214: CS = F2S / G2S
215: * Make sure abs(FF) = 1
216: * Do complex/real division explicitly with 2 real divisions
217: IF( ABS1( F ).GT.ONE ) THEN
218: D = DLAPY2( DBLE( F ), DIMAG( F ) )
219: FF = DCMPLX( DBLE( F ) / D, DIMAG( F ) / D )
220: ELSE
221: DR = SAFMX2*DBLE( F )
222: DI = SAFMX2*DIMAG( F )
223: D = DLAPY2( DR, DI )
224: FF = DCMPLX( DR / D, DI / D )
225: END IF
226: SN = FF*DCMPLX( DBLE( GS ) / G2S, -DIMAG( GS ) / G2S )
227: R = CS*F + SN*G
228: ELSE
229: *
230: * This is the most common case.
231: * Neither F2 nor F2/G2 are less than SAFMIN
232: * F2S cannot overflow, and it is accurate
233: *
234: F2S = SQRT( ONE+G2 / F2 )
235: * Do the F2S(real)*FS(complex) multiply with two real multiplies
236: R = DCMPLX( F2S*DBLE( FS ), F2S*DIMAG( FS ) )
237: CS = ONE / F2S
238: D = F2 + G2
239: * Do complex/real division explicitly with two real divisions
240: SN = DCMPLX( DBLE( R ) / D, DIMAG( R ) / D )
241: SN = SN*DCONJG( GS )
242: IF( COUNT.NE.0 ) THEN
243: IF( COUNT.GT.0 ) THEN
244: DO 30 I = 1, COUNT
245: R = R*SAFMX2
246: 30 CONTINUE
247: ELSE
248: DO 40 I = 1, -COUNT
249: R = R*SAFMN2
250: 40 CONTINUE
251: END IF
252: END IF
253: END IF
254: RETURN
255: *
256: * End of ZLARTG
257: *
258: END
CVSweb interface <joel.bertrand@systella.fr>