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>