Blame view

fvn_fnlib/r9admp.f 11.3 KB
38581db0c   daniau   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
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
        subroutine r9admp (x, ampl, phi)
  c july 1980 edition.  w. fullerton, bell labs.  revised nov 1983.
  c
  c evaluate the derivative of the airy function modulus and
  c phase for x .le. -1.0.
  c
        dimension an22cs(33), an21cs(24), an20cs(16), aph2cs(32),
       1  aph1cs(22), aph0cs(15)
        external csevl, inits, r1mach
  c
  c series for an22 on the interval -1.00000e+00 to -1.25000e-01
  c                                        with weighted error   3.30e-17
  c                                         log weighted error  16.48
  c                               significant figures required  14.95
  c                                    decimal places required  17.24
  c
        data an22cs(  1) /     0.0537418629 629794329 e0/
        data an22cs(  2) /    -0.0126661435 859883193 e0/
        data an22cs(  3) /    -0.0011924334 106593007 e0/
        data an22cs(  4) /    -0.0002032327 627275655 e0/
        data an22cs(  5) /    -0.0000446468 963075164 e0/
        data an22cs(  6) /    -0.0000113359 036053123 e0/
        data an22cs(  7) /    -0.0000031641 352378546 e0/
        data an22cs(  8) /    -0.0000009446 708886149 e0/
        data an22cs(  9) /    -0.0000002966 562236472 e0/
        data an22cs( 10) /    -0.0000000969 118892024 e0/
        data an22cs( 11) /    -0.0000000326 822538653 e0/
        data an22cs( 12) /    -0.0000000113 144618964 e0/
        data an22cs( 13) /    -0.0000000040 042691002 e0/
        data an22cs( 14) /    -0.0000000014 440333684 e0/
        data an22cs( 15) /    -0.0000000005 292853746 e0/
        data an22cs( 16) /    -0.0000000001 967763374 e0/
        data an22cs( 17) /    -0.0000000000 740800096 e0/
        data an22cs( 18) /    -0.0000000000 282016314 e0/
        data an22cs( 19) /    -0.0000000000 108440066 e0/
        data an22cs( 20) /    -0.0000000000 042074801 e0/
        data an22cs( 21) /    -0.0000000000 016459150 e0/
        data an22cs( 22) /    -0.0000000000 006486827 e0/
        data an22cs( 23) /    -0.0000000000 002574095 e0/
        data an22cs( 24) /    -0.0000000000 001027889 e0/
        data an22cs( 25) /    -0.0000000000 000412846 e0/
        data an22cs( 26) /    -0.0000000000 000166711 e0/
        data an22cs( 27) /    -0.0000000000 000067657 e0/
        data an22cs( 28) /    -0.0000000000 000027585 e0/
        data an22cs( 29) /    -0.0000000000 000011296 e0/
        data an22cs( 30) /    -0.0000000000 000004645 e0/
        data an22cs( 31) /    -0.0000000000 000001917 e0/
        data an22cs( 32) /    -0.0000000000 000000794 e0/
        data an22cs( 33) /    -0.0000000000 000000330 e0/
  c
  c series for an21 on the interval -1.25000e-01 to -1.56250e-02
  c                                        with weighted error   3.43e-17
  c                                         log weighted error  16.47
  c                               significant figures required  14.48
  c                                    decimal places required  17.16
  c
        data an21cs(  1) /     0.0198313155 263169394 e0/
        data an21cs(  2) /    -0.0029376249 067087533 e0/
        data an21cs(  3) /    -0.0001136260 695958196 e0/
        data an21cs(  4) /    -0.0000100554 451087156 e0/
        data an21cs(  5) /    -0.0000013048 787116563 e0/
        data an21cs(  6) /    -0.0000002123 881993151 e0/
        data an21cs(  7) /    -0.0000000402 270833384 e0/
        data an21cs(  8) /    -0.0000000084 996745953 e0/
        data an21cs(  9) /    -0.0000000019 514839426 e0/
        data an21cs( 10) /    -0.0000000004 783865344 e0/
        data an21cs( 11) /    -0.0000000001 236733992 e0/
        data an21cs( 12) /    -0.0000000000 334137486 e0/
        data an21cs( 13) /    -0.0000000000 093702824 e0/
        data an21cs( 14) /    -0.0000000000 027130128 e0/
        data an21cs( 15) /    -0.0000000000 008075954 e0/
        data an21cs( 16) /    -0.0000000000 002463214 e0/
        data an21cs( 17) /    -0.0000000000 000767656 e0/
        data an21cs( 18) /    -0.0000000000 000243883 e0/
        data an21cs( 19) /    -0.0000000000 000078831 e0/
        data an21cs( 20) /    -0.0000000000 000025882 e0/
        data an21cs( 21) /    -0.0000000000 000008619 e0/
        data an21cs( 22) /    -0.0000000000 000002908 e0/
        data an21cs( 23) /    -0.0000000000 000000993 e0/
        data an21cs( 24) /    -0.0000000000 000000343 e0/
  c
  c series for an20 on the interval -1.56250e-02 to  0.00000e+00
  c                                        with weighted error   4.41e-17
  c                                         log weighted error  16.36
  c                               significant figures required  14.16
  c                                    decimal places required  16.96
  c
        data an20cs(  1) /     0.0126732217 145738027 e0/
        data an20cs(  2) /    -0.0005212847 072615621 e0/
        data an20cs(  3) /    -0.0000052672 111140370 e0/
        data an20cs(  4) /    -0.0000001628 202185026 e0/
        data an20cs(  5) /    -0.0000000090 991442687 e0/
        data an20cs(  6) /    -0.0000000007 438647126 e0/
        data an20cs(  7) /    -0.0000000000 795494752 e0/
        data an20cs(  8) /    -0.0000000000 104050944 e0/
        data an20cs(  9) /    -0.0000000000 015932426 e0/
        data an20cs( 10) /    -0.0000000000 002770648 e0/
        data an20cs( 11) /    -0.0000000000 000535343 e0/
        data an20cs( 12) /    -0.0000000000 000113062 e0/
        data an20cs( 13) /    -0.0000000000 000025772 e0/
        data an20cs( 14) /    -0.0000000000 000006278 e0/
        data an20cs( 15) /    -0.0000000000 000001621 e0/
        data an20cs( 16) /    -0.0000000000 000000441 e0/
  c
  c series for aph2 on the interval -1.00000e+00 to -1.25000e-01
  c                                        with weighted error   2.94e-17
  c                                         log weighted error  16.53
  c                               significant figures required  15.58
  c                                    decimal places required  17.28
  c
        data aph2cs(  1) /    -0.2057088719 781465107 e0/
        data aph2cs(  2) /     0.0422196961 357771922 e0/
        data aph2cs(  3) /     0.0020482560 511207275 e0/
        data aph2cs(  4) /     0.0002607800 735165006 e0/
        data aph2cs(  5) /     0.0000474824 268004729 e0/
        data aph2cs(  6) /     0.0000105102 756431612 e0/
        data aph2cs(  7) /     0.0000026353 534014668 e0/
        data aph2cs(  8) /     0.0000007208 824863499 e0/
        data aph2cs(  9) /     0.0000002103 236664473 e0/
        data aph2cs( 10) /     0.0000000644 975634555 e0/
        data aph2cs( 11) /     0.0000000205 802377264 e0/
        data aph2cs( 12) /     0.0000000067 836273921 e0/
        data aph2cs( 13) /     0.0000000022 974015284 e0/
        data aph2cs( 14) /     0.0000000007 961306765 e0/
        data aph2cs( 15) /     0.0000000002 813860610 e0/
        data aph2cs( 16) /     0.0000000001 011749057 e0/
        data aph2cs( 17) /     0.0000000000 369306738 e0/
        data aph2cs( 18) /     0.0000000000 136615066 e0/
        data aph2cs( 19) /     0.0000000000 051142751 e0/
        data aph2cs( 20) /     0.0000000000 019351689 e0/
        data aph2cs( 21) /     0.0000000000 007393607 e0/
        data aph2cs( 22) /     0.0000000000 002849792 e0/
        data aph2cs( 23) /     0.0000000000 001107281 e0/
        data aph2cs( 24) /     0.0000000000 000433412 e0/
        data aph2cs( 25) /     0.0000000000 000170801 e0/
        data aph2cs( 26) /     0.0000000000 000067733 e0/
        data aph2cs( 27) /     0.0000000000 000027017 e0/
        data aph2cs( 28) /     0.0000000000 000010835 e0/
        data aph2cs( 29) /     0.0000000000 000004367 e0/
        data aph2cs( 30) /     0.0000000000 000001769 e0/
        data aph2cs( 31) /     0.0000000000 000000719 e0/
        data aph2cs( 32) /     0.0000000000 000000294 e0/
  c
  c series for aph1 on the interval -1.25000e-01 to -1.56250e-02
  c                                        with weighted error   6.38e-17
  c                                         log weighted error  16.20
  c                               significant figures required  14.91
  c                                    decimal places required  16.87
  c
        data aph1cs(  1) /    -0.1024172908 077571694 e0/
        data aph1cs(  2) /     0.0071697275 146591248 e0/
        data aph1cs(  3) /     0.0001209959 363122329 e0/
        data aph1cs(  4) /     0.0000073361 512841220 e0/
        data aph1cs(  5) /     0.0000007535 382954272 e0/
        data aph1cs(  6) /     0.0000001041 478171741 e0/
        data aph1cs(  7) /     0.0000000174 358728519 e0/
        data aph1cs(  8) /     0.0000000033 399795033 e0/
        data aph1cs(  9) /     0.0000000007 073075174 e0/
        data aph1cs( 10) /     0.0000000001 619187515 e0/
        data aph1cs( 11) /     0.0000000000 394539982 e0/
        data aph1cs( 12) /     0.0000000000 101192282 e0/
        data aph1cs( 13) /     0.0000000000 027092778 e0/
        data aph1cs( 14) /     0.0000000000 007523806 e0/
        data aph1cs( 15) /     0.0000000000 002156369 e0/
        data aph1cs( 16) /     0.0000000000 000635283 e0/
        data aph1cs( 17) /     0.0000000000 000191757 e0/
        data aph1cs( 18) /     0.0000000000 000059143 e0/
        data aph1cs( 19) /     0.0000000000 000018597 e0/
        data aph1cs( 20) /     0.0000000000 000005950 e0/
        data aph1cs( 21) /     0.0000000000 000001934 e0/
        data aph1cs( 22) /     0.0000000000 000000638 e0/
  c
  c series for aph0 on the interval -1.56250e-02 to  0.00000e+00
  c                                        with weighted error   2.29e-17
  c                                         log weighted error  16.64
  c                               significant figures required  15.27
  c                                    decimal places required  17.23
  c
        data aph0cs(  1) /    -0.0855849241 130933257 e0/
        data aph0cs(  2) /     0.0011214378 867065261 e0/
        data aph0cs(  3) /     0.0000042721 029353664 e0/
        data aph0cs(  4) /     0.0000000817 607381483 e0/
        data aph0cs(  5) /     0.0000000033 907645000 e0/
        data aph0cs(  6) /     0.0000000002 253264423 e0/
        data aph0cs(  7) /     0.0000000000 206284209 e0/
        data aph0cs(  8) /     0.0000000000 023858763 e0/
        data aph0cs(  9) /     0.0000000000 003301618 e0/
        data aph0cs( 10) /     0.0000000000 000527010 e0/
        data aph0cs( 11) /     0.0000000000 000094555 e0/
        data aph0cs( 12) /     0.0000000000 000018709 e0/
        data aph0cs( 13) /     0.0000000000 000004024 e0/
        data aph0cs( 14) /     0.0000000000 000000930 e0/
        data aph0cs( 15) /     0.0000000000 000000229 e0/
  c
        data pi34 / 2.3561944 901923449e0 /
  c     data pi34 / 2.3561944 9019234492 8846982537 4596271631 3d0  /
        data nan20, nan21, nan22, naph0, naph1, naph2 / 6*0 /
        data xsml / 0.0 /
  c
        if (nan20.ne.0) go to 10
        eta = 0.1*r1mach(3)
        nan20 = inits (an20cs, 16, eta)
        nan21 = inits (an21cs, 24, eta)
        nan22 = inits (an22cs, 33, eta)
        naph0 = inits (aph0cs, 15, eta)
        naph1 = inits (aph1cs, 22, eta)
        naph2 = inits (aph2cs, 32, eta)
  c
        xsml = -(128.0/r1mach(3))**0.3333
  c
   10   if (x.ge.(-4.0)) go to 20
        z = 1.0
        if (x.gt.xsml) z = 128.0/x**3 + 1.0
        ampl = 0.3125 + csevl (z, an20cs, nan20)
        phi = -0.625 + csevl (z, aph0cs, naph0)
        go to 40
  c
   20   if (x.ge.(-2.0)) go to 30
        z = (128.0/x**3 + 9.0) / 7.0
        ampl = 0.3125 + csevl (z, an21cs, nan21)
        phi = -0.625 + csevl (z, aph1cs, naph1)
        go to 40
  c
   30   if (x.gt.(-1.0)) call seteru (
       1  25hr9admp  x must be le -1.0, 25, 1, 2)
  c
        z = (16.0/x**3 + 9.0) / 7.0
        ampl = 0.3125 + csevl (z, an22cs, nan22)
        phi = -0.625 + csevl (z, aph2cs, naph2)
  c
   40   sqrtx = sqrt (-x)
        ampl = sqrt (ampl*sqrtx)
        phi = pi34 - x*sqrtx*phi
  c
        return
        end