Annotation of coherent/b/lib/libm/j0.c, revision 1.1.1.1

1.1       root        1: /*
                      2:  * Evaluate the Bessel function of the first kind for the zeroth order.
                      3:  */
                      4: 
                      5: #include <math.h>
                      6: 
                      7: /*
                      8:  * (Hart 5848, 18.54)
                      9:  */
                     10: static readonly double j0sntab[] ={
                     11:         0.1208181340866561224763662419e+12,
                     12:        -0.2956513002312076810191727211e+11,
                     13:         0.1729413174598080383355729444e+10,
                     14:        -0.4281611621547871420502838045e+08,
                     15:         0.5645169313685735094277826749e+06,
                     16:        -0.4471963251278787165486324342e+04,
                     17:         0.2281027164345610253338043760e+02,
                     18:        -0.7777570245675629906097285039e-01,
                     19:         0.1792464784997734953753734861e-03,
                     20:        -0.2735011670747987792661294323e-06,
                     21:         0.2553996162031530552738418047e-09,
                     22:        -0.1135416951138795305302383379e-12
                     23: };
                     24: static readonly double j0smtab[] ={
                     25:         0.1208181340866561225104607422e+12,
                     26:         0.6394034985432622416780183619e+09,
                     27:         0.1480704129894421521840387092e+07,
                     28:         0.1806405145147135549477896097e+04,
                     29:         0.1000000000000000000000000000e+01
                     30: };
                     31: 
                     32: /*
                     33:  * (Hart 6547, 17.03)
                     34:  */
                     35: static readonly double j0pntab[] ={
                     36:         0.2180773647883051611610017444e+07,
                     37:         0.3034316360847526998223619545e+07,
                     38:         0.1103544421040585180513992902e+07,
                     39:         0.1125251525566438051490397793e+06,
                     40:         0.2239903669775096469121786611e+04
                     41: };
                     42: static readonly double j0pmtab[] ={
                     43:         0.2180773647883051631694575061e+07,
                     44:         0.3036712230333721250883918177e+07,
                     45:         0.1106820941229570783867456372e+07,
                     46:         0.1136627571261390604805029875e+06,
                     47:         0.2340314010639454134522071924e+04,
                     48:         0.1000000000000000000000000000e+01
                     49: };
                     50: 
                     51: /*
                     52:  * (Hart 6947, 17.23)
                     53:  */
                     54: static readonly double j0qntab[] ={
                     55:        -0.1766545623380246464429878870e+04,
                     56:        -0.2872031612145666435287793370e+04,
                     57:        -0.1259882860132553867070968820e+04,
                     58:        -0.1626342106227059314963763000e+03,
                     59:        -0.4423758485693335344324390000e+01
                     60: };
                     61: static readonly double j0qmtab[] ={
                     62:         0.1130589198963358159265207160e+06,
                     63:         0.1848451085035102526452988821e+06,
                     64:         0.8227466098014465706873825955e+05,
                     65:         0.1108580583675148668209607374e+05,
                     66:         0.3565514005857633096065669500e+03,
                     67:         0.1000000000000000000000000000e+01
                     68: };
                     69: 
                     70: double
                     71: j0(x)
                     72: double x;
                     73: {
                     74:        double r, p, q, z, xn;
                     75: 
                     76:        if (x < 0.0)
                     77:                x = -x;
                     78:        if (x <= 8.0) {
                     79:                r = x*x;
                     80:                r = _pol(r, j0sntab, 12)/_pol(r, j0smtab, 5);
                     81:        } else {
                     82:                /* N.B. misprint in Hart 1968 edition, corrected 1978 reprint */
                     83:                z = 8.0 / x;
                     84:                xn = x - PI/4.0;
                     85:                r = z*z;
                     86:                p = _pol(r, j0pntab, 5) / _pol(r, j0pmtab, 6);
                     87:                q = z*_pol(r, j0qntab, 5) / _pol(r, j0qmtab, 6);
                     88:                r = sqrt(2.0/(PI*x)) * (p*cos(xn) - q*sin(xn));
                     89:        }
                     90:        return (r);
                     91: }

unix.superglobalmegacorp.com

This archive runs on limited infrastructure. Preserving old code on modern bandwidth. Automated agents are requested to crawl responsibly.