Annotation of rpl/lapack/lapack/zlartg.f, revision 1.17

1.11      bertrand    1: *> \brief \b ZLARTG generates a plane rotation with real cosine and complex sine.
1.8       bertrand    2: *
                      3: *  =========== DOCUMENTATION ===========
                      4: *
1.16      bertrand    5: * Online html documentation available at
                      6: *            http://www.netlib.org/lapack/explore-html/
1.8       bertrand    7: *
                      8: *> \htmlonly
1.16      bertrand    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">
1.8       bertrand   15: *> [TXT]</a>
1.16      bertrand   16: *> \endhtmlonly
1.8       bertrand   17: *
                     18: *  Definition:
                     19: *  ===========
                     20: *
                     21: *       SUBROUTINE ZLARTG( F, G, CS, SN, R )
1.16      bertrand   22: *
1.8       bertrand   23: *       .. Scalar Arguments ..
                     24: *       DOUBLE PRECISION   CS
                     25: *       COMPLEX*16         F, G, R, SN
                     26: *       ..
1.16      bertrand   27: *
1.8       bertrand   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: *
1.16      bertrand   83: *> \author Univ. of Tennessee
                     84: *> \author Univ. of California Berkeley
                     85: *> \author Univ. of Colorado Denver
                     86: *> \author NAG Ltd.
