Blame view
fvn_fnlib/d9b0mp.f
14 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 |
subroutine d9b0mp (x, ampl, theta)
c july 1977 edition. w. fullerton, c3, los alamos scientific lab.
c
c evaluate the modulus and phase for the bessel j0 and y0 functions.
c
double precision x, ampl, theta, bm0cs(37), bt02cs(39),
1 bm02cs(40), bth0cs(44), xmax, pi4, z, d1mach, dcsevl,
2 dsqrt
external d1mach, dcsevl, initds
c
c series for bm0 on the interval 1.56250e-02 to 6.25000e-02
c with weighted error 4.40e-32
c log weighted error 31.36
c significant figures required 30.02
c decimal places required 32.14
c
data bm0 cs( 1) / +.9211656246 8277427125 7376773018 2 d-1 /
data bm0 cs( 2) / -.1050590997 2719051024 8071637175 5 d-2 /
data bm0 cs( 3) / +.1470159840 7687597540 5639285095 2 d-4 /
data bm0 cs( 4) / -.5058557606 0385542233 4792932770 2 d-6 /
data bm0 cs( 5) / +.2787254538 6324441766 3035613788 1 d-7 /
data bm0 cs( 6) / -.2062363611 7809148026 1884101897 3 d-8 /
data bm0 cs( 7) / +.1870214313 1388796751 3817259626 1 d-9 /
data bm0 cs( 8) / -.1969330971 1356362002 4173077782 5 d-10 /
data bm0 cs( 9) / +.2325973793 9992754440 1250881805 2 d-11 /
data bm0 cs( 10) / -.3009520344 9382502728 5122473448 2 d-12 /
data bm0 cs( 11) / +.4194521333 8506691814 7120676864 6 d-13 /
data bm0 cs( 12) / -.6219449312 1884458259 7326742956 4 d-14 /
data bm0 cs( 13) / +.9718260411 3360684696 0176588526 9 d-15 /
data bm0 cs( 14) / -.1588478585 7010752073 6663596693 7 d-15 /
data bm0 cs( 15) / +.2700072193 6713088900 8621732445 8 d-16 /
data bm0 cs( 16) / -.4750092365 2340089924 7750478677 3 d-17 /
data bm0 cs( 17) / +.8615128162 6043708731 9170374656 0 d-18 /
data bm0 cs( 18) / -.1605608686 9561448157 4560270335 9 d-18 /
data bm0 cs( 19) / +.3066513987 3144829751 8853980159 9 d-19 /
data bm0 cs( 20) / -.5987764223 1939564306 9650561706 6 d-20 /
data bm0 cs( 21) / +.1192971253 7482483064 8906984106 6 d-20 /
data bm0 cs( 22) / -.2420969142 0448054894 8468258133 3 d-21 /
data bm0 cs( 23) / +.4996751760 5106164533 7100287999 9 d-22 /
data bm0 cs( 24) / -.1047493639 3511585100 9504051199 9 d-22 /
data bm0 cs( 25) / +.2227786843 7974681010 4818346666 6 d-23 /
data bm0 cs( 26) / -.4801813239 3981628623 7054293333 3 d-24 /
data bm0 cs( 27) / +.1047962723 4709599564 7699626666 6 d-24 /
data bm0 cs( 28) / -.2313858165 6786153251 0126080000 0 d-25 /
data bm0 cs( 29) / +.5164823088 4626742116 3519999999 9 d-26 /
data bm0 cs( 30) / -.1164691191 8500653895 2540159999 9 d-26 /
data bm0 cs( 31) / +.2651788486 0433192829 5833600000 0 d-27 /
data bm0 cs( 32) / -.6092559503 8257284976 9130666666 6 d-28 /
data bm0 cs( 33) / +.1411804686 1442593080 3882666666 6 d-28 /
data bm0 cs( 34) / -.3298094961 2317372457 5061333333 3 d-29 /
data bm0 cs( 35) / +.7763931143 0740650317 1413333333 3 d-30 /
data bm0 cs( 36) / -.1841031343 6614584784 2133333333 3 d-30 /
data bm0 cs( 37) / +.4395880138 5943107371 0079999999 9 d-31 /
c
c series for bth0 on the interval 0. to 1.56250e-02
c with weighted error 2.66e-32
c log weighted error 31.57
c significant figures required 30.67
c decimal places required 32.40
c
data bth0cs( 1) / -.2490178086 2128936717 7097937899 67 d+0 /
data bth0cs( 2) / +.4855029960 9623749241 0486155354 85 d-3 /
data bth0cs( 3) / -.5451183734 5017204950 6562735635 05 d-5 /
data bth0cs( 4) / +.1355867305 9405964054 3774459299 03 d-6 /
data bth0cs( 5) / -.5569139890 2227626227 5832184149 20 d-8 /
data bth0cs( 6) / +.3260903182 4994335304 0042057194 68 d-9 /
data bth0cs( 7) / -.2491880786 2461341125 2379038779 93 d-10 /
data bth0cs( 8) / +.2344937742 0882520554 3524135648 91 d-11 /
data bth0cs( 9) / -.2609653444 4310387762 1775747661 36 d-12 /
data bth0cs( 10) / +.3335314042 0097395105 8699550149 23 d-13 /
data bth0cs( 11) / -.4789000044 0572684646 7507705574 09 d-14 /
data bth0cs( 12) / +.7595617843 6192215972 6425685452 48 d-15 /
data bth0cs( 13) / -.1313155601 6891440382 7733974876 33 d-15 /
data bth0cs( 14) / +.2448361834 5240857495 4268207383 55 d-16 /
data bth0cs( 15) / -.4880572981 0618777683 2567619183 31 d-17 /
data bth0cs( 16) / +.1032728502 9786316149 2237563612 04 d-17 /
data bth0cs( 17) / -.2305763381 5057217157 0047445270 25 d-18 /
data bth0cs( 18) / +.5404444300 1892693993 0171084837 65 d-19 /
data bth0cs( 19) / -.1324069519 4366572724 1550328823 85 d-19 /
data bth0cs( 20) / +.3378079562 1371970203 4247921247 22 d-20 /
data bth0cs( 21) / -.8945762915 7111779003 0269262922 99 d-21 /
data bth0cs( 22) / +.2451990688 9219317090 8999086514 05 d-21 /
data bth0cs( 23) / -.6938842287 6866318680 1399331576 57 d-22 /
data bth0cs( 24) / +.2022827871 4890138392 9463033377 91 d-22 /
data bth0cs( 25) / -.6062850000 2335483105 7941953717 64 d-23 /
data bth0cs( 26) / +.1864974896 4037635381 8237883962 70 d-23 /
data bth0cs( 27) / -.5878373238 4849894560 2450365308 67 d-24 /
data bth0cs( 28) / +.1895859144 7999563485 5311795035 13 d-24 /
data bth0cs( 29) / -.6248197937 2258858959 2916207285 65 d-25 /
data bth0cs( 30) / +.2101790168 4551024686 6386335290 74 d-25 /
data bth0cs( 31) / -.7208430093 5209253690 8139339924 46 d-26 /
data bth0cs( 32) / +.2518136389 2474240867 1564059767 46 d-26 /
data bth0cs( 33) / -.8951804225 8785778806 1439459536 43 d-27 /
data bth0cs( 34) / +.3235723747 9762298533 2562358685 87 d-27 /
data bth0cs( 35) / -.1188301051 9855353657 0471441137 96 d-27 /
data bth0cs( 36) / +.4430628690 7358104820 5792319417 31 d-28 /
data bth0cs( 37) / -.1676100964 8834829495 7920101356 81 d-28 /
data bth0cs( 38) / +.6429294692 1207466972 5323939660 88 d-29 /
data bth0cs( 39) / -.2499226116 6978652421 2072136827 63 d-29 /
data bth0cs( 40) / +.9839979429 9521955672 8282603553 18 d-30 /
data bth0cs( 41) / -.3922037524 2408016397 9891316261 58 d-30 /
data bth0cs( 42) / +.1581810703 0056522138 5906188456 92 d-30 /
data bth0cs( 43) / -.6452550614 4890715944 3440983654 26 d-31 /
data bth0cs( 44) / +.2661111136 9199356137 1770183463 67 d-31 /
c
c series for bm02 on the interval 0. to 1.56250e-02
c with weighted error 4.72e-32
c log weighted error 31.33
c significant figures required 30.00
c decimal places required 32.13
c
data bm02cs( 1) / +.9500415145 2283813693 3086133556 0 d-1 /
data bm02cs( 2) / -.3801864682 3656709917 4808156685 1 d-3 /
data bm02cs( 3) / +.2258339301 0314811929 5182992722 4 d-5 /
data bm02cs( 4) / -.3895725802 3722287647 3062141260 5 d-7 /
data bm02cs( 5) / +.1246886416 5120816979 3099052972 5 d-8 /
data bm02cs( 6) / -.6065949022 1025037798 0383505838 7 d-10 /
data bm02cs( 7) / +.4008461651 4217469910 1527597104 5 d-11 /
data bm02cs( 8) / -.3350998183 3980942184 6729879457 4 d-12 /
data bm02cs( 9) / +.3377119716 5174173670 6326434199 6 d-13 /
data bm02cs( 10) / -.3964585901 6350127005 6935629582 3 d-14 /
data bm02cs( 11) / +.5286111503 8838572173 8793974473 5 d-15 /
data bm02cs( 12) / -.7852519083 4508523136 5464024349 3 d-16 /
data bm02cs( 13) / +.1280300573 3866822010 1163407344 9 d-16 /
data bm02cs( 14) / -.2263996296 3914297762 8709924488 4 d-17 /
data bm02cs( 15) / +.4300496929 6567903886 4641029047 7 d-18 /
data bm02cs( 16) / -.8705749805 1325870797 4753545145 5 d-19 /
data bm02cs( 17) / +.1865862713 9620951411 8144277205 0 d-19 /
data bm02cs( 18) / -.4210482486 0930654573 4508697230 1 d-20 /
data bm02cs( 19) / +.9956676964 2284009915 8162741784 2 d-21 /
data bm02cs( 20) / -.2457357442 8053133596 0592147854 7 d-21 /
data bm02cs( 21) / +.6307692160 7620315680 8735370705 9 d-22 /
data bm02cs( 22) / -.1678773691 4407401426 9333117238 8 d-22 /
data bm02cs( 23) / +.4620259064 6739044337 7087813608 7 d-23 /
data bm02cs( 24) / -.1311782266 8603087322 3769340249 6 d-23 /
data bm02cs( 25) / +.3834087564 1163028277 4792244027 6 d-24 /
data bm02cs( 26) / -.1151459324 0777412710 7261329357 6 d-24 /
data bm02cs( 27) / +.3547210007 5233385230 7697134521 3 d-25 /
data bm02cs( 28) / -.1119218385 8150046462 6435594217 6 d-25 /
data bm02cs( 29) / +.3611879427 6298378316 9840499425 7 d-26 /
data bm02cs( 30) / -.1190687765 9133331500 9264176246 3 d-26 /
data bm02cs( 31) / +.4005094059 4039681318 0247644953 6 d-27 /
data bm02cs( 32) / -.1373169422 4522123905 9519391601 7 d-27 /
data bm02cs( 33) / +.4794199088 7425315859 9649152643 7 d-28 /
data bm02cs( 34) / -.1702965627 6241095840 0699447645 2 d-28 /
data bm02cs( 35) / +.6149512428 9363300715 0357516132 4 d-29 /
data bm02cs( 36) / -.2255766896 5818283499 4430023724 2 d-29 /
data bm02cs( 37) / +.8399707509 2942994860 6165835320 0 d-30 /
data bm02cs( 38) / -.3172997595 5626023555 6742393615 2 d-30 /
data bm02cs( 39) / +.1215205298 8812985545 8333302651 4 d-30 /
data bm02cs( 40) / -.4715852749 7544386930 1321056804 5 d-31 /
c
c series for bt02 on the interval 1.56250e-02 to 6.25000e-02
c with weighted error 2.99e-32
c log weighted error 31.52
c significant figures required 30.61
c decimal places required 32.32
c
data bt02cs( 1) / -.2454829521 3424597462 0504672493 24 d+0 /
data bt02cs( 2) / +.1254412103 9084615780 7853317782 99 d-2 /
data bt02cs( 3) / -.3125395041 4871522854 9734467095 71 d-4 /
data bt02cs( 4) / +.1470977824 9940831164 4534269693 14 d-5 /
data bt02cs( 5) / -.9954348893 7950033643 4688503511 58 d-7 /
data bt02cs( 6) / +.8549316673 3203041247 5787113977 51 d-8 /
data bt02cs( 7) / -.8698975952 6554334557 9855121791 92 d-9 /
data bt02cs( 8) / +.1005209953 3559791084 5401010821 53 d-9 /
data bt02cs( 9) / -.1282823060 1708892903 4836236855 44 d-10 /
data bt02cs( 10) / +.1773170078 1805131705 6557504510 23 d-11 /
data bt02cs( 11) / -.2617457456 9485577488 6362841809 25 d-12 /
data bt02cs( 12) / +.4082835138 9972059621 9664812211 03 d-13 /
data bt02cs( 13) / -.6675166823 9742720054 6067495542 61 d-14 /
data bt02cs( 14) / +.1136576139 3071629448 3924695499 51 d-14 /
data bt02cs( 15) / -.2005118962 0647160250 5592664121 17 d-15 /
data bt02cs( 16) / +.3649797879 4766269635 7205914641 06 d-16 /
data bt02cs( 17) / -.6830963756 4582303169 3558437888 00 d-17 /
data bt02cs( 18) / +.1310758314 5670756620 0571042679 46 d-17 /
data bt02cs( 19) / -.2572336310 1850607778 7571306495 99 d-18 /
data bt02cs( 20) / +.5152165744 1863959925 2677809493 33 d-19 /
data bt02cs( 21) / -.1051301756 3758802637 9407414613 33 d-19 /
data bt02cs( 22) / +.2182038199 1194813847 3010845013 33 d-20 /
data bt02cs( 23) / -.4600470121 0362160577 2259054933 33 d-21 /
data bt02cs( 24) / +.9840700692 5466818520 9536511999 99 d-22 /
data bt02cs( 25) / -.2133403803 5728375844 7359863466 66 d-22 /
data bt02cs( 26) / +.4683103642 3973365296 0662869333 33 d-23 /
data bt02cs( 27) / -.1040021369 1985747236 5133823999 99 d-23 /
data bt02cs( 28) / +.2334910567 7301510051 7777408000 00 d-24 /
data bt02cs( 29) / -.5295682532 3318615788 0497493333 33 d-25 /
data bt02cs( 30) / +.1212634195 2959756829 1962879999 99 d-25 /
data bt02cs( 31) / -.2801889708 2289428760 2756266666 66 d-26 /
data bt02cs( 32) / +.6529267898 7012873342 5937066666 66 d-27 /
data bt02cs( 33) / -.1533798006 1873346427 8357333333 33 d-27 /
data bt02cs( 34) / +.3630588430 6364536682 3594666666 66 d-28 /
data bt02cs( 35) / -.8656075571 3629122479 1722666666 66 d-29 /
data bt02cs( 36) / +.2077990997 2536284571 2383999999 99 d-29 /
data bt02cs( 37) / -.5021117022 1417221674 3253333333 33 d-30 /
data bt02cs( 38) / +.1220836027 9441714184 1919999999 99 d-30 /
data bt02cs( 39) / -.2986005626 7039913454 2506666666 66 d-31 /
c
data pi4 / 0.7853981633 9744830961 5660845819 876 d0 /
data nbm0, nbt02, nbm02, nbth0, xmax / 4*0, 0.d0 /
c
if (nbm0.ne.0) go to 10
eta = 0.1*sngl(d1mach(3))
nbm0 = initds (bm0cs, 37, eta)
nbt02 = initds (bt02cs, 39, eta)
nbm02 = initds (bm02cs, 40, eta)
nbth0 = initds (bth0cs, 44, eta)
c
xmax = 1.0d0/d1mach(4)
c
10 if (x.lt.4.d0) call seteru (22hd9b0mp x must be ge 4, 22, 1, 2)
c
if (x.gt.8.d0) go to 20
z = (128.d0/(x*x) - 5.d0)/3.d0
ampl = (.75d0 + dcsevl (z, bm0cs, nbm0))/dsqrt(x)
theta = x - pi4 + dcsevl (z, bt02cs, nbt02)/x
return
c
20 if (x.gt.xmax) call seteru (
1 37hd9b0mp no precision because x is big, 37, 2, 2)
c
z = 128.d0/(x*x) - 1.d0
ampl = (.75d0 + dcsevl (z, bm02cs, nbm02))/dsqrt(x)
theta = x - pi4 + dcsevl (z, bth0cs, nbth0)/x
return
c
end
|