Blame view
fvn_fnlib/bie.f
8.42 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 |
function bie (x)
c july 1977 edition. w. fullerton, c3, los alamos scientific lab.
c
c evaluate bi(x) for x .le. 0 and bi(x)*exp(-zeta) where
c zeta = 2/3 * x**(3/2) for x .ge. 0.0
c
dimension bifcs(9), bigcs(8), bif2cs(10), big2cs(10), bipcs(24),
1 bip2cs(29)
external csevl, inits, r1mach
c
c series for bif on the interval -1.00000d+00 to 1.00000d+00
c with weighted error 1.88e-19
c log weighted error 18.72
c significant figures required 17.74
c decimal places required 19.20
c
data bif cs( 1) / -.0167302164 7198664948e0 /
data bif cs( 2) / .1025233583 424944561e0 /
data bif cs( 3) / .0017083092 5073815165e0 /
data bif cs( 4) / .0000118625 4546774468e0 /
data bif cs( 5) / .0000000449 3290701779e0 /
data bif cs( 6) / .0000000001 0698207143e0 /
data bif cs( 7) / .0000000000 0017480643e0 /
data bif cs( 8) / .0000000000 0000020810e0 /
data bif cs( 9) / .0000000000 0000000018e0 /
c
c series for big on the interval -1.00000d+00 to 1.00000d+00
c with weighted error 2.61e-17
c log weighted error 16.58
c significant figures required 15.17
c decimal places required 17.03
c
data big cs( 1) / .0224662232 4857452e0 /
data big cs( 2) / .0373647754 5301955e0 /
data big cs( 3) / .0004447621 8957212e0 /
data big cs( 4) / .0000024708 0756363e0 /
data big cs( 5) / .0000000079 1913533e0 /
data big cs( 6) / .0000000000 1649807e0 /
data big cs( 7) / .0000000000 0002411e0 /
data big cs( 8) / .0000000000 0000002e0 /
c
c series for bif2 on the interval 1.00000d+00 to 8.00000d+00
c with weighted error 1.11e-17
c log weighted error 16.95
c approx significant figures required 16.5
c decimal places required 17.45
c
data bif2cs( 1) / 0.0998457269 3816041e0 /
data bif2cs( 2) / .4786249778 63005538e0 /
data bif2cs( 3) / .0251552119 604330118e0 /
data bif2cs( 4) / .0005820693 885232645e0 /
data bif2cs( 5) / .0000074997 659644377e0 /
data bif2cs( 6) / .0000000613 460287034e0 /
data bif2cs( 7) / .0000000003 462753885e0 /
data bif2cs( 8) / .0000000000 014288910e0 /
data bif2cs( 9) / .0000000000 000044962e0 /
data bif2cs(10) / .0000000000 000000111e0 /
c
c series for big2 on the interval 1.00000d+00 to 8.00000d+00
c with weighted error 1.19e-18
c log weighted error 17.92
c approx significant figures required 17.2
c decimal places required 18.42
c
data big2cs( 1) / .0333056621 45514340e0 /
data big2cs( 2) / .1613092151 23197068e0 /
data big2cs( 3) / .0063190073 096134286e0 /
data big2cs( 4) / .0001187904 568162517e0 /
data big2cs( 5) / .0000013045 345886200e0 /
data big2cs( 6) / .0000000093 741259955e0 /
data big2cs( 7) / .0000000000 474580188e0 /
data big2cs( 8) / .0000000000 001783107e0 /
data big2cs( 9) / .0000000000 000005167e0 /
data big2cs(10) / .0000000000 000000011e0 /
c
c series for bip on the interval 1.25000d-01 to 3.53553d-01
c with weighted error 1.91e-17
c log weighted error 16.72
c significant figures required 15.35
c decimal places required 17.41
c
data bip cs( 1) / -.0832204747 7943447e0 /
data bip cs( 2) / .0114611892 7371174e0 /
data bip cs( 3) / .0004289644 0718911e0 /
data bip cs( 4) / -.0001490663 9379950e0 /
data bip cs( 5) / -.0000130765 9726787e0 /
data bip cs( 6) / .0000063275 9839610e0 /
data bip cs( 7) / -.0000004222 6696982e0 /
data bip cs( 8) / -.0000001914 7186298e0 /
data bip cs( 9) / .0000000645 3106284e0 /
data bip cs(10) / -.0000000078 4485467e0 /
data bip cs(11) / -.0000000009 6077216e0 /
data bip cs(12) / .0000000007 0004713e0 /
data bip cs(13) / -.0000000001 7731789e0 /
data bip cs(14) / .0000000000 2272089e0 /
data bip cs(15) / .0000000000 0165404e0 /
data bip cs(16) / -.0000000000 0185171e0 /
data bip cs(17) / .0000000000 0059576e0 /
data bip cs(18) / -.0000000000 0012194e0 /
data bip cs(19) / .0000000000 0001334e0 /
data bip cs(20) / .0000000000 0000172e0 /
data bip cs(21) / -.0000000000 0000145e0 /
data bip cs(22) / .0000000000 0000049e0 /
data bip cs(23) / -.0000000000 0000011e0 /
data bip cs(24) / .0000000000 0000001e0 /
c
c series for bip2 on the interval 0. to 1.25000d-01
c with weighted error 1.05e-18
c log weighted error 17.98
c significant figures required 16.74
c decimal places required 18.71
c
data bip2cs( 1) / -.1135967375 85988679e0 /
data bip2cs( 2) / .0041381473 947881595e0 /
data bip2cs( 3) / .0001353470 622119332e0 /
data bip2cs( 4) / .0000104273 166530153e0 /
data bip2cs( 5) / .0000013474 954767849e0 /
data bip2cs( 6) / .0000001696 537405438e0 /
data bip2cs( 7) / -.0000000100 965008656e0 /
data bip2cs( 8) / -.0000000167 291194937e0 /
data bip2cs( 9) / -.0000000045 815364485e0 /
data bip2cs(10) / .0000000003 736681366e0 /
data bip2cs(11) / .0000000005 766930320e0 /
data bip2cs(12) / .0000000000 621812650e0 /
data bip2cs(13) / -.0000000000 632941202e0 /
data bip2cs(14) / -.0000000000 149150479e0 /
data bip2cs(15) / .0000000000 078896213e0 /
data bip2cs(16) / .0000000000 024960513e0 /
data bip2cs(17) / -.0000000000 012130075e0 /
data bip2cs(18) / -.0000000000 003740493e0 /
data bip2cs(19) / .0000000000 002237727e0 /
data bip2cs(20) / .0000000000 000474902e0 /
data bip2cs(21) / -.0000000000 000452616e0 /
data bip2cs(22) / -.0000000000 000030172e0 /
data bip2cs(23) / .0000000000 000091058e0 /
data bip2cs(24) / -.0000000000 000009814e0 /
data bip2cs(25) / -.0000000000 000016429e0 /
data bip2cs(26) / .0000000000 000005533e0 /
data bip2cs(27) / .0000000000 000002175e0 /
data bip2cs(28) / -.0000000000 000001737e0 /
data bip2cs(29) / -.0000000000 000000010e0 /
c
data atr / 8.750690570 8484345 e0 /
c atr = 16.0/(sqrt(8.)-1.) and btr = -(sqrt(8.)+1.)/(sqrt(8.)-1.)
data btr / -2.093836321 356054 e0 /
c
data nbif, nbig, nbif2, nbig2, nbip, nbip2 / 6*0 /
data x3sml, x32sml, xbig / 3*0.0 /
c
if (nbif.ne.0) go to 10
eta = 0.1*r1mach(3)
nbif = inits (bifcs, 9, eta)
nbig = inits (bigcs, 8, eta)
nbif2 = inits (bif2cs, 10, eta)
nbig2 = inits (big2cs, 10, eta)
nbip = inits (bipcs , 24, eta)
nbip2 = inits (bip2cs, 29, eta)
c
x3sml = eta**0.3333
x32sml = 1.3104*x3sml**2
xbig = r1mach(2)**0.6666
c
10 if (x.ge.(-1.0)) go to 20
call r9aimp (x, xm, theta)
bie = xm * sin(theta)
return
c
20 if (x.gt.1.0) go to 30
z = 0.0
if (abs(x).gt.x3sml) z = x**3
bie = 0.625 + csevl (z, bifcs, nbif) + x*(0.4375 +
1 csevl (z, bigcs, nbig))
if (x.gt.x32sml) bie = bie * exp(-2.0*x*sqrt(x)/3.0)
return
c
30 if (x.gt.2.0) go to 40
z = (2.0*x**3 - 9.0) / 7.0
bie = exp(-2.0*x*sqrt(x)/3.0) * (1.125 + csevl (z, bif2cs, nbif2)
1 + x*(0.625 + csevl (z, big2cs, nbig2)) )
return
c
40 if (x.gt.4.0) go to 50
sqrtx = sqrt(x)
z = atr/(x*sqrtx) + btr
bie = (0.625 + csevl (z, bipcs, nbip)) / sqrt(sqrtx)
return
c
50 sqrtx = sqrt(x)
z = -1.0
if (x.lt.xbig) z = 16.0/(x*sqrtx) - 1.0
bie = (0.625 + csevl (z, bip2cs, nbip2))/sqrt(sqrtx)
return
c
end
|