drotmg.c revision acf83bd421067e7a2b828c5fe61594a4fda5e9a5
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; 107L10: 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.. */ 116L20: 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; 134L30: 135 dflag = zero; 136 *dd1 /= du; 137 *dd2 /= du; 138 *dx1 *= du; 139/* GO SCALE-CHECK.. */ 140 goto L100; 141L40: 142 if (! (dq2 < zero)) { 143 goto L50; 144 } 145/* GO ZERO-H-D-AND-DX1.. */ 146 goto L60; 147L50: 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.. */ 159L60: 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.. */ 172L70: 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; 184L80: 185 dh21 = -one; 186 dh12 = one; 187 dflag = -one; 188L90: 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 */ 196L100: 197L110: 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; 208L120: 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; 216L130: 217L140: 218 if (! (*dd1 >= gamsq)) { 219 goto L160; 220 } 221 igo = 1; 222 igo_fmt = fmt_150; 223/* FIX-H.. */ 224 goto L70; 225L150: 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; 233L160: 234L170: 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; 245L180: 246/* Computing 2nd power */ 247 d__1 = gam; 248 *dd2 *= d__1 * d__1; 249 dh21 /= gam; 250 dh22 /= gam; 251 goto L170; 252L190: 253L200: 254 if (! (abs(*dd2) >= gamsq)) { 255 goto L220; 256 } 257 igo = 3; 258 igo_fmt = fmt_210; 259/* FIX-H.. */ 260 goto L70; 261L210: 262/* Computing 2nd power */ 263 d__1 = gam; 264 *dd2 /= d__1 * d__1; 265 dh21 *= gam; 266 dh22 *= gam; 267 goto L200; 268L220: 269 if (dflag < 0.) { 270 goto L250; 271 } else if (dflag == 0) { 272 goto L230; 273 } else { 274 goto L240; 275 } 276L230: 277 dparam[3] = dh21; 278 dparam[4] = dh12; 279 goto L260; 280L240: 281 dparam[2] = dh11; 282 dparam[5] = dh22; 283 goto L260; 284L250: 285 dparam[2] = dh11; 286 dparam[3] = dh21; 287 dparam[4] = dh12; 288 dparam[5] = dh22; 289L260: 290 dparam[1] = dflag; 291 return 0; 292} /* drotmg_ */ 293 294