dai.f
3.66 KB
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