1.8       bertrand   87: *
1.16      bertrand   88: *> \date December 2016
1.8       bertrand   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.16      bertrand  106: *  -- LAPACK auxiliary routine (version 3.7.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.16      bertrand  109: *     December 2016
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
1.13      bertrand  133:       LOGICAL            DISNAN
                    134:       EXTERNAL           DLAMCH, DLAPY2, DISNAN
1.1       bertrand  135: *     ..
                    136: *     .. Intrinsic Functions ..
                    137:       INTRINSIC          ABS, DBLE, DCMPLX, DCONJG, DIMAG, INT, LOG,
                    138:      $                   MAX, SQRT
                    139: *     ..
                    140: *     .. Statement Functions ..
                    141:       DOUBLE PRECISION   ABS1, ABSSQ
                    142: *     ..
                    143: *     .. Statement Function definitions ..
                    144:       ABS1( FF ) = MAX( ABS( DBLE( FF ) ), ABS( DIMAG( FF ) ) )
                    145:       ABSSQ( FF ) = DBLE( FF )**2 + DIMAG( FF )**2
                    146: *     ..
                    147: *     .. Executable Statements ..
                    148: *
1.13      bertrand  149:       SAFMIN = DLAMCH( 'S' )
                    150:       EPS = DLAMCH( 'E' )
                    151:       SAFMN2 = DLAMCH( 'B' )**INT( LOG( SAFMIN / EPS ) /
                    152:      $         LOG( DLAMCH( 'B' ) ) / TWO )
                    153:       SAFMX2 = ONE / SAFMN2
1.1       bertrand  154:       SCALE = MAX( ABS1( F ), ABS1( G ) )
                    155:       FS = F
                    156:       GS = G
                    157:       COUNT = 0
                    158:       IF( SCALE.GE.SAFMX2 ) THEN
                    159:    10    CONTINUE
                    160:          COUNT = COUNT + 1
                    161:          FS = FS*SAFMN2
                    162:          GS = GS*SAFMN2
                    163:          SCALE = SCALE*SAFMN2
                    164:          IF( SCALE.GE.SAFMX2 )
                    165:      $      GO TO 10
                    166:       ELSE IF( SCALE.LE.SAFMN2 ) THEN
1.13      bertrand  167:          IF( G.EQ.CZERO.OR.DISNAN( ABS( G ) ) ) THEN
1.1       bertrand  168:             CS = ONE
                    169:             SN = CZERO
                    170:             R = F
                    171:             RETURN
                    172:          END IF
                    173:    20    CONTINUE
                    174:          COUNT = COUNT - 1
                    175:          FS = FS*SAFMX2
                    176:          GS = GS*SAFMX2
                    177:          SCALE = SCALE*SAFMX2
                    178:          IF( SCALE.LE.SAFMN2 )
                    179:      $      GO TO 20
                    180:       END IF
                    181:       F2 = ABSSQ( FS )
                    182:       G2 = ABSSQ( GS )
                    183:       IF( F2.LE.MAX( G2, ONE )*SAFMIN ) THEN
                    184: *
                    185: *        This is a rare case: F is very small.
                    186: *
                    187:          IF( F.EQ.CZERO ) THEN
                    188:             CS = ZERO
                    189:             R = DLAPY2( DBLE( G ), DIMAG( G ) )
                    190: *           Do complex/real division explicitly with two real divisions
                    191:             D = DLAPY2( DBLE( GS ), DIMAG( GS ) )
                    192:             SN = DCMPLX( DBLE( GS ) / D, -DIMAG( GS ) / D )
                    193:             RETURN
                    194:          END IF
                    195:          F2S = DLAPY2( DBLE( FS ), DIMAG( FS ) )
                    196: *        G2 and G2S are accurate
                    197: *        G2 is at least SAFMIN, and G2S is at least SAFMN2
                    198:          G2S = SQRT( G2 )
                    199: *        Error in CS from underflow in F2S is at most
                    200: *        UNFL / SAFMN2 .lt. sqrt(UNFL*EPS) .lt. EPS
                    201: *        If MAX(G2,ONE)=G2, then F2 .lt. G2*SAFMIN,
                    202: *        and so CS .lt. sqrt(SAFMIN)
                    203: *        If MAX(G2,ONE)=ONE, then F2 .lt. SAFMIN
                    204: *        and so CS .lt. sqrt(SAFMIN)/SAFMN2 = sqrt(EPS)
                    205: *        Therefore, CS = F2S/G2S / sqrt( 1 + (F2S/G2S)**2 ) = F2S/G2S
                    206:          CS = F2S / G2S
                    207: *        Make sure abs(FF) = 1
                    208: *        Do complex/real division explicitly with 2 real divisions
                    209:          IF( ABS1( F ).GT.ONE ) THEN
                    210:             D = DLAPY2( DBLE( F ), DIMAG( F ) )
                    211:             FF = DCMPLX( DBLE( F ) / D, DIMAG( F ) / D )
                    212:          ELSE
                    213:             DR = SAFMX2*DBLE( F )
                    214:             DI = SAFMX2*DIMAG( F )
                    215:             D = DLAPY2( DR, DI )
                    216:             FF = DCMPLX( DR / D, DI / D )
                    217:          END IF
                    218:          SN = FF*DCMPLX( DBLE( GS ) / G2S, -DIMAG( GS ) / G2S )
                    219:          R = CS*F + SN*G
                    220:       ELSE
                    221: *
                    222: *        This is the most common case.
                    223: *        Neither F2 nor F2/G2 are less than SAFMIN
                    224: *        F2S cannot overflow, and it is accurate
                    225: *
                    226:          F2S = SQRT( ONE+G2 / F2 )
                    227: *        Do the F2S(real)*FS(complex) multiply with two real multiplies
                    228:          R = DCMPLX( F2S*DBLE( FS ), F2S*DIMAG( FS ) )
                    229:          CS = ONE / F2S
                    230:          D = F2 + G2
                    231: *        Do complex/real division explicitly with two real divisions
                    232:          SN = DCMPLX( DBLE( R ) / D, DIMAG( R ) / D )
                    233:          SN = SN*DCONJG( GS )
                    234:          IF( COUNT.NE.0 ) THEN
                    235:             IF( COUNT.GT.0 ) THEN
                    236:                DO 30 I = 1, COUNT
                    237:                   R = R*SAFMX2
                    238:    30          CONTINUE
                    239:             ELSE
                    240:                DO 40 I = 1, -COUNT
                    241:                   R = R*SAFMN2
                    242:    40          CONTINUE
                    243:             END IF
                    244:          END IF
                    245:       END IF
                    246:       RETURN
                    247: *
                    248: *     End of ZLARTG
                    249: *
                    250:       END

CVSweb interface <joel.bertrand@systella.fr>