![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b CLARTG 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 *> \htmlonly 00009 *> Download CLARTG + dependencies 00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/clartg.f"> 00011 *> [TGZ]</a> 00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/clartg.f"> 00013 *> [ZIP]</a> 00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/clartg.f"> 00015 *> [TXT]</a> 00016 *> \endhtmlonly 00017 * 00018 * Definition: 00019 * =========== 00020 * 00021 * SUBROUTINE CLARTG( F, G, CS, SN, R ) 00022 * 00023 * .. Scalar Arguments .. 00024 * REAL CS 00025 * COMPLEX F, G, R, SN 00026 * .. 00027 * 00028 * 00029 *> \par Purpose: 00030 * ============= 00031 *> 00032 *> \verbatim 00033 *> 00034 *> CLARTG generates a plane rotation so that 00035 *> 00036 *> [ CS SN ] [ F ] [ R ] 00037 *> [ __ ] . [ ] = [ ] where CS**2 + |SN|**2 = 1. 00038 *> [ -SN CS ] [ G ] [ 0 ] 00039 *> 00040 *> This is a faster version of the BLAS1 routine CROTG, except for 00041 *> the following differences: 00042 *> F and G are unchanged on return. 00043 *> If G=0, then CS=1 and SN=0. 00044 *> If F=0, then CS=0 and SN is chosen so that R is real. 00045 *> \endverbatim 00046 * 00047 * Arguments: 00048 * ========== 00049 * 00050 *> \param[in] F 00051 *> \verbatim 00052 *> F is COMPLEX 00053 *> The first component of vector to be rotated. 00054 *> \endverbatim 00055 *> 00056 *> \param[in] G 00057 *> \verbatim 00058 *> G is COMPLEX 00059 *> The second component of vector to be rotated. 00060 *> \endverbatim 00061 *> 00062 *> \param[out] CS 00063 *> \verbatim 00064 *> CS is REAL 00065 *> The cosine of the rotation. 00066 *> \endverbatim 00067 *> 00068 *> \param[out] SN 00069 *> \verbatim 00070 *> SN is COMPLEX 00071 *> The sine of the rotation. 00072 *> \endverbatim 00073 *> 00074 *> \param[out] R 00075 *> \verbatim 00076 *> R is COMPLEX 00077 *> The nonzero component of the rotated vector. 00078 *> \endverbatim 00079 * 00080 * Authors: 00081 * ======== 00082 * 00083 *> \author Univ. of Tennessee 00084 *> \author Univ. of California Berkeley 00085 *> \author Univ. of Colorado Denver 00086 *> \author NAG Ltd. 00087 * 00088 *> \date November 2011 00089 * 00090 *> \ingroup complexOTHERauxiliary 00091 * 00092 *> \par Further Details: 00093 * ===================== 00094 *> 00095 *> \verbatim 00096 *> 00097 *> 3-5-96 - Modified with a new algorithm by W. Kahan and J. Demmel 00098 *> 00099 *> This version has a few statements commented out for thread safety 00100 *> (machine parameters are computed on each entry). 10 feb 03, SJH. 00101 *> \endverbatim 00102 *> 00103 * ===================================================================== 00104 SUBROUTINE CLARTG( F, G, CS, SN, R ) 00105 * 00106 * -- LAPACK auxiliary routine (version 3.4.0) -- 00107 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00108 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00109 * November 2011 00110 * 00111 * .. Scalar Arguments .. 00112 REAL CS 00113 COMPLEX F, G, R, SN 00114 * .. 00115 * 00116 * ===================================================================== 00117 * 00118 * .. Parameters .. 00119 REAL TWO, ONE, ZERO 00120 PARAMETER ( TWO = 2.0E+0, ONE = 1.0E+0, ZERO = 0.0E+0 ) 00121 COMPLEX CZERO 00122 PARAMETER ( CZERO = ( 0.0E+0, 0.0E+0 ) ) 00123 * .. 00124 * .. Local Scalars .. 00125 * LOGICAL FIRST 00126 INTEGER COUNT, I 00127 REAL D, DI, DR, EPS, F2, F2S, G2, G2S, SAFMIN, 00128 $ SAFMN2, SAFMX2, SCALE 00129 COMPLEX FF, FS, GS 00130 * .. 00131 * .. External Functions .. 00132 REAL SLAMCH, SLAPY2 00133 EXTERNAL SLAMCH, SLAPY2 00134 * .. 00135 * .. Intrinsic Functions .. 00136 INTRINSIC ABS, AIMAG, CMPLX, CONJG, INT, LOG, MAX, REAL, 00137 $ SQRT 00138 * .. 00139 * .. Statement Functions .. 00140 REAL ABS1, ABSSQ 00141 * .. 00142 * .. Save statement .. 00143 * SAVE FIRST, SAFMX2, SAFMIN, SAFMN2 00144 * .. 00145 * .. Data statements .. 00146 * DATA FIRST / .TRUE. / 00147 * .. 00148 * .. Statement Function definitions .. 00149 ABS1( FF ) = MAX( ABS( REAL( FF ) ), ABS( AIMAG( FF ) ) ) 00150 ABSSQ( FF ) = REAL( FF )**2 + AIMAG( FF )**2 00151 * .. 00152 * .. Executable Statements .. 00153 * 00154 * IF( FIRST ) THEN 00155 SAFMIN = SLAMCH( 'S' ) 00156 EPS = SLAMCH( 'E' ) 00157 SAFMN2 = SLAMCH( 'B' )**INT( LOG( SAFMIN / EPS ) / 00158 $ LOG( SLAMCH( 'B' ) ) / TWO ) 00159 SAFMX2 = ONE / SAFMN2 00160 * FIRST = .FALSE. 00161 * END IF 00162 SCALE = MAX( ABS1( F ), ABS1( G ) ) 00163 FS = F 00164 GS = G 00165 COUNT = 0 00166 IF( SCALE.GE.SAFMX2 ) THEN 00167 10 CONTINUE 00168 COUNT = COUNT + 1 00169 FS = FS*SAFMN2 00170 GS = GS*SAFMN2 00171 SCALE = SCALE*SAFMN2 00172 IF( SCALE.GE.SAFMX2 ) 00173 $ GO TO 10 00174 ELSE IF( SCALE.LE.SAFMN2 ) THEN 00175 IF( G.EQ.CZERO ) THEN 00176 CS = ONE 00177 SN = CZERO 00178 R = F 00179 RETURN 00180 END IF 00181 20 CONTINUE 00182 COUNT = COUNT - 1 00183 FS = FS*SAFMX2 00184 GS = GS*SAFMX2 00185 SCALE = SCALE*SAFMX2 00186 IF( SCALE.LE.SAFMN2 ) 00187 $ GO TO 20 00188 END IF 00189 F2 = ABSSQ( FS ) 00190 G2 = ABSSQ( GS ) 00191 IF( F2.LE.MAX( G2, ONE )*SAFMIN ) THEN 00192 * 00193 * This is a rare case: F is very small. 00194 * 00195 IF( F.EQ.CZERO ) THEN 00196 CS = ZERO 00197 R = SLAPY2( REAL( G ), AIMAG( G ) ) 00198 * Do complex/real division explicitly with two real divisions 00199 D = SLAPY2( REAL( GS ), AIMAG( GS ) ) 00200 SN = CMPLX( REAL( GS ) / D, -AIMAG( GS ) / D ) 00201 RETURN 00202 END IF 00203 F2S = SLAPY2( REAL( FS ), AIMAG( FS ) ) 00204 * G2 and G2S are accurate 00205 * G2 is at least SAFMIN, and G2S is at least SAFMN2 00206 G2S = SQRT( G2 ) 00207 * Error in CS from underflow in F2S is at most 00208 * UNFL / SAFMN2 .lt. sqrt(UNFL*EPS) .lt. EPS 00209 * If MAX(G2,ONE)=G2, then F2 .lt. G2*SAFMIN, 00210 * and so CS .lt. sqrt(SAFMIN) 00211 * If MAX(G2,ONE)=ONE, then F2 .lt. SAFMIN 00212 * and so CS .lt. sqrt(SAFMIN)/SAFMN2 = sqrt(EPS) 00213 * Therefore, CS = F2S/G2S / sqrt( 1 + (F2S/G2S)**2 ) = F2S/G2S 00214 CS = F2S / G2S 00215 * Make sure abs(FF) = 1 00216 * Do complex/real division explicitly with 2 real divisions 00217 IF( ABS1( F ).GT.ONE ) THEN 00218 D = SLAPY2( REAL( F ), AIMAG( F ) ) 00219 FF = CMPLX( REAL( F ) / D, AIMAG( F ) / D ) 00220 ELSE 00221 DR = SAFMX2*REAL( F ) 00222 DI = SAFMX2*AIMAG( F ) 00223 D = SLAPY2( DR, DI ) 00224 FF = CMPLX( DR / D, DI / D ) 00225 END IF 00226 SN = FF*CMPLX( REAL( GS ) / G2S, -AIMAG( GS ) / G2S ) 00227 R = CS*F + SN*G 00228 ELSE 00229 * 00230 * This is the most common case. 00231 * Neither F2 nor F2/G2 are less than SAFMIN 00232 * F2S cannot overflow, and it is accurate 00233 * 00234 F2S = SQRT( ONE+G2 / F2 ) 00235 * Do the F2S(real)*FS(complex) multiply with two real multiplies 00236 R = CMPLX( F2S*REAL( FS ), F2S*AIMAG( FS ) ) 00237 CS = ONE / F2S 00238 D = F2 + G2 00239 * Do complex/real division explicitly with two real divisions 00240 SN = CMPLX( REAL( R ) / D, AIMAG( R ) / D ) 00241 SN = SN*CONJG( GS ) 00242 IF( COUNT.NE.0 ) THEN 00243 IF( COUNT.GT.0 ) THEN 00244 DO 30 I = 1, COUNT 00245 R = R*SAFMX2 00246 30 CONTINUE 00247 ELSE 00248 DO 40 I = 1, -COUNT 00249 R = R*SAFMN2 00250 40 CONTINUE 00251 END IF 00252 END IF 00253 END IF 00254 RETURN 00255 * 00256 * End of CLARTG 00257 * 00258 END