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