SUBROUTINE SLARTG( F, G, CS, SN, R ) * * -- LAPACK auxiliary routine (version 3.0) -- * Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., * Courant Institute, Argonne National Lab, and Rice University * Jan 17, 2001 * * .. Scalar Arguments .. REAL CS, F, G, R, SN * .. * * Purpose * ======= * * SLARTG generate a plane rotation so that * * [ CS SN ] . [ F ] = [ R ] where CS**2 + SN**2 = 1. * [ -SN CS ] [ G ] [ 0 ] * * This is a slower, more accurate version of the BLAS1 routine SROTG, * with the following other differences: * F and G are unchanged on return. * If F=0 and G=0, then CS=1, SN=0, and R=0. * If F .ne. 0 and G=0, then CS=1, SN=0, and R=F. * If F=0 and G .ne. 0, then CS=0, SN=sign(G), and R=abs(G). * If F .ne. 0 and (G .ne. 0), then * CS = abs(F)/sqrt(F**2 + G**2) * SN = sign(F)*G/sqrt(F**2 + G**2) * R = sign(F)*sqrt(F**2 + G**2) * * If a NaN occurs in the input, then R and possibly also * CS and SN will be a NaN. * * If an infinity occurs in the input, then R and possibly also * CS and SN will be infinite or NaNs. * * The complex routine CLARTG returns the same * CS and SN on complex inputs (F,0) and (G,0). * * Arguments * ========= * * F (input) REAL * The first component of vector to be rotated. * * G (input) REAL * The second component of vector to be rotated. * * CS (output) REAL * The cosine of the rotation. * * SN (output) REAL * The sine of the rotation. * * R (output) REAL * The nonzero component of the rotated vector. * * ===================================================================== * * .. Parameters .. REAL ZERO PARAMETER ( ZERO = 0.0E0 ) REAL ONE PARAMETER ( ONE = 1.0E0 ) REAL TWO PARAMETER ( TWO = 2.0E0 ) * .. * .. Local Scalars .. LOGICAL FIRST INTEGER COUNT, I, MAXCNT REAL EPS, F1, G1, SAFMIN, SAFMN2, SAFMX2, SCALE REAL SAFMN, SAFMX, SCL, ESFMN2, RR * .. * .. External Functions .. REAL SLAMCH EXTERNAL SLAMCH * .. * .. Intrinsic Functions .. INTRINSIC ABS, INT, LOG, MAX, SQRT, SIGN * .. * .. Save statement .. SAVE FIRST, EPS, SAFMX2, SAFMIN, SAFMN2, SAFMN SAVE SAFMX * .. * .. Data statements .. DATA FIRST / .TRUE. / * .. * .. Executable Statements .. * IF( FIRST ) THEN * * On first call to SLARTG, compute * SAFMN2 = sqrt(SAFMIN/EPS) rounded down to the nearest power * of the floating point radix * This means that scaling by multiplication by SAFMN2 and its * reciprocal SAFMX2 cause no roundoff error * FIRST = .FALSE. SAFMIN = SLAMCH( 'S' ) EPS = SLAMCH( 'E' ) ESFMN2 = INT( LOG( SAFMIN / EPS ) / LOG( SLAMCH( 'B' ) ) / TWO) SAFMN2 = SLAMCH( 'B' )**ESFMN2 SAFMN = SAFMN2**2 SAFMX2 = ONE / SAFMN2 SAFMX = SAFMX2**2 MAXCNT = INT( -MAX( SLAMCH('L'), SLAMCH('N')+ONE-SLAMCH('M') ) + /ESFMN2 - ONE - ONE ) END IF IF( G.EQ.ZERO ) THEN * * Includes the case F=G=0 * CS = ONE SN = ZERO R = F ELSE IF( F.EQ.ZERO ) THEN * * G must be nonzero * CS = ZERO SN = SIGN( ONE, G ) R = ABS(G) ELSE * * Both F and G must be nonzero * F1 = F G1 = G SCALE = MAX( ABS( F1 ), ABS( G1 ) ) COUNT = 0 IF( SCALE.GE.SAFMX2 ) THEN * * Handle case where F1**2 + G1**2 might overflow * SCL = SAFMX2 * COUNT = COUNT + 1 F1 = F1*SAFMN2 G1 = G1*SAFMN2 SCALE = SCALE*SAFMN2 IF( SCALE.LT.SAFMX2 ) GOTO 100 * COUNT = COUNT + 1 F1 = F1*SAFMN2 G1 = G1*SAFMN2 SCALE = SCALE*SAFMN2 IF( SCALE.LT.SAFMX2 ) GOTO 100 * 10 CONTINUE COUNT = COUNT + 1 F1 = F1*SAFMN2 G1 = G1*SAFMN2 SCALE = SCALE*SAFMN2 IF( SCALE.GT.SAFMX2 .AND. COUNT.LE.MAXCNT ) GOTO 10 ELSE IF( SCALE.LE.SAFMN2 ) THEN * * Handle case where F1**2 + G1**2 might underflow * SCL = SAFMN2 * COUNT = COUNT + 1 F1 = F1*SAFMX2 G1 = G1*SAFMX2 SCALE = SCALE*SAFMX2 IF( SCALE.GT.SAFMN2 ) GOTO 100 * COUNT = COUNT + 1 F1 = F1*SAFMX2 G1 = G1*SAFMX2 SCALE = SCALE*SAFMX2 IF( SCALE.GT.SAFMN2 ) GOTO 100 * 20 CONTINUE COUNT = COUNT + 1 F1 = F1*SAFMX2 G1 = G1*SAFMX2 SCALE = SCALE*SAFMX2 IF( SCALE.LT.SAFMN2 .AND. COUNT.LE.MAXCNT ) GOTO 20 ENDIF 100 CONTINUE R = SQRT( F1**2+G1**2 ) RR = ONE/R CS = ABS(F1) * RR SN = G1 * RR IF (F .LT. ZERO) THEN R = -R SN = -SN ENDIF DO 40 I = 1, COUNT R = R*SCL 40 CONTINUE ENDIF RETURN * * End of SLARTG * END