Blame view
fvn_fnlib/dai.f
3.66 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 |
double precision function dai (x) c july 1977 edition. w. fullerton, c3, los alamos scientific lab. double precision x, aifcs(13), aigcs(13), theta, xm, xmax, x3sml, 1 z, d1mach, dcsevl, daie, dcos, dexp, dlog, dsqrt external d1mach, daie, dcsevl, initds c c series for aif on the interval -1.00000e+00 to 1.00000e+00 c with weighted error 8.37e-33 c log weighted error 32.08 c significant figures required 30.87 c decimal places required 32.63 c data aif cs( 1) / -.3797135849 6669997496 1970894694 14 d-1 / data aif cs( 2) / +.5919188853 7263638574 3197280137 77 d-1 / data aif cs( 3) / +.9862928057 7279975365 6038910440 60 d-3 / data aif cs( 4) / +.6848843819 0765667554 8548301824 12 d-5 / data aif cs( 5) / +.2594202596 2194713019 4892790814 03 d-7 / data aif cs( 6) / +.6176612774 0813750329 4457496972 36 d-10 / data aif cs( 7) / +.1009245417 2466117901 4295562246 01 d-12 / data aif cs( 8) / +.1201479251 1179938141 2880332253 33 d-15 / data aif cs( 9) / +.1088294558 8716991878 5252954666 66 d-18 / data aif cs( 10) / +.7751377219 6684887039 2384000000 00 d-22 / data aif cs( 11) / +.4454811203 7175638391 4666666666 66 d-25 / data aif cs( 12) / +.2109284523 1692343466 6666666666 66 d-28 / data aif cs( 13) / +.8370173591 0741333333 3333333333 33 d-32 / c c series for aig on the interval -1.00000e+00 to 1.00000e+00 c with weighted error 7.47e-34 c log weighted error 33.13 c significant figures required 31.50 c decimal places required 33.68 c data aig cs( 1) / +.1815236558 1161273011 5562099578 64 d-1 / data aig cs( 2) / +.2157256316 6010755534 0306388199 68 d-1 / data aig cs( 3) / +.2567835698 7483249659 0524280901 33 d-3 / data aig cs( 4) / +.1426521411 9792403898 8294969217 21 d-5 / data aig cs( 5) / +.4572114920 0180426070 4340975581 91 d-8 / data aig cs( 6) / +.9525170843 5647098607 3922788405 92 d-11 / data aig cs( 7) / +.1392563460 5771399051 1504206861 90 d-13 / data aig cs( 8) / +.1507099914 2762379592 3069911386 66 d-16 / data aig cs( 9) / +.1255914831 2567778822 7032053333 33 d-19 / data aig cs( 10) / +.8306307377 0821340343 8293333333 33 d-23 / data aig cs( 11) / +.4465753849 3718567445 3333333333 33 d-26 / data aig cs( 12) / +.1990085503 4518869333 3333333333 33 d-29 / data aig cs( 13) / +.7470288525 6533333333 3333333333 33 d-33 / c data naif, naig, x3sml, xmax / 2*0, 2*0.d0 / c if (naif.ne.0) go to 10 naif = initds (aifcs, 13, 0.1*sngl(d1mach(3))) naig = initds (aigcs, 13, 0.1*sngl(d1mach(3))) c x3sml = d1mach(3)**0.3334d0 xmax = (-1.5d0*dlog(d1mach(1)))**0.6667d0 xmax = xmax - xmax*dlog(xmax)/(4.0d0*xmax*dsqrt(xmax)+1.0d0) 1 - 0.01d0 c 10 if (x.ge.(-1.d0)) go to 20 call d9aimp (x, xm, theta) dai = xm * dcos(theta) return c 20 if (x.gt.1.0d0) go to 30 z = 0.0d0 if (dabs(x).gt.x3sml) z = x**3 dai = 0.375d0 + (dcsevl (z, aifcs, naif) - x*(0.25d0 + 1 dcsevl (z, aigcs, naig)) ) return c 30 if (x.gt.xmax) go to 40 dai = daie(x) * dexp(-2.0d0*x*dsqrt(x)/3.0d0) return c 40 dai = 0.0d0 call seteru (30hdai x so big ai underflows, 30, 1, 0) return c end |