[BACK]Return to drotmg.f CVS log [TXT] [DIR] Up to [local] / rpl / lapack / blas

Diff for /rpl/lapack/blas/drotmg.f between versions 1.4 and 1.17

version 1.4, 2010/08/07 13:18:05 version 1.17, 2023/08/07 08:38:43
Line 1 Line 1
       SUBROUTINE DROTMG(DD1,DD2,DX1,DY1,DPARAM)  *> \brief \b DROTMG
 *     .. Scalar Arguments ..  
       DOUBLE PRECISION DD1,DD2,DX1,DY1  
 *     ..  
 *     .. Array Arguments ..  
       DOUBLE PRECISION DPARAM(5)  
 *     ..  
 *  
 *  Purpose  
 *  =======  
 *  *
 *     CONSTRUCT THE MODIFIED GIVENS TRANSFORMATION MATRIX H WHICH ZEROS  *  =========== DOCUMENTATION ===========
 *     THE SECOND COMPONENT OF THE 2-VECTOR  (DSQRT(DD1)*DX1,DSQRT(DD2)*  
 *     DY2)**T.  
 *     WITH DPARAM(1)=DFLAG, H HAS ONE OF THE FOLLOWING FORMS..  
 *  *
 *     DFLAG=-1.D0     DFLAG=0.D0        DFLAG=1.D0     DFLAG=-2.D0  * Online html documentation available at
   *            http://www.netlib.org/lapack/explore-html/
 *  *
 *       (DH11  DH12)    (1.D0  DH12)    (DH11  1.D0)    (1.D0  0.D0)  *  Definition:
 *     H=(          )    (          )    (          )    (          )  *  ===========
 *       (DH21  DH22),   (DH21  1.D0),   (-1.D0 DH22),   (0.D0  1.D0).  *
 *     LOCATIONS 2-4 OF DPARAM CONTAIN DH11, DH21, DH12, AND DH22  *       SUBROUTINE DROTMG(DD1,DD2,DX1,DY1,DPARAM)
 *     RESPECTIVELY. (VALUES OF 1.D0, -1.D0, OR 0.D0 IMPLIED BY THE  *
 *     VALUE OF DPARAM(1) ARE NOT STORED IN DPARAM.)  *       .. Scalar Arguments ..
   *       DOUBLE PRECISION DD1,DD2,DX1,DY1
   *       ..
   *       .. Array Arguments ..
   *       DOUBLE PRECISION DPARAM(5)
   *       ..
   *
   *
   *> \par Purpose:
   *  =============
   *>
   *> \verbatim
   *>
   *>    CONSTRUCT THE MODIFIED GIVENS TRANSFORMATION MATRIX H WHICH ZEROS
   *>    THE SECOND COMPONENT OF THE 2-VECTOR  (DSQRT(DD1)*DX1,DSQRT(DD2)*>    DY2)**T.
   *>    WITH DPARAM(1)=DFLAG, H HAS ONE OF THE FOLLOWING FORMS..
   *>
   *>    DFLAG=-1.D0     DFLAG=0.D0        DFLAG=1.D0     DFLAG=-2.D0
   *>
   *>      (DH11  DH12)    (1.D0  DH12)    (DH11  1.D0)    (1.D0  0.D0)
   *>    H=(          )    (          )    (          )    (          )
   *>      (DH21  DH22),   (DH21  1.D0),   (-1.D0 DH22),   (0.D0  1.D0).
   *>    LOCATIONS 2-4 OF DPARAM CONTAIN DH11, DH21, DH12, AND DH22
   *>    RESPECTIVELY. (VALUES OF 1.D0, -1.D0, OR 0.D0 IMPLIED BY THE
   *>    VALUE OF DPARAM(1) ARE NOT STORED IN DPARAM.)
   *>
   *>    THE VALUES OF GAMSQ AND RGAMSQ SET IN THE DATA STATEMENT MAY BE
   *>    INEXACT.  THIS IS OK AS THEY ARE ONLY USED FOR TESTING THE SIZE
   *>    OF DD1 AND DD2.  ALL ACTUAL SCALING OF DATA IS DONE USING GAM.
   *>
   *> \endverbatim
   *
   *  Arguments:
   *  ==========
   *
   *> \param[in,out] DD1
   *> \verbatim
   *>          DD1 is DOUBLE PRECISION
   *> \endverbatim
   *>
   *> \param[in,out] DD2
   *> \verbatim
   *>          DD2 is DOUBLE PRECISION
   *> \endverbatim
   *>
   *> \param[in,out] DX1
   *> \verbatim
   *>          DX1 is DOUBLE PRECISION
   *> \endverbatim
   *>
   *> \param[in] DY1
   *> \verbatim
   *>          DY1 is DOUBLE PRECISION
   *> \endverbatim
   *>
   *> \param[out] DPARAM
   *> \verbatim
   *>          DPARAM is DOUBLE PRECISION array, dimension (5)
   *>     DPARAM(1)=DFLAG
   *>     DPARAM(2)=DH11
   *>     DPARAM(3)=DH21
   *>     DPARAM(4)=DH12
   *>     DPARAM(5)=DH22
   *> \endverbatim
   *
   *  Authors:
   *  ========
   *
   *> \author Univ. of Tennessee
   *> \author Univ. of California Berkeley
   *> \author Univ. of Colorado Denver
   *> \author NAG Ltd.
 *  *
 *     THE VALUES OF GAMSQ AND RGAMSQ SET IN THE DATA STATEMENT MAY BE  *> \ingroup double_blas_level1
 *     INEXACT.  THIS IS OK AS THEY ARE ONLY USED FOR TESTING THE SIZE  
 *     OF DD1 AND DD2.  ALL ACTUAL SCALING OF DATA IS DONE USING GAM.  
 *  *
   *  =====================================================================
         SUBROUTINE DROTMG(DD1,DD2,DX1,DY1,DPARAM)
 *  *
 *  Arguments  *  -- Reference BLAS level1 routine --
 *  =========  *  -- Reference BLAS is a software package provided by Univ. of Tennessee,    --
 *  *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 *  DD1    (input/output) DOUBLE PRECISION  
 *  
 *  DD2    (input/output) DOUBLE PRECISION  
 *  
 *  DX1    (input/output) DOUBLE PRECISION  
 *  
 *  DY1    (input) DOUBLE PRECISION  
 *  *
 *  DPARAM (input/output)  DOUBLE PRECISION array, dimension 5  *     .. Scalar Arguments ..
 *     DPARAM(1)=DFLAG        DOUBLE PRECISION DD1,DD2,DX1,DY1
 *     DPARAM(2)=DH11  *     ..
 *     DPARAM(3)=DH21  *     .. Array Arguments ..
 *     DPARAM(4)=DH12        DOUBLE PRECISION DPARAM(5)
 *     DPARAM(5)=DH22  *     ..
 *  *
 *  =====================================================================  *  =====================================================================
 *  *
 *     .. Local Scalars ..  *     .. Local Scalars ..
       DOUBLE PRECISION DFLAG,DH11,DH12,DH21,DH22,DP1,DP2,DQ1,DQ2,DTEMP,        DOUBLE PRECISION DFLAG,DH11,DH12,DH21,DH22,DP1,DP2,DQ1,DQ2,DTEMP,
      +                 DU,GAM,GAMSQ,ONE,RGAMSQ,TWO,ZERO       $                 DU,GAM,GAMSQ,ONE,RGAMSQ,TWO,ZERO
       INTEGER IGO  
 *     ..  *     ..
 *     .. Intrinsic Functions ..  *     .. Intrinsic Functions ..
       INTRINSIC DABS        INTRINSIC DABS
