source: sasmodels/sasmodels/special.py @ 859567e

core_shell_microgelsmagnetic_modelticket-1257-vesicle-productticket_1156ticket_1265_superballticket_822_more_unit_tests
Last change on this file since 859567e was 859567e, checked in by Paul Kienzle <pkienzle@…>, 7 years ago

explore quadrature options for bcc at low q

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