File:  [local] / rpl / lapack / lapack / zlartg.f
Revision 1.8: download - view: text, annotated - select for diffs - revision graph
Mon Nov 21 20:43:17 2011 UTC (12 years, 5 months ago) by bertrand
Branches: MAIN
CVS tags: HEAD
Mise à jour de Lapack.

    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: *  =====================================================================
  104:       SUBROUTINE ZLARTG( F, G, CS, SN, R )
  105: *
  106: *  -- LAPACK auxiliary routine (version 3.4.0) --
  107: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
  108: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
  109: *     November 2011
  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>