function aide (x) c july 1980 edition. w. fullerton, bell labs. c c evaluate the derivative of the airy function for x .le. 0 c and evaluate aid(x)*exp(zeta) for x .ge. 0 where c zeta = 2/3 * x**(3/2) c dimension aifcs(8), aigcs(9), aip1cs(25), aip2cs(15) external csevl, inits, r1mach c c series for aif on the interval -1.00000e+00 to 1.00000e+00 c with weighted error 5.22e-18 c log weighted error 17.28 c significant figures required 16.01 c decimal places required 17.73 c data aifcs( 1) / 0.1052746122 6531408809 e0/ data aifcs( 2) / 0.0118361362 8152997844 e0/ data aifcs( 3) / 0.0001232810 4173225664 e0/ data aifcs( 4) / 0.0000006226 1225638140 e0/ data aifcs( 5) / 0.0000000018 5298887844 e0/ data aifcs( 6) / 0.0000000000 0363328873 e0/ data aifcs( 7) / 0.0000000000 0000504622 e0/ data aifcs( 8) / 0.0000000000 0000000522 e0/ c c series for aig on the interval -1.00000e+00 to 1.00000e+00 c with weighted error 3.14e-19 c log weighted error 18.50 c significant figures required 17.44 c decimal places required 18.98 c data aigcs( 1) / 0.0212338781 5091866685 2e0/ data aigcs( 2) / 0.0863159303 3521440675 2e0/ data aigcs( 3) / 0.0017975947 2038323135 8e0/ data aigcs( 4) / 0.0000142654 9987555069 3e0/ data aigcs( 5) / 0.0000000594 3799528368 3e0/ data aigcs( 6) / 0.0000000001 5240336647 9e0/ data aigcs( 7) / 0.0000000000 0026458766 0e0/ data aigcs( 8) / 0.0000000000 0000033156 2e0/ data aigcs( 9) / 0.0000000000 0000000031 4e0/ c c series for aip2 on the interval 0.00000e+00 to 1.25000e-01 c with weighted error 2.15e-17 c log weighted error 16.67 c significant figures required 14.27 c decimal places required 17.26 c data aip2cs( 1) / 0.0065457691 989713757 e0/ data aip2cs( 2) / 0.0023833724 120774592 e0/ data aip2cs( 3) / -0.0000430700 770220586 e0/ data aip2cs( 4) / 0.0000015629 125858629 e0/ data aip2cs( 5) / -0.0000000815 417186163 e0/ data aip2cs( 6) / 0.0000000054 103738057 e0/ data aip2cs( 7) / -0.0000000004 284130883 e0/ data aip2cs( 8) / 0.0000000000 389497963 e0/ data aip2cs( 9) / -0.0000000000 039623161 e0/ data aip2cs( 10) / 0.0000000000 004428184 e0/ data aip2cs( 11) / -0.0000000000 000536297 e0/ data aip2cs( 12) / 0.0000000000 000069650 e0/ data aip2cs( 13) / -0.0000000000 000009620 e0/ data aip2cs( 14) / 0.0000000000 000001403 e0/ data aip2cs( 15) / -0.0000000000 000000215 e0/ c c series for aip1 on the interval 1.25000e-01 to 1.00000e+00 c with weighted error 2.60e-17 c log weighted error 16.58 c significant figures required 14.91 c decimal places required 17.28 c data aip1cs( 1) / 0.0358865097 808301538 e0/ data aip1cs( 2) / 0.0114668575 627764899 e0/ data aip1cs( 3) / -0.0007592073 583861400 e0/ data aip1cs( 4) / 0.0000869517 610893841 e0/ data aip1cs( 5) / -0.0000128237 294298592 e0/ data aip1cs( 6) / 0.0000022062 695681038 e0/ data aip1cs( 7) / -0.0000004222 295185921 e0/ data aip1cs( 8) / 0.0000000874 686415726 e0/ data aip1cs( 9) / -0.0000000192 773588418 e0/ data aip1cs( 10) / 0.0000000044 668460054 e0/ data aip1cs( 11) / -0.0000000010 790108052 e0/ data aip1cs( 12) / 0.0000000002 700029447 e0/ data aip1cs( 13) / -0.0000000000 696480108 e0/ data aip1cs( 14) / 0.0000000000 184489907 e0/ data aip1cs( 15) / -0.0000000000 050027817 e0/ data aip1cs( 16) / 0.0000000000 013852243 e0/ data aip1cs( 17) / -0.0000000000 003908218 e0/ data aip1cs( 18) / 0.0000000000 001121536 e0/ data aip1cs( 19) / -0.0000000000 000326862 e0/ data aip1cs( 20) / 0.0000000000 000096619 e0/ data aip1cs( 21) / -0.0000000000 000028935 e0/ data aip1cs( 22) / 0.0000000000 000008770 e0/ data aip1cs( 23) / -0.0000000000 000002688 e0/ data aip1cs( 24) / 0.0000000000 000000832 e0/ data aip1cs( 25) / -0.0000000000 000000260 e0/ c data naif, naig, naip1, naip2 / 4*0 / data x2sml, x3sml, x32sml, xbig / 4*0.0 / c if (naif.ne.0) go to 10 eta = 0.1*r1mach(3) naif = inits (aifcs, 8, eta) naig = inits (aigcs, 9, eta) naip1 = inits (aip1cs, 25, eta) naip2 = inits (aip2cs, 15, 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) aide = xn * cos(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 aide = (x2*(0.125 + csevl (x3, aifcs, naif)) - 1 csevl (x3, aigcs, naig)) - 0.25 if (x.gt.x32sml) aide = aide * exp (2.0*x*sqrt(x)/3.0) return c 30 if (x.gt.4.0) go to 40 sqrtx = sqrt(x) z = (16.0/(x*sqrtx) - 9.0)/7.0 aide = (-0.28125 - csevl (z, aip1cs, naip1)) * sqrt(sqrtx) return c 40 sqrtx = sqrt(x) z = -1.0 if (x.lt.xbig) z = 16.0/(x*sqrtx) - 1.0 aide = (-0.28125 - csevl (z, aip2cs, naip2)) * sqrt(sqrtx) return c end