1 /* drotmg.f -- translated by f2c (version 20100827). 2 You must link the resulting object file with libf2c: 3 on Microsoft Windows system, link with libf2c.lib; 4 on Linux or Unix systems, link with .../path/to/libf2c.a -lm 5 or, if you install libf2c.a in a standard place, with -lf2c -lm 6 -- in that order, at the end of the command line, as in 7 cc *.o -lf2c -lm 8 Source for libf2c is in /netlib/f2c/libf2c.zip, e.g., 9 10 http://www.netlib.org/f2c/libf2c.zip 11 */ 12 13 #include "datatypes.h" 14 15 /* Subroutine */ int drotmg_(doublereal *dd1, doublereal *dd2, doublereal * 16 dx1, doublereal *dy1, doublereal *dparam) 17 { 18 /* Initialized data */ 19 20 static doublereal zero = 0.; 21 static doublereal one = 1.; 22 static doublereal two = 2.; 23 static doublereal gam = 4096.; 24 static doublereal gamsq = 16777216.; 25 static doublereal rgamsq = 5.9604645e-8; 26 27 /* Format strings */ 28 static char fmt_120[] = ""; 29 static char fmt_150[] = ""; 30 static char fmt_180[] = ""; 31 static char fmt_210[] = ""; 32 33 /* System generated locals */ 34 doublereal d__1; 35 36 /* Local variables */ 37 doublereal du, dp1, dp2, dq1, dq2, dh11, dh12, dh21, dh22; 38 integer igo; 39 doublereal dflag, dtemp; 40 41 /* Assigned format variables */ 42 static char *igo_fmt; 43 44 /* .. Scalar Arguments .. */ 45 /* .. */ 46 /* .. Array Arguments .. */ 47 /* .. */ 48 49 /* Purpose */ 50 /* ======= */ 51 52 /* CONSTRUCT THE MODIFIED GIVENS TRANSFORMATION MATRIX H WHICH ZEROS */ 53 /* THE SECOND COMPONENT OF THE 2-VECTOR (DSQRT(DD1)*DX1,DSQRT(DD2)* */ 54 /* DY2)**T. */ 55 /* WITH DPARAM(1)=DFLAG, H HAS ONE OF THE FOLLOWING FORMS.. */ 56 57 /* DFLAG=-1.D0 DFLAG=0.D0 DFLAG=1.D0 DFLAG=-2.D0 */ 58 59 /* (DH11 DH12) (1.D0 DH12) (DH11 1.D0) (1.D0 0.D0) */ 60 /* H=( ) ( ) ( ) ( ) */ 61 /* (DH21 DH22), (DH21 1.D0), (-1.D0 DH22), (0.D0 1.D0). */ 62 /* LOCATIONS 2-4 OF DPARAM CONTAIN DH11, DH21, DH12, AND DH22 */ 63 /* RESPECTIVELY. (VALUES OF 1.D0, -1.D0, OR 0.D0 IMPLIED BY THE */ 64 /* VALUE OF DPARAM(1) ARE NOT STORED IN DPARAM.) */ 65 66 /* THE VALUES OF GAMSQ AND RGAMSQ SET IN THE DATA STATEMENT MAY BE */ 67 /* INEXACT. THIS IS OK AS THEY ARE ONLY USED FOR TESTING THE SIZE */ 68 /* OF DD1 AND DD2. ALL ACTUAL SCALING OF DATA IS DONE USING GAM. */ 69 70 71 /* Arguments */ 72 /* ========= */ 73 74 /* DD1 (input/output) DOUBLE PRECISION */ 75 76 /* DD2 (input/output) DOUBLE PRECISION */ 77 78 /* DX1 (input/output) DOUBLE PRECISION */ 79 80 /* DY1 (input) DOUBLE PRECISION */ 81 82 /* DPARAM (input/output) DOUBLE PRECISION array, dimension 5 */ 83 /* DPARAM(1)=DFLAG */ 84 /* DPARAM(2)=DH11 */ 85 /* DPARAM(3)=DH21 */ 86 /* DPARAM(4)=DH12 */ 87 /* DPARAM(5)=DH22 */ 88 89 /* ===================================================================== */ 90 91 /* .. Local Scalars .. */ 92 /* .. */ 93 /* .. Intrinsic Functions .. */ 94 /* .. */ 95 /* .. Data statements .. */ 96 97 /* Parameter adjustments */ 98 --dparam; 99 100 /* Function Body */ 101 /* .. */ 102 if (! (*dd1 < zero)) { 103 goto L10; 104 } 105 /* GO ZERO-H-D-AND-DX1.. */ 106 goto L60; 107 L10: 108 /* CASE-DD1-NONNEGATIVE */ 109 dp2 = *dd2 * *dy1; 110 if (! (dp2 == zero)) { 111 goto L20; 112 } 113 dflag = -two; 114 goto L260; 115 /* REGULAR-CASE.. */ 116 L20: 117 dp1 = *dd1 * *dx1; 118 dq2 = dp2 * *dy1; 119 dq1 = dp1 * *dx1; 120 121 if (! (abs(dq1) > abs(dq2))) { 122 goto L40; 123 } 124 dh21 = -(*dy1) / *dx1; 125 dh12 = dp2 / dp1; 126 127 du = one - dh12 * dh21; 128 129 if (! (du <= zero)) { 130 goto L30; 131 } 132 /* GO ZERO-H-D-AND-DX1.. */ 133 goto L60; 134 L30: 135 dflag = zero; 136 *dd1 /= du; 137 *dd2 /= du; 138 *dx1 *= du; 139 /* GO SCALE-CHECK.. */ 140 goto L100; 141 L40: 142 if (! (dq2 < zero)) { 143 goto L50; 144 } 145 /* GO ZERO-H-D-AND-DX1.. */ 146 goto L60; 147 L50: 148 dflag = one; 149 dh11 = dp1 / dp2; 150 dh22 = *dx1 / *dy1; 151 du = one + dh11 * dh22; 152 dtemp = *dd2 / du; 153 *dd2 = *dd1 / du; 154 *dd1 = dtemp; 155 *dx1 = *dy1 * du; 156 /* GO SCALE-CHECK */ 157 goto L100; 158 /* PROCEDURE..ZERO-H-D-AND-DX1.. */ 159 L60: 160 dflag = -one; 161 dh11 = zero; 162 dh12 = zero; 163 dh21 = zero; 164 dh22 = zero; 165 166 *dd1 = zero; 167 *dd2 = zero; 168 *dx1 = zero; 169 /* RETURN.. */ 170 goto L220; 171 /* PROCEDURE..FIX-H.. */ 172 L70: 173 if (! (dflag >= zero)) { 174 goto L90; 175 } 176 177 if (! (dflag == zero)) { 178 goto L80; 179 } 180 dh11 = one; 181 dh22 = one; 182 dflag = -one; 183 goto L90; 184 L80: 185 dh21 = -one; 186 dh12 = one; 187 dflag = -one; 188 L90: 189 switch (igo) { 190 case 0: goto L120; 191 case 1: goto L150; 192 case 2: goto L180; 193 case 3: goto L210; 194 } 195 /* PROCEDURE..SCALE-CHECK */ 196 L100: 197 L110: 198 if (! (*dd1 <= rgamsq)) { 199 goto L130; 200 } 201 if (*dd1 == zero) { 202 goto L160; 203 } 204 igo = 0; 205 igo_fmt = fmt_120; 206 /* FIX-H.. */ 207 goto L70; 208 L120: 209 /* Computing 2nd power */ 210 d__1 = gam; 211 *dd1 *= d__1 * d__1; 212 *dx1 /= gam; 213 dh11 /= gam; 214 dh12 /= gam; 215 goto L110; 216 L130: 217 L140: 218 if (! (*dd1 >= gamsq)) { 219 goto L160; 220 } 221 igo = 1; 222 igo_fmt = fmt_150; 223 /* FIX-H.. */ 224 goto L70; 225 L150: 226 /* Computing 2nd power */ 227 d__1 = gam; 228 *dd1 /= d__1 * d__1; 229 *dx1 *= gam; 230 dh11 *= gam; 231 dh12 *= gam; 232 goto L140; 233 L160: 234 L170: 235 if (! (abs(*dd2) <= rgamsq)) { 236 goto L190; 237 } 238 if (*dd2 == zero) { 239 goto L220; 240 } 241 igo = 2; 242 igo_fmt = fmt_180; 243 /* FIX-H.. */ 244 goto L70; 245 L180: 246 /* Computing 2nd power */ 247 d__1 = gam; 248 *dd2 *= d__1 * d__1; 249 dh21 /= gam; 250 dh22 /= gam; 251 goto L170; 252 L190: 253 L200: 254 if (! (abs(*dd2) >= gamsq)) { 255 goto L220; 256 } 257 igo = 3; 258 igo_fmt = fmt_210; 259 /* FIX-H.. */ 260 goto L70; 261 L210: 262 /* Computing 2nd power */ 263 d__1 = gam; 264 *dd2 /= d__1 * d__1; 265 dh21 *= gam; 266 dh22 *= gam; 267 goto L200; 268 L220: 269 if (dflag < 0.) { 270 goto L250; 271 } else if (dflag == 0) { 272 goto L230; 273 } else { 274 goto L240; 275 } 276 L230: 277 dparam[3] = dh21; 278 dparam[4] = dh12; 279 goto L260; 280 L240: 281 dparam[2] = dh11; 282 dparam[5] = dh22; 283 goto L260; 284 L250: 285 dparam[2] = dh11; 286 dparam[3] = dh21; 287 dparam[4] = dh12; 288 dparam[5] = dh22; 289 L260: 290 dparam[1] = dflag; 291 return 0; 292 } /* drotmg_ */ 293 294