Blame view
fvn_fnlib/bide.f
8.86 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 |
function bide (x)
c july 1980 edition. w. fullerton, bell labs.
c
c evaluate the dervative of the airy function bi(x) for x.lt.0
c and evaluate bid(x) * exp(-zeta) for x .ge. 0 where
c zeta = 2/3 * x**(3/2)
c
dimension bifcs(8), bigcs(9), bif2cs(10), big2cs(10), bip1cs(24),
1 bip2cs(29)
external csevl, inits, r1mach
c
c series for bif on the interval -1.00000e+00 to 1.00000e+00
c with weighted error 9.05e-18
c log weighted error 17.04
c significant figures required 15.83
c decimal places required 17.49
c
data bifcs( 1) / 0.1153536790 828570243 e0/
data bifcs( 2) / 0.0205007894 049192875 e0/
data bifcs( 3) / 0.0002135290 278902876 e0/
data bifcs( 4) / 0.0000010783 960614677 e0/
data bifcs( 5) / 0.0000000032 094708833 e0/
data bifcs( 6) / 0.0000000000 062930407 e0/
data bifcs( 7) / 0.0000000000 000087403 e0/
data bifcs( 8) / 0.0000000000 000000090 e0/
c
c series for big on the interval -1.00000e+00 to 1.00000e+00
c with weighted error 5.44e-19
c log weighted error 18.26
c significant figures required 17.46
c decimal places required 18.74
c
data bigcs( 1) / -0.0971964404 1644353739 0e0/
data bigcs( 2) / 0.1495035768 4316706657 1e0/
data bigcs( 3) / 0.0031135253 8712132604 2e0/
data bigcs( 4) / 0.0000247085 7057982129 7e0/
data bigcs( 5) / 0.0000001029 4962773137 9e0/
data bigcs( 6) / 0.0000000002 6397037398 7e0/
data bigcs( 7) / 0.0000000000 0045827927 1e0/
data bigcs( 8) / 0.0000000000 0000057428 3e0/
data bigcs( 9) / 0.0000000000 0000000054 4e0/
c
c series for bif2 on the interval 1.00000e+00 to 8.00000e+00
c with weighted error 3.82e-19
c log weighted error 18.42
c significant figures required 17.68
c decimal places required 18.92
c
data bif2cs( 1) / 0.3234939876 0352203352 1e0/
data bif2cs( 2) / 0.0862978715 3556355913 9e0/
data bif2cs( 3) / 0.0029940255 5265539742 6e0/
data bif2cs( 4) / 0.0000514305 2836466163 7e0/
data bif2cs( 5) / 0.0000005258 4025003681 1e0/
data bif2cs( 6) / 0.0000000035 6175137395 8e0/
data bif2cs( 7) / 0.0000000000 1714686400 7e0/
data bif2cs( 8) / 0.0000000000 0006166352 0e0/
data bif2cs( 9) / 0.0000000000 0000017191 1e0/
data bif2cs( 10) / 0.0000000000 0000000038 2e0/
c
c series for big2 on the interval 1.00000e+00 to 8.00000e+00
c with weighted error 3.35e-17
c log weighted error 16.48
c significant figures required 16.52
c decimal places required 16.98
c
data big2cs( 1) / 1.6062999463 621294578 e0/
data big2cs( 2) / 0.7449088819 876088652 e0/
data big2cs( 3) / 0.0470138738 610277380 e0/
data big2cs( 4) / 0.0012284422 062548239 e0/
data big2cs( 5) / 0.0000173222 412256624 e0/
data big2cs( 6) / 0.0000001521 901652368 e0/
data big2cs( 7) / 0.0000000009 113560249 e0/
data big2cs( 8) / 0.0000000000 039547918 e0/
data big2cs( 9) / 0.0000000000 000130017 e0/
data big2cs( 10) / 0.0000000000 000000335 e0/
c
c series for bip2 on the interval 0.00000e+00 to 1.25000e-01
c with weighted error 2.07e-18
c log weighted error 17.69
c significant figures required 16.51
c decimal places required 18.42
c
data bip2cs( 1) / -0.1326970544 3526630495 e0/
data bip2cs( 2) / -0.0056844362 6045977481 e0/
data bip2cs( 3) / -0.0001564360 1119611610 e0/
data bip2cs( 4) / -0.0000113673 7203679562 e0/
data bip2cs( 5) / -0.0000014346 4350991284 e0/
data bip2cs( 6) / -0.0000001809 8531185164 e0/
data bip2cs( 7) / 0.0000000092 6177343611 e0/
data bip2cs( 8) / 0.0000000171 0005490721 e0/
data bip2cs( 9) / 0.0000000047 6698163504 e0/
data bip2cs( 10) / -0.0000000003 5195022023 e0/
data bip2cs( 11) / -0.0000000005 8890614316 e0/
data bip2cs( 12) / -0.0000000000 6678499608 e0/
data bip2cs( 13) / 0.0000000000 6395565102 e0/
data bip2cs( 14) / 0.0000000000 1554529427 e0/
data bip2cs( 15) / -0.0000000000 0792397000 e0/
data bip2cs( 16) / -0.0000000000 0258326243 e0/
data bip2cs( 17) / 0.0000000000 0121655048 e0/
data bip2cs( 18) / 0.0000000000 0038707207 e0/
data bip2cs( 19) / -0.0000000000 0022487045 e0/
data bip2cs( 20) / -0.0000000000 0004953477 e0/
data bip2cs( 21) / 0.0000000000 0004563782 e0/
data bip2cs( 22) / 0.0000000000 0000332998 e0/
data bip2cs( 23) / -0.0000000000 0000921750 e0/
data bip2cs( 24) / 0.0000000000 0000094157 e0/
data bip2cs( 25) / 0.0000000000 0000167154 e0/
data bip2cs( 26) / -0.0000000000 0000055134 e0/
data bip2cs( 27) / -0.0000000000 0000022369 e0/
data bip2cs( 28) / 0.0000000000 0000017487 e0/
data bip2cs( 29) / 0.0000000000 0000000207 e0/
c
c series for bip1 on the interval 1.25000e-01 to 3.53553e-01
c with weighted error 1.86e-17
c log weighted error 16.73
c significant figures required 15.67
c decimal places required 17.42
c
data bip1cs( 1) / -0.1729187351 079553719 e0/
data bip1cs( 2) / -0.0149358492 984694364 e0/
data bip1cs( 3) / -0.0005471104 951678566 e0/
data bip1cs( 4) / 0.0001537966 292958408 e0/
data bip1cs( 5) / 0.0000154353 476192179 e0/
data bip1cs( 6) / -0.0000065434 113851906 e0/
data bip1cs( 7) / 0.0000003728 082407879 e0/
data bip1cs( 8) / 0.0000002072 078388189 e0/
data bip1cs( 9) / -0.0000000658 173336470 e0/
data bip1cs( 10) / 0.0000000074 926746354 e0/
data bip1cs( 11) / 0.0000000011 101336884 e0/
data bip1cs( 12) / -0.0000000007 265140553 e0/
data bip1cs( 13) / 0.0000000001 782723560 e0/
data bip1cs( 14) / -0.0000000000 217346352 e0/
data bip1cs( 15) / -0.0000000000 020302035 e0/
data bip1cs( 16) / 0.0000000000 019311827 e0/
data bip1cs( 17) / -0.0000000000 006044953 e0/
data bip1cs( 18) / 0.0000000000 001209450 e0/
data bip1cs( 19) / -0.0000000000 000125109 e0/
data bip1cs( 20) / -0.0000000000 000019917 e0/
data bip1cs( 21) / 0.0000000000 000015154 e0/
data bip1cs( 22) / -0.0000000000 000004977 e0/
data bip1cs( 23) / 0.0000000000 000001155 e0/
data bip1cs( 24) / -0.0000000000 000000186 e0/
c
data atr / 8.750690570 8484345 e0 /
c atr = 16./(sqrt(8)-1) and btr = -(sqrt(8)+1)/(sqrt(8)-1)
data btr /-2.093836321 3560543 e0 /
c
data nbif, nbig, nbif2, nbig2, nbip1, nbip2 / 6*0 /
data x2sml, x3sml, x32sml, xbig / 4*0.0 /
c
if (nbif.ne.0) go to 10
eta = 0.1*r1mach(3)
nbif = inits (bifcs, 8, eta)
nbig = inits (bigcs, 9, eta)
nbif2 = inits (bif2cs, 10, eta)
nbig2 = inits (big2cs, 10, eta)
nbip1 = inits (bip1cs, 24, eta)
nbip2 = inits (bip2cs, 29, eta)
c
x2sml = sqrt (eta)
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 r9admp (x, xn, phi)
bide = xn * sin (phi)
return
c
20 if (x.gt.1.0) go to 30
x3 = 0.0
if (abs(x).gt.x3sml) x3 = x**3
x2 = 0.0
if (abs(x).gt.x2sml) x2 = x*x
bide = x2 * (csevl(x3, bifcs, nbif) + 0.25) +
1 csevl (x3, bigcs, nbig) + 0.5
if (x.gt.x32sml) bide = bide * 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
bide = exp (-2.0*x*sqrt(x)/3.0) * (x*x * (0.25 +
1 csevl (z, bif2cs, nbif2)) + 0.5 + csevl (z, big2cs, nbig2))
return
c
40 if (x.gt.4.0) go to 50
sqrtx = sqrt(x)
z = atr/(x*sqrtx) + btr
bide = (0.625 + csevl (z, bip1cs, nbip1)) * sqrt(sqrtx)
return
c
50 sqrtx = sqrt(x)
z = -1.0
if (x.lt.xbig) z = 16.0/(x*sqrtx) - 1.0
bide = (0.625 + csevl (z, bip2cs, nbip2)) * sqrt(sqrtx)
return
c
end
|