File:  [local] / rpl / lapack / blas / drotmg.f
Revision 1.8: download - view: text, annotated - select for diffs - revision graph
Fri Jul 22 07:38:01 2011 UTC (12 years, 9 months ago) by bertrand
Branches: MAIN
CVS tags: rpl-4_1_3, rpl-4_1_2, rpl-4_1_1, HEAD
En route vers la 4.4.1.

    1:       SUBROUTINE DROTMG(DD1,DD2,DX1,DY1,DPARAM)
    2: *     .. Scalar Arguments ..
    3:       DOUBLE PRECISION DD1,DD2,DX1,DY1
    4: *     ..
    5: *     .. Array Arguments ..
    6:       DOUBLE PRECISION DPARAM(5)
    7: *     ..
    8: *
    9: *  Purpose
   10: *  =======
   11: *
   12: *     CONSTRUCT THE MODIFIED GIVENS TRANSFORMATION MATRIX H WHICH ZEROS
   13: *     THE SECOND COMPONENT OF THE 2-VECTOR  (DSQRT(DD1)*DX1,DSQRT(DD2)*
   14: *     DY2)**T.
   15: *     WITH DPARAM(1)=DFLAG, H HAS ONE OF THE FOLLOWING FORMS..
   16: *
   17: *     DFLAG=-1.D0     DFLAG=0.D0        DFLAG=1.D0     DFLAG=-2.D0
   18: *
   19: *       (DH11  DH12)    (1.D0  DH12)    (DH11  1.D0)    (1.D0  0.D0)
   20: *     H=(          )    (          )    (          )    (          )
   21: *       (DH21  DH22),   (DH21  1.D0),   (-1.D0 DH22),   (0.D0  1.D0).
   22: *     LOCATIONS 2-4 OF DPARAM CONTAIN DH11, DH21, DH12, AND DH22
   23: *     RESPECTIVELY. (VALUES OF 1.D0, -1.D0, OR 0.D0 IMPLIED BY THE
   24: *     VALUE OF DPARAM(1) ARE NOT STORED IN DPARAM.)
   25: *
   26: *     THE VALUES OF GAMSQ AND RGAMSQ SET IN THE DATA STATEMENT MAY BE
   27: *     INEXACT.  THIS IS OK AS THEY ARE ONLY USED FOR TESTING THE SIZE
   28: *     OF DD1 AND DD2.  ALL ACTUAL SCALING OF DATA IS DONE USING GAM.
   29: *
   30: *
   31: *  Arguments
   32: *  =========
   33: *
   34: *  DD1    (input/output) DOUBLE PRECISION
   35: *
   36: *  DD2    (input/output) DOUBLE PRECISION
   37: *
   38: *  DX1    (input/output) DOUBLE PRECISION
   39: *
   40: *  DY1    (input) DOUBLE PRECISION
   41: *
   42: *  DPARAM (input/output)  DOUBLE PRECISION array, dimension 5
   43: *     DPARAM(1)=DFLAG
   44: *     DPARAM(2)=DH11
   45: *     DPARAM(3)=DH21
   46: *     DPARAM(4)=DH12
   47: *     DPARAM(5)=DH22
   48: *
   49: *  =====================================================================
   50: *
   51: *     .. Local Scalars ..
   52:       DOUBLE PRECISION DFLAG,DH11,DH12,DH21,DH22,DP1,DP2,DQ1,DQ2,DTEMP,
   53:      $                 DU,GAM,GAMSQ,ONE,RGAMSQ,TWO,ZERO
   54: *     ..
   55: *     .. Intrinsic Functions ..
   56:       INTRINSIC DABS
   57: *     ..
   58: *     .. Data statements ..
   59: *
   60:       DATA ZERO,ONE,TWO/0.D0,1.D0,2.D0/
   61:       DATA GAM,GAMSQ,RGAMSQ/4096.D0,16777216.D0,5.9604645D-8/
   62: *     ..
   63: 
   64:       IF (DD1.LT.ZERO) THEN
   65: *        GO ZERO-H-D-AND-DX1..
   66:          DFLAG = -ONE
   67:          DH11 = ZERO
   68:          DH12 = ZERO
   69:          DH21 = ZERO
   70:          DH22 = ZERO
   71: *
   72:          DD1 = ZERO
   73:          DD2 = ZERO
   74:          DX1 = ZERO
   75:       ELSE
   76: *        CASE-DD1-NONNEGATIVE
   77:          DP2 = DD2*DY1
   78:          IF (DP2.EQ.ZERO) THEN
   79:             DFLAG = -TWO
   80:             DPARAM(1) = DFLAG
   81:             RETURN
   82:          END IF 
   83: *        REGULAR-CASE..
   84:          DP1 = DD1*DX1
   85:          DQ2 = DP2*DY1
   86:          DQ1 = DP1*DX1
   87: *
   88:          IF (DABS(DQ1).GT.DABS(DQ2)) THEN
   89:             DH21 = -DY1/DX1
   90:             DH12 = DP2/DP1
   91: *
   92:             DU = ONE - DH12*DH21
   93: *
   94:            IF (DU.GT.ZERO) THEN
   95:              DFLAG = ZERO
   96:              DD1 = DD1/DU
   97:              DD2 = DD2/DU
   98:              DX1 = DX1*DU
   99:            END IF
  100:          ELSE
  101: 
  102:             IF (DQ2.LT.ZERO) THEN
  103: *              GO ZERO-H-D-AND-DX1..
  104:                DFLAG = -ONE
  105:                DH11 = ZERO
  106:                DH12 = ZERO
  107:                DH21 = ZERO
  108:                DH22 = ZERO
  109: *
  110:                DD1 = ZERO
  111:                DD2 = ZERO
  112:                DX1 = ZERO
  113:             ELSE
  114:                DFLAG = ONE
  115:                DH11 = DP1/DP2
  116:                DH22 = DX1/DY1
  117:                DU = ONE + DH11*DH22
  118:                DTEMP = DD2/DU
  119:                DD2 = DD1/DU
  120:                DD1 = DTEMP
  121:                DX1 = DY1*DU
  122:             END IF
  123:          END IF
  124: 
  125: *     PROCEDURE..SCALE-CHECK
  126:          IF (DD1.NE.ZERO) THEN
  127:             DO WHILE ((DD1.LE.RGAMSQ) .OR. (DD1.GE.GAMSQ))
  128:                IF (DFLAG.EQ.ZERO) THEN
  129:                   DH11 = ONE
  130:                   DH22 = ONE
  131:                   DFLAG = -ONE
  132:                ELSE
  133:                   DH21 = -ONE
  134:                   DH12 = ONE
  135:                   DFLAG = -ONE
  136:                END IF
  137:                IF (DD1.LE.RGAMSQ) THEN
  138:                   DD1 = DD1*GAM**2
  139:                   DX1 = DX1/GAM
  140:                   DH11 = DH11/GAM
  141:                   DH12 = DH12/GAM
  142:                ELSE
  143:                   DD1 = DD1/GAM**2
  144:                   DX1 = DX1*GAM
  145:                   DH11 = DH11*GAM
  146:                   DH12 = DH12*GAM
  147:                END IF
  148:             ENDDO
  149:          END IF
  150:   
  151:          IF (DD2.NE.ZERO) THEN
  152:             DO WHILE ( (DABS(DD2).LE.RGAMSQ) .OR. (DABS(DD2).GE.GAMSQ) )
  153:                IF (DFLAG.EQ.ZERO) THEN
  154:                   DH11 = ONE
  155:                   DH22 = ONE
  156:                   DFLAG = -ONE
  157:                ELSE
  158:                   DH21 = -ONE
  159:                   DH12 = ONE
  160:                   DFLAG = -ONE
  161:                END IF
  162:                IF (DABS(DD2).LE.RGAMSQ) THEN
  163:                   DD2 = DD2*GAM**2
  164:                   DH21 = DH21/GAM
  165:                   DH22 = DH22/GAM
  166:                ELSE
  167:                   DD2 = DD2/GAM**2
  168:                   DH21 = DH21*GAM
  169:                   DH22 = DH22*GAM
  170:                END IF      
  171:             END DO
  172:          END IF
  173:      
  174:       END IF
  175: 
  176:       IF (DFLAG.LT.ZERO) THEN
  177:          DPARAM(2) = DH11
  178:          DPARAM(3) = DH21
  179:          DPARAM(4) = DH12
  180:          DPARAM(5) = DH22
  181:       ELSE IF (DFLAG.EQ.ZERO) THEN
  182:          DPARAM(3) = DH21
  183:          DPARAM(4) = DH12 
  184:       ELSE
  185:          DPARAM(2) = DH11
  186:          DPARAM(5) = DH22
  187:       END IF
  188: 
  189:   260 CONTINUE
  190:       DPARAM(1) = DFLAG
  191:       RETURN
  192:       END
  193:       
  194:      
  195:      
  196:      

CVSweb interface <joel.bertrand@systella.fr>