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

1.1     ! bertrand    1:       SUBROUTINE ZLARTG( F, G, CS, SN, R )
        !             2: *
        !             3: *  -- LAPACK auxiliary routine (version 3.2) --
        !             4: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
        !             5: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
        !             6: *     November 2006
        !             7: *
        !             8: *     .. Scalar Arguments ..
        !             9:       DOUBLE PRECISION   CS
        !            10:       COMPLEX*16         F, G, R, SN
        !            11: *     ..
        !            12: *
        !            13: *  Purpose
        !            14: *  =======
        !            15: *
        !            16: *  ZLARTG generates a plane rotation so that
        !            17: *
        !            18: *     [  CS  SN  ]     [ F ]     [ R ]
        !            19: *     [  __      ]  .  [   ]  =  [   ]   where CS**2 + |SN|**2 = 1.
        !            20: *     [ -SN  CS  ]     [ G ]     [ 0 ]
        !            21: *
        !            22: *  This is a faster version of the BLAS1 routine ZROTG, except for
        !            23: *  the following differences:
        !            24: *     F and G are unchanged on return.
        !            25: *     If G=0, then CS=1 and SN=0.
        !            26: *     If F=0, then CS=0 and SN is chosen so that R is real.
        !            27: *
        !            28: *  Arguments
        !            29: *  =========
        !            30: *
        !            31: *  F       (input) COMPLEX*16
        !            32: *          The first component of vector to be rotated.
        !            33: *
        !            34: *  G       (input) COMPLEX*16
        !            35: *          The second component of vector to be rotated.
        !            36: *
        !            37: *  CS      (output) DOUBLE PRECISION
        !            38: *          The cosine of the rotation.
        !            39: *
        !            40: *  SN      (output) COMPLEX*16
        !            41: *          The sine of the rotation.
        !            42: *
        !            43: *  R       (output) COMPLEX*16
        !            44: *          The nonzero component of the rotated vector.
        !            45: *
        !            46: *  Further Details
        !            47: *  ======= =======
        !            48: *
        !            49: *  3-5-96 - Modified with a new algorithm by W. Kahan and J. Demmel
        !            50: *
        !            51: *  This version has a few statements commented out for thread safety
        !            52: *  (machine parameters are computed on each entry). 10 feb 03, SJH.
        !            53: *
        !            54: *  =====================================================================
        !            55: *
        !            56: *     .. Parameters ..
        !            57:       DOUBLE PRECISION   TWO, ONE, ZERO
        !            58:       PARAMETER          ( TWO = 2.0D+0, ONE = 1.0D+0, ZERO = 0.0D+0 )
        !            59:       COMPLEX*16         CZERO
        !            60:       PARAMETER          ( CZERO = ( 0.0D+0, 0.0D+0 ) )
        !            61: *     ..
        !            62: *     .. Local Scalars ..
        !            63: *     LOGICAL            FIRST
        !            64:       INTEGER            COUNT, I
        !            65:       DOUBLE PRECISION   D, DI, DR, EPS, F2, F2S, G2, G2S, SAFMIN,
        !            66:      $                   SAFMN2, SAFMX2, SCALE
        !            67:       COMPLEX*16         FF, FS, GS
        !            68: *     ..
        !            69: *     .. External Functions ..
        !            70:       DOUBLE PRECISION   DLAMCH, DLAPY2
        !            71:       EXTERNAL           DLAMCH, DLAPY2
        !            72: *     ..
        !            73: *     .. Intrinsic Functions ..
        !            74:       INTRINSIC          ABS, DBLE, DCMPLX, DCONJG, DIMAG, INT, LOG,
        !            75:      $                   MAX, SQRT
        !            76: *     ..
        !            77: *     .. Statement Functions ..
        !            78:       DOUBLE PRECISION   ABS1, ABSSQ
        !            79: *     ..
        !            80: *     .. Save statement ..
        !            81: *     SAVE               FIRST, SAFMX2, SAFMIN, SAFMN2
        !            82: *     ..
        !            83: *     .. Data statements ..
        !            84: *     DATA               FIRST / .TRUE. /
        !            85: *     ..
        !            86: *     .. Statement Function definitions ..
        !            87:       ABS1( FF ) = MAX( ABS( DBLE( FF ) ), ABS( DIMAG( FF ) ) )
        !            88:       ABSSQ( FF ) = DBLE( FF )**2 + DIMAG( FF )**2
        !            89: *     ..
        !            90: *     .. Executable Statements ..
        !            91: *
        !            92: *     IF( FIRST ) THEN
        !            93:          SAFMIN = DLAMCH( 'S' )
        !            94:          EPS = DLAMCH( 'E' )
        !            95:          SAFMN2 = DLAMCH( 'B' )**INT( LOG( SAFMIN / EPS ) /
        !            96:      $            LOG( DLAMCH( 'B' ) ) / TWO )
        !            97:          SAFMX2 = ONE / SAFMN2
        !            98: *        FIRST = .FALSE.
        !            99: *     END IF
        !           100:       SCALE = MAX( ABS1( F ), ABS1( G ) )
        !           101:       FS = F
        !           102:       GS = G
        !           103:       COUNT = 0
        !           104:       IF( SCALE.GE.SAFMX2 ) THEN
        !           105:    10    CONTINUE
        !           106:          COUNT = COUNT + 1
        !           107:          FS = FS*SAFMN2
        !           108:          GS = GS*SAFMN2
        !           109:          SCALE = SCALE*SAFMN2
        !           110:          IF( SCALE.GE.SAFMX2 )
        !           111:      $      GO TO 10
        !           112:       ELSE IF( SCALE.LE.SAFMN2 ) THEN
        !           113:          IF( G.EQ.CZERO ) THEN
        !           114:             CS = ONE
        !           115:             SN = CZERO
        !           116:             R = F
        !           117:             RETURN
        !           118:          END IF
        !           119:    20    CONTINUE
        !           120:          COUNT = COUNT - 1
        !           121:          FS = FS*SAFMX2
        !           122:          GS = GS*SAFMX2
        !           123:          SCALE = SCALE*SAFMX2
        !           124:          IF( SCALE.LE.SAFMN2 )
        !           125:      $      GO TO 20
        !           126:       END IF
        !           127:       F2 = ABSSQ( FS )
        !           128:       G2 = ABSSQ( GS )
        !           129:       IF( F2.LE.MAX( G2, ONE )*SAFMIN ) THEN
        !           130: *
        !           131: *        This is a rare case: F is very small.
        !           132: *
        !           133:          IF( F.EQ.CZERO ) THEN
        !           134:             CS = ZERO
        !           135:             R = DLAPY2( DBLE( G ), DIMAG( G ) )
        !           136: *           Do complex/real division explicitly with two real divisions
        !           137:             D = DLAPY2( DBLE( GS ), DIMAG( GS ) )
        !           138:             SN = DCMPLX( DBLE( GS ) / D, -DIMAG( GS ) / D )
        !           139:             RETURN
        !           140:          END IF
        !           141:          F2S = DLAPY2( DBLE( FS ), DIMAG( FS ) )
        !           142: *        G2 and G2S are accurate
        !           143: *        G2 is at least SAFMIN, and G2S is at least SAFMN2
        !           144:          G2S = SQRT( G2 )
        !           145: *        Error in CS from underflow in F2S is at most
        !           146: *        UNFL / SAFMN2 .lt. sqrt(UNFL*EPS) .lt. EPS
        !           147: *        If MAX(G2,ONE)=G2, then F2 .lt. G2*SAFMIN,
        !           148: *        and so CS .lt. sqrt(SAFMIN)
        !           149: *        If MAX(G2,ONE)=ONE, then F2 .lt. SAFMIN
        !           150: *        and so CS .lt. sqrt(SAFMIN)/SAFMN2 = sqrt(EPS)
        !           151: *        Therefore, CS = F2S/G2S / sqrt( 1 + (F2S/G2S)**2 ) = F2S/G2S
        !           152:          CS = F2S / G2S
        !           153: *        Make sure abs(FF) = 1
        !           154: *        Do complex/real division explicitly with 2 real divisions
        !           155:          IF( ABS1( F ).GT.ONE ) THEN
        !           156:             D = DLAPY2( DBLE( F ), DIMAG( F ) )
        !           157:             FF = DCMPLX( DBLE( F ) / D, DIMAG( F ) / D )
        !           158:          ELSE
        !           159:             DR = SAFMX2*DBLE( F )
        !           160:             DI = SAFMX2*DIMAG( F )
        !           161:             D = DLAPY2( DR, DI )
        !           162:             FF = DCMPLX( DR / D, DI / D )
        !           163:          END IF
        !           164:          SN = FF*DCMPLX( DBLE( GS ) / G2S, -DIMAG( GS ) / G2S )
        !           165:          R = CS*F + SN*G
        !           166:       ELSE
        !           167: *
        !           168: *        This is the most common case.
        !           169: *        Neither F2 nor F2/G2 are less than SAFMIN
        !           170: *        F2S cannot overflow, and it is accurate
        !           171: *
        !           172:          F2S = SQRT( ONE+G2 / F2 )
        !           173: *        Do the F2S(real)*FS(complex) multiply with two real multiplies
        !           174:          R = DCMPLX( F2S*DBLE( FS ), F2S*DIMAG( FS ) )
        !           175:          CS = ONE / F2S
        !           176:          D = F2 + G2
        !           177: *        Do complex/real division explicitly with two real divisions
        !           178:          SN = DCMPLX( DBLE( R ) / D, DIMAG( R ) / D )
        !           179:          SN = SN*DCONJG( GS )
        !           180:          IF( COUNT.NE.0 ) THEN
        !           181:             IF( COUNT.GT.0 ) THEN
        !           182:                DO 30 I = 1, COUNT
        !           183:                   R = R*SAFMX2
        !           184:    30          CONTINUE
        !           185:             ELSE
        !           186:                DO 40 I = 1, -COUNT
        !           187:                   R = R*SAFMN2
        !           188:    40          CONTINUE
        !           189:             END IF
        !           190:          END IF
        !           191:       END IF
        !           192:       RETURN
        !           193: *
        !           194: *     End of ZLARTG
        !           195: *
        !           196:       END

CVSweb interface <joel.bertrand@systella.fr>