Blame view
fvn_fnlib/r9sifg.f
8.16 KB
38581db0c git-svn-id: https... |
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 |
subroutine r9sifg (x, f, g) c december 1980 edition. w. fullerton, bell labs dimension f1cs(20), f2cs(29), g1cs(21), g2cs(34) external csevl, inits, r1mach c c series for f1 on the interval 2.00000e-02 to 6.25000e-02 c with weighted error 2.82e-17 c log weighted error 16.55 c significant figures required 15.36 c decimal places required 17.20 c data f1 cs( 1) / -0.1191081969 051363610 e0/ data f1 cs( 2) / -0.0247823144 996236248 e0/ data f1 cs( 3) / 0.0011910281 453357821 e0/ data f1 cs( 4) / -0.0000927027 714388562 e0/ data f1 cs( 5) / 0.0000093373 141568271 e0/ data f1 cs( 6) / -0.0000011058 287820557 e0/ data f1 cs( 7) / 0.0000001464 772071460 e0/ data f1 cs( 8) / -0.0000000210 694496288 e0/ data f1 cs( 9) / 0.0000000032 293492367 e0/ data f1 cs( 10) / -0.0000000005 206529618 e0/ data f1 cs( 11) / 0.0000000000 874878885 e0/ data f1 cs( 12) / -0.0000000000 152176187 e0/ data f1 cs( 13) / 0.0000000000 027257192 e0/ data f1 cs( 14) / -0.0000000000 005007053 e0/ data f1 cs( 15) / 0.0000000000 000940241 e0/ data f1 cs( 16) / -0.0000000000 000180014 e0/ data f1 cs( 17) / 0.0000000000 000035063 e0/ data f1 cs( 18) / -0.0000000000 000006935 e0/ data f1 cs( 19) / 0.0000000000 000001391 e0/ data f1 cs( 20) / -0.0000000000 000000282 e0/ c c series for f2 on the interval 0.00000e+00 to 2.00000e-02 c with weighted error 4.32e-17 c log weighted error 16.36 c significant figures required 14.75 c decimal places required 17.10 c data f2 cs( 1) / -0.0348409253 897013234 e0/ data f2 cs( 2) / -0.0166842205 677959686 e0/ data f2 cs( 3) / 0.0006752901 241237738 e0/ data f2 cs( 4) / -0.0000535066 622544701 e0/ data f2 cs( 5) / 0.0000062693 421779007 e0/ data f2 cs( 6) / -0.0000009526 638801991 e0/ data f2 cs( 7) / 0.0000001745 629224251 e0/ data f2 cs( 8) / -0.0000000368 795403065 e0/ data f2 cs( 9) / 0.0000000087 202677705 e0/ data f2 cs( 10) / -0.0000000022 601970392 e0/ data f2 cs( 11) / 0.0000000006 324624977 e0/ data f2 cs( 12) / -0.0000000001 888911889 e0/ data f2 cs( 13) / 0.0000000000 596774674 e0/ data f2 cs( 14) / -0.0000000000 198044313 e0/ data f2 cs( 15) / 0.0000000000 068641396 e0/ data f2 cs( 16) / -0.0000000000 024731020 e0/ data f2 cs( 17) / 0.0000000000 009226360 e0/ data f2 cs( 18) / -0.0000000000 003552364 e0/ data f2 cs( 19) / 0.0000000000 001407606 e0/ data f2 cs( 20) / -0.0000000000 000572623 e0/ data f2 cs( 21) / 0.0000000000 000238654 e0/ data f2 cs( 22) / -0.0000000000 000101714 e0/ data f2 cs( 23) / 0.0000000000 000044259 e0/ data f2 cs( 24) / -0.0000000000 000019634 e0/ data f2 cs( 25) / 0.0000000000 000008868 e0/ data f2 cs( 26) / -0.0000000000 000004074 e0/ data f2 cs( 27) / 0.0000000000 000001901 e0/ data f2 cs( 28) / -0.0000000000 000000900 e0/ data f2 cs( 29) / 0.0000000000 000000432 e0/ c c series for g1 on the interval 2.00000e-02 to 6.25000e-02 c with weighted error 5.48e-17 c log weighted error 16.26 c significant figures required 15.47 c decimal places required 16.92 c data g1 cs( 1) / -0.3040578798 253495954 e0/ data g1 cs( 2) / -0.0566890984 597120588 e0/ data g1 cs( 3) / 0.0039046158 173275644 e0/ data g1 cs( 4) / -0.0003746075 959202261 e0/ data g1 cs( 5) / 0.0000435431 556559844 e0/ data g1 cs( 6) / -0.0000057417 294453025 e0/ data g1 cs( 7) / 0.0000008282 552104503 e0/ data g1 cs( 8) / -0.0000001278 245892595 e0/ data g1 cs( 9) / 0.0000000207 978352949 e0/ data g1 cs( 10) / -0.0000000035 313205922 e0/ data g1 cs( 11) / 0.0000000006 210824236 e0/ data g1 cs( 12) / -0.0000000001 125215474 e0/ data g1 cs( 13) / 0.0000000000 209088918 e0/ data g1 cs( 14) / -0.0000000000 039715832 e0/ data g1 cs( 15) / 0.0000000000 007690431 e0/ data g1 cs( 16) / -0.0000000000 001514697 e0/ data g1 cs( 17) / 0.0000000000 000302892 e0/ data g1 cs( 18) / -0.0000000000 000061400 e0/ data g1 cs( 19) / 0.0000000000 000012601 e0/ data g1 cs( 20) / -0.0000000000 000002615 e0/ data g1 cs( 21) / 0.0000000000 000000548 e0/ c c series for g2 on the interval 0.00000e+00 to 2.00000e-02 c with weighted error 5.01e-17 c log weighted error 16.30 c significant figures required 15.12 c decimal places required 17.07 c data g2 cs( 1) / -0.0967329367 532432218 e0/ data g2 cs( 2) / -0.0452077907 957459871 e0/ data g2 cs( 3) / 0.0028190005 352706523 e0/ data g2 cs( 4) / -0.0002899167 740759160 e0/ data g2 cs( 5) / 0.0000407444 664601121 e0/ data g2 cs( 6) / -0.0000071056 382192354 e0/ data g2 cs( 7) / 0.0000014534 723163019 e0/ data g2 cs( 8) / -0.0000003364 116512503 e0/ data g2 cs( 9) / 0.0000000859 774367886 e0/ data g2 cs( 10) / -0.0000000238 437656302 e0/ data g2 cs( 11) / 0.0000000070 831906340 e0/ data g2 cs( 12) / -0.0000000022 318068154 e0/ data g2 cs( 13) / 0.0000000007 401087359 e0/ data g2 cs( 14) / -0.0000000002 567171162 e0/ data g2 cs( 15) / 0.0000000000 926707021 e0/ data g2 cs( 16) / -0.0000000000 346693311 e0/ data g2 cs( 17) / 0.0000000000 133950573 e0/ data g2 cs( 18) / -0.0000000000 053290754 e0/ data g2 cs( 19) / 0.0000000000 021775312 e0/ data g2 cs( 20) / -0.0000000000 009118621 e0/ data g2 cs( 21) / 0.0000000000 003905864 e0/ data g2 cs( 22) / -0.0000000000 001708459 e0/ data g2 cs( 23) / 0.0000000000 000762015 e0/ data g2 cs( 24) / -0.0000000000 000346151 e0/ data g2 cs( 25) / 0.0000000000 000159996 e0/ data g2 cs( 26) / -0.0000000000 000075213 e0/ data g2 cs( 27) / 0.0000000000 000035970 e0/ data g2 cs( 28) / -0.0000000000 000017530 e0/ data g2 cs( 29) / 0.0000000000 000008738 e0/ data g2 cs( 30) / -0.0000000000 000004487 e0/ data g2 cs( 31) / 0.0000000000 000002397 e0/ data g2 cs( 32) / -0.0000000000 000001347 e0/ data g2 cs( 33) / 0.0000000000 000000801 e0/ data g2 cs( 34) / -0.0000000000 000000501 e0/ c data nf1, nf2, ng1, ng2 / 4*0 / data xbnd, xbig, xmaxf, xmaxg / 4*0.0 / c if (nf1.ne.0) go to 10 tol = 0.1*r1mach(3) nf1 = inits (f1cs, 20, tol) nf2 = inits (f2cs, 29, tol) ng1 = inits (g1cs, 21, tol) ng2 = inits (g2cs, 34, tol) c xbig = sqrt (1.0/r1mach(3)) xmaxf = exp (amin1(-alog(r1mach(1)), alog(r1mach(2))) - 0.01) xmaxg = 1.0/sqrt(r1mach(1)) xbnd = sqrt(50.0) c 10 if (x.lt.4.0) call seteru ( 1 34hr9sifg approxs invalid for x lt 4, 34, 1, 2) c if (x.gt.xbnd) go to 20 f = (1.0 + csevl ((1.0/x**2-0.04125)/0.02125, f1cs, nf1))/x g = (1.0 + csevl ((1.0/x**2-0.04125)/0.02125, g1cs, ng1))/x**2 return c 20 if (x.gt.xbig) go to 30 f = (1.0 + csevl (100./x**2-1., f2cs, nf2))/x g = (1.0 + csevl (100./x**2-1., g2cs, ng2))/x**2 return c 30 f = 0.0 if (x.lt.xmaxf) f = 1.0/x g = 0.0 if (x.lt.xmaxg) g = 1.0/x**2 return c end |