Blame view
fvn_fnlib/r9admp.f
11.3 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 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 |
subroutine r9admp (x, ampl, phi) c july 1980 edition. w. fullerton, bell labs. revised nov 1983. c c evaluate the derivative of the airy function modulus and c phase for x .le. -1.0. c dimension an22cs(33), an21cs(24), an20cs(16), aph2cs(32), 1 aph1cs(22), aph0cs(15) external csevl, inits, r1mach c c series for an22 on the interval -1.00000e+00 to -1.25000e-01 c with weighted error 3.30e-17 c log weighted error 16.48 c significant figures required 14.95 c decimal places required 17.24 c data an22cs( 1) / 0.0537418629 629794329 e0/ data an22cs( 2) / -0.0126661435 859883193 e0/ data an22cs( 3) / -0.0011924334 106593007 e0/ data an22cs( 4) / -0.0002032327 627275655 e0/ data an22cs( 5) / -0.0000446468 963075164 e0/ data an22cs( 6) / -0.0000113359 036053123 e0/ data an22cs( 7) / -0.0000031641 352378546 e0/ data an22cs( 8) / -0.0000009446 708886149 e0/ data an22cs( 9) / -0.0000002966 562236472 e0/ data an22cs( 10) / -0.0000000969 118892024 e0/ data an22cs( 11) / -0.0000000326 822538653 e0/ data an22cs( 12) / -0.0000000113 144618964 e0/ data an22cs( 13) / -0.0000000040 042691002 e0/ data an22cs( 14) / -0.0000000014 440333684 e0/ data an22cs( 15) / -0.0000000005 292853746 e0/ data an22cs( 16) / -0.0000000001 967763374 e0/ data an22cs( 17) / -0.0000000000 740800096 e0/ data an22cs( 18) / -0.0000000000 282016314 e0/ data an22cs( 19) / -0.0000000000 108440066 e0/ data an22cs( 20) / -0.0000000000 042074801 e0/ data an22cs( 21) / -0.0000000000 016459150 e0/ data an22cs( 22) / -0.0000000000 006486827 e0/ data an22cs( 23) / -0.0000000000 002574095 e0/ data an22cs( 24) / -0.0000000000 001027889 e0/ data an22cs( 25) / -0.0000000000 000412846 e0/ data an22cs( 26) / -0.0000000000 000166711 e0/ data an22cs( 27) / -0.0000000000 000067657 e0/ data an22cs( 28) / -0.0000000000 000027585 e0/ data an22cs( 29) / -0.0000000000 000011296 e0/ data an22cs( 30) / -0.0000000000 000004645 e0/ data an22cs( 31) / -0.0000000000 000001917 e0/ data an22cs( 32) / -0.0000000000 000000794 e0/ data an22cs( 33) / -0.0000000000 000000330 e0/ c c series for an21 on the interval -1.25000e-01 to -1.56250e-02 c with weighted error 3.43e-17 c log weighted error 16.47 c significant figures required 14.48 c decimal places required 17.16 c data an21cs( 1) / 0.0198313155 263169394 e0/ data an21cs( 2) / -0.0029376249 067087533 e0/ data an21cs( 3) / -0.0001136260 695958196 e0/ data an21cs( 4) / -0.0000100554 451087156 e0/ data an21cs( 5) / -0.0000013048 787116563 e0/ data an21cs( 6) / -0.0000002123 881993151 e0/ data an21cs( 7) / -0.0000000402 270833384 e0/ data an21cs( 8) / -0.0000000084 996745953 e0/ data an21cs( 9) / -0.0000000019 514839426 e0/ data an21cs( 10) / -0.0000000004 783865344 e0/ data an21cs( 11) / -0.0000000001 236733992 e0/ data an21cs( 12) / -0.0000000000 334137486 e0/ data an21cs( 13) / -0.0000000000 093702824 e0/ data an21cs( 14) / -0.0000000000 027130128 e0/ data an21cs( 15) / -0.0000000000 008075954 e0/ data an21cs( 16) / -0.0000000000 002463214 e0/ data an21cs( 17) / -0.0000000000 000767656 e0/ data an21cs( 18) / -0.0000000000 000243883 e0/ data an21cs( 19) / -0.0000000000 000078831 e0/ data an21cs( 20) / -0.0000000000 000025882 e0/ data an21cs( 21) / -0.0000000000 000008619 e0/ data an21cs( 22) / -0.0000000000 000002908 e0/ data an21cs( 23) / -0.0000000000 000000993 e0/ data an21cs( 24) / -0.0000000000 000000343 e0/ c c series for an20 on the interval -1.56250e-02 to 0.00000e+00 c with weighted error 4.41e-17 c log weighted error 16.36 c significant figures required 14.16 c decimal places required 16.96 c data an20cs( 1) / 0.0126732217 145738027 e0/ data an20cs( 2) / -0.0005212847 072615621 e0/ data an20cs( 3) / -0.0000052672 111140370 e0/ data an20cs( 4) / -0.0000001628 202185026 e0/ data an20cs( 5) / -0.0000000090 991442687 e0/ data an20cs( 6) / -0.0000000007 438647126 e0/ data an20cs( 7) / -0.0000000000 795494752 e0/ data an20cs( 8) / -0.0000000000 104050944 e0/ data an20cs( 9) / -0.0000000000 015932426 e0/ data an20cs( 10) / -0.0000000000 002770648 e0/ data an20cs( 11) / -0.0000000000 000535343 e0/ data an20cs( 12) / -0.0000000000 000113062 e0/ data an20cs( 13) / -0.0000000000 000025772 e0/ data an20cs( 14) / -0.0000000000 000006278 e0/ data an20cs( 15) / -0.0000000000 000001621 e0/ data an20cs( 16) / -0.0000000000 000000441 e0/ c c series for aph2 on the interval -1.00000e+00 to -1.25000e-01 c with weighted error 2.94e-17 c log weighted error 16.53 c significant figures required 15.58 c decimal places required 17.28 c data aph2cs( 1) / -0.2057088719 781465107 e0/ data aph2cs( 2) / 0.0422196961 357771922 e0/ data aph2cs( 3) / 0.0020482560 511207275 e0/ data aph2cs( 4) / 0.0002607800 735165006 e0/ data aph2cs( 5) / 0.0000474824 268004729 e0/ data aph2cs( 6) / 0.0000105102 756431612 e0/ data aph2cs( 7) / 0.0000026353 534014668 e0/ data aph2cs( 8) / 0.0000007208 824863499 e0/ data aph2cs( 9) / 0.0000002103 236664473 e0/ data aph2cs( 10) / 0.0000000644 975634555 e0/ data aph2cs( 11) / 0.0000000205 802377264 e0/ data aph2cs( 12) / 0.0000000067 836273921 e0/ data aph2cs( 13) / 0.0000000022 974015284 e0/ data aph2cs( 14) / 0.0000000007 961306765 e0/ data aph2cs( 15) / 0.0000000002 813860610 e0/ data aph2cs( 16) / 0.0000000001 011749057 e0/ data aph2cs( 17) / 0.0000000000 369306738 e0/ data aph2cs( 18) / 0.0000000000 136615066 e0/ data aph2cs( 19) / 0.0000000000 051142751 e0/ data aph2cs( 20) / 0.0000000000 019351689 e0/ data aph2cs( 21) / 0.0000000000 007393607 e0/ data aph2cs( 22) / 0.0000000000 002849792 e0/ data aph2cs( 23) / 0.0000000000 001107281 e0/ data aph2cs( 24) / 0.0000000000 000433412 e0/ data aph2cs( 25) / 0.0000000000 000170801 e0/ data aph2cs( 26) / 0.0000000000 000067733 e0/ data aph2cs( 27) / 0.0000000000 000027017 e0/ data aph2cs( 28) / 0.0000000000 000010835 e0/ data aph2cs( 29) / 0.0000000000 000004367 e0/ data aph2cs( 30) / 0.0000000000 000001769 e0/ data aph2cs( 31) / 0.0000000000 000000719 e0/ data aph2cs( 32) / 0.0000000000 000000294 e0/ c c series for aph1 on the interval -1.25000e-01 to -1.56250e-02 c with weighted error 6.38e-17 c log weighted error 16.20 c significant figures required 14.91 c decimal places required 16.87 c data aph1cs( 1) / -0.1024172908 077571694 e0/ data aph1cs( 2) / 0.0071697275 146591248 e0/ data aph1cs( 3) / 0.0001209959 363122329 e0/ data aph1cs( 4) / 0.0000073361 512841220 e0/ data aph1cs( 5) / 0.0000007535 382954272 e0/ data aph1cs( 6) / 0.0000001041 478171741 e0/ data aph1cs( 7) / 0.0000000174 358728519 e0/ data aph1cs( 8) / 0.0000000033 399795033 e0/ data aph1cs( 9) / 0.0000000007 073075174 e0/ data aph1cs( 10) / 0.0000000001 619187515 e0/ data aph1cs( 11) / 0.0000000000 394539982 e0/ data aph1cs( 12) / 0.0000000000 101192282 e0/ data aph1cs( 13) / 0.0000000000 027092778 e0/ data aph1cs( 14) / 0.0000000000 007523806 e0/ data aph1cs( 15) / 0.0000000000 002156369 e0/ data aph1cs( 16) / 0.0000000000 000635283 e0/ data aph1cs( 17) / 0.0000000000 000191757 e0/ data aph1cs( 18) / 0.0000000000 000059143 e0/ data aph1cs( 19) / 0.0000000000 000018597 e0/ data aph1cs( 20) / 0.0000000000 000005950 e0/ data aph1cs( 21) / 0.0000000000 000001934 e0/ data aph1cs( 22) / 0.0000000000 000000638 e0/ c c series for aph0 on the interval -1.56250e-02 to 0.00000e+00 c with weighted error 2.29e-17 c log weighted error 16.64 c significant figures required 15.27 c decimal places required 17.23 c data aph0cs( 1) / -0.0855849241 130933257 e0/ data aph0cs( 2) / 0.0011214378 867065261 e0/ data aph0cs( 3) / 0.0000042721 029353664 e0/ data aph0cs( 4) / 0.0000000817 607381483 e0/ data aph0cs( 5) / 0.0000000033 907645000 e0/ data aph0cs( 6) / 0.0000000002 253264423 e0/ data aph0cs( 7) / 0.0000000000 206284209 e0/ data aph0cs( 8) / 0.0000000000 023858763 e0/ data aph0cs( 9) / 0.0000000000 003301618 e0/ data aph0cs( 10) / 0.0000000000 000527010 e0/ data aph0cs( 11) / 0.0000000000 000094555 e0/ data aph0cs( 12) / 0.0000000000 000018709 e0/ data aph0cs( 13) / 0.0000000000 000004024 e0/ data aph0cs( 14) / 0.0000000000 000000930 e0/ data aph0cs( 15) / 0.0000000000 000000229 e0/ c data pi34 / 2.3561944 901923449e0 / c data pi34 / 2.3561944 9019234492 8846982537 4596271631 3d0 / data nan20, nan21, nan22, naph0, naph1, naph2 / 6*0 / data xsml / 0.0 / c if (nan20.ne.0) go to 10 eta = 0.1*r1mach(3) nan20 = inits (an20cs, 16, eta) nan21 = inits (an21cs, 24, eta) nan22 = inits (an22cs, 33, eta) naph0 = inits (aph0cs, 15, eta) naph1 = inits (aph1cs, 22, eta) naph2 = inits (aph2cs, 32, eta) c xsml = -(128.0/r1mach(3))**0.3333 c 10 if (x.ge.(-4.0)) go to 20 z = 1.0 if (x.gt.xsml) z = 128.0/x**3 + 1.0 ampl = 0.3125 + csevl (z, an20cs, nan20) phi = -0.625 + csevl (z, aph0cs, naph0) go to 40 c 20 if (x.ge.(-2.0)) go to 30 z = (128.0/x**3 + 9.0) / 7.0 ampl = 0.3125 + csevl (z, an21cs, nan21) phi = -0.625 + csevl (z, aph1cs, naph1) go to 40 c 30 if (x.gt.(-1.0)) call seteru ( 1 25hr9admp x must be le -1.0, 25, 1, 2) c z = (16.0/x**3 + 9.0) / 7.0 ampl = 0.3125 + csevl (z, an22cs, nan22) phi = -0.625 + csevl (z, aph2cs, naph2) c 40 sqrtx = sqrt (-x) ampl = sqrt (ampl*sqrtx) phi = pi34 - x*sqrtx*phi c return end |