Line 62 Line 114
       DATA GAM,GAMSQ,RGAMSQ/4096.D0,16777216.D0,5.9604645D-8/        DATA GAM,GAMSQ,RGAMSQ/4096.D0,16777216.D0,5.9604645D-8/
 *     ..  *     ..
   
       IF (.NOT.DD1.LT.ZERO) GO TO 10        IF (DD1.LT.ZERO) THEN
 *       GO ZERO-H-D-AND-DX1..  *        GO ZERO-H-D-AND-DX1..
       GO TO 60           DFLAG = -ONE
    10 CONTINUE           DH11 = ZERO
 *     CASE-DD1-NONNEGATIVE           DH12 = ZERO
       DP2 = DD2*DY1           DH21 = ZERO
       IF (.NOT.DP2.EQ.ZERO) GO TO 20           DH22 = ZERO
       DFLAG = -TWO  *
       GO TO 260           DD1 = ZERO
    20 CONTINUE           DD2 = ZERO
 *     REGULAR-CASE..           DX1 = ZERO
       DP1 = DD1*DX1        ELSE
       DQ2 = DP2*DY1  *        CASE-DD1-NONNEGATIVE
       DQ1 = DP1*DX1           DP2 = DD2*DY1
 *           IF (DP2.EQ.ZERO) THEN
       IF (.NOT.DABS(DQ1).GT.DABS(DQ2)) GO TO 40              DFLAG = -TWO
       DH21 = -DY1/DX1              DPARAM(1) = DFLAG
       DH12 = DP2/DP1              RETURN
 *           END IF
       DU = ONE - DH12*DH21  *        REGULAR-CASE..
 *           DP1 = DD1*DX1
       IF (.NOT.DU.LE.ZERO) GO TO 30           DQ2 = DP2*DY1
 *         GO ZERO-H-D-AND-DX1..           DQ1 = DP1*DX1
       GO TO 60  *
    30 CONTINUE           IF (DABS(DQ1).GT.DABS(DQ2)) THEN
       DFLAG = ZERO              DH21 = -DY1/DX1
       DD1 = DD1/DU              DH12 = DP2/DP1
       DD2 = DD2/DU  *
       DX1 = DX1*DU              DU = ONE - DH12*DH21
 *         GO SCALE-CHECK..  *
       GO TO 100             IF (DU.GT.ZERO) THEN
    40 CONTINUE               DFLAG = ZERO
       IF (.NOT.DQ2.LT.ZERO) GO TO 50               DD1 = DD1/DU
 *         GO ZERO-H-D-AND-DX1..               DD2 = DD2/DU
       GO TO 60               DX1 = DX1*DU
    50 CONTINUE             ELSE
       DFLAG = ONE  *            This code path if here for safety. We do not expect this
       DH11 = DP1/DP2  *            condition to ever hold except in edge cases with rounding
       DH22 = DX1/DY1  *            errors. See DOI: 10.1145/355841.355847
       DU = ONE + DH11*DH22               DFLAG = -ONE
       DTEMP = DD2/DU               DH11 = ZERO
       DD2 = DD1/DU               DH12 = ZERO
       DD1 = DTEMP               DH21 = ZERO
       DX1 = DY1*DU               DH22 = ZERO
 *         GO SCALE-CHECK  *
       GO TO 100               DD1 = ZERO
    60 CONTINUE               DD2 = ZERO
 *     PROCEDURE..ZERO-H-D-AND-DX1..               DX1 = ZERO
       DFLAG = -ONE             END IF
       DH11 = ZERO           ELSE
       DH12 = ZERO  
       DH21 = ZERO              IF (DQ2.LT.ZERO) THEN
       DH22 = ZERO  *              GO ZERO-H-D-AND-DX1..
 *                 DFLAG = -ONE
       DD1 = ZERO                 DH11 = ZERO
       DD2 = ZERO                 DH12 = ZERO
       DX1 = ZERO                 DH21 = ZERO
 *         RETURN..                 DH22 = ZERO
       GO TO 220  *
    70 CONTINUE                 DD1 = ZERO
 *     PROCEDURE..FIX-H..                 DD2 = ZERO
       IF (.NOT.DFLAG.GE.ZERO) GO TO 90                 DX1 = ZERO
 *              ELSE
       IF (.NOT.DFLAG.EQ.ZERO) GO TO 80                 DFLAG = ONE
       DH11 = ONE                 DH11 = DP1/DP2
       DH22 = ONE                 DH22 = DX1/DY1
       DFLAG = -ONE                 DU = ONE + DH11*DH22
       GO TO 90                 DTEMP = DD2/DU
    80 CONTINUE                 DD2 = DD1/DU
       DH21 = -ONE                 DD1 = DTEMP
       DH12 = ONE                 DX1 = DY1*DU
       DFLAG = -ONE              END IF
    90 CONTINUE           END IF
       GO TO (150,180,210) IGO  
       GO TO 120  
   100 CONTINUE  
 *     PROCEDURE..SCALE-CHECK  *     PROCEDURE..SCALE-CHECK
   110 CONTINUE           IF (DD1.NE.ZERO) THEN
       IF (.NOT.DD1.LE.RGAMSQ) GO TO 130              DO WHILE ((DD1.LE.RGAMSQ) .OR. (DD1.GE.GAMSQ))
       IF (DD1.EQ.ZERO) GO TO 160                 IF (DFLAG.EQ.ZERO) THEN
       IGO = 0                    DH11 = ONE
 *              FIX-H..                    DH22 = ONE
       GO TO 70                    DFLAG = -ONE
   120 CONTINUE                 ELSE
       DD1 = DD1*GAM**2                    DH21 = -ONE
       DX1 = DX1/GAM                    DH12 = ONE
       DH11 = DH11/GAM                    DFLAG = -ONE
       DH12 = DH12/GAM                 END IF
       GO TO 110                 IF (DD1.LE.RGAMSQ) THEN
   130 CONTINUE                    DD1 = DD1*GAM**2
   140 CONTINUE                    DX1 = DX1/GAM
       IF (.NOT.DD1.GE.GAMSQ) GO TO 160                    DH11 = DH11/GAM
       IGO = 1                    DH12 = DH12/GAM
 *              FIX-H..                 ELSE
       GO TO 70                    DD1 = DD1/GAM**2
   150 CONTINUE                    DX1 = DX1*GAM
       DD1 = DD1/GAM**2                    DH11 = DH11*GAM
       DX1 = DX1*GAM                    DH12 = DH12*GAM
       DH11 = DH11*GAM                 END IF
       DH12 = DH12*GAM              ENDDO
       GO TO 140           END IF
   160 CONTINUE  
   170 CONTINUE           IF (DD2.NE.ZERO) THEN
       IF (.NOT.DABS(DD2).LE.RGAMSQ) GO TO 190              DO WHILE ( (DABS(DD2).LE.RGAMSQ) .OR. (DABS(DD2).GE.GAMSQ) )
       IF (DD2.EQ.ZERO) GO TO 220                 IF (DFLAG.EQ.ZERO) THEN
       IGO = 2                    DH11 = ONE
 *              FIX-H..                    DH22 = ONE
       GO TO 70                    DFLAG = -ONE
   180 CONTINUE                 ELSE
       DD2 = DD2*GAM**2                    DH21 = -ONE
       DH21 = DH21/GAM                    DH12 = ONE
       DH22 = DH22/GAM                    DFLAG = -ONE
       GO TO 170                 END IF
   190 CONTINUE                 IF (DABS(DD2).LE.RGAMSQ) THEN
   200 CONTINUE                    DD2 = DD2*GAM**2
       IF (.NOT.DABS(DD2).GE.GAMSQ) GO TO 220                    DH21 = DH21/GAM
       IGO = 3                    DH22 = DH22/GAM
 *              FIX-H..                 ELSE
       GO TO 70                    DD2 = DD2/GAM**2
   210 CONTINUE                    DH21 = DH21*GAM
       DD2 = DD2/GAM**2                    DH22 = DH22*GAM
       DH21 = DH21*GAM                 END IF
       DH22 = DH22*GAM              END DO
       GO TO 200           END IF
   220 CONTINUE  
         END IF
   
       IF (DFLAG.LT.ZERO) THEN        IF (DFLAG.LT.ZERO) THEN
          GO TO 250           DPARAM(2) = DH11
            DPARAM(3) = DH21
            DPARAM(4) = DH12
            DPARAM(5) = DH22
       ELSE IF (DFLAG.EQ.ZERO) THEN        ELSE IF (DFLAG.EQ.ZERO) THEN
          GO TO 230            DPARAM(3) = DH21
            DPARAM(4) = DH12
       ELSE        ELSE
          GO TO 240           DPARAM(2) = DH11
            DPARAM(5) = DH22
       END IF        END IF
   230 CONTINUE  
       DPARAM(3) = DH21  
       DPARAM(4) = DH12  
       GO TO 260  
   240 CONTINUE  
       DPARAM(2) = DH11  
       DPARAM(5) = DH22  
       GO TO 260  
   250 CONTINUE  
       DPARAM(2) = DH11  
       DPARAM(3) = DH21  
       DPARAM(4) = DH12  
       DPARAM(5) = DH22  
   260 CONTINUE  
       DPARAM(1) = DFLAG        DPARAM(1) = DFLAG
       RETURN        RETURN
   *
   *     End of DROTMG
   *
       END        END
Removed from v.1.4  
changed lines
  Added in v.1.17

CVSweb interface <joel.bertrand@systella.fr>