|
Packit |
67cb25 |
/* specfunc/bessel_zero.c
|
|
Packit |
67cb25 |
*
|
|
Packit |
67cb25 |
* Copyright (C) 1996, 1997, 1998, 1999, 2000 Gerard Jungman
|
|
Packit |
67cb25 |
*
|
|
Packit |
67cb25 |
* This program is free software; you can redistribute it and/or modify
|
|
Packit |
67cb25 |
* it under the terms of the GNU General Public License as published by
|
|
Packit |
67cb25 |
* the Free Software Foundation; either version 3 of the License, or (at
|
|
Packit |
67cb25 |
* your option) any later version.
|
|
Packit |
67cb25 |
*
|
|
Packit |
67cb25 |
* This program is distributed in the hope that it will be useful, but
|
|
Packit |
67cb25 |
* WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
Packit |
67cb25 |
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
|
|
Packit |
67cb25 |
* General Public License for more details.
|
|
Packit |
67cb25 |
*
|
|
Packit |
67cb25 |
* You should have received a copy of the GNU General Public License
|
|
Packit |
67cb25 |
* along with this program; if not, write to the Free Software
|
|
Packit |
67cb25 |
* Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
|
|
Packit |
67cb25 |
*/
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Author: G. Jungman */
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
#include <config.h>
|
|
Packit |
67cb25 |
#include <gsl/gsl_math.h>
|
|
Packit |
67cb25 |
#include <gsl/gsl_errno.h>
|
|
Packit |
67cb25 |
#include <gsl/gsl_sf_airy.h>
|
|
Packit |
67cb25 |
#include <gsl/gsl_sf_pow_int.h>
|
|
Packit |
67cb25 |
#include <gsl/gsl_sf_bessel.h>
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
#include "error.h"
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
#include "bessel_olver.h"
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* For Chebyshev expansions of the roots as functions of nu,
|
|
Packit |
67cb25 |
* see [G. Nemeth, Mathematical Approximation of Special Functions].
|
|
Packit |
67cb25 |
* This gives the fits for all nu and s <= 10.
|
|
Packit |
67cb25 |
* I made the fits for other values of s myself [GJ].
|
|
Packit |
67cb25 |
*/
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,1} = c_k T_k*(nu/2), nu <= 2 */
|
|
Packit |
67cb25 |
static const double coef_jnu1_a[] = {
|
|
Packit |
67cb25 |
3.801775243633476,
|
|
Packit |
67cb25 |
1.360704737511120,
|
|
Packit |
67cb25 |
-0.030707710261106,
|
|
Packit |
67cb25 |
0.004526823746202,
|
|
Packit |
67cb25 |
-0.000808682832134,
|
|
Packit |
67cb25 |
0.000159218792489,
|
|
Packit |
67cb25 |
-0.000033225189761,
|
|
Packit |
67cb25 |
0.000007205599763,
|
|
Packit |
67cb25 |
-0.000001606110397,
|
|
Packit |
67cb25 |
0.000000365439424,
|
|
Packit |
67cb25 |
-0.000000084498039,
|
|
Packit |
67cb25 |
0.000000019793815,
|
|
Packit |
67cb25 |
-0.000000004687054,
|
|
Packit |
67cb25 |
0.000000001120052,
|
|
Packit |
67cb25 |
-0.000000000269767,
|
|
Packit |
67cb25 |
0.000000000065420,
|
|
Packit |
67cb25 |
-0.000000000015961,
|
|
Packit |
67cb25 |
0.000000000003914,
|
|
Packit |
67cb25 |
-0.000000000000965,
|
|
Packit |
67cb25 |
0.000000000000239,
|
|
Packit |
67cb25 |
-0.000000000000059,
|
|
Packit |
67cb25 |
0.000000000000015,
|
|
Packit |
67cb25 |
-0.000000000000004,
|
|
Packit |
67cb25 |
0.000000000000001
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,1} = nu c_k T_k*((2/nu)^(2/3)), nu >= 2 */
|
|
Packit |
67cb25 |
static const double coef_jnu1_b[] = {
|
|
Packit |
67cb25 |
1.735063412537096,
|
|
Packit |
67cb25 |
0.784478100951978,
|
|
Packit |
67cb25 |
0.048881473180370,
|
|
Packit |
67cb25 |
-0.000578279783021,
|
|
Packit |
67cb25 |
-0.000038984957864,
|
|
Packit |
67cb25 |
0.000005758297879,
|
|
Packit |
67cb25 |
-0.000000327583229,
|
|
Packit |
67cb25 |
-0.000000003853878,
|
|
Packit |
67cb25 |
0.000000002284653,
|
|
Packit |
67cb25 |
-0.000000000153079,
|
|
Packit |
67cb25 |
-0.000000000000895,
|
|
Packit |
67cb25 |
0.000000000000283,
|
|
Packit |
67cb25 |
0.000000000000043,
|
|
Packit |
67cb25 |
0.000000000000010,
|
|
Packit |
67cb25 |
-0.000000000000003
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,2} = c_k T_k*(nu/2), nu <= 2 */
|
|
Packit |
67cb25 |
static const double coef_jnu2_a[] = {
|
|
Packit |
67cb25 |
6.992370244046161,
|
|
Packit |
67cb25 |
1.446379282056534,
|
|
Packit |
67cb25 |
-0.023458616207293,
|
|
Packit |
67cb25 |
0.002172149448700,
|
|
Packit |
67cb25 |
-0.000246262775620,
|
|
Packit |
67cb25 |
0.000030990180959,
|
|
Packit |
67cb25 |
-0.000004154183047,
|
|
Packit |
67cb25 |
0.000000580766328,
|
|
Packit |
67cb25 |
-0.000000083648175,
|
|
Packit |
67cb25 |
0.000000012317355,
|
|
Packit |
67cb25 |
-0.000000001844887,
|
|
Packit |
67cb25 |
0.000000000280076,
|
|
Packit |
67cb25 |
-0.000000000042986,
|
|
Packit |
67cb25 |
0.000000000006658,
|
|
Packit |
67cb25 |
-0.000000000001039,
|
|
Packit |
67cb25 |
0.000000000000163,
|
|
Packit |
67cb25 |
-0.000000000000026,
|
|
Packit |
67cb25 |
0.000000000000004,
|
|
Packit |
67cb25 |
-0.000000000000001
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,2} = nu c_k T_k*((2/nu)^(2/3)), nu >= 2 */
|
|
Packit |
67cb25 |
static const double coef_jnu2_b[] = {
|
|
Packit |
67cb25 |
2.465611864263400,
|
|
Packit |
67cb25 |
1.607952988471069,
|
|
Packit |
67cb25 |
0.138758034431497,
|
|
Packit |
67cb25 |
-0.003687791182054,
|
|
Packit |
67cb25 |
-0.000051276007868,
|
|
Packit |
67cb25 |
0.000045113570749,
|
|
Packit |
67cb25 |
-0.000007579172152,
|
|
Packit |
67cb25 |
0.000000736469208,
|
|
Packit |
67cb25 |
-0.000000011118527,
|
|
Packit |
67cb25 |
-0.000000011919884,
|
|
Packit |
67cb25 |
0.000000002696788,
|
|
Packit |
67cb25 |
-0.000000000314488,
|
|
Packit |
67cb25 |
0.000000000008124,
|
|
Packit |
67cb25 |
0.000000000005211,
|
|
Packit |
67cb25 |
-0.000000000001292,
|
|
Packit |
67cb25 |
0.000000000000158,
|
|
Packit |
67cb25 |
-0.000000000000004,
|
|
Packit |
67cb25 |
-0.000000000000003,
|
|
Packit |
67cb25 |
0.000000000000001
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,3} = c_k T_k*(nu/3), nu <= 3 */
|
|
Packit |
67cb25 |
static const double coef_jnu3_a[] = {
|
|
Packit |
67cb25 |
10.869647065239236,
|
|
Packit |
67cb25 |
2.177524286141710,
|
|
Packit |
67cb25 |
-0.034822817125293,
|
|
Packit |
67cb25 |
0.003167249102413,
|
|
Packit |
67cb25 |
-0.000353960349344,
|
|
Packit |
67cb25 |
0.000044039086085,
|
|
Packit |
67cb25 |
-0.000005851380981,
|
|
Packit |
67cb25 |
0.000000812575483,
|
|
Packit |
67cb25 |
-0.000000116463617,
|
|
Packit |
67cb25 |
0.000000017091246,
|
|
Packit |
67cb25 |
-0.000000002554376,
|
|
Packit |
67cb25 |
0.000000000387335,
|
|
Packit |
67cb25 |
-0.000000000059428,
|
|
Packit |
67cb25 |
0.000000000009207,
|
|
Packit |
67cb25 |
-0.000000000001438,
|
|
Packit |
67cb25 |
0.000000000000226,
|
|
Packit |
67cb25 |
-0.000000000000036,
|
|
Packit |
67cb25 |
0.000000000000006,
|
|
Packit |
67cb25 |
-0.000000000000001
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,3} = nu c_k T_k*((3/nu)^(2/3)), nu >= 3 */
|
|
Packit |
67cb25 |
static const double coef_jnu3_b[] = {
|
|
Packit |
67cb25 |
2.522816775173244,
|
|
Packit |
67cb25 |
1.673199424973720,
|
|
Packit |
67cb25 |
0.146431617506314,
|
|
Packit |
67cb25 |
-0.004049001763912,
|
|
Packit |
67cb25 |
-0.000039517767244,
|
|
Packit |
67cb25 |
0.000048781729288,
|
|
Packit |
67cb25 |
-0.000008729705695,
|
|
Packit |
67cb25 |
0.000000928737310,
|
|
Packit |
67cb25 |
-0.000000028388244,
|
|
Packit |
67cb25 |
-0.000000012927432,
|
|
Packit |
67cb25 |
0.000000003441008,
|
|
Packit |
67cb25 |
-0.000000000471695,
|
|
Packit |
67cb25 |
0.000000000025590,
|
|
Packit |
67cb25 |
0.000000000005502,
|
|
Packit |
67cb25 |
-0.000000000001881,
|
|
Packit |
67cb25 |
0.000000000000295,
|
|
Packit |
67cb25 |
-0.000000000000020,
|
|
Packit |
67cb25 |
-0.000000000000003,
|
|
Packit |
67cb25 |
0.000000000000001
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,4} = c_k T_k*(nu/4), nu <= 4 */
|
|
Packit |
67cb25 |
static const double coef_jnu4_a[] = {
|
|
Packit |
67cb25 |
14.750310252773009,
|
|
Packit |
67cb25 |
2.908010932941708,
|
|
Packit |
67cb25 |
-0.046093293420315,
|
|
Packit |
67cb25 |
0.004147172321412,
|
|
Packit |
67cb25 |
-0.000459092310473,
|
|
Packit |
67cb25 |
0.000056646951906,
|
|
Packit |
67cb25 |
-0.000007472351546,
|
|
Packit |
67cb25 |
0.000001031210065,
|
|
Packit |
67cb25 |
-0.000000147008137,
|
|
Packit |
67cb25 |
0.000000021475218,
|
|
Packit |
67cb25 |
-0.000000003197208,
|
|
Packit |
67cb25 |
0.000000000483249,
|
|
Packit |
67cb25 |
-0.000000000073946,
|
|
Packit |
67cb25 |
0.000000000011431,
|
|
Packit |
67cb25 |
-0.000000000001782,
|
|
Packit |
67cb25 |
0.000000000000280,
|
|
Packit |
67cb25 |
-0.000000000000044,
|
|
Packit |
67cb25 |
0.000000000000007,
|
|
Packit |
67cb25 |
-0.000000000000001
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,4} = nu c_k T_k*((4/nu)^(2/3)), nu >= 4 */
|
|
Packit |
67cb25 |
static const double coef_jnu4_b[] = {
|
|
Packit |
67cb25 |
2.551681323117914,
|
|
Packit |
67cb25 |
1.706177978336572,
|
|
Packit |
67cb25 |
0.150357658406131,
|
|
Packit |
67cb25 |
-0.004234001378590,
|
|
Packit |
67cb25 |
-0.000033854229898,
|
|
Packit |
67cb25 |
0.000050763551485,
|
|
Packit |
67cb25 |
-0.000009337464057,
|
|
Packit |
67cb25 |
0.000001029717834,
|
|
Packit |
67cb25 |
-0.000000037474196,
|
|
Packit |
67cb25 |
-0.000000013450153,
|
|
Packit |
67cb25 |
0.000000003836180,
|
|
Packit |
67cb25 |
-0.000000000557404,
|
|
Packit |
67cb25 |
0.000000000035748,
|
|
Packit |
67cb25 |
0.000000000005487,
|
|
Packit |
67cb25 |
-0.000000000002187,
|
|
Packit |
67cb25 |
0.000000000000374,
|
|
Packit |
67cb25 |
-0.000000000000031,
|
|
Packit |
67cb25 |
-0.000000000000003,
|
|
Packit |
67cb25 |
0.000000000000001
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,5} = c_k T_k*(nu/5), nu <= 5 */
|
|
Packit |
67cb25 |
static const double coef_jnu5_a[] = {
|
|
Packit |
67cb25 |
18.632261081028211,
|
|
Packit |
67cb25 |
3.638249012596966,
|
|
Packit |
67cb25 |
-0.057329705998828,
|
|
Packit |
67cb25 |
0.005121709126820,
|
|
Packit |
67cb25 |
-0.000563325259487,
|
|
Packit |
67cb25 |
0.000069100826174,
|
|
Packit |
67cb25 |
-0.000009066603030,
|
|
Packit |
67cb25 |
0.000001245181383,
|
|
Packit |
67cb25 |
-0.000000176737282,
|
|
Packit |
67cb25 |
0.000000025716695,
|
|
Packit |
67cb25 |
-0.000000003815184,
|
|
Packit |
67cb25 |
0.000000000574839,
|
|
Packit |
67cb25 |
-0.000000000087715,
|
|
Packit |
67cb25 |
0.000000000013526,
|
|
Packit |
67cb25 |
-0.000000000002104,
|
|
Packit |
67cb25 |
0.000000000000330,
|
|
Packit |
67cb25 |
-0.000000000000052,
|
|
Packit |
67cb25 |
0.000000000000008,
|
|
Packit |
67cb25 |
-0.000000000000001
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,5} = nu c_k T_k*((5/nu)^(2/3)), nu >= 5 */
|
|
Packit |
67cb25 |
/* FIXME: There is something wrong with this fit, in about the
|
|
Packit |
67cb25 |
* 9th or 10th decimal place.
|
|
Packit |
67cb25 |
*/
|
|
Packit |
67cb25 |
static const double coef_jnu5_b[] = {
|
|
Packit |
67cb25 |
2.569079487591442,
|
|
Packit |
67cb25 |
1.726073360882134,
|
|
Packit |
67cb25 |
0.152740776809531,
|
|
Packit |
67cb25 |
-0.004346449660148,
|
|
Packit |
67cb25 |
-0.000030512461856,
|
|
Packit |
67cb25 |
0.000052000821080,
|
|
Packit |
67cb25 |
-0.000009713343981,
|
|
Packit |
67cb25 |
0.000001091997863,
|
|
Packit |
67cb25 |
-0.000000043061707,
|
|
Packit |
67cb25 |
-0.000000013779413,
|
|
Packit |
67cb25 |
0.000000004082870,
|
|
Packit |
67cb25 |
-0.000000000611259,
|
|
Packit |
67cb25 |
0.000000000042242,
|
|
Packit |
67cb25 |
0.000000000005448,
|
|
Packit |
67cb25 |
-0.000000000002377,
|
|
Packit |
67cb25 |
0.000000000000424,
|
|
Packit |
67cb25 |
-0.000000000000038,
|
|
Packit |
67cb25 |
-0.000000000000002,
|
|
Packit |
67cb25 |
0.000000000000002
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,6} = c_k T_k*(nu/6), nu <= 6 */
|
|
Packit |
67cb25 |
static const double coef_jnu6_a[] = {
|
|
Packit |
67cb25 |
22.514836143374042,
|
|
Packit |
67cb25 |
4.368367257557198,
|
|
Packit |
67cb25 |
-0.068550155285562,
|
|
Packit |
67cb25 |
0.006093776505822,
|
|
Packit |
67cb25 |
-0.000667152784957,
|
|
Packit |
67cb25 |
0.000081486022398,
|
|
Packit |
67cb25 |
-0.000010649011647,
|
|
Packit |
67cb25 |
0.000001457089679,
|
|
Packit |
67cb25 |
-0.000000206105082,
|
|
Packit |
67cb25 |
0.000000029894724,
|
|
Packit |
67cb25 |
-0.000000004422012,
|
|
Packit |
67cb25 |
0.000000000664471,
|
|
Packit |
67cb25 |
-0.000000000101140,
|
|
Packit |
67cb25 |
0.000000000015561,
|
|
Packit |
67cb25 |
-0.000000000002416,
|
|
Packit |
67cb25 |
0.000000000000378,
|
|
Packit |
67cb25 |
-0.000000000000060,
|
|
Packit |
67cb25 |
0.000000000000009,
|
|
Packit |
67cb25 |
-0.000000000000002
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,6} = nu c_k T_k*((6/nu)^(2/3)), nu >= 6 */
|
|
Packit |
67cb25 |
static const double coef_jnu6_b[] = {
|
|
Packit |
67cb25 |
2.580710285494837,
|
|
Packit |
67cb25 |
1.739380728566154,
|
|
Packit |
67cb25 |
0.154340696401691,
|
|
Packit |
67cb25 |
-0.004422028860168,
|
|
Packit |
67cb25 |
-0.000028305272624,
|
|
Packit |
67cb25 |
0.000052845975269,
|
|
Packit |
67cb25 |
-0.000009968794373,
|
|
Packit |
67cb25 |
0.000001134252926,
|
|
Packit |
67cb25 |
-0.000000046841241,
|
|
Packit |
67cb25 |
-0.000000014007555,
|
|
Packit |
67cb25 |
0.000000004251816,
|
|
Packit |
67cb25 |
-0.000000000648213,
|
|
Packit |
67cb25 |
0.000000000046728,
|
|
Packit |
67cb25 |
0.000000000005414,
|
|
Packit |
67cb25 |
-0.000000000002508,
|
|
Packit |
67cb25 |
0.000000000000459,
|
|
Packit |
67cb25 |
-0.000000000000043,
|
|
Packit |
67cb25 |
-0.000000000000002,
|
|
Packit |
67cb25 |
0.000000000000002
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,7} = c_k T_k*(nu/7), nu <= 7 */
|
|
Packit |
67cb25 |
static const double coef_jnu7_a[] = {
|
|
Packit |
67cb25 |
26.397760539730869,
|
|
Packit |
67cb25 |
5.098418721711790,
|
|
Packit |
67cb25 |
-0.079761896398948,
|
|
Packit |
67cb25 |
0.007064521280487,
|
|
Packit |
67cb25 |
-0.000770766522482,
|
|
Packit |
67cb25 |
0.000093835449636,
|
|
Packit |
67cb25 |
-0.000012225308542,
|
|
Packit |
67cb25 |
0.000001667939800,
|
|
Packit |
67cb25 |
-0.000000235288157,
|
|
Packit |
67cb25 |
0.000000034040347,
|
|
Packit |
67cb25 |
-0.000000005023142,
|
|
Packit |
67cb25 |
0.000000000753101,
|
|
Packit |
67cb25 |
-0.000000000114389,
|
|
Packit |
67cb25 |
0.000000000017564,
|
|
Packit |
67cb25 |
-0.000000000002722,
|
|
Packit |
67cb25 |
0.000000000000425,
|
|
Packit |
67cb25 |
-0.000000000000067,
|
|
Packit |
67cb25 |
0.000000000000011,
|
|
Packit |
67cb25 |
-0.000000000000002
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,7} = nu c_k T_k*((7/nu)^(2/3)), nu >= 7 */
|
|
Packit |
67cb25 |
static const double coef_jnu7_b[] = {
|
|
Packit |
67cb25 |
2.589033335856773,
|
|
Packit |
67cb25 |
1.748907007612678,
|
|
Packit |
67cb25 |
0.155488900387653,
|
|
Packit |
67cb25 |
-0.004476317805688,
|
|
Packit |
67cb25 |
-0.000026737952924,
|
|
Packit |
67cb25 |
0.000053459680946,
|
|
Packit |
67cb25 |
-0.000010153699240,
|
|
Packit |
67cb25 |
0.000001164804272,
|
|
Packit |
67cb25 |
-0.000000049566917,
|
|
Packit |
67cb25 |
-0.000000014175403,
|
|
Packit |
67cb25 |
0.000000004374840,
|
|
Packit |
67cb25 |
-0.000000000675135,
|
|
Packit |
67cb25 |
0.000000000050004,
|
|
Packit |
67cb25 |
0.000000000005387,
|
|
Packit |
67cb25 |
-0.000000000002603,
|
|
Packit |
67cb25 |
0.000000000000485,
|
|
Packit |
67cb25 |
-0.000000000000047,
|
|
Packit |
67cb25 |
-0.000000000000002,
|
|
Packit |
67cb25 |
0.000000000000002
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,8} = c_k T_k*(nu/8), nu <= 8 */
|
|
Packit |
67cb25 |
static const double coef_jnu8_a[] = {
|
|
Packit |
67cb25 |
30.280900001606662,
|
|
Packit |
67cb25 |
5.828429205461221,
|
|
Packit |
67cb25 |
-0.090968381181069,
|
|
Packit |
67cb25 |
0.008034479731033,
|
|
Packit |
67cb25 |
-0.000874254899080,
|
|
Packit |
67cb25 |
0.000106164151611,
|
|
Packit |
67cb25 |
-0.000013798098749,
|
|
Packit |
67cb25 |
0.000001878187386,
|
|
Packit |
67cb25 |
-0.000000264366627,
|
|
Packit |
67cb25 |
0.000000038167685,
|
|
Packit |
67cb25 |
-0.000000005621060,
|
|
Packit |
67cb25 |
0.000000000841165,
|
|
Packit |
67cb25 |
-0.000000000127538,
|
|
Packit |
67cb25 |
0.000000000019550,
|
|
Packit |
67cb25 |
-0.000000000003025,
|
|
Packit |
67cb25 |
0.000000000000472,
|
|
Packit |
67cb25 |
-0.000000000000074,
|
|
Packit |
67cb25 |
0.000000000000012,
|
|
Packit |
67cb25 |
-0.000000000000002
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,8} = nu c_k T_k*((8/nu)^(2/3)), nu >= 8 */
|
|
Packit |
67cb25 |
static const double coef_jnu8_b[] = {
|
|
Packit |
67cb25 |
2.595283877150078,
|
|
Packit |
67cb25 |
1.756063044986928,
|
|
Packit |
67cb25 |
0.156352972371030,
|
|
Packit |
67cb25 |
-0.004517201896761,
|
|
Packit |
67cb25 |
-0.000025567187878,
|
|
Packit |
67cb25 |
0.000053925472558,
|
|
Packit |
67cb25 |
-0.000010293734486,
|
|
Packit |
67cb25 |
0.000001187923085,
|
|
Packit |
67cb25 |
-0.000000051625122,
|
|
Packit |
67cb25 |
-0.000000014304212,
|
|
Packit |
67cb25 |
0.000000004468450,
|
|
Packit |
67cb25 |
-0.000000000695620,
|
|
Packit |
67cb25 |
0.000000000052500,
|
|
Packit |
67cb25 |
0.000000000005367,
|
|
Packit |
67cb25 |
-0.000000000002676,
|
|
Packit |
67cb25 |
0.000000000000505,
|
|
Packit |
67cb25 |
-0.000000000000050,
|
|
Packit |
67cb25 |
-0.000000000000002,
|
|
Packit |
67cb25 |
0.000000000000002
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,9} = c_k T_k*(nu/9), nu <= 9 */
|
|
Packit |
67cb25 |
static const double coef_jnu9_a[] = {
|
|
Packit |
67cb25 |
34.164181213238386,
|
|
Packit |
67cb25 |
6.558412747925228,
|
|
Packit |
67cb25 |
-0.102171455365016,
|
|
Packit |
67cb25 |
0.009003934361201,
|
|
Packit |
67cb25 |
-0.000977663914535,
|
|
Packit |
67cb25 |
0.000118479876579,
|
|
Packit |
67cb25 |
-0.000015368714220,
|
|
Packit |
67cb25 |
0.000002088064285,
|
|
Packit |
67cb25 |
-0.000000293381154,
|
|
Packit |
67cb25 |
0.000000042283900,
|
|
Packit |
67cb25 |
-0.000000006217033,
|
|
Packit |
67cb25 |
0.000000000928887,
|
|
Packit |
67cb25 |
-0.000000000140627,
|
|
Packit |
67cb25 |
0.000000000021526,
|
|
Packit |
67cb25 |
-0.000000000003326,
|
|
Packit |
67cb25 |
0.000000000000518,
|
|
Packit |
67cb25 |
-0.000000000000081,
|
|
Packit |
67cb25 |
0.000000000000013,
|
|
Packit |
67cb25 |
-0.000000000000002
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,9} = nu c_k T_k*((9/nu)^(2/3)), nu >= 9 */
|
|
Packit |
67cb25 |
static const double coef_jnu9_b[] = {
|
|
Packit |
67cb25 |
2.600150240905079,
|
|
Packit |
67cb25 |
1.761635491694032,
|
|
Packit |
67cb25 |
0.157026743724010,
|
|
Packit |
67cb25 |
-0.004549100368716,
|
|
Packit |
67cb25 |
-0.000024659248617,
|
|
Packit |
67cb25 |
0.000054291035068,
|
|
Packit |
67cb25 |
-0.000010403464334,
|
|
Packit |
67cb25 |
0.000001206027524,
|
|
Packit |
67cb25 |
-0.000000053234089,
|
|
Packit |
67cb25 |
-0.000000014406241,
|
|
Packit |
67cb25 |
0.000000004542078,
|
|
Packit |
67cb25 |
-0.000000000711728,
|
|
Packit |
67cb25 |
0.000000000054464,
|
|
Packit |
67cb25 |
0.000000000005350,
|
|
Packit |
67cb25 |
-0.000000000002733,
|
|
Packit |
67cb25 |
0.000000000000521,
|
|
Packit |
67cb25 |
-0.000000000000052,
|
|
Packit |
67cb25 |
-0.000000000000002,
|
|
Packit |
67cb25 |
0.000000000000002
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,10} = c_k T_k*(nu/10), nu <= 10 */
|
|
Packit |
67cb25 |
static const double coef_jnu10_a[] = {
|
|
Packit |
67cb25 |
38.047560766184647,
|
|
Packit |
67cb25 |
7.288377637926008,
|
|
Packit |
67cb25 |
-0.113372193277897,
|
|
Packit |
67cb25 |
0.009973047509098,
|
|
Packit |
67cb25 |
-0.001081019701335,
|
|
Packit |
67cb25 |
0.000130786983847,
|
|
Packit |
67cb25 |
-0.000016937898538,
|
|
Packit |
67cb25 |
0.000002297699179,
|
|
Packit |
67cb25 |
-0.000000322354218,
|
|
Packit |
67cb25 |
0.000000046392941,
|
|
Packit |
67cb25 |
-0.000000006811759,
|
|
Packit |
67cb25 |
0.000000001016395,
|
|
Packit |
67cb25 |
-0.000000000153677,
|
|
Packit |
67cb25 |
0.000000000023486,
|
|
Packit |
67cb25 |
-0.000000000003616,
|
|
Packit |
67cb25 |
0.000000000000561,
|
|
Packit |
67cb25 |
-0.000000000000095,
|
|
Packit |
67cb25 |
0.000000000000027,
|
|
Packit |
67cb25 |
-0.000000000000013,
|
|
Packit |
67cb25 |
0.000000000000005
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,10} = nu c_k T_k*((10/nu)^(2/3)), nu >= 10 */
|
|
Packit |
67cb25 |
static const double coef_jnu10_b[] = {
|
|
Packit |
67cb25 |
2.604046346867949,
|
|
Packit |
67cb25 |
1.766097596481182,
|
|
Packit |
67cb25 |
0.157566834446511,
|
|
Packit |
67cb25 |
-0.004574682244089,
|
|
Packit |
67cb25 |
-0.000023934500688,
|
|
Packit |
67cb25 |
0.000054585558231,
|
|
Packit |
67cb25 |
-0.000010491765415,
|
|
Packit |
67cb25 |
0.000001220589364,
|
|
Packit |
67cb25 |
-0.000000054526331,
|
|
Packit |
67cb25 |
-0.000000014489078,
|
|
Packit |
67cb25 |
0.000000004601510,
|
|
Packit |
67cb25 |
-0.000000000724727,
|
|
Packit |
67cb25 |
0.000000000056049,
|
|
Packit |
67cb25 |
0.000000000005337,
|
|
Packit |
67cb25 |
-0.000000000002779,
|
|
Packit |
67cb25 |
0.000000000000533,
|
|
Packit |
67cb25 |
-0.000000000000054,
|
|
Packit |
67cb25 |
-0.000000000000002,
|
|
Packit |
67cb25 |
0.000000000000002
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,11} = c_k T_k*(nu/22), nu <= 22 */
|
|
Packit |
67cb25 |
static const double coef_jnu11_a[] = {
|
|
Packit |
67cb25 |
49.5054081076848637,
|
|
Packit |
67cb25 |
15.33692279367165101,
|
|
Packit |
67cb25 |
-0.33677234163517130,
|
|
Packit |
67cb25 |
0.04623235772920729,
|
|
Packit |
67cb25 |
-0.00781084960665093,
|
|
Packit |
67cb25 |
0.00147217395434708,
|
|
Packit |
67cb25 |
-0.00029695043846867,
|
|
Packit |
67cb25 |
0.00006273356860235,
|
|
Packit |
67cb25 |
-0.00001370575125628,
|
|
Packit |
67cb25 |
3.07171282012e-6,
|
|
Packit |
67cb25 |
-7.0235041249e-7,
|
|
Packit |
67cb25 |
1.6320559339e-7,
|
|
Packit |
67cb25 |
-3.843117306e-8,
|
|
Packit |
67cb25 |
9.15083800e-9,
|
|
Packit |
67cb25 |
-2.19957642e-9,
|
|
Packit |
67cb25 |
5.3301703e-10,
|
|
Packit |
67cb25 |
-1.3007541e-10,
|
|
Packit |
67cb25 |
3.193827e-11,
|
|
Packit |
67cb25 |
-7.88605e-12,
|
|
Packit |
67cb25 |
1.95918e-12,
|
|
Packit |
67cb25 |
-4.9020e-13,
|
|
Packit |
67cb25 |
1.2207e-13,
|
|
Packit |
67cb25 |
-2.820e-14,
|
|
Packit |
67cb25 |
5.25e-15,
|
|
Packit |
67cb25 |
-1.88e-15,
|
|
Packit |
67cb25 |
2.80e-15,
|
|
Packit |
67cb25 |
-2.45e-15
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,12} = c_k T_k*(nu/24), nu <= 24 */
|
|
Packit |
67cb25 |
static const double coef_jnu12_a[] = {
|
|
Packit |
67cb25 |
54.0787833216641519,
|
|
Packit |
67cb25 |
16.7336367772863598,
|
|
Packit |
67cb25 |
-0.36718411124537953,
|
|
Packit |
67cb25 |
0.05035523375053820,
|
|
Packit |
67cb25 |
-0.00849884978867533,
|
|
Packit |
67cb25 |
0.00160027692813434,
|
|
Packit |
67cb25 |
-0.00032248114889921,
|
|
Packit |
67cb25 |
0.00006806354127199,
|
|
Packit |
67cb25 |
-0.00001485665901339,
|
|
Packit |
67cb25 |
3.32668783672e-6,
|
|
Packit |
67cb25 |
-7.5998952729e-7,
|
|
Packit |
67cb25 |
1.7644939709e-7,
|
|
Packit |
67cb25 |
-4.151538210e-8,
|
|
Packit |
67cb25 |
9.87722772e-9,
|
|
Packit |
67cb25 |
-2.37230133e-9,
|
|
Packit |
67cb25 |
5.7442875e-10,
|
|
Packit |
67cb25 |
-1.4007767e-10,
|
|
Packit |
67cb25 |
3.437166e-11,
|
|
Packit |
67cb25 |
-8.48215e-12,
|
|
Packit |
67cb25 |
2.10554e-12,
|
|
Packit |
67cb25 |
-5.2623e-13,
|
|
Packit |
67cb25 |
1.3189e-13,
|
|
Packit |
67cb25 |
-3.175e-14,
|
|
Packit |
67cb25 |
5.73e-15,
|
|
Packit |
67cb25 |
5.6e-16,
|
|
Packit |
67cb25 |
-8.7e-16,
|
|
Packit |
67cb25 |
-6.5e-16
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,13} = c_k T_k*(nu/26), nu <= 26 */
|
|
Packit |
67cb25 |
static const double coef_jnu13_a[] = {
|
|
Packit |
67cb25 |
58.6521941921708890,
|
|
Packit |
67cb25 |
18.1303398137970284,
|
|
Packit |
67cb25 |
-0.39759381380126650,
|
|
Packit |
67cb25 |
0.05447765240465494,
|
|
Packit |
67cb25 |
-0.00918674227679980,
|
|
Packit |
67cb25 |
0.00172835361420579,
|
|
Packit |
67cb25 |
-0.00034800528297612,
|
|
Packit |
67cb25 |
0.00007339183835188,
|
|
Packit |
67cb25 |
-0.00001600713368099,
|
|
Packit |
67cb25 |
3.58154960392e-6,
|
|
Packit |
67cb25 |
-8.1759873497e-7,
|
|
Packit |
67cb25 |
1.8968523220e-7,
|
|
Packit |
67cb25 |
-4.459745253e-8,
|
|
Packit |
67cb25 |
1.060304419e-8,
|
|
Packit |
67cb25 |
-2.54487624e-9,
|
|
Packit |
67cb25 |
6.1580214e-10,
|
|
Packit |
67cb25 |
-1.5006751e-10,
|
|
Packit |
67cb25 |
3.679707e-11,
|
|
Packit |
67cb25 |
-9.07159e-12,
|
|
Packit |
67cb25 |
2.24713e-12,
|
|
Packit |
67cb25 |
-5.5943e-13,
|
|
Packit |
67cb25 |
1.4069e-13,
|
|
Packit |
67cb25 |
-3.679e-14,
|
|
Packit |
67cb25 |
1.119e-14,
|
|
Packit |
67cb25 |
-4.99e-15,
|
|
Packit |
67cb25 |
3.43e-15,
|
|
Packit |
67cb25 |
-2.85e-15,
|
|
Packit |
67cb25 |
2.3e-15,
|
|
Packit |
67cb25 |
-1.7e-15,
|
|
Packit |
67cb25 |
8.7e-16
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,14} = c_k T_k*(nu/28), nu <= 28 */
|
|
Packit |
67cb25 |
static const double coef_jnu14_a[] = {
|
|
Packit |
67cb25 |
63.2256329577315566,
|
|
Packit |
67cb25 |
19.5270342832914901,
|
|
Packit |
67cb25 |
-0.42800190567884337,
|
|
Packit |
67cb25 |
0.05859971627729398,
|
|
Packit |
67cb25 |
-0.00987455163523582,
|
|
Packit |
67cb25 |
0.00185641011402081,
|
|
Packit |
67cb25 |
-0.00037352439419968,
|
|
Packit |
67cb25 |
0.00007871886257265,
|
|
Packit |
67cb25 |
-0.00001715728110045,
|
|
Packit |
67cb25 |
3.83632624437e-6,
|
|
Packit |
67cb25 |
-8.7518558668e-7,
|
|
Packit |
67cb25 |
2.0291515353e-7,
|
|
Packit |
67cb25 |
-4.767795233e-8,
|
|
Packit |
67cb25 |
1.132844415e-8,
|
|
Packit |
67cb25 |
-2.71734219e-9,
|
|
Packit |
67cb25 |
6.5714886e-10,
|
|
Packit |
67cb25 |
-1.6005342e-10,
|
|
Packit |
67cb25 |
3.922557e-11,
|
|
Packit |
67cb25 |
-9.66637e-12,
|
|
Packit |
67cb25 |
2.39379e-12,
|
|
Packit |
67cb25 |
-5.9541e-13,
|
|
Packit |
67cb25 |
1.4868e-13,
|
|
Packit |
67cb25 |
-3.726e-14,
|
|
Packit |
67cb25 |
9.37e-15,
|
|
Packit |
67cb25 |
-2.36e-15,
|
|
Packit |
67cb25 |
6.0e-16
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,15} = c_k T_k*(nu/30), nu <= 30 */
|
|
Packit |
67cb25 |
static const double coef_jnu15_a[] = {
|
|
Packit |
67cb25 |
67.7990939565631635,
|
|
Packit |
67cb25 |
20.9237219226859859,
|
|
Packit |
67cb25 |
-0.45840871823085836,
|
|
Packit |
67cb25 |
0.06272149946755639,
|
|
Packit |
67cb25 |
-0.01056229551143042,
|
|
Packit |
67cb25 |
0.00198445078693100,
|
|
Packit |
67cb25 |
-0.00039903958650729,
|
|
Packit |
67cb25 |
0.00008404489865469,
|
|
Packit |
67cb25 |
-0.00001830717574922,
|
|
Packit |
67cb25 |
4.09103745566e-6,
|
|
Packit |
67cb25 |
-9.3275533309e-7,
|
|
Packit |
67cb25 |
2.1614056403e-7,
|
|
Packit |
67cb25 |
-5.075725222e-8,
|
|
Packit |
67cb25 |
1.205352081e-8,
|
|
Packit |
67cb25 |
-2.88971837e-9,
|
|
Packit |
67cb25 |
6.9846848e-10,
|
|
Packit |
67cb25 |
-1.7002946e-10,
|
|
Packit |
67cb25 |
4.164941e-11,
|
|
Packit |
67cb25 |
-1.025859e-11,
|
|
Packit |
67cb25 |
2.53921e-12,
|
|
Packit |
67cb25 |
-6.3128e-13,
|
|
Packit |
67cb25 |
1.5757e-13,
|
|
Packit |
67cb25 |
-3.947e-14,
|
|
Packit |
67cb25 |
9.92e-15,
|
|
Packit |
67cb25 |
-2.50e-15,
|
|
Packit |
67cb25 |
6.3e-16
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,16} = c_k T_k*(nu/32), nu <= 32 */
|
|
Packit |
67cb25 |
static const double coef_jnu16_a[] = {
|
|
Packit |
67cb25 |
72.3725729616724770,
|
|
Packit |
67cb25 |
22.32040402918608585,
|
|
Packit |
67cb25 |
-0.48881449782358690,
|
|
Packit |
67cb25 |
0.06684305681828766,
|
|
Packit |
67cb25 |
-0.01124998690363398,
|
|
Packit |
67cb25 |
0.00211247882775445,
|
|
Packit |
67cb25 |
-0.00042455166484632,
|
|
Packit |
67cb25 |
0.00008937015316346,
|
|
Packit |
67cb25 |
-0.00001945687139551,
|
|
Packit |
67cb25 |
4.34569739281e-6,
|
|
Packit |
67cb25 |
-9.9031173548e-7,
|
|
Packit |
67cb25 |
2.2936247195e-7,
|
|
Packit |
67cb25 |
-5.383562595e-8,
|
|
Packit |
67cb25 |
1.277835103e-8,
|
|
Packit |
67cb25 |
-3.06202860e-9,
|
|
Packit |
67cb25 |
7.3977037e-10,
|
|
Packit |
67cb25 |
-1.8000071e-10,
|
|
Packit |
67cb25 |
4.407196e-11,
|
|
Packit |
67cb25 |
-1.085046e-11,
|
|
Packit |
67cb25 |
2.68453e-12,
|
|
Packit |
67cb25 |
-6.6712e-13,
|
|
Packit |
67cb25 |
1.6644e-13,
|
|
Packit |
67cb25 |
-4.168e-14,
|
|
Packit |
67cb25 |
1.047e-14,
|
|
Packit |
67cb25 |
-2.64e-15,
|
|
Packit |
67cb25 |
6.7e-16
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,17} = c_k T_k*(nu/34), nu <= 34 */
|
|
Packit |
67cb25 |
static const double coef_jnu17_a[] = {
|
|
Packit |
67cb25 |
76.9460667535209549,
|
|
Packit |
67cb25 |
23.71708159112252670,
|
|
Packit |
67cb25 |
-0.51921943142405352,
|
|
Packit |
67cb25 |
0.07096442978067622,
|
|
Packit |
67cb25 |
-0.01193763559341369,
|
|
Packit |
67cb25 |
0.00224049662974902,
|
|
Packit |
67cb25 |
-0.00045006122941781,
|
|
Packit |
67cb25 |
0.00009469477941684,
|
|
Packit |
67cb25 |
-0.00002060640777107,
|
|
Packit |
67cb25 |
4.60031647195e-6,
|
|
Packit |
67cb25 |
-1.04785755046e-6,
|
|
Packit |
67cb25 |
2.4258161247e-7,
|
|
Packit |
67cb25 |
-5.691327087e-8,
|
|
Packit |
67cb25 |
1.350298805e-8,
|
|
Packit |
67cb25 |
-3.23428733e-9,
|
|
Packit |
67cb25 |
7.8105847e-10,
|
|
Packit |
67cb25 |
-1.8996825e-10,
|
|
Packit |
67cb25 |
4.649350e-11,
|
|
Packit |
67cb25 |
-1.144205e-11,
|
|
Packit |
67cb25 |
2.82979e-12,
|
|
Packit |
67cb25 |
-7.0294e-13,
|
|
Packit |
67cb25 |
1.7531e-13,
|
|
Packit |
67cb25 |
-4.388e-14,
|
|
Packit |
67cb25 |
1.102e-14,
|
|
Packit |
67cb25 |
-2.78e-15,
|
|
Packit |
67cb25 |
7.0e-16
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,18} = c_k T_k*(nu/36), nu <= 36 */
|
|
Packit |
67cb25 |
static const double coef_jnu18_a[] = {
|
|
Packit |
67cb25 |
81.5195728368096659,
|
|
Packit |
67cb25 |
25.11375537470259305,
|
|
Packit |
67cb25 |
-0.54962366347317668,
|
|
Packit |
67cb25 |
0.07508565026117689,
|
|
Packit |
67cb25 |
-0.01262524908033818,
|
|
Packit |
67cb25 |
0.00236850602019778,
|
|
Packit |
67cb25 |
-0.00047556873651929,
|
|
Packit |
67cb25 |
0.00010001889347161,
|
|
Packit |
67cb25 |
-0.00002175581482429,
|
|
Packit |
67cb25 |
4.85490251239e-6,
|
|
Packit |
67cb25 |
-1.10539483940e-6,
|
|
Packit |
67cb25 |
2.5579853343e-7,
|
|
Packit |
67cb25 |
-5.999033352e-8,
|
|
Packit |
67cb25 |
1.422747129e-8,
|
|
Packit |
67cb25 |
-3.40650521e-9,
|
|
Packit |
67cb25 |
8.2233565e-10,
|
|
Packit |
67cb25 |
-1.9993286e-10,
|
|
Packit |
67cb25 |
4.891426e-11,
|
|
Packit |
67cb25 |
-1.203343e-11,
|
|
Packit |
67cb25 |
2.97498e-12,
|
|
Packit |
67cb25 |
-7.3875e-13,
|
|
Packit |
67cb25 |
1.8418e-13,
|
|
Packit |
67cb25 |
-4.608e-14,
|
|
Packit |
67cb25 |
1.157e-14,
|
|
Packit |
67cb25 |
-2.91e-15,
|
|
Packit |
67cb25 |
7.4e-16
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,19} = c_k T_k*(nu/38), nu <= 38 */
|
|
Packit |
67cb25 |
static const double coef_jnu19_a[] = {
|
|
Packit |
67cb25 |
86.0930892477047512,
|
|
Packit |
67cb25 |
26.51042598308271729,
|
|
Packit |
67cb25 |
-0.58002730731948358,
|
|
Packit |
67cb25 |
0.07920674321589394,
|
|
Packit |
67cb25 |
-0.01331283320930301,
|
|
Packit |
67cb25 |
0.00249650841778073,
|
|
Packit |
67cb25 |
-0.00050107453900793,
|
|
Packit |
67cb25 |
0.00010534258471335,
|
|
Packit |
67cb25 |
-0.00002290511552874,
|
|
Packit |
67cb25 |
5.10946148897e-6,
|
|
Packit |
67cb25 |
-1.16292517157e-6,
|
|
Packit |
67cb25 |
2.6901365037e-7,
|
|
Packit |
67cb25 |
-6.306692473e-8,
|
|
Packit |
67cb25 |
1.495183048e-8,
|
|
Packit |
67cb25 |
-3.57869025e-9,
|
|
Packit |
67cb25 |
8.6360410e-10,
|
|
Packit |
67cb25 |
-2.0989514e-10,
|
|
Packit |
67cb25 |
5.133439e-11,
|
|
Packit |
67cb25 |
-1.262465e-11,
|
|
Packit |
67cb25 |
3.12013e-12,
|
|
Packit |
67cb25 |
-7.7455e-13,
|
|
Packit |
67cb25 |
1.9304e-13,
|
|
Packit |
67cb25 |
-4.829e-14,
|
|
Packit |
67cb25 |
1.212e-14,
|
|
Packit |
67cb25 |
-3.05e-15,
|
|
Packit |
67cb25 |
7.7e-16
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Chebyshev expansion: j_{nu,20} = c_k T_k*(nu/40), nu <= 40 */
|
|
Packit |
67cb25 |
static const double coef_jnu20_a[] = {
|
|
Packit |
67cb25 |
90.6666144195163770,
|
|
Packit |
67cb25 |
27.9070938975436823,
|
|
Packit |
67cb25 |
-0.61043045315390591,
|
|
Packit |
67cb25 |
0.08332772844325554,
|
|
Packit |
67cb25 |
-0.01400039260208282,
|
|
Packit |
67cb25 |
0.00262450494035660,
|
|
Packit |
67cb25 |
-0.00052657891389470,
|
|
Packit |
67cb25 |
0.00011066592304919,
|
|
Packit |
67cb25 |
-0.00002405432778364,
|
|
Packit |
67cb25 |
5.36399803946e-6,
|
|
Packit |
67cb25 |
-1.22044976064e-6,
|
|
Packit |
67cb25 |
2.8222728362e-7,
|
|
Packit |
67cb25 |
-6.614312964e-8,
|
|
Packit |
67cb25 |
1.567608839e-8,
|
|
Packit |
67cb25 |
-3.75084856e-9,
|
|
Packit |
67cb25 |
9.0486546e-10,
|
|
Packit |
67cb25 |
-2.1985553e-10,
|
|
Packit |
67cb25 |
5.375401e-11,
|
|
Packit |
67cb25 |
-1.321572e-11,
|
|
Packit |
67cb25 |
3.26524e-12,
|
|
Packit |
67cb25 |
-8.1033e-13,
|
|
Packit |
67cb25 |
2.0190e-13,
|
|
Packit |
67cb25 |
-5.049e-14,
|
|
Packit |
67cb25 |
1.267e-14,
|
|
Packit |
67cb25 |
-3.19e-15,
|
|
Packit |
67cb25 |
8.0e-16,
|
|
Packit |
67cb25 |
-2.0e-16
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
static const double * coef_jnu_a[] = {
|
|
Packit |
67cb25 |
0,
|
|
Packit |
67cb25 |
coef_jnu1_a,
|
|
Packit |
67cb25 |
coef_jnu2_a,
|
|
Packit |
67cb25 |
coef_jnu3_a,
|
|
Packit |
67cb25 |
coef_jnu4_a,
|
|
Packit |
67cb25 |
coef_jnu5_a,
|
|
Packit |
67cb25 |
coef_jnu6_a,
|
|
Packit |
67cb25 |
coef_jnu7_a,
|
|
Packit |
67cb25 |
coef_jnu8_a,
|
|
Packit |
67cb25 |
coef_jnu9_a,
|
|
Packit |
67cb25 |
coef_jnu10_a,
|
|
Packit |
67cb25 |
coef_jnu11_a,
|
|
Packit |
67cb25 |
coef_jnu12_a,
|
|
Packit |
67cb25 |
coef_jnu13_a,
|
|
Packit |
67cb25 |
coef_jnu14_a,
|
|
Packit |
67cb25 |
coef_jnu15_a,
|
|
Packit |
67cb25 |
coef_jnu16_a,
|
|
Packit |
67cb25 |
coef_jnu17_a,
|
|
Packit |
67cb25 |
coef_jnu18_a,
|
|
Packit |
67cb25 |
coef_jnu19_a,
|
|
Packit |
67cb25 |
coef_jnu20_a
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
static const size_t size_jnu_a[] = {
|
|
Packit |
67cb25 |
0,
|
|
Packit |
67cb25 |
sizeof(coef_jnu1_a)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu2_a)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu3_a)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu4_a)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu5_a)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu6_a)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu7_a)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu8_a)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu9_a)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu10_a)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu11_a)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu12_a)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu13_a)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu14_a)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu15_a)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu16_a)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu17_a)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu18_a)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu19_a)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu20_a)/sizeof(double)
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
static const double * coef_jnu_b[] = {
|
|
Packit |
67cb25 |
0,
|
|
Packit |
67cb25 |
coef_jnu1_b,
|
|
Packit |
67cb25 |
coef_jnu2_b,
|
|
Packit |
67cb25 |
coef_jnu3_b,
|
|
Packit |
67cb25 |
coef_jnu4_b,
|
|
Packit |
67cb25 |
coef_jnu5_b,
|
|
Packit |
67cb25 |
coef_jnu6_b,
|
|
Packit |
67cb25 |
coef_jnu7_b,
|
|
Packit |
67cb25 |
coef_jnu8_b,
|
|
Packit |
67cb25 |
coef_jnu9_b,
|
|
Packit |
67cb25 |
coef_jnu10_b
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
static const size_t size_jnu_b[] = {
|
|
Packit |
67cb25 |
0,
|
|
Packit |
67cb25 |
sizeof(coef_jnu1_b)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu2_b)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu3_b)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu4_b)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu5_b)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu6_b)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu7_b)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu8_b)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu9_b)/sizeof(double),
|
|
Packit |
67cb25 |
sizeof(coef_jnu10_b)/sizeof(double)
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Evaluate Clenshaw recurrence for
|
|
Packit |
67cb25 |
* a T* Chebyshev series.
|
|
Packit |
67cb25 |
* sizeof(c) = N+1
|
|
Packit |
67cb25 |
*/
|
|
Packit |
67cb25 |
static double
|
|
Packit |
67cb25 |
clenshaw(const double * c, int N, double u)
|
|
Packit |
67cb25 |
{
|
|
Packit |
67cb25 |
double B_np1 = 0.0;
|
|
Packit |
67cb25 |
double B_n = c[N];
|
|
Packit |
67cb25 |
double B_nm1;
|
|
Packit |
67cb25 |
int n;
|
|
Packit |
67cb25 |
for(n=N; n>0; n--) {
|
|
Packit |
67cb25 |
B_nm1 = 2.0*(2.0*u-1.0) * B_n - B_np1 + c[n-1];
|
|
Packit |
67cb25 |
B_np1 = B_n;
|
|
Packit |
67cb25 |
B_n = B_nm1;
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
return B_n - (2.0*u-1.0)*B_np1;
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* correction terms to leading McMahon expansion
|
|
Packit |
67cb25 |
* [Abramowitz+Stegun 9.5.12]
|
|
Packit |
67cb25 |
* [Olver, Royal Society Math. Tables, v. 7]
|
|
Packit |
67cb25 |
* We factor out a beta, so that this is a multiplicative
|
|
Packit |
67cb25 |
* correction:
|
|
Packit |
67cb25 |
* j_{nu,s} = beta(s,nu) * mcmahon_correction(nu, beta(s,nu))
|
|
Packit |
67cb25 |
* macmahon_correction --> 1 as s --> Inf
|
|
Packit |
67cb25 |
*/
|
|
Packit |
67cb25 |
static double
|
|
Packit |
67cb25 |
mcmahon_correction(const double mu, const double beta)
|
|
Packit |
67cb25 |
{
|
|
Packit |
67cb25 |
const double eb = 8.0*beta;
|
|
Packit |
67cb25 |
const double ebsq = eb*eb;
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
if(mu < GSL_DBL_EPSILON) {
|
|
Packit |
67cb25 |
/* Prevent division by zero below. */
|
|
Packit |
67cb25 |
const double term1 = 1.0/ebsq;
|
|
Packit |
67cb25 |
const double term2 = -4.0*31.0/(3*ebsq*ebsq);
|
|
Packit |
67cb25 |
const double term3 = 32.0*3779.0/(15.0*ebsq*ebsq*ebsq);
|
|
Packit |
67cb25 |
const double term4 = -64.0*6277237.0/(105.0*ebsq*ebsq*ebsq*ebsq);
|
|
Packit |
67cb25 |
const double term5 = 512.0*2092163573.0/(315.0*ebsq*ebsq*ebsq*ebsq*ebsq);
|
|
Packit |
67cb25 |
return 1.0 + 8.0*(term1 + term2 + term3 + term4 + term5);
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
else {
|
|
Packit |
67cb25 |
/* Here we do things in terms of 1/mu, which
|
|
Packit |
67cb25 |
* is purely to prevent overflow in the very
|
|
Packit |
67cb25 |
* unlikely case that mu is really big.
|
|
Packit |
67cb25 |
*/
|
|
Packit |
67cb25 |
const double mi = 1.0/mu;
|
|
Packit |
67cb25 |
const double r = mu/ebsq;
|
|
Packit |
67cb25 |
const double n2 = 4.0/3.0 * (7.0 - 31.0*mi);
|
|
Packit |
67cb25 |
const double n3 = 32.0/15.0 * (83.0 + (-982.0 + 3779.0*mi)*mi);
|
|
Packit |
67cb25 |
const double n4 = 64.0/105.0 * (6949.0 + (-153855.0 + (1585743.0 - 6277237.0*mi)*mi)*mi);
|
|
Packit |
67cb25 |
const double n5 = 512.0/315.0 * (70197.0 + (-2479316.0 + (48010494.0 + (-512062548.0 + 2092163573.0*mi)*mi)*mi)*mi);
|
|
Packit |
67cb25 |
const double n6 = 2048.0/3465.0 * (5592657.0 + (-287149133.0 + (8903961290.0 + (-179289628602.0 + (1982611456181.0 - 8249725736393.0*mi)*mi)*mi)*mi)*mi);
|
|
Packit |
67cb25 |
const double term1 = (1.0 - mi) * r;
|
|
Packit |
67cb25 |
const double term2 = term1 * n2 * r;
|
|
Packit |
67cb25 |
const double term3 = term1 * n3 * r*r;
|
|
Packit |
67cb25 |
const double term4 = term1 * n4 * r*r*r;
|
|
Packit |
67cb25 |
const double term5 = term1 * n5 * r*r*r*r;
|
|
Packit |
67cb25 |
const double term6 = term1 * n6 * r*r*r*r*r;
|
|
Packit |
67cb25 |
return 1.0 - 8.0*(term1 + term2 + term3 + term4 + term5 + term6);
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* Assumes z >= 1.0 */
|
|
Packit |
67cb25 |
static double
|
|
Packit |
67cb25 |
olver_b0(double z, double minus_zeta)
|
|
Packit |
67cb25 |
{
|
|
Packit |
67cb25 |
if(z < 1.02) {
|
|
Packit |
67cb25 |
const double a = 1.0-z;
|
|
Packit |
67cb25 |
const double c0 = 0.0179988721413553309252458658183;
|
|
Packit |
67cb25 |
const double c1 = 0.0111992982212877614645974276203;
|
|
Packit |
67cb25 |
const double c2 = 0.0059404069786014304317781160605;
|
|
Packit |
67cb25 |
const double c3 = 0.0028676724516390040844556450173;
|
|
Packit |
67cb25 |
const double c4 = 0.0012339189052567271708525111185;
|
|
Packit |
67cb25 |
const double c5 = 0.0004169250674535178764734660248;
|
|
Packit |
67cb25 |
const double c6 = 0.0000330173385085949806952777365;
|
|
Packit |
67cb25 |
const double c7 = -0.0001318076238578203009990106425;
|
|
Packit |
67cb25 |
const double c8 = -0.0001906870370050847239813945647;
|
|
Packit |
67cb25 |
return c0 + a*(c1 + a*(c2 + a*(c3 + a*(c4 + a*(c5 + a*(c6 + a*(c7 + a*c8)))))));
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
else {
|
|
Packit |
67cb25 |
const double abs_zeta = minus_zeta;
|
|
Packit |
67cb25 |
const double t = 1.0/(z*sqrt(1.0 - 1.0/(z*z)));
|
|
Packit |
67cb25 |
return -5.0/(48.0*abs_zeta*abs_zeta) + t*(3.0 + 5.0*t*t)/(24.0*sqrt(abs_zeta));
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
inline
|
|
Packit |
67cb25 |
static double
|
|
Packit |
67cb25 |
olver_f1(double z, double minus_zeta)
|
|
Packit |
67cb25 |
{
|
|
Packit |
67cb25 |
const double b0 = olver_b0(z, minus_zeta);
|
|
Packit |
67cb25 |
const double h2 = sqrt(4.0*minus_zeta/(z*z-1.0)); /* FIXME */
|
|
Packit |
67cb25 |
return 0.5 * z * h2 * b0;
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
int
|
|
Packit |
67cb25 |
gsl_sf_bessel_zero_J0_e(unsigned int s, gsl_sf_result * result)
|
|
Packit |
67cb25 |
{
|
|
Packit |
67cb25 |
/* CHECK_POINTER(result) */
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
if(s == 0){
|
|
Packit |
67cb25 |
result->val = 0.0;
|
|
Packit |
67cb25 |
result->err = 0.0;
|
|
Packit |
67cb25 |
GSL_ERROR ("error", GSL_EINVAL);
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
else {
|
|
Packit |
67cb25 |
/* See [F. Lether, J. Comp. Appl .Math. 67, 167 (1996)]. */
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
static const double P[] = { 1567450796.0/12539606369.0,
|
|
Packit |
67cb25 |
8903660.0/2365861.0,
|
|
Packit |
67cb25 |
10747040.0/536751.0,
|
|
Packit |
67cb25 |
17590991.0/1696654.0
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
static const double Q[] = { 1.0,
|
|
Packit |
67cb25 |
29354255.0/954518.0,
|
|
Packit |
67cb25 |
76900001.0/431847.0,
|
|
Packit |
67cb25 |
67237052.0/442411.0
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
const double beta = (s - 0.25) * M_PI;
|
|
Packit |
67cb25 |
const double bi2 = 1.0/(beta*beta);
|
|
Packit |
67cb25 |
const double R33num = P[0] + bi2 * (P[1] + bi2 * (P[2] + P[3] * bi2));
|
|
Packit |
67cb25 |
const double R33den = Q[0] + bi2 * (Q[1] + bi2 * (Q[2] + Q[3] * bi2));
|
|
Packit |
67cb25 |
const double R33 = R33num/R33den;
|
|
Packit |
67cb25 |
result->val = beta + R33/beta;
|
|
Packit |
67cb25 |
result->err = fabs(3.0e-15 * result->val);
|
|
Packit |
67cb25 |
return GSL_SUCCESS;
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
int
|
|
Packit |
67cb25 |
gsl_sf_bessel_zero_J1_e(unsigned int s, gsl_sf_result * result)
|
|
Packit |
67cb25 |
{
|
|
Packit |
67cb25 |
/* CHECK_POINTER(result) */
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
if(s == 0) {
|
|
Packit |
67cb25 |
result->val = 0.0;
|
|
Packit |
67cb25 |
result->err = 0.0;
|
|
Packit |
67cb25 |
return GSL_SUCCESS;
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
else {
|
|
Packit |
67cb25 |
/* See [M. Branders et al., J. Comp. Phys. 42, 403 (1981)]. */
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
static const double a[] = { -0.362804405737084,
|
|
Packit |
67cb25 |
0.120341279038597,
|
|
Packit |
67cb25 |
0.439454547101171e-01,
|
|
Packit |
67cb25 |
0.159340088474713e-02
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
static const double b[] = { 1.0,
|
|
Packit |
67cb25 |
-0.325641790801361,
|
|
Packit |
67cb25 |
-0.117453445968927,
|
|
Packit |
67cb25 |
-0.424906902601794e-02
|
|
Packit |
67cb25 |
};
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
const double beta = (s + 0.25) * M_PI;
|
|
Packit |
67cb25 |
const double bi2 = 1.0/(beta*beta);
|
|
Packit |
67cb25 |
const double Rnum = a[3] + bi2 * (a[2] + bi2 * (a[1] + bi2 * a[0]));
|
|
Packit |
67cb25 |
const double Rden = b[3] + bi2 * (b[2] + bi2 * (b[1] + bi2 * b[0]));
|
|
Packit |
67cb25 |
const double R = Rnum/Rden;
|
|
Packit |
67cb25 |
result->val = beta * (1.0 + R*bi2);
|
|
Packit |
67cb25 |
result->err = fabs(2.0e-14 * result->val);
|
|
Packit |
67cb25 |
return GSL_SUCCESS;
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
int
|
|
Packit |
67cb25 |
gsl_sf_bessel_zero_Jnu_e(double nu, unsigned int s, gsl_sf_result * result)
|
|
Packit |
67cb25 |
{
|
|
Packit |
67cb25 |
/* CHECK_POINTER(result) */
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
if(nu <= -1.0) {
|
|
Packit |
67cb25 |
DOMAIN_ERROR(result);
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
else if(s == 0) {
|
|
Packit |
67cb25 |
result->val = 0.0;
|
|
Packit |
67cb25 |
result->err = 0.0;
|
|
Packit |
67cb25 |
if (nu == 0.0) {
|
|
Packit |
67cb25 |
GSL_ERROR ("no zero-th root for nu = 0.0", GSL_EINVAL);
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
return GSL_SUCCESS;
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
else if(nu < 0.0) {
|
|
Packit |
67cb25 |
/* This can be done, I'm just lazy now. */
|
|
Packit |
67cb25 |
result->val = 0.0;
|
|
Packit |
67cb25 |
result->err = 0.0;
|
|
Packit |
67cb25 |
GSL_ERROR("unimplemented", GSL_EUNIMPL);
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
else if(s == 1) {
|
|
Packit |
67cb25 |
/* Chebyshev fits for the first positive zero.
|
|
Packit |
67cb25 |
* For some reason Nemeth made this different from the others.
|
|
Packit |
67cb25 |
*/
|
|
Packit |
67cb25 |
if(nu < 2.0) {
|
|
Packit |
67cb25 |
const double * c = coef_jnu_a[s];
|
|
Packit |
67cb25 |
const size_t L = size_jnu_a[s];
|
|
Packit |
67cb25 |
const double arg = nu/2.0;
|
|
Packit |
67cb25 |
const double chb = clenshaw(c, L-1, arg);
|
|
Packit |
67cb25 |
result->val = chb;
|
|
Packit |
67cb25 |
result->err = 2.0e-15 * result->val;
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
else {
|
|
Packit |
67cb25 |
const double * c = coef_jnu_b[s];
|
|
Packit |
67cb25 |
const size_t L = size_jnu_b[s];
|
|
Packit |
67cb25 |
const double arg = pow(2.0/nu, 2.0/3.0);
|
|
Packit |
67cb25 |
const double chb = clenshaw(c, L-1, arg);
|
|
Packit |
67cb25 |
result->val = nu * chb;
|
|
Packit |
67cb25 |
result->err = 2.0e-15 * result->val;
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
return GSL_SUCCESS;
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
else if(s <= 10) {
|
|
Packit |
67cb25 |
/* Chebyshev fits for the first 10 positive zeros. */
|
|
Packit |
67cb25 |
if(nu < s) {
|
|
Packit |
67cb25 |
const double * c = coef_jnu_a[s];
|
|
Packit |
67cb25 |
const size_t L = size_jnu_a[s];
|
|
Packit |
67cb25 |
const double arg = nu/s;
|
|
Packit |
67cb25 |
const double chb = clenshaw(c, L-1, arg);
|
|
Packit |
67cb25 |
result->val = chb;
|
|
Packit |
67cb25 |
result->err = 2.0e-15 * result->val;
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
else {
|
|
Packit |
67cb25 |
const double * c = coef_jnu_b[s];
|
|
Packit |
67cb25 |
const size_t L = size_jnu_b[s];
|
|
Packit |
67cb25 |
const double arg = pow(s/nu, 2.0/3.0);
|
|
Packit |
67cb25 |
const double chb = clenshaw(c, L-1, arg);
|
|
Packit |
67cb25 |
result->val = nu * chb;
|
|
Packit |
67cb25 |
result->err = 2.0e-15 * result->val;
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/* FIXME: truth in advertising for the screwed up
|
|
Packit |
67cb25 |
* s = 5 fit. Need to fix that.
|
|
Packit |
67cb25 |
*/
|
|
Packit |
67cb25 |
if(s == 5) {
|
|
Packit |
67cb25 |
result->err *= 5.0e+06;
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
return GSL_SUCCESS;
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
else if(s > 0.5*nu && s <= 20) {
|
|
Packit |
67cb25 |
/* Chebyshev fits for 10 < s <= 20. */
|
|
Packit |
67cb25 |
const double * c = coef_jnu_a[s];
|
|
Packit |
67cb25 |
const size_t L = size_jnu_a[s];
|
|
Packit |
67cb25 |
const double arg = nu/(2.0*s);
|
|
Packit |
67cb25 |
const double chb = clenshaw(c, L-1, arg);
|
|
Packit |
67cb25 |
result->val = chb;
|
|
Packit |
67cb25 |
result->err = 4.0e-15 * chb;
|
|
Packit |
67cb25 |
return GSL_SUCCESS;
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
else if(s > 2.0 * nu) {
|
|
Packit |
67cb25 |
/* McMahon expansion if s is large compared to nu. */
|
|
Packit |
67cb25 |
const double beta = (s + 0.5*nu - 0.25) * M_PI;
|
|
Packit |
67cb25 |
const double mc = mcmahon_correction(4.0*nu*nu, beta);
|
|
Packit |
67cb25 |
gsl_sf_result rat12;
|
|
Packit |
67cb25 |
gsl_sf_pow_int_e(nu/beta, 14, &rat12);
|
|
Packit |
67cb25 |
result->val = beta * mc;
|
|
Packit |
67cb25 |
result->err = 4.0 * fabs(beta) * rat12.val;
|
|
Packit |
67cb25 |
result->err += 4.0 * fabs(GSL_DBL_EPSILON * result->val);
|
|
Packit |
67cb25 |
return GSL_SUCCESS;
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
else {
|
|
Packit |
67cb25 |
/* Olver uniform asymptotic. */
|
|
Packit |
67cb25 |
gsl_sf_result as;
|
|
Packit |
67cb25 |
const int stat_as = gsl_sf_airy_zero_Ai_e(s, &as);
|
|
Packit |
67cb25 |
const double minus_zeta = -pow(nu,-2.0/3.0) * as.val;
|
|
Packit |
67cb25 |
const double z = gsl_sf_bessel_Olver_zofmzeta(minus_zeta);
|
|
Packit |
67cb25 |
const double f1 = olver_f1(z, minus_zeta);
|
|
Packit |
67cb25 |
result->val = nu * (z + f1/(nu*nu));
|
|
Packit |
67cb25 |
result->err = 0.001/(nu*nu*nu);
|
|
Packit |
67cb25 |
result->err += 2.0 * GSL_DBL_EPSILON * fabs(result->val);
|
|
Packit |
67cb25 |
return stat_as;
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
/*-*-*-*-*-*-*-*-*-* Functions w/ Natural Prototypes *-*-*-*-*-*-*-*-*-*-*/
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
#include "eval.h"
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
double gsl_sf_bessel_zero_J0(unsigned int s)
|
|
Packit |
67cb25 |
{
|
|
Packit |
67cb25 |
EVAL_RESULT(gsl_sf_bessel_zero_J0_e(s, &result));
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
double gsl_sf_bessel_zero_J1(unsigned int s)
|
|
Packit |
67cb25 |
{
|
|
Packit |
67cb25 |
EVAL_RESULT(gsl_sf_bessel_zero_J1_e(s, &result));
|
|
Packit |
67cb25 |
}
|
|
Packit |
67cb25 |
|
|
Packit |
67cb25 |
double gsl_sf_bessel_zero_Jnu(double nu, unsigned int s)
|
|
Packit |
67cb25 |
{
|
|
Packit |
67cb25 |
EVAL_RESULT(gsl_sf_bessel_zero_Jnu_e(nu, s, &result));
|
|
Packit |
67cb25 |
}
|