1 | r""" |
---|
2 | Special Functions |
---|
3 | ................. |
---|
4 | |
---|
5 | This following standard C99 math functions are available: |
---|
6 | |
---|
7 | M_PI, M_PI_2, M_PI_4, M_SQRT1_2, M_E: |
---|
8 | $\pi$, $\pi/2$, $\pi/4$, $1/\sqrt{2}$ and Euler's constant $e$ |
---|
9 | |
---|
10 | exp, log, pow(x,y), expm1, log1p, sqrt, cbrt: |
---|
11 | Power functions $e^x$, $\ln x$, $x^y$, $e^x - 1$, $\ln 1 + x$, |
---|
12 | $\sqrt{x}$, $\sqrt[3]{x}$. The functions expm1(x) and log1p(x) |
---|
13 | are accurate across all $x$, including $x$ very close to zero. |
---|
14 | |
---|
15 | sin, cos, tan, asin, acos, atan: |
---|
16 | Trigonometry functions and inverses, operating on radians. |
---|
17 | |
---|
18 | sinh, cosh, tanh, asinh, acosh, atanh: |
---|
19 | Hyperbolic trigonometry functions. |
---|
20 | |
---|
21 | atan2(y,x): |
---|
22 | Angle from the $x$\ -axis to the point $(x,y)$, which is equal to |
---|
23 | $\tan^{-1}(y/x)$ corrected for quadrant. That is, if $x$ and $y$ are |
---|
24 | both negative, then atan2(y,x) returns a value in quadrant III where |
---|
25 | atan(y/x) would return a value in quadrant I. Similarly for |
---|
26 | quadrants II and IV when $x$ and $y$ have opposite sign. |
---|
27 | |
---|
28 | fabs(x), fmin(x,y), fmax(x,y), trunc, rint: |
---|
29 | Floating point functions. rint(x) returns the nearest integer. |
---|
30 | |
---|
31 | NAN: |
---|
32 | NaN, Not a Number, $0/0$. Use isnan(x) to test for NaN. Note that |
---|
33 | you cannot use :code:`x == NAN` to test for NaN values since that |
---|
34 | will always return false. NAN does not equal NAN! The alternative, |
---|
35 | :code:`x != x` may fail if the compiler optimizes the test away. |
---|
36 | |
---|
37 | INFINITY: |
---|
38 | $\infty, 1/0$. Use isinf(x) to test for infinity, or isfinite(x) |
---|
39 | to test for finite and not NaN. |
---|
40 | |
---|
41 | erf, erfc, tgamma, lgamma: **do not use** |
---|
42 | Special functions that should be part of the standard, but are missing |
---|
43 | or inaccurate on some platforms. Use sas_erf, sas_erfc and sas_gamma |
---|
44 | instead (see below). Note: lgamma(x) has not yet been tested. |
---|
45 | |
---|
46 | Some non-standard constants and functions are also provided: |
---|
47 | |
---|
48 | M_PI_180, M_4PI_3: |
---|
49 | $\frac{\pi}{180}$, $\frac{4\pi}{3}$ |
---|
50 | |
---|
51 | SINCOS(x, s, c): |
---|
52 | Macro which sets s=sin(x) and c=cos(x). The variables *c* and *s* |
---|
53 | must be declared first. |
---|
54 | |
---|
55 | square(x): |
---|
56 | $x^2$ |
---|
57 | |
---|
58 | cube(x): |
---|
59 | $x^3$ |
---|
60 | |
---|
61 | sas_sinx_x(x): |
---|
62 | $\sin(x)/x$, with limit $\sin(0)/0 = 1$. |
---|
63 | |
---|
64 | powr(x, y): |
---|
65 | $x^y$ for $x \ge 0$; this is faster than general $x^y$ on some GPUs. |
---|
66 | |
---|
67 | pown(x, n): |
---|
68 | $x^n$ for $n$ integer; this is faster than general $x^n$ on some GPUs. |
---|
69 | |
---|
70 | FLOAT_SIZE: |
---|
71 | The number of bytes in a floating point value. Even though all |
---|
72 | variables are declared double, they may be converted to single |
---|
73 | precision float before running. If your algorithm depends on |
---|
74 | precision (which is not uncommon for numerical algorithms), use |
---|
75 | the following:: |
---|
76 | |
---|
77 | #if FLOAT_SIZE>4 |
---|
78 | ... code for double precision ... |
---|
79 | #else |
---|
80 | ... code for single precision ... |
---|
81 | #endif |
---|
82 | |
---|
83 | SAS_DOUBLE: |
---|
84 | A replacement for :code:`double` so that the declared variable will |
---|
85 | stay double precision; this should generally not be used since some |
---|
86 | graphics cards do not support double precision. There is no provision |
---|
87 | for forcing a constant to stay double precision. |
---|
88 | |
---|
89 | The following special functions and scattering calculations are defined. |
---|
90 | These functions have been tuned to be fast and numerically stable down |
---|
91 | to $q=0$ even in single precision. In some cases they work around bugs |
---|
92 | which appear on some platforms but not others, so use them where needed. |
---|
93 | Add the files listed in :code:`source = ["lib/file.c", ...]` to your *model.py* |
---|
94 | file in the order given, otherwise these functions will not be available. |
---|
95 | |
---|
96 | polevl(x, c, n): |
---|
97 | Polynomial evaluation $p(x) = \sum_{i=0}^n c_i x^i$ using Horner's |
---|
98 | method so it is faster and more accurate. |
---|
99 | |
---|
100 | $c = \{c_n, c_{n-1}, \ldots, c_0 \}$ is the table of coefficients, |
---|
101 | sorted from highest to lowest. |
---|
102 | |
---|
103 | p1evl(x, c, n): |
---|
104 | Evaluate normalized polynomial $p(x) = x^n + \sum_{i=0}^{n-1} c_i x^i$ |
---|
105 | using Horner's method so it is faster and more accurate. |
---|
106 | |
---|
107 | $c = \{c_{n-1}, c_{n-2} \ldots, c_0 \}$ is the table of coefficients, |
---|
108 | sorted from highest to lowest. |
---|
109 | |
---|
110 | sas_gamma(x): |
---|
111 | Gamma function $\text{sas_gamma}(x) = \Gamma(x)$. |
---|
112 | |
---|
113 | The standard math function, tgamma(x) is unstable for $x < 1$ |
---|
114 | on some platforms. |
---|
115 | |
---|
116 | sas_gammaln(x): |
---|
117 | log gamma function sas_gammaln\ $(x) = \log \Gamma(|x|)$. |
---|
118 | |
---|
119 | The standard math function, lgamma(x), is incorrect for single |
---|
120 | precision on some platforms. |
---|
121 | |
---|
122 | sas_gammainc(a, x), sas_gammaincc(a, x): |
---|
123 | Incomplete gamma function |
---|
124 | sas_gammainc\ $(a, x) = \int_0^x t^{a-1}e^{-t}\,dt / \Gamma(a)$ |
---|
125 | and complementary incomplete gamma function |
---|
126 | sas_gammaincc\ $(a, x) = \int_x^\infty t^{a-1}e^{-t}\,dt / \Gamma(a)$ |
---|
127 | |
---|
128 | sas_erf(x), sas_erfc(x): |
---|
129 | Error function |
---|
130 | $\text{sas_erf}(x) = \frac{2}{\sqrt\pi}\int_0^x e^{-t^2}\,dt$ |
---|
131 | and complementary error function |
---|
132 | $\text{sas_erfc}(x) = \frac{2}{\sqrt\pi}\int_x^{\infty} e^{-t^2}\,dt$. |
---|
133 | |
---|
134 | The standard math functions erf(x) and erfc(x) are slower and broken |
---|
135 | on some platforms. |
---|
136 | |
---|
137 | sas_J0(x): |
---|
138 | Bessel function of the first kind $\text{sas_J0}(x)=J_0(x)$ where |
---|
139 | $J_0(x) = \frac{1}{\pi}\int_0^\pi \cos(x\sin(\tau))\,d\tau$. |
---|
140 | |
---|
141 | The standard math function j0(x) is not available on all platforms. |
---|
142 | |
---|
143 | sas_J1(x): |
---|
144 | Bessel function of the first kind $\text{sas_J1}(x)=J_1(x)$ where |
---|
145 | $J_1(x) = \frac{1}{\pi}\int_0^\pi \cos(\tau - x\sin(\tau))\,d\tau$. |
---|
146 | |
---|
147 | The standard math function j1(x) is not available on all platforms. |
---|
148 | |
---|
149 | sas_JN(n, x): |
---|
150 | Bessel function of the first kind and integer order $n$: |
---|
151 | $\text{sas_JN}(n, x)=J_n(x)$ where |
---|
152 | $J_n(x) = \frac{1}{\pi}\int_0^\pi \cos(n\tau - x\sin(\tau))\,d\tau$. |
---|
153 | If $n$ = 0 or 1, it uses sas_J0(x) or sas_J1(x), respectively. |
---|
154 | |
---|
155 | The standard math function jn(n, x) is not available on all platforms. |
---|
156 | |
---|
157 | sas_Si(x): |
---|
158 | Sine integral $\text{Si}(x) = \int_0^x \tfrac{\sin t}{t}\,dt$. |
---|
159 | |
---|
160 | This function uses Taylor series for small and large arguments: |
---|
161 | |
---|
162 | For large arguments, |
---|
163 | |
---|
164 | .. math:: |
---|
165 | |
---|
166 | \text{Si}(x) \sim \frac{\pi}{2} |
---|
167 | - \frac{\cos(x)}{x} |
---|
168 | \left(1 - \frac{2!}{x^2} + \frac{4!}{x^4} - \frac{6!}{x^6} \right) |
---|
169 | - \frac{\sin(x)}{x} |
---|
170 | \left(\frac{1}{x} - \frac{3!}{x^3} + \frac{5!}{x^5} - \frac{7!}{x^7}\right) |
---|
171 | |
---|
172 | For small arguments, |
---|
173 | |
---|
174 | .. math:: |
---|
175 | |
---|
176 | \text{Si}(x) \sim x |
---|
177 | - \frac{x^3}{3\times 3!} + \frac{x^5}{5 \times 5!} - \frac{x^7}{7 \times 7!} |
---|
178 | + \frac{x^9}{9\times 9!} - \frac{x^{11}}{11\times 11!} |
---|
179 | |
---|
180 | sas_3j1x_x(x): |
---|
181 | Spherical Bessel form |
---|
182 | $\text{sph_j1c}(x) = 3 j_1(x)/x = 3 (\sin(x) - x \cos(x))/x^3$, |
---|
183 | with a limiting value of 1 at $x=0$, where $j_1(x)$ is the spherical |
---|
184 | Bessel function of the first kind and first order. |
---|
185 | |
---|
186 | This function uses a Taylor series for small $x$ for numerical accuracy. |
---|
187 | |
---|
188 | |
---|
189 | sas_2J1x_x(x): |
---|
190 | Bessel form $\text{sas_J1c}(x) = 2 J_1(x)/x$, with a limiting value |
---|
191 | of 1 at $x=0$, where $J_1(x)$ is the Bessel function of first kind |
---|
192 | and first order. |
---|
193 | |
---|
194 | |
---|
195 | gauss76.n, gauss76.z[i], gauss76.w[i]: |
---|
196 | Points $z_i$ and weights $w_i$ for 76-point Gaussian quadrature, respectively, |
---|
197 | computing $\int_{-1}^1 f(z)\,dz \approx \sum_{i=1}^{76} w_i\,f(z_i)$. |
---|
198 | When translating the model to C, include 'lib/gauss76.c' in the source |
---|
199 | and use :code:`GAUSS_N`, :code:`GAUSS_Z`, and :code:`GAUSS_W`. |
---|
200 | |
---|
201 | Similar arrays are available in :code:`gauss20` for 20-point quadrature |
---|
202 | and :code:`gauss150.c` for 150-point quadrature. By using |
---|
203 | :code:`import gauss76 as gauss` it is easy to change the number of |
---|
204 | points in the integration. |
---|
205 | """ |
---|
206 | # pylint: disable=unused-import |
---|
207 | |
---|
208 | import numpy as np |
---|
209 | |
---|
210 | # Functions to add to our standard set |
---|
211 | from numpy import degrees, radians |
---|
212 | |
---|
213 | # C99 standard math library functions |
---|
214 | from numpy import exp, log, power as pow, expm1, log1p, sqrt, cbrt |
---|
215 | from numpy import sin, cos, tan, arcsin as asin, arccos as acos, arctan as atan |
---|
216 | from numpy import sinh, cosh, tanh, arcsinh as asinh, arccosh as acosh, arctanh as atanh |
---|
217 | from numpy import arctan2 as atan2 |
---|
218 | from numpy import fabs, fmin, fmax, trunc, rint |
---|
219 | from numpy import pi, nan, inf |
---|
220 | from scipy.special import gamma as sas_gamma |
---|
221 | from scipy.special import gammaln as sas_gammaln |
---|
222 | from scipy.special import gammainc as sas_gammainc |
---|
223 | from scipy.special import gammaincc as sas_gammaincc |
---|
224 | from scipy.special import erf as sas_erf |
---|
225 | from scipy.special import erfc as sas_erfc |
---|
226 | from scipy.special import j0 as sas_J0 |
---|
227 | from scipy.special import j1 as sas_J1 |
---|
228 | from scipy.special import jn as sas_JN |
---|
229 | |
---|
230 | # erf, erfc, tgamma, lgamma **do not use** |
---|
231 | |
---|
232 | # C99 standard math constants |
---|
233 | M_PI, M_PI_2, M_PI_4, M_SQRT1_2, M_E = np.pi, np.pi/2, np.pi/4, np.sqrt(0.5), np.e |
---|
234 | NAN = nan |
---|
235 | INFINITY = inf |
---|
236 | |
---|
237 | # non-standard constants |
---|
238 | M_PI_180, M_4PI_3 = M_PI/180, 4*M_PI/3 |
---|
239 | |
---|
240 | # can't do SINCOS in python; use "s, c = SINCOS(x)" instead |
---|
241 | def SINCOS(x): |
---|
242 | """return sin(x), cos(x)""" |
---|
243 | return sin(x), cos(x) |
---|
244 | sincos = SINCOS |
---|
245 | |
---|
246 | def square(x): |
---|
247 | """return x^2""" |
---|
248 | return x*x |
---|
249 | |
---|
250 | def cube(x): |
---|
251 | """return x^3""" |
---|
252 | return x*x*x |
---|
253 | |
---|
254 | def sas_sinx_x(x): |
---|
255 | """return sin(x)/x""" |
---|
256 | from numpy import sinc as _sinc |
---|
257 | return _sinc(x/M_PI) |
---|
258 | |
---|
259 | def powr(x, y): |
---|
260 | """return x^y for x>0""" |
---|
261 | return x**y |
---|
262 | def pown(x, n): |
---|
263 | """return x^n for n integer""" |
---|
264 | return x**n |
---|
265 | |
---|
266 | FLOAT_SIZE = 8 |
---|
267 | |
---|
268 | def polevl(x, c, n): |
---|
269 | """return p(x) for polynomial p of degree n-1 with coefficients c""" |
---|
270 | return np.polyval(c[:n], x) |
---|
271 | |
---|
272 | def p1evl(x, c, n): |
---|
273 | """return x^n + p(x) for polynomial p of degree n-1 with coefficients c""" |
---|
274 | return np.polyval(np.hstack(([1.], c))[:n], x) |
---|
275 | |
---|
276 | def sas_Si(x): |
---|
277 | """return Si(x)""" |
---|
278 | from scipy.special import sici |
---|
279 | return sici(x)[0] |
---|
280 | |
---|
281 | def sas_j1(x): |
---|
282 | """return j1(x)""" |
---|
283 | if np.isscalar(x): |
---|
284 | retvalue = (sin(x) - x*cos(x))/x**2 if x != 0. else 0. |
---|
285 | else: |
---|
286 | with np.errstate(all='ignore'): |
---|
287 | retvalue = (sin(x) - x*cos(x))/x**2 |
---|
288 | retvalue[x == 0.] = 0. |
---|
289 | return retvalue |
---|
290 | |
---|
291 | def sas_3j1x_x(x): |
---|
292 | """return 3*j1(x)/x""" |
---|
293 | if np.isscalar(x): |
---|
294 | retvalue = 3*(sin(x) - x*cos(x))/x**3 if x != 0. else 1. |
---|
295 | else: |
---|
296 | with np.errstate(all='ignore'): |
---|
297 | retvalue = 3*(sin(x) - x*cos(x))/x**3 |
---|
298 | retvalue[x == 0.] = 1. |
---|
299 | return retvalue |
---|
300 | |
---|
301 | def sas_2J1x_x(x): |
---|
302 | """return 2*J1(x)/x""" |
---|
303 | if np.isscalar(x): |
---|
304 | retvalue = 2*sas_J1(x)/x if x != 0 else 1. |
---|
305 | else: |
---|
306 | with np.errstate(all='ignore'): |
---|
307 | retvalue = 2*sas_J1(x)/x |
---|
308 | retvalue[x == 0] = 1. |
---|
309 | return retvalue |
---|
310 | |
---|
311 | |
---|
312 | # Gaussians |
---|
313 | class Gauss: |
---|
314 | def __init__(self, w, z): |
---|
315 | self.n = len(w) |
---|
316 | self.w = w |
---|
317 | self.z = z |
---|
318 | |
---|
319 | gauss20 = Gauss( |
---|
320 | w=np.array([ |
---|
321 | .0176140071391521, |
---|
322 | .0406014298003869, |
---|
323 | .0626720483341091, |
---|
324 | .0832767415767047, |
---|
325 | .10193011981724, |
---|
326 | .118194531961518, |
---|
327 | .131688638449177, |
---|
328 | .142096109318382, |
---|
329 | .149172986472604, |
---|
330 | .152753387130726, |
---|
331 | .152753387130726, |
---|
332 | .149172986472604, |
---|
333 | .142096109318382, |
---|
334 | .131688638449177, |
---|
335 | .118194531961518, |
---|
336 | .10193011981724, |
---|
337 | .0832767415767047, |
---|
338 | .0626720483341091, |
---|
339 | .0406014298003869, |
---|
340 | .0176140071391521 |
---|
341 | ]), |
---|
342 | z=np.array([ |
---|
343 | -.993128599185095, |
---|
344 | -.963971927277914, |
---|
345 | -.912234428251326, |
---|
346 | -.839116971822219, |
---|
347 | -.746331906460151, |
---|
348 | -.636053680726515, |
---|
349 | -.510867001950827, |
---|
350 | -.37370608871542, |
---|
351 | -.227785851141645, |
---|
352 | -.076526521133497, |
---|
353 | .0765265211334973, |
---|
354 | .227785851141645, |
---|
355 | .37370608871542, |
---|
356 | .510867001950827, |
---|
357 | .636053680726515, |
---|
358 | .746331906460151, |
---|
359 | .839116971822219, |
---|
360 | .912234428251326, |
---|
361 | .963971927277914, |
---|
362 | .993128599185095 |
---|
363 | ]) |
---|
364 | ) |
---|
365 | |
---|
366 | gauss76 = Gauss( |
---|
367 | w=np.array([ |
---|
368 | .00126779163408536, #0 |
---|
369 | .00294910295364247, |
---|
370 | .00462793522803742, |
---|
371 | .00629918049732845, |
---|
372 | .00795984747723973, |
---|
373 | .00960710541471375, |
---|
374 | .0112381685696677, |
---|
375 | .0128502838475101, |
---|
376 | .0144407317482767, |
---|
377 | .0160068299122486, |
---|
378 | .0175459372914742, #10 |
---|
379 | .0190554584671906, |
---|
380 | .020532847967908, |
---|
381 | .0219756145344162, |
---|
382 | .0233813253070112, |
---|
383 | .0247476099206597, |
---|
384 | .026072164497986, |
---|
385 | .0273527555318275, |
---|
386 | .028587223650054, |
---|
387 | .029773487255905, |
---|
388 | .0309095460374916, #20 |
---|
389 | .0319934843404216, |
---|
390 | .0330234743977917, |
---|
391 | .0339977794120564, |
---|
392 | .0349147564835508, |
---|
393 | .0357728593807139, |
---|
394 | .0365706411473296, |
---|
395 | .0373067565423816, |
---|
396 | .0379799643084053, |
---|
397 | .0385891292645067, |
---|
398 | .0391332242205184, #30 |
---|
399 | .0396113317090621, |
---|
400 | .0400226455325968, |
---|
401 | .040366472122844, |
---|
402 | .0406422317102947, |
---|
403 | .0408494593018285, |
---|
404 | .040987805464794, |
---|
405 | .0410570369162294, |
---|
406 | .0410570369162294, |
---|
407 | .040987805464794, |
---|
408 | .0408494593018285, #40 |
---|
409 | .0406422317102947, |
---|
410 | .040366472122844, |
---|
411 | .0400226455325968, |
---|
412 | .0396113317090621, |
---|
413 | .0391332242205184, |
---|
414 | .0385891292645067, |
---|
415 | .0379799643084053, |
---|
416 | .0373067565423816, |
---|
417 | .0365706411473296, |
---|
418 | .0357728593807139, #50 |
---|
419 | .0349147564835508, |
---|
420 | .0339977794120564, |
---|
421 | .0330234743977917, |
---|
422 | .0319934843404216, |
---|
423 | .0309095460374916, |
---|
424 | .029773487255905, |
---|
425 | .028587223650054, |
---|
426 | .0273527555318275, |
---|
427 | .026072164497986, |
---|
428 | .0247476099206597, #60 |
---|
429 | .0233813253070112, |
---|
430 | .0219756145344162, |
---|
431 | .020532847967908, |
---|
432 | .0190554584671906, |
---|
433 | .0175459372914742, |
---|
434 | .0160068299122486, |
---|
435 | .0144407317482767, |
---|
436 | .0128502838475101, |
---|
437 | .0112381685696677, |
---|
438 | .00960710541471375, #70 |
---|
439 | .00795984747723973, |
---|
440 | .00629918049732845, |
---|
441 | .00462793522803742, |
---|
442 | .00294910295364247, |
---|
443 | .00126779163408536 #75 (indexed from 0) |
---|
444 | ]), |
---|
445 | z=np.array([ |
---|
446 | -.999505948362153, #0 |
---|
447 | -.997397786355355, |
---|
448 | -.993608772723527, |
---|
449 | -.988144453359837, |
---|
450 | -.981013938975656, |
---|
451 | -.972229228520377, |
---|
452 | -.961805126758768, |
---|
453 | -.949759207710896, |
---|
454 | -.936111781934811, |
---|
455 | -.92088586125215, |
---|
456 | -.904107119545567, #10 |
---|
457 | -.885803849292083, |
---|
458 | -.866006913771982, |
---|
459 | -.844749694983342, |
---|
460 | -.822068037328975, |
---|
461 | -.7980001871612, |
---|
462 | -.77258672828181, |
---|
463 | -.74587051350361, |
---|
464 | -.717896592387704, |
---|
465 | -.688712135277641, |
---|
466 | -.658366353758143, #20 |
---|
467 | -.626910417672267, |
---|
468 | -.594397368836793, |
---|
469 | -.560882031601237, |
---|
470 | -.526420920401243, |
---|
471 | -.491072144462194, |
---|
472 | -.454895309813726, |
---|
473 | -.417951418780327, |
---|
474 | -.380302767117504, |
---|
475 | -.342012838966962, |
---|
476 | -.303146199807908, #30 |
---|
477 | -.263768387584994, |
---|
478 | -.223945802196474, |
---|
479 | -.183745593528914, |
---|
480 | -.143235548227268, |
---|
481 | -.102483975391227, |
---|
482 | -.0615595913906112, |
---|
483 | -.0205314039939986, |
---|
484 | .0205314039939986, |
---|
485 | .0615595913906112, |
---|
486 | .102483975391227, #40 |
---|
487 | .143235548227268, |
---|
488 | .183745593528914, |
---|
489 | .223945802196474, |
---|
490 | .263768387584994, |
---|
491 | .303146199807908, |
---|
492 | .342012838966962, |
---|
493 | .380302767117504, |
---|
494 | .417951418780327, |
---|
495 | .454895309813726, |
---|
496 | .491072144462194, #50 |
---|
497 | .526420920401243, |
---|
498 | .560882031601237, |
---|
499 | .594397368836793, |
---|
500 | .626910417672267, |
---|
501 | .658366353758143, |
---|
502 | .688712135277641, |
---|
503 | .717896592387704, |
---|
504 | .74587051350361, |
---|
505 | .77258672828181, |
---|
506 | .7980001871612, #60 |
---|
507 | .822068037328975, |
---|
508 | .844749694983342, |
---|
509 | .866006913771982, |
---|
510 | .885803849292083, |
---|
511 | .904107119545567, |
---|
512 | .92088586125215, |
---|
513 | .936111781934811, |
---|
514 | .949759207710896, |
---|
515 | .961805126758768, |
---|
516 | .972229228520377, #70 |
---|
517 | .981013938975656, |
---|
518 | .988144453359837, |
---|
519 | .993608772723527, |
---|
520 | .997397786355355, |
---|
521 | .999505948362153 #75 |
---|
522 | ]) |
---|
523 | ) |
---|
524 | |
---|
525 | gauss150 = Gauss( |
---|
526 | z=np.array([ |
---|
527 | -0.9998723404457334, |
---|
528 | -0.9993274305065947, |
---|
529 | -0.9983473449340834, |
---|
530 | -0.9969322929775997, |
---|
531 | -0.9950828645255290, |
---|
532 | -0.9927998590434373, |
---|
533 | -0.9900842691660192, |
---|
534 | -0.9869372772712794, |
---|
535 | -0.9833602541697529, |
---|
536 | -0.9793547582425894, |
---|
537 | -0.9749225346595943, |
---|
538 | -0.9700655145738374, |
---|
539 | -0.9647858142586956, |
---|
540 | -0.9590857341746905, |
---|
541 | -0.9529677579610971, |
---|
542 | -0.9464345513503147, |
---|
543 | -0.9394889610042837, |
---|
544 | -0.9321340132728527, |
---|
545 | -0.9243729128743136, |
---|
546 | -0.9162090414984952, |
---|
547 | -0.9076459563329236, |
---|
548 | -0.8986873885126239, |
---|
549 | -0.8893372414942055, |
---|
550 | -0.8795995893549102, |
---|
551 | -0.8694786750173527, |
---|
552 | -0.8589789084007133, |
---|
553 | -0.8481048644991847, |
---|
554 | -0.8368612813885015, |
---|
555 | -0.8252530581614230, |
---|
556 | -0.8132852527930605, |
---|
557 | -0.8009630799369827, |
---|
558 | -0.7882919086530552, |
---|
559 | -0.7752772600680049, |
---|
560 | -0.7619248049697269, |
---|
561 | -0.7482403613363824, |
---|
562 | -0.7342298918013638, |
---|
563 | -0.7198995010552305, |
---|
564 | -0.7052554331857488, |
---|
565 | -0.6903040689571928, |
---|
566 | -0.6750519230300931, |
---|
567 | -0.6595056411226444, |
---|
568 | -0.6436719971150083, |
---|
569 | -0.6275578900977726, |
---|
570 | -0.6111703413658551, |
---|
571 | -0.5945164913591590, |
---|
572 | -0.5776035965513142, |
---|
573 | -0.5604390262878617, |
---|
574 | -0.5430302595752546, |
---|
575 | -0.5253848818220803, |
---|
576 | -0.5075105815339176, |
---|
577 | -0.4894151469632753, |
---|
578 | -0.4711064627160663, |
---|
579 | -0.4525925063160997, |
---|
580 | -0.4338813447290861, |
---|
581 | -0.4149811308476706, |
---|
582 | -0.3959000999390257, |
---|
583 | -0.3766465660565522, |
---|
584 | -0.3572289184172501, |
---|
585 | -0.3376556177463400, |
---|
586 | -0.3179351925907259, |
---|
587 | -0.2980762356029071, |
---|
588 | -0.2780873997969574, |
---|
589 | -0.2579773947782034, |
---|
590 | -0.2377549829482451, |
---|
591 | -0.2174289756869712, |
---|
592 | -0.1970082295132342, |
---|
593 | -0.1765016422258567, |
---|
594 | -0.1559181490266516, |
---|
595 | -0.1352667186271445, |
---|
596 | -0.1145563493406956, |
---|
597 | -0.0937960651617229, |
---|
598 | -0.0729949118337358, |
---|
599 | -0.0521619529078925, |
---|
600 | -0.0313062657937972, |
---|
601 | -0.0104369378042598, |
---|
602 | 0.0104369378042598, |
---|
603 | 0.0313062657937972, |
---|
604 | 0.0521619529078925, |
---|
605 | 0.0729949118337358, |
---|
606 | 0.0937960651617229, |
---|
607 | 0.1145563493406956, |
---|
608 | 0.1352667186271445, |
---|
609 | 0.1559181490266516, |
---|
610 | 0.1765016422258567, |
---|
611 | 0.1970082295132342, |
---|
612 | 0.2174289756869712, |
---|
613 | 0.2377549829482451, |
---|
614 | 0.2579773947782034, |
---|
615 | 0.2780873997969574, |
---|
616 | 0.2980762356029071, |
---|
617 | 0.3179351925907259, |
---|
618 | 0.3376556177463400, |
---|
619 | 0.3572289184172501, |
---|
620 | 0.3766465660565522, |
---|
621 | 0.3959000999390257, |
---|
622 | 0.4149811308476706, |
---|
623 | 0.4338813447290861, |
---|
624 | 0.4525925063160997, |
---|
625 | 0.4711064627160663, |
---|
626 | 0.4894151469632753, |
---|
627 | 0.5075105815339176, |
---|
628 | 0.5253848818220803, |
---|
629 | 0.5430302595752546, |
---|
630 | 0.5604390262878617, |
---|
631 | 0.5776035965513142, |
---|
632 | 0.5945164913591590, |
---|
633 | 0.6111703413658551, |
---|
634 | 0.6275578900977726, |
---|
635 | 0.6436719971150083, |
---|
636 | 0.6595056411226444, |
---|
637 | 0.6750519230300931, |
---|
638 | 0.6903040689571928, |
---|
639 | 0.7052554331857488, |
---|
640 | 0.7198995010552305, |
---|
641 | 0.7342298918013638, |
---|
642 | 0.7482403613363824, |
---|
643 | 0.7619248049697269, |
---|
644 | 0.7752772600680049, |
---|
645 | 0.7882919086530552, |
---|
646 | 0.8009630799369827, |
---|
647 | 0.8132852527930605, |
---|
648 | 0.8252530581614230, |
---|
649 | 0.8368612813885015, |
---|
650 | 0.8481048644991847, |
---|
651 | 0.8589789084007133, |
---|
652 | 0.8694786750173527, |
---|
653 | 0.8795995893549102, |
---|
654 | 0.8893372414942055, |
---|
655 | 0.8986873885126239, |
---|
656 | 0.9076459563329236, |
---|
657 | 0.9162090414984952, |
---|
658 | 0.9243729128743136, |
---|
659 | 0.9321340132728527, |
---|
660 | 0.9394889610042837, |
---|
661 | 0.9464345513503147, |
---|
662 | 0.9529677579610971, |
---|
663 | 0.9590857341746905, |
---|
664 | 0.9647858142586956, |
---|
665 | 0.9700655145738374, |
---|
666 | 0.9749225346595943, |
---|
667 | 0.9793547582425894, |
---|
668 | 0.9833602541697529, |
---|
669 | 0.9869372772712794, |
---|
670 | 0.9900842691660192, |
---|
671 | 0.9927998590434373, |
---|
672 | 0.9950828645255290, |
---|
673 | 0.9969322929775997, |
---|
674 | 0.9983473449340834, |
---|
675 | 0.9993274305065947, |
---|
676 | 0.9998723404457334 |
---|
677 | ]), |
---|
678 | w=np.array([ |
---|
679 | 0.0003276086705538, |
---|
680 | 0.0007624720924706, |
---|
681 | 0.0011976474864367, |
---|
682 | 0.0016323569986067, |
---|
683 | 0.0020663664924131, |
---|
684 | 0.0024994789888943, |
---|
685 | 0.0029315036836558, |
---|
686 | 0.0033622516236779, |
---|
687 | 0.0037915348363451, |
---|
688 | 0.0042191661429919, |
---|
689 | 0.0046449591497966, |
---|
690 | 0.0050687282939456, |
---|
691 | 0.0054902889094487, |
---|
692 | 0.0059094573005900, |
---|
693 | 0.0063260508184704, |
---|
694 | 0.0067398879387430, |
---|
695 | 0.0071507883396855, |
---|
696 | 0.0075585729801782, |
---|
697 | 0.0079630641773633, |
---|
698 | 0.0083640856838475, |
---|
699 | 0.0087614627643580, |
---|
700 | 0.0091550222717888, |
---|
701 | 0.0095445927225849, |
---|
702 | 0.0099300043714212, |
---|
703 | 0.0103110892851360, |
---|
704 | 0.0106876814158841, |
---|
705 | 0.0110596166734735, |
---|
706 | 0.0114267329968529, |
---|
707 | 0.0117888704247183, |
---|
708 | 0.0121458711652067, |
---|
709 | 0.0124975796646449, |
---|
710 | 0.0128438426753249, |
---|
711 | 0.0131845093222756, |
---|
712 | 0.0135194311690004, |
---|
713 | 0.0138484622795371, |
---|
714 | 0.0141714592928592, |
---|
715 | 0.0144882814685445, |
---|
716 | 0.0147987907597169, |
---|
717 | 0.0151028518701744, |
---|
718 | 0.0154003323133401, |
---|
719 | 0.0156911024699895, |
---|
720 | 0.0159750356447283, |
---|
721 | 0.0162520081211971, |
---|
722 | 0.0165218992159766, |
---|
723 | 0.0167845913311726, |
---|
724 | 0.0170399700056559, |
---|
725 | 0.0172879239649355, |
---|
726 | 0.0175283451696437, |
---|
727 | 0.0177611288626114, |
---|
728 | 0.0179861736145128, |
---|
729 | 0.0182033813680609, |
---|
730 | 0.0184126574807331, |
---|
731 | 0.0186139107660094, |
---|
732 | 0.0188070535331042, |
---|
733 | 0.0189920016251754, |
---|
734 | 0.0191686744559934, |
---|
735 | 0.0193369950450545, |
---|
736 | 0.0194968900511231, |
---|
737 | 0.0196482898041878, |
---|
738 | 0.0197911283358190, |
---|
739 | 0.0199253434079123, |
---|
740 | 0.0200508765398072, |
---|
741 | 0.0201676730337687, |
---|
742 | 0.0202756819988200, |
---|
743 | 0.0203748563729175, |
---|
744 | 0.0204651529434560, |
---|
745 | 0.0205465323660984, |
---|
746 | 0.0206189591819181, |
---|
747 | 0.0206824018328499, |
---|
748 | 0.0207368326754401, |
---|
749 | 0.0207822279928917, |
---|
750 | 0.0208185680053983, |
---|
751 | 0.0208458368787627, |
---|
752 | 0.0208640227312962, |
---|
753 | 0.0208731176389954, |
---|
754 | 0.0208731176389954, |
---|
755 | 0.0208640227312962, |
---|
756 | 0.0208458368787627, |
---|
757 | 0.0208185680053983, |
---|
758 | 0.0207822279928917, |
---|
759 | 0.0207368326754401, |
---|
760 | 0.0206824018328499, |
---|
761 | 0.0206189591819181, |
---|
762 | 0.0205465323660984, |
---|
763 | 0.0204651529434560, |
---|
764 | 0.0203748563729175, |
---|
765 | 0.0202756819988200, |
---|
766 | 0.0201676730337687, |
---|
767 | 0.0200508765398072, |
---|
768 | 0.0199253434079123, |
---|
769 | 0.0197911283358190, |
---|
770 | 0.0196482898041878, |
---|
771 | 0.0194968900511231, |
---|
772 | 0.0193369950450545, |
---|
773 | 0.0191686744559934, |
---|
774 | 0.0189920016251754, |
---|
775 | 0.0188070535331042, |
---|
776 | 0.0186139107660094, |
---|
777 | 0.0184126574807331, |
---|
778 | 0.0182033813680609, |
---|
779 | 0.0179861736145128, |
---|
780 | 0.0177611288626114, |
---|
781 | 0.0175283451696437, |
---|
782 | 0.0172879239649355, |
---|
783 | 0.0170399700056559, |
---|
784 | 0.0167845913311726, |
---|
785 | 0.0165218992159766, |
---|
786 | 0.0162520081211971, |
---|
787 | 0.0159750356447283, |
---|
788 | 0.0156911024699895, |
---|
789 | 0.0154003323133401, |
---|
790 | 0.0151028518701744, |
---|
791 | 0.0147987907597169, |
---|
792 | 0.0144882814685445, |
---|
793 | 0.0141714592928592, |
---|
794 | 0.0138484622795371, |
---|
795 | 0.0135194311690004, |
---|
796 | 0.0131845093222756, |
---|
797 | 0.0128438426753249, |
---|
798 | 0.0124975796646449, |
---|
799 | 0.0121458711652067, |
---|
800 | 0.0117888704247183, |
---|
801 | 0.0114267329968529, |
---|
802 | 0.0110596166734735, |
---|
803 | 0.0106876814158841, |
---|
804 | 0.0103110892851360, |
---|
805 | 0.0099300043714212, |
---|
806 | 0.0095445927225849, |
---|
807 | 0.0091550222717888, |
---|
808 | 0.0087614627643580, |
---|
809 | 0.0083640856838475, |
---|
810 | 0.0079630641773633, |
---|
811 | 0.0075585729801782, |
---|
812 | 0.0071507883396855, |
---|
813 | 0.0067398879387430, |
---|
814 | 0.0063260508184704, |
---|
815 | 0.0059094573005900, |
---|
816 | 0.0054902889094487, |
---|
817 | 0.0050687282939456, |
---|
818 | 0.0046449591497966, |
---|
819 | 0.0042191661429919, |
---|
820 | 0.0037915348363451, |
---|
821 | 0.0033622516236779, |
---|
822 | 0.0029315036836558, |
---|
823 | 0.0024994789888943, |
---|
824 | 0.0020663664924131, |
---|
825 | 0.0016323569986067, |
---|
826 | 0.0011976474864367, |
---|
827 | 0.0007624720924706, |
---|
828 | 0.0003276086705538 |
---|
829 | ]) |
---|
830 | ) |
---|