46 #include <gsl/gsl_sf.h> 50 #define PI 3.14159265358979323846264338328 54 int compare(
const std::string &
name,
double v1,
double v2,
double scale = 2.0) {
57 std::cout << std::setw(50) << std::left << name <<
":\t";
62 double delta = v2 - v1;
64 if (delta < 0 ) delta = - delta;
65 if (v1 == 0 || v2 == 0) {
76 if ( delta/d > eps && delta > eps )
81 std::cout <<
" OK" << std::endl;
83 std::cout <<
" FAILED " << std::endl;
84 int pr = std::cout.precision (18);
85 std::cout <<
"\nDiscrepancy in " << name <<
"() :\n " << v1 <<
" != " << v2 <<
" discr = " << int(delta/d/eps)
86 <<
" (Allowed discrepancy is " << eps <<
")\n\n";
87 std::cout.precision (pr);
116 std::cout.precision(20);
132 const int inc_gamma_scale = 100;
134 const int inc_gamma_scale = 200;
167 iret |=
compare(
"conf_hyperg(1, 1.5, 1) ",
conf_hyperg(1, 1.5, 1), 2.0300784692787049755);
171 iret |=
compare(
"cyl_bessel_j(0.75, 1.0) ",
cyl_bessel_j(0.75, 1.0), 0.5586524932048917478, 16);
175 iret |=
compare(
"cyl_neumann(0.75, 1.0) ",
cyl_neumann(0.75, 1.0), -0.6218694174429746383 );
177 iret |=
compare(
"ellint_1(0.50, PI/3.0) ",
ellint_1(0.50,
PI/3.0), 1.0895506700518854093);
179 iret |=
compare(
"ellint_2(0.50, PI/3.0) ",
ellint_2(0.50,
PI/3.0), 1.0075555551444720293);
181 iret |=
compare(
"ellint_3(-0.50, 0.5, PI/3.0) ",
ellint_3(-0.50, 0.5,
PI/3.0), 0.9570574331323584890);
183 iret |=
compare(
"expint(1.0) ",
expint(1.0), 1.8951178163559367555);
187 iret |=
compare(
"hyperg(8, -8, 1, 0.5) ",
hyperg(8, -8, 1, 0.5), 0.13671875);
191 iret |=
compare(
"legendre(10, -0.5) ",
legendre(10, -0.5), -0.1882286071777345);
193 iret |=
compare(
"riemann_zeta(-0.5) ",
riemann_zeta(-0.5), -0.207886224977354566017307, 16);
195 iret |=
compare(
"sph_bessel(1, 10.0) ",
sph_bessel(1, 10.0), 0.07846694179875154709000);
201 iret |=
compare(
"airy_Ai(-0.5) ",
airy_Ai(-0.5), 0.475728091610539583);
203 iret |=
compare(
"airy_Bi(0.5) ",
airy_Bi(0.5), 0.854277043103155553);
218 std::cout <<
"\n\nError: Special Functions Test FAILED !!!!!" << std::endl;
226 #ifdef CHECK_WITH_GSL 230 iret = gsl_sf_ellint_P_e(
PI/2.0, 0.5, -0.5, GSL_PREC_DOUBLE, &r);
231 std::cout <<
"comp_ellint_3(0.50, 0.5) : " << r.val <<
" err: " << r.err <<
" iret: " << iret << std::endl;
233 iret = gsl_sf_ellint_P_e(
PI/3.0, 0.5, 0.5, GSL_PREC_DOUBLE, &r);
234 std::cout <<
"ellint_3(0.50, 0.5, PI/3.0) : " << r.val <<
" err: " << r.err <<
" iret: " << iret << std::endl;
236 iret = gsl_sf_zeta_e(-0.5, &r);
237 std::cout <<
"riemann_zeta(-0.5) : " << r.val <<
" err: " << r.err <<
" iret: " << iret << std::endl;
double erf(double x)
Error function encountered in integrating the normal distribution.
double airy_zero_Bi(unsigned int s)
Calculates the zeroes of the Airy function Bi.
double sph_legendre(unsigned l, unsigned m, double theta)
Computes the spherical (normalized) associated Legendre polynomials, or spherical harmonic without az...
double cyl_bessel_j(double nu, double x)
Calculates the (cylindrical) Bessel functions of the first kind (also called regular (cylindrical) Be...
double inc_beta(double x, double a, double b)
Calculates the normalized (regularized) incomplete beta function.
double airy_zero_Ai_deriv(unsigned int s)
Calculates the zeroes of the derivative of the Airy function Ai.
double ellint_2(double k, double phi)
Calculates the complete elliptic integral of the second kind.
double tgamma(double x)
The gamma function is defined to be the extension of the factorial to real numbers.
double legendre(unsigned l, double x)
Calculates the Legendre polynomials.
double beta(double x, double y)
Calculates the beta function.
double erfc(double x)
Complementary error function.
double inc_gamma_c(double a, double x)
Calculates the normalized (regularized) upper incomplete gamma function (upper integral) ...
double comp_ellint_1(double k)
Calculates the complete elliptic integral of the first kind.
double sph_bessel(unsigned n, double x)
Calculates the spherical Bessel functions of the first kind (also called regular spherical Bessel fun...
double conf_hyperg(double a, double b, double z)
Calculates the confluent hypergeometric functions of the first kind.
double hyperg(double a, double b, double c, double x)
Calculates Gauss' hypergeometric function.
double airy_zero_Bi_deriv(unsigned int s)
Calculates the zeroes of the derivative of the Airy function Bi.
double expint(double x)
Calculates the exponential integral.
double ellint_3(double n, double k, double phi)
Calculates the complete elliptic integral of the third kind.
double riemann_zeta(double x)
Calculates the Riemann zeta function.
double comp_ellint_2(double k)
Calculates the complete elliptic integral of the second kind.
double assoc_laguerre(unsigned n, double m, double x)
Computes the generalized Laguerre polynomials for .
double assoc_legendre(unsigned l, unsigned m, double x)
Computes the associated Legendre polynomials.
double airy_zero_Ai(unsigned int s)
Calculates the zeroes of the Airy function Ai.
double airy_Bi(double x)
Calculates the Airy function Bi.
double laguerre(unsigned n, double x)
Calculates the Laguerre polynomials.
double airy_Bi_deriv(double x)
Calculates the derivative of the Airy function Bi.
double comp_ellint_3(double n, double k)
Calculates the complete elliptic integral of the third kind.
double airy_Ai_deriv(double x)
Calculates the derivative of the Airy function Ai.
double ellint_1(double k, double phi)
Calculates the incomplete elliptic integral of the first kind.
int compare(const std::string &name, double v1, double v2, double scale=2.0)
double cyl_neumann(double nu, double x)
Calculates the (cylindrical) Bessel functions of the second kind (also called irregular (cylindrical)...
#define PI
Test file for the special functions implemented in MathMore.
double sph_neumann(unsigned n, double x)
Calculates the spherical Bessel functions of the second kind (also called irregular spherical Bessel ...
double cyl_bessel_i(double nu, double x)
Calculates the modified Bessel function of the first kind (also called regular modified (cylindrical)...
double lgamma(double x)
Calculates the logarithm of the gamma function.
double airy_Ai(double x)
Calculates the Airy function Ai.
double inc_gamma(double a, double x)
Calculates the normalized (regularized) lower incomplete gamma function (lower integral) ...
double cyl_bessel_k(double nu, double x)
Calculates the modified Bessel functions of the second kind (also called irregular modified (cylindri...