File indexing completed on 2024-05-12 03:47:53
0001 /* 0002 File : nsl_sf_basic.c 0003 Project : LabPlot 0004 Description : NSL special basic functions 0005 -------------------------------------------------------------------- 0006 SPDX-FileCopyrightText: 2018-2022 Stefan Gerlach <stefan.gerlach@uni.kn> 0007 SPDX-License-Identifier: GPL-2.0-or-later 0008 */ 0009 0010 #include "nsl_sf_basic.h" 0011 #include <gsl/gsl_math.h> 0012 #include <gsl/gsl_randist.h> 0013 #include <gsl/gsl_sf.h> 0014 #include <math.h> 0015 #include <stdlib.h> 0016 #ifdef HAVE_LIBCERF 0017 #include <cerf.h> 0018 #elif !defined(_MSC_VER) 0019 #include "Faddeeva.h" 0020 #endif 0021 0022 double nsl_sf_dummy(double x) { 0023 return 0 * x; // "use" parameter 0024 } 0025 0026 double nsl_sf_dummy2(double p, double x) { 0027 return 0 * p * x; // "use" parameter 0028 } 0029 0030 /* stdlib.h */ 0031 double nsl_sf_rand(void) { 0032 return rand(); 0033 } 0034 #if defined(HAVE_RANDOM_FUNCTION) 0035 double nsl_sf_random(void) { 0036 return random(); 0037 } 0038 double nsl_sf_drand(void) { 0039 return random() / (double)RAND_MAX; 0040 } 0041 #else 0042 double nsl_sf_random(void) { 0043 return rand(); 0044 } 0045 double nsl_sf_drand(void) { 0046 return rand() / (double)RAND_MAX; 0047 } 0048 #endif 0049 0050 double nsl_sf_sgn(double x) { 0051 #ifndef _WIN32 0052 return copysign(1.0, x); 0053 #else 0054 if (x == 0) 0055 return 0; 0056 else 0057 return GSL_SIGN(x); 0058 #endif 0059 } 0060 0061 double nsl_sf_theta(double x) { 0062 if (x >= 0) 0063 return 1; 0064 else 0065 return 0; 0066 } 0067 0068 double nsl_sf_exp10(double x) { 0069 #ifdef __linux__ 0070 return exp10(x); 0071 #else 0072 return pow(10., x); 0073 #endif 0074 } 0075 0076 /* non-uniform random number generation 0077 * using a generator chosen by the environment variable GSL_RNG_TYPE 0078 */ 0079 #define SETUP_GSL_RNG \ 0080 gsl_rng_env_setup(); \ 0081 const gsl_rng_type* T = gsl_rng_default; \ 0082 gsl_rng* r = gsl_rng_alloc(T); \ 0083 gsl_rng_set(r, rand()); /*seed*/ 0084 0085 double nsl_sf_ran_gaussian(double sigma) { 0086 SETUP_GSL_RNG 0087 return gsl_ran_gaussian_ziggurat(r, sigma); 0088 } 0089 double nsl_sf_ran_exponential(double mu) { 0090 SETUP_GSL_RNG 0091 return gsl_ran_exponential(r, mu); 0092 } 0093 double nsl_sf_ran_laplace(double a) { 0094 SETUP_GSL_RNG 0095 return gsl_ran_laplace(r, a); 0096 } 0097 double nsl_sf_ran_cauchy(double a) { 0098 SETUP_GSL_RNG 0099 return gsl_ran_cauchy(r, a); 0100 } 0101 double nsl_sf_ran_rayleigh(double sigma) { 0102 SETUP_GSL_RNG 0103 return gsl_ran_rayleigh(r, sigma); 0104 } 0105 double nsl_sf_ran_landau(void) { 0106 SETUP_GSL_RNG 0107 return gsl_ran_landau(r); 0108 } 0109 double nsl_sf_ran_levy(double c, double alpha) { 0110 SETUP_GSL_RNG 0111 return gsl_ran_levy(r, c, alpha); 0112 } 0113 double nsl_sf_ran_gamma(double a, double b) { 0114 SETUP_GSL_RNG 0115 return gsl_ran_gamma(r, a, b); 0116 } 0117 double nsl_sf_ran_flat(double a, double b) { 0118 SETUP_GSL_RNG 0119 return gsl_ran_flat(r, a, b); 0120 } 0121 double nsl_sf_ran_lognormal(double zeta, double sigma) { 0122 SETUP_GSL_RNG 0123 return gsl_ran_lognormal(r, zeta, sigma); 0124 } 0125 double nsl_sf_ran_chisq(double nu) { 0126 SETUP_GSL_RNG 0127 return gsl_ran_chisq(r, nu); 0128 } 0129 double nsl_sf_ran_tdist(double nu) { 0130 SETUP_GSL_RNG 0131 return gsl_ran_tdist(r, nu); 0132 } 0133 double nsl_sf_ran_logistic(double a) { 0134 SETUP_GSL_RNG 0135 return gsl_ran_logistic(r, a); 0136 } 0137 0138 double nsl_sf_ran_poisson(double mu) { 0139 SETUP_GSL_RNG 0140 return (double)gsl_ran_poisson(r, mu); 0141 } 0142 double nsl_sf_ran_bernoulli(double p) { 0143 SETUP_GSL_RNG 0144 return (double)gsl_ran_bernoulli(r, p); 0145 } 0146 double nsl_sf_ran_binomial(double p, double n) { 0147 SETUP_GSL_RNG 0148 return (double)gsl_ran_binomial(r, p, (unsigned int)round(n)); 0149 } 0150 0151 /* 0152 * source: https://stackoverflow.com/questions/11376288/fast-computing-of-log2-for-64-bit-integers 0153 * source: https://graphics.stanford.edu/~seander/bithacks.html#IntegerLogLookup 0154 */ 0155 int nsl_sf_log2_int(unsigned int x) { 0156 #ifdef _MSC_VER 0157 return nsl_sf_log2_int2(x); 0158 #else 0159 return nsl_sf_log2_int0(x); 0160 #endif 0161 } 0162 int nsl_sf_log2_int0(unsigned int x) { 0163 #ifdef _MSC_VER 0164 return 0; /* not implemented yet */ 0165 #else 0166 return (int)(8 * sizeof(unsigned int) - __builtin_clz(x) - 1); 0167 #endif 0168 } 0169 int nsl_sf_log2_int2(int x) { 0170 const signed char LogTable256[256] = {-1, 0, 1, 1, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3, 3, 3, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 0171 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 5, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 0172 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 0173 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 6, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 0174 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 0175 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 0176 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7}; 0177 0178 unsigned int r; // r will be lg(v) 0179 unsigned int t, tt; // temporaries 0180 if ((tt = x >> 16)) 0181 r = (t = tt >> 8) ? 24 + LogTable256[t] : 16 + LogTable256[tt]; 0182 else 0183 r = (t = x >> 8) ? 8 + LogTable256[t] : LogTable256[x]; 0184 0185 return r; 0186 } 0187 int nsl_sf_log2_int3(uint64_t value) { 0188 const int tab64[64] = {63, 0, 58, 1, 59, 47, 53, 2, 60, 39, 48, 27, 54, 33, 42, 3, 61, 51, 37, 40, 49, 18, 28, 20, 55, 30, 34, 11, 43, 14, 22, 4, 0189 62, 57, 46, 52, 38, 26, 32, 41, 50, 36, 17, 19, 29, 10, 13, 21, 56, 45, 25, 31, 35, 16, 9, 12, 44, 24, 15, 8, 23, 7, 6, 5}; 0190 0191 value |= value >> 1; 0192 value |= value >> 2; 0193 value |= value >> 4; 0194 value |= value >> 8; 0195 value |= value >> 16; 0196 value |= value >> 32; 0197 0198 return tab64[((uint64_t)((value - (value >> 1)) * 0x07EDD5E59A4E28C2)) >> 58]; 0199 } 0200 int nsl_sf_log2p1_int(int x) { 0201 // fastest method 0202 return nsl_sf_log2_int(x) + 1; 0203 // TODO: why is this so slow 0204 // return (int)log2(x) + 1; 0205 } 0206 int nsl_sf_log2_longlong(unsigned long long x) { 0207 #ifdef _MSC_VER 0208 return 0; /* not implemented yet */ 0209 #else 0210 return (int)(8 * sizeof(unsigned long long) - __builtin_clzll(x) - 1); 0211 #endif 0212 } 0213 0214 double nsl_sf_sec(double x) { 0215 return 1. / cos(x); 0216 } 0217 double nsl_sf_csc(double x) { 0218 return 1. / sin(x); 0219 } 0220 double nsl_sf_cot(double x) { 0221 return 1. / tan(x); 0222 } 0223 double nsl_sf_asec(double x) { 0224 return acos(1. / x); 0225 } 0226 double nsl_sf_acsc(double x) { 0227 return asin(1. / x); 0228 } 0229 double nsl_sf_acot(double x) { 0230 if (x > 0) 0231 return atan(1. / x); 0232 else 0233 return atan(1. / x) + M_PI; 0234 } 0235 double nsl_sf_sech(double x) { 0236 return 1. / cosh(x); 0237 } 0238 double nsl_sf_csch(double x) { 0239 return 1. / sinh(x); 0240 } 0241 double nsl_sf_coth(double x) { 0242 return 1. / tanh(x); 0243 } 0244 double nsl_sf_asech(double x) { 0245 return gsl_acosh(1. / x); 0246 } 0247 double nsl_sf_acsch(double x) { 0248 return gsl_asinh(1. / x); 0249 } 0250 double nsl_sf_acoth(double x) { 0251 return gsl_atanh(1. / x); 0252 } 0253 0254 double nsl_sf_harmonic(double x) { 0255 // check if x is a negative integer 0256 if (x < 0 && !gsl_fcmp(round(x) - x, 0., 1.e-16)) 0257 return GSL_POSINF; 0258 0259 return gsl_sf_psi(x + 1) + M_EULER; 0260 } 0261 0262 /* error functions and related */ 0263 double nsl_sf_erfcx(double x) { 0264 #ifdef HAVE_LIBCERF 0265 return erfcx(x); 0266 #elif defined(_MSC_VER) 0267 return 0.; // not supported yet 0268 #else 0269 return Faddeeva_erfcx_re(x); 0270 #endif 0271 } 0272 0273 double nsl_sf_erfi(double x) { 0274 #ifdef HAVE_LIBCERF 0275 return erfi(x); 0276 #elif defined(_MSC_VER) 0277 return 0.; // not supported yet 0278 #else 0279 return Faddeeva_erfi_re(x); 0280 #endif 0281 } 0282 0283 double nsl_sf_im_w_of_x(double x) { 0284 #ifdef HAVE_LIBCERF 0285 return im_w_of_x(x); 0286 #elif defined(_MSC_VER) 0287 return 0.; // not supported yet 0288 #else 0289 return Faddeeva_w_im(x); 0290 #endif 0291 } 0292 0293 #if !defined(_MSC_VER) 0294 double nsl_sf_im_w_of_z(COMPLEX z) { 0295 #ifdef HAVE_LIBCERF 0296 return cimag(w_of_z(z)); 0297 #else 0298 return cimag(Faddeeva_w(z, 0)); 0299 #endif 0300 } 0301 #endif 0302 0303 double nsl_sf_dawson(double x) { 0304 #ifdef HAVE_LIBCERF 0305 return dawson(x); 0306 #elif defined(_MSC_VER) 0307 return 0.; // not supported yet 0308 #else 0309 return Faddeeva_Dawson_re(x); 0310 #endif 0311 } 0312 0313 double nsl_sf_voigt(double x, double sigma, double gamma) { 0314 #ifdef HAVE_LIBCERF 0315 return voigt(x, sigma, gamma); 0316 #elif defined(_MSC_VER) 0317 return 0.; // not supported yet 0318 #else 0319 COMPLEX z = (x + I * gamma) / (M_SQRT2 * sigma); 0320 return creal(Faddeeva_w(z, 0)) / (M_SQRT2 * M_SQRTPI * sigma); 0321 #endif 0322 } 0323 0324 double nsl_sf_pseudovoigt(double x, double eta, double sigma, double gamma) { 0325 if (sigma == 0 || gamma == 0) 0326 return 0; 0327 // TODO: what if eta < 0 or > 1? 0328 0329 return (1. - eta) * gsl_ran_gaussian_pdf(x, sigma) + eta * gsl_ran_cauchy_pdf(x, gamma); 0330 } 0331 0332 double nsl_sf_pseudovoigt1(double x, double eta, double w) { 0333 // 2w - FWHM, sigma_G = w/sqrt(2ln(2)) 0334 return nsl_sf_pseudovoigt(x, eta, w / sqrt(2. * M_LN2), w); 0335 } 0336 0337 /* wrapper for GSL functions with integer parameters */ 0338 #define MODE GSL_PREC_DOUBLE 0339 /* mathematical functions */ 0340 double nsl_sf_ldexp(double x, double expo) { 0341 return gsl_ldexp(x, (int)round(expo)); 0342 } 0343 double nsl_sf_powint(double x, double n) { 0344 return gsl_sf_pow_int(x, (int)round(n)); 0345 } 0346 /* Airy functions */ 0347 double nsl_sf_airy_Ai(double x) { 0348 return gsl_sf_airy_Ai(x, MODE); 0349 } 0350 double nsl_sf_airy_Bi(double x) { 0351 return gsl_sf_airy_Bi(x, MODE); 0352 } 0353 double nsl_sf_airy_Ais(double x) { 0354 return gsl_sf_airy_Ai_scaled(x, MODE); 0355 } 0356 double nsl_sf_airy_Bis(double x) { 0357 return gsl_sf_airy_Bi_scaled(x, MODE); 0358 } 0359 double nsl_sf_airy_Aid(double x) { 0360 return gsl_sf_airy_Ai_deriv(x, MODE); 0361 } 0362 double nsl_sf_airy_Bid(double x) { 0363 return gsl_sf_airy_Bi_deriv(x, MODE); 0364 } 0365 double nsl_sf_airy_Aids(double x) { 0366 return gsl_sf_airy_Ai_deriv_scaled(x, MODE); 0367 } 0368 double nsl_sf_airy_Bids(double x) { 0369 return gsl_sf_airy_Bi_deriv_scaled(x, MODE); 0370 } 0371 double nsl_sf_airy_0_Ai(double s) { 0372 return gsl_sf_airy_zero_Ai((unsigned int)round(s)); 0373 } 0374 double nsl_sf_airy_0_Bi(double s) { 0375 return gsl_sf_airy_zero_Bi((unsigned int)round(s)); 0376 } 0377 double nsl_sf_airy_0_Aid(double s) { 0378 return gsl_sf_airy_zero_Ai_deriv((unsigned int)round(s)); 0379 } 0380 double nsl_sf_airy_0_Bid(double s) { 0381 return gsl_sf_airy_zero_Bi_deriv((unsigned int)round(s)); 0382 } 0383 /* Bessel functions */ 0384 double nsl_sf_bessel_Jn(double n, double x) { 0385 return gsl_sf_bessel_Jn((int)round(n), x); 0386 } 0387 double nsl_sf_bessel_Yn(double n, double x) { 0388 return gsl_sf_bessel_Yn((int)round(n), x); 0389 } 0390 double nsl_sf_bessel_In(double n, double x) { 0391 return gsl_sf_bessel_In((int)round(n), x); 0392 } 0393 double nsl_sf_bessel_Ins(double n, double x) { 0394 return gsl_sf_bessel_In_scaled((int)round(n), x); 0395 } 0396 double nsl_sf_bessel_Kn(double n, double x) { 0397 return gsl_sf_bessel_Kn((int)round(n), x); 0398 } 0399 double nsl_sf_bessel_Kns(double n, double x) { 0400 return gsl_sf_bessel_Kn_scaled((int)round(n), x); 0401 } 0402 double nsl_sf_bessel_jl(double l, double x) { 0403 return gsl_sf_bessel_jl((int)round(l), x); 0404 } 0405 double nsl_sf_bessel_yl(double l, double x) { 0406 return gsl_sf_bessel_yl((int)round(l), x); 0407 } 0408 double nsl_sf_bessel_ils(double l, double x) { 0409 return gsl_sf_bessel_il_scaled((int)round(l), x); 0410 } 0411 double nsl_sf_bessel_kls(double l, double x) { 0412 return gsl_sf_bessel_kl_scaled((int)round(l), x); 0413 } 0414 double nsl_sf_bessel_0_J0(double s) { 0415 return gsl_sf_bessel_zero_J0((unsigned int)round(s)); 0416 } 0417 double nsl_sf_bessel_0_J1(double s) { 0418 return gsl_sf_bessel_zero_J1((unsigned int)round(s)); 0419 } 0420 double nsl_sf_bessel_0_Jnu(double nu, double s) { 0421 return gsl_sf_bessel_zero_Jnu(nu, (unsigned int)round(s)); 0422 } 0423 0424 double nsl_sf_hydrogenicR(double n, double l, double z, double r) { 0425 return gsl_sf_hydrogenicR((int)round(n), (int)round(l), z, r); 0426 } 0427 /* elliptic integrals */ 0428 double nsl_sf_ellint_Kc(double x) { 0429 return gsl_sf_ellint_Kcomp(x, MODE); 0430 } 0431 double nsl_sf_ellint_Ec(double x) { 0432 return gsl_sf_ellint_Ecomp(x, MODE); 0433 } 0434 double nsl_sf_ellint_Pc(double x, double n) { 0435 return gsl_sf_ellint_Pcomp(x, n, MODE); 0436 } 0437 double nsl_sf_ellint_F(double phi, double k) { 0438 return gsl_sf_ellint_F(phi, k, MODE); 0439 } 0440 double nsl_sf_ellint_E(double phi, double k) { 0441 return gsl_sf_ellint_E(phi, k, MODE); 0442 } 0443 double nsl_sf_ellint_P(double phi, double k, double n) { 0444 return gsl_sf_ellint_P(phi, k, n, MODE); 0445 } 0446 double nsl_sf_ellint_D(double phi, double k) { 0447 #if GSL_MAJOR_VERSION >= 2 0448 return gsl_sf_ellint_D(phi, k, MODE); 0449 #else 0450 return gsl_sf_ellint_D(phi, k, 0.0, MODE); 0451 #endif 0452 } 0453 double nsl_sf_ellint_RC(double x, double y) { 0454 return gsl_sf_ellint_RC(x, y, MODE); 0455 } 0456 double nsl_sf_ellint_RD(double x, double y, double z) { 0457 return gsl_sf_ellint_RD(x, y, z, MODE); 0458 } 0459 double nsl_sf_ellint_RF(double x, double y, double z) { 0460 return gsl_sf_ellint_RF(x, y, z, MODE); 0461 } 0462 double nsl_sf_ellint_RJ(double x, double y, double z, double p) { 0463 return gsl_sf_ellint_RJ(x, y, z, p, MODE); 0464 } 0465 0466 double nsl_sf_exprel_n(double n, double x) { 0467 return gsl_sf_exprel_n((int)round(n), x); 0468 } 0469 double nsl_sf_fermi_dirac_int(double j, double x) { 0470 return gsl_sf_fermi_dirac_int((int)round(j), x); 0471 } 0472 /* Gamma */ 0473 double nsl_sf_fact(double n) { 0474 return gsl_sf_fact((unsigned int)round(n)); 0475 } 0476 double nsl_sf_doublefact(double n) { 0477 return gsl_sf_doublefact((unsigned int)round(n)); 0478 } 0479 double nsl_sf_lnfact(double n) { 0480 return gsl_sf_lnfact((unsigned int)round(n)); 0481 } 0482 double nsl_sf_lndoublefact(double n) { 0483 return gsl_sf_lndoublefact((unsigned int)round(n)); 0484 } 0485 double nsl_sf_choose(double n, double m) { 0486 return gsl_sf_choose((unsigned int)round(n), (unsigned int)round(m)); 0487 } 0488 double nsl_sf_lnchoose(double n, double m) { 0489 return gsl_sf_lnchoose((unsigned int)round(n), (unsigned int)round(m)); 0490 } 0491 double nsl_sf_taylorcoeff(double n, double x) { 0492 return gsl_sf_taylorcoeff((int)round(n), x); 0493 } 0494 0495 double nsl_sf_gegenpoly_n(double n, double l, double x) { 0496 return gsl_sf_gegenpoly_n((int)round(n), l, x); 0497 } 0498 0499 /* Hermite polynomials and functions since GSL 2.4 */ 0500 #if (GSL_MAJOR_VERSION > 2) || (GSL_MAJOR_VERSION == 2) && (GSL_MINOR_VERSION >= 4) 0501 double nsl_sf_hermite_prob(double n, double x) { 0502 return gsl_sf_hermite_prob((int)round(n), x); 0503 } 0504 double nsl_sf_hermite_func(double n, double x) { 0505 return gsl_sf_hermite_func((int)round(n), x); 0506 } 0507 double nsl_sf_hermite_func_der(double m, double n, double x) { 0508 return gsl_sf_hermite_func_der((int)round(m), (int)round(n), x); 0509 } 0510 #if (GSL_MAJOR_VERSION == 2) && (GSL_MINOR_VERSION < 6) 0511 double nsl_sf_hermite(double n, double x) { 0512 return gsl_sf_hermite_phys((int)round(n), x); 0513 } 0514 double nsl_sf_hermite_deriv(double m, double n, double x) { 0515 return gsl_sf_hermite_phys_der((int)round(m), (int)round(n), x); 0516 } 0517 double nsl_sf_hermite_prob_deriv(double m, double n, double x) { 0518 return gsl_sf_hermite_prob_der((int)round(m), (int)round(n), x); 0519 } 0520 #else /* renamed and added func_fast in GSL 2.6 */ 0521 double nsl_sf_hermite(double n, double x) { 0522 return gsl_sf_hermite((int)round(n), x); 0523 } 0524 double nsl_sf_hermite_deriv(double m, double n, double x) { 0525 return gsl_sf_hermite_deriv((int)round(m), (int)round(n), x); 0526 } 0527 double nsl_sf_hermite_prob_deriv(double m, double n, double x) { 0528 return gsl_sf_hermite_prob_deriv((int)round(m), (int)round(n), x); 0529 } 0530 double nsl_sf_hermite_func_fast(double n, double x) { 0531 return gsl_sf_hermite_func_fast((int)round(n), x); 0532 } 0533 #endif 0534 #endif 0535 0536 double nsl_sf_hyperg_1F1i(double m, double n, double x) { 0537 return gsl_sf_hyperg_1F1_int((int)round(m), (int)round(n), x); 0538 } 0539 double nsl_sf_hyperg_Ui(double m, double n, double x) { 0540 return gsl_sf_hyperg_U_int((int)round(m), (int)round(n), x); 0541 } 0542 double nsl_sf_laguerre_n(double n, double a, double x) { 0543 return gsl_sf_laguerre_n((int)round(n), a, x); 0544 } 0545 0546 double nsl_sf_legendre_Pl(double l, double x) { 0547 return gsl_sf_legendre_Pl((int)round(l), x); 0548 } 0549 double nsl_sf_legendre_Ql(double l, double x) { 0550 return gsl_sf_legendre_Ql((int)round(l), x); 0551 } 0552 double nsl_sf_legendre_Plm(double l, double m, double x) { 0553 return gsl_sf_legendre_Plm((int)round(l), (int)round(m), x); 0554 } 0555 double nsl_sf_legendre_sphPlm(double l, double m, double x) { 0556 return gsl_sf_legendre_sphPlm((int)round(l), (int)round(m), x); 0557 } 0558 double nsl_sf_conicalP_sphreg(double l, double L, double x) { 0559 return gsl_sf_conicalP_sph_reg((int)round(l), L, x); 0560 } 0561 double nsl_sf_conicalP_cylreg(double m, double l, double x) { 0562 return gsl_sf_conicalP_sph_reg((int)round(m), l, x); 0563 } 0564 double nsl_sf_legendre_H3d(double l, double L, double e) { 0565 return gsl_sf_legendre_H3d((int)round(l), L, e); 0566 } 0567 0568 #if (GSL_MAJOR_VERSION >= 2) 0569 double nsl_sf_mathieu_a(double n, double q) { 0570 return gsl_sf_mathieu_a((int)round(n), q); 0571 } 0572 double nsl_sf_mathieu_b(double n, double q) { 0573 return gsl_sf_mathieu_b((int)round(n), q); 0574 } 0575 double nsl_sf_mathieu_ce(double n, double q, double x) { 0576 return gsl_sf_mathieu_ce((int)round(n), q, x); 0577 } 0578 double nsl_sf_mathieu_se(double n, double q, double x) { 0579 return gsl_sf_mathieu_se((int)round(n), q, x); 0580 } 0581 double nsl_sf_mathieu_Mc(double j, double n, double q, double x) { 0582 return gsl_sf_mathieu_Mc((int)round(j), (int)round(n), q, x); 0583 } 0584 double nsl_sf_mathieu_Ms(double j, double n, double q, double x) { 0585 return gsl_sf_mathieu_Ms((int)round(j), (int)round(n), q, x); 0586 } 0587 #endif 0588 0589 double nsl_sf_psiint(double n) { 0590 return gsl_sf_psi_int((int)round(n)); 0591 } 0592 double nsl_sf_psi1int(double n) { 0593 return gsl_sf_psi_1_int((int)round(n)); 0594 } 0595 double nsl_sf_psin(double n, double x) { 0596 return gsl_sf_psi_n((int)round(n), x); 0597 } 0598 0599 double nsl_sf_zetaint(double n) { 0600 return gsl_sf_zeta_int((int)round(n)); 0601 } 0602 double nsl_sf_zetam1int(double n) { 0603 return gsl_sf_zetam1_int((int)round(n)); 0604 } 0605 double nsl_sf_etaint(double n) { 0606 return gsl_sf_eta_int((int)round(n)); 0607 } 0608 0609 /* random number distributions */ 0610 double nsl_sf_poisson(double k, double m) { 0611 return gsl_ran_poisson_pdf((unsigned int)round(k), m); 0612 } 0613 double nsl_sf_bernoulli(double k, double p) { 0614 return gsl_ran_bernoulli_pdf((unsigned int)round(k), p); 0615 } 0616 double nsl_sf_binomial(double k, double p, double n) { 0617 return gsl_ran_binomial_pdf((unsigned int)round(k), p, (unsigned int)round(n)); 0618 } 0619 double nsl_sf_negative_binomial(double k, double p, double n) { 0620 return gsl_ran_negative_binomial_pdf((unsigned int)round(k), p, n); 0621 } 0622 double nsl_sf_pascal(double k, double p, double n) { 0623 return gsl_ran_pascal_pdf((unsigned int)round(k), p, (unsigned int)round(n)); 0624 } 0625 double nsl_sf_geometric(double k, double p) { 0626 return gsl_ran_geometric_pdf((unsigned int)round(k), p); 0627 } 0628 double nsl_sf_hypergeometric(double k, double n1, double n2, double t) { 0629 return gsl_ran_hypergeometric_pdf((unsigned int)round(k), (unsigned int)round(n1), (unsigned int)round(n2), (unsigned int)round(t)); 0630 } 0631 double nsl_sf_logarithmic(double k, double p) { 0632 return gsl_ran_logarithmic_pdf((unsigned int)round(k), p); 0633 }