double precision function dbid (x) c july 1980 edition. w. fullerton, bell labs. c c evaluate the derivative of the airy function bi(x). c double precision x, bifcs(13), bigcs(13), bif2cs(15), big2cs(16), 1 x2sml, x3sml, xmax, x3, x2, z, xn, phi, d1mach, dbide, dcsevl, 2 dexp, dlog, dsin, dsqrt external d1mach, dbide, dcsevl, initds c 1 sqrt c c series for bif on the interval -1.00000e+00 to 1.00000e+00 c with weighted error 3.82e-34 c log weighted error 33.42 c significant figures required 32.21 c decimal places required 33.98 c data bif cs( 1) / 0.1153536790 8285702426 7474446284 908879d0/ data bif cs( 2) / 0.0205007894 0491928753 0357789445 940252d0/ data bif cs( 3) / 0.0002135290 2789028758 1892679619 451158d0/ data bif cs( 4) / 0.0000010783 9606146768 3042209155 523569d0/ data bif cs( 5) / 0.0000000032 0947088332 0666783353 670420d0/ data bif cs( 6) / 0.0000000000 0629304067 1833540390 213316d0/ data bif cs( 7) / 0.0000000000 0000874030 4300063083 340121d0/ data bif cs( 8) / 0.0000000000 0000000904 7915683496 049529d0/ data bif cs( 9) / 0.0000000000 0000000000 7249923164 709251d0/ data bif cs( 10) / 0.0000000000 0000000000 0004629576 649604d0/ data bif cs( 11) / 0.0000000000 0000000000 0000002411 236436d0/ data bif cs( 12) / 0.0000000000 0000000000 0000000001 043825d0/ data bif cs( 13) / 0.0000000000 0000000000 0000000000 000382d0/ c c series for big on the interval -1.00000e+00 to 1.00000e+00 c with weighted error 4.79e-32 c log weighted error 31.32 c significant figures required 30.52 c decimal places required 31.88 c data big cs( 1) / -0.0971964404 1644353738 9779097460 6802d0/ data big cs( 2) / 0.1495035768 4316706657 1084344532 6264d0/ data big cs( 3) / 0.0031135253 8712132604 1941917683 9631d0/ data big cs( 4) / 0.0000247085 7057982129 6777702192 0569d0/ data big cs( 5) / 0.0000001029 4962773137 8608198732 4295d0/ data big cs( 6) / 0.0000000002 6397037398 6943289267 6778d0/ data big cs( 7) / 0.0000000000 0045827927 0780320660 8181d0/ data big cs( 8) / 0.0000000000 0000057428 2974089344 7321d0/ data big cs( 9) / 0.0000000000 0000000054 3827538523 8549d0/ data big cs( 10) / 0.0000000000 0000000000 0402834726 7083d0/ data big cs( 11) / 0.0000000000 0000000000 0000239782 3826d0/ data big cs( 12) / 0.0000000000 0000000000 0000000117 1956d0/ data big cs( 13) / 0.0000000000 0000000000 0000000000 0479d0/ c c series for bif2 on the interval 1.00000e+00 to 8.00000e+00 c with weighted error 1.38e-33 c log weighted error 32.86 c significant figures required 32.12 c decimal places required 33.45 c data bif2cs( 1) / 0.3234939876 0352203352 1191935962 66015d0/ data bif2cs( 2) / 0.0862978715 3556355913 8888353238 11100d0/ data bif2cs( 3) / 0.0029940255 5265539742 6138210507 27155d0/ data bif2cs( 4) / 0.0000514305 2836466163 7204643169 50821d0/ data bif2cs( 5) / 0.0000005258 4025003681 1460260330 98613d0/ data bif2cs( 6) / 0.0000000035 6175137395 7700281027 30600d0/ data bif2cs( 7) / 0.0000000000 1714686400 7145848305 18308d0/ data bif2cs( 8) / 0.0000000000 0006166351 9692325554 06693d0/ data bif2cs( 9) / 0.0000000000 0000017191 0821543159 85806d0/ data bif2cs( 10) / 0.0000000000 0000000038 2368895188 03943d0/ data bif2cs( 11) / 0.0000000000 0000000000 0694241736 24884d0/ data bif2cs( 12) / 0.0000000000 0000000000 0001048339 32510d0/ data bif2cs( 13) / 0.0000000000 0000000000 0000001337 21972d0/ data bif2cs( 14) / 0.0000000000 0000000000 0000000001 45986d0/ data bif2cs( 15) / 0.0000000000 0000000000 0000000000 00138d0/ c c series for big2 on the interval 1.00000e+00 to 8.00000e+00 c with weighted error 1.91e-34 c log weighted error 33.72 c significant figures required 33.76 c decimal places required 34.32 c data big2cs( 1) / 1.6062999463 6212945775 9284537862 622883d0/ data big2cs( 2) / 0.7449088819 8760886520 1476685194 753972d0/ data big2cs( 3) / 0.0470138738 6102773796 4095177635 353019d0/ data big2cs( 4) / 0.0012284422 0625482390 7016188785 848091d0/ data big2cs( 5) / 0.0000173222 4122566236 2670987355 613727d0/ data big2cs( 6) / 0.0000001521 9016523680 1893711508 366563d0/ data big2cs( 7) / 0.0000000009 1135602491 1957704145 528786d0/ data big2cs( 8) / 0.0000000000 0395479184 2356644201 722554d0/ data big2cs( 9) / 0.0000000000 0001300173 7033862320 007309d0/ data big2cs( 10) / 0.0000000000 0000003349 3506858269 079763d0/ data big2cs( 11) / 0.0000000000 0000000006 9419094403 694057d0/ data big2cs( 12) / 0.0000000000 0000000000 0118248256 604581d0/ data big2cs( 13) / 0.0000000000 0000000000 0000168462 493472d0/ data big2cs( 14) / 0.0000000000 0000000000 0000000203 684674d0/ data big2cs( 15) / 0.0000000000 0000000000 0000000000 211619d0/ data big2cs( 16) / 0.0000000000 0000000000 0000000000 000191d0/ c data nbif, nbig, nbif2, nbig2 / 4*0 / data x2sml, x3sml, xmax / 3*0.0d0 / c if (nbif.ne.0) go to 10 eta = 0.1*sngl(d1mach(3)) nbif = initds (bifcs, 13, eta) nbig = initds (bigcs, 13, eta) nbif2 = initds (bif2cs, 15, eta) nbig2 = initds (big2cs, 16, eta) c x2sml = sqrt (eta) x3sml = eta**0.3333 xmax = (1.5d0*dlog(d1mach(2)))**0.6666d0 c 10 if (x.ge.(-1.0d0)) go to 20 call d9admp (x, xn, phi) dbid = xn * dsin (phi) return c 20 if (x.gt.1.0d0) go to 30 x3 = 0.0d0 if (dabs(x).gt.x3sml) x3 = x**3 x2 = 0.0d0 if (dabs(x).gt.x2sml) x2 = x*x dbid = x2*(dcsevl (x3, bifcs, nbif) + 0.25d0) + 1 dcsevl (x3, bigcs, nbig) + 0.5d0 return c 30 if (x.gt.2.0d0) go to 40 z = (2.0d0*x**3 - 9.0d0) / 7.0d0 dbid = x*x*(dcsevl (z, bif2cs, nbif2) + 0.25d0) + 1 dcsevl (z, big2cs, nbig2) + 0.5d0 return c 40 if (x.gt.xmax) call seteru ( 1 36hdbid x so big that dbid overflows, 36, 1, 2) c dbid = dbide(x) * dexp (2.0d0*x*dsqrt(x)/3.0d0) return c end