166 #include "../../Dynamics/Utils/Faddeeva.hh" 174 #define Inf numeric_limits<double>::infinity() 175 #define NaN numeric_limits<double>::quiet_NaN() 177 typedef complex<double>
cmplx;
180 #define cexp(z) exp(z) 181 #define creal(z) real(z) 182 #define cimag(z) imag(z) 183 #define cpolar(r, t) polar(r, t) 185 #define C(a, b) cmplx(a, b) 187 #define FADDEEVA(name) Faddeeva::name 188 #define FADDEEVA_RE(name) Faddeeva::name 191 #if (__cplusplus < 201103L) && (!defined(HAVE_ISNAN) || !defined(HAVE_ISINF)) 192 static inline bool my_isnan(
double x) {
return x != x; }
193 #define isnan my_isnan 194 static inline bool my_isinf(
double x) {
return 1 / x == 0.; }
195 #define isinf my_isinf 196 #elif (__cplusplus >= 201103L) 198 #define isnan std::isnan 199 #define isinf std::isinf 203 #if defined(_WIN32) || defined(__WIN32__) 204 #define copysign _copysign // of course MS had to be different 205 #elif defined(GNULIB_NAMESPACE) // we are using using gnulib <cmath> 206 #define copysign GNULIB_NAMESPACE::copysign 207 #elif (__cplusplus < 201103L) && !defined(HAVE_COPYSIGN) && \ 208 !defined(__linux__) && !(defined(__APPLE__) && defined(__MACH__)) && \ 210 static inline double my_copysign(
double x,
double y) {
211 return x < 0 != y < 0 ? -x : x;
213 #define copysign my_copysign 222 #ifdef GNULIB_NAMESPACE 223 #define floor GNULIB_NAMESPACE::floor 226 #else // !__cplusplus, i.e. pure C (requires C99 features) 228 #include "Faddeeva.h" 230 #define _GNU_SOURCE // enable GNU libc NAN extension if possible 237 #define FADDEEVA(name) Faddeeva_##name 238 #define FADDEEVA_RE(name) Faddeeva_##name##_re 251 #if !defined(CMPLX) && \ 252 (__GNUC__ > 4 || (__GNUC__ == 4 && __GNUC_MINOR__ >= 7)) && \ 253 !(defined(__ICC) || defined(__INTEL_COMPILER)) 254 #define CMPLX(a, b) __builtin_complex((double)(a), (double)(b)) 258 #define C(a, b) CMPLX(a, b) 259 #define Inf INFINITY // C99 infinity 260 #ifdef NAN // GNU libc extension 263 #define NaN (0. / 0.) // NaN 266 #define C(a, b) ((a) + I * (b)) 267 #define Inf (1. / 0.) 268 #define NaN (0. / 0.) 272 if (r == 0.0 && !isnan(t))
275 return C(r * cos(t), r * sin(t));
278 #endif // !__cplusplus, i.e. pure C (requires C99 features) 285 return FADDEEVA(
w)(C(-cimag(z), creal(z)), relerr);
289 double FADDEEVA_RE(
erf)(
double x) {
290 #if !defined(__cplusplus) 292 #elif (__cplusplus >= 201103L) || defined(HAVE_ERF) 297 return (x >= 0 ? 1.0 : -1.0);
302 return 1.0 - exp(mx2) * FADDEEVA_RE(
erfcx)(x);
306 return exp(mx2) * FADDEEVA_RE(
erfcx)(-x) - 1.0;
312 return x * (1.1283791670955125739 +
313 mx2 * (0.37612638903183752464 +
314 mx2 * (0.11283791670955125739 +
315 mx2 * (0.026866170645131251760 +
316 mx2 * 0.0052239776254421878422))));
322 double x = creal(z), y = cimag(z);
325 return C(FADDEEVA_RE(
erf)(x),
332 y * y > 720 ? (y > 0 ?
Inf : -
Inf)
333 : exp(y * y) * FADDEEVA(
w_im)(y));
335 double mRe_z2 = (y - x) * (x + y);
336 double mIm_z2 = -2 * x * y;
338 return (x >= 0 ? 1.0 : -1.0);
347 else if (fabs(mIm_z2) < 5e-3 && x < 5e-3)
352 return 1.0 - exp(mRe_z2) * (C(cos(mIm_z2), sin(mIm_z2)) *
353 FADDEEVA(
w)(C(-y, x), relerr));
358 else if (fabs(mIm_z2) < 5e-3 && x > -5e-3)
361 return C(
NaN, y == 0 ? 0 :
NaN);
365 (C(cos(mIm_z2), sin(mIm_z2)) * FADDEEVA(
w)(C(y, -x), relerr)) -
372 cmplx mz2 = C(mRe_z2, mIm_z2);
373 return z * (1.1283791670955125739 +
374 mz2 * (0.37612638903183752464 +
375 mz2 * (0.11283791670955125739 +
376 mz2 * (0.026866170645131251760 +
377 mz2 * 0.0052239776254421878422))));
390 double x2 = x * x, y2 = y * y;
391 double expy2 = exp(y2);
394 (1.1283791670955125739 -
395 x2 * (0.37612638903183752464 + 0.75225277806367504925 * y2) +
397 (0.11283791670955125739 +
398 y2 * (0.45135166683820502956 + 0.15045055561273500986 * y2))),
399 expy2 * (FADDEEVA(
w_im)(y) - x2 * y *
400 (1.1283791670955125739 -
401 x2 * (0.56418958354775628695 +
402 0.37612638903183752464 * y2))));
408 cmplx e = FADDEEVA(
erf)(C(-cimag(z), creal(z)), relerr);
409 return C(cimag(e), -creal(e));
413 double FADDEEVA_RE(
erfi)(
double x) {
414 return x * x > 720 ? (x > 0 ?
Inf : -
Inf) : exp(x * x) * FADDEEVA(
w_im)(x);
418 double FADDEEVA_RE(
erfc)(
double x) {
419 #if !defined(__cplusplus) 421 #elif (__cplusplus >= 201103L) || defined(HAVE_ERFC) 425 return (x >= 0 ? 0.0 : 2.0);
426 return x >= 0 ? exp(-x * x) * FADDEEVA_RE(
erfcx)(x)
427 : 2. - exp(-x * x) * FADDEEVA_RE(
erfcx)(-x);
433 double x = creal(z), y = cimag(z);
440 y * y > 720 ? (y > 0 ? -
Inf :
Inf)
441 : -exp(y * y) * FADDEEVA(
w_im)(y));
444 return C(x >= 0 ? 0.0 : 2.0,
446 return C(x >= 0 ? exp(-x * x) * FADDEEVA_RE(
erfcx)(x)
447 : 2. - exp(-x * x) * FADDEEVA_RE(
erfcx)(-x),
451 double mRe_z2 = (y - x) * (x + y);
452 double mIm_z2 = -2 * x * y;
454 return (x >= 0 ? 0.0 : 2.0);
457 return cexp(C(mRe_z2, mIm_z2)) * FADDEEVA(
w)(C(-y, x), relerr);
459 return 2.0 - cexp(C(mRe_z2, mIm_z2)) * FADDEEVA(
w)(C(y, -x), relerr);
464 const double spi2 = 0.8862269254527580136490837416705725913990;
465 return spi2 * FADDEEVA(
w_im)(x);
470 const double spi2 = 0.8862269254527580136490837416705725913990;
471 double x = creal(z), y = cimag(z);
475 return C(spi2 * FADDEEVA(
w_im)(x),
481 y * (1. + y2 * (0.6666666666666666666666666666666666666667 +
482 y2 * 0.26666666666666666666666666666666666667)));
485 spi2 * (y >= 0 ? exp(y2) - FADDEEVA_RE(
erfcx)(y)
486 : FADDEEVA_RE(
erfcx)(-y) - exp(y2)));
489 double mRe_z2 = (y - x) * (x + y);
490 double mIm_z2 = -2 * x * y;
491 cmplx mz2 = C(mRe_z2, mIm_z2);
500 else if (fabs(mIm_z2) < 5e-3)
501 goto taylor_realaxis;
503 cmplx res = cexp(mz2) - FADDEEVA(
w)(z, relerr);
504 return spi2 * C(-cimag(res), creal(res));
509 else if (fabs(mIm_z2) < 5e-3)
510 goto taylor_realaxis;
512 return C(x == 0 ? 0 :
NaN,
NaN);
513 cmplx res = FADDEEVA(
w)(-z, relerr) - cexp(mz2);
514 return spi2 * C(-cimag(res), creal(res));
520 return z * (1. + mz2 * (0.6666666666666666666666666666666666666667 +
521 mz2 * 0.2666666666666666666666666666666666666667));
561 double xy2 = (x * y) * (x * y);
563 (0.5 + y2 * (0.5 + 0.25 * y2 - 0.16666666666666666667 * xy2)) / x,
566 y2 * (-0.66666666666666666667 + 0.13333333333333333333 * xy2 -
567 0.26666666666666666667 * y2)) /
570 return (1. / (-15 + x2 * (90 + x2 * (-60 + 8 * x2)))) *
571 C(x * (33 + x2 * (-28 + 4 * x2) + y2 * (18 - 4 * x2 + 4 * y2)),
572 y * (-15 + x2 * (24 - 4 * x2) + y2 * (4 * x2 - 10 - 4 * y2)));
574 double D = spi2 * FADDEEVA(
w_im)(x);
577 D + y2 * (D + x - 2 * D * x2) +
579 (D * (0.5 - x2 * (2 - 0.66666666666666666667 * x2)) +
580 x * (0.83333333333333333333 - 0.33333333333333333333 * x2)),
582 y2 * 0.66666666666666666667 * (1 - x2 - D * x * (3 - 2 * x2)) +
584 (0.26666666666666666667 -
585 x2 * (0.6 - 0.13333333333333333333 * x2) -
587 (1 - x2 * (1.3333333333333333333 -
588 0.26666666666666666667 * x2)))));
597 static inline double sinc(
double x,
double sinx) {
598 return fabs(x) < 1e-4 ? 1 - (0.1666666666666666666667) * x * x : sinx / x;
603 return x * (1 + (x * x) * (0.1666666666666666666667 +
604 0.00833333333333333333333 * (x * x)));
607 static inline double sqr(
double x) {
return x * x; }
612 7.64405281671221563e-01,
613 3.41424527166548425e-01,
614 8.91072646929412548e-02,
615 1.35887299055460086e-02,
616 1.21085455253437481e-03,
617 6.30452613933449404e-05,
618 1.91805156577114683e-06,
619 3.40969447714832381e-08,
620 3.54175089099469393e-10,
621 2.14965079583260682e-12,
622 7.62368911833724354e-15,
623 1.57982797110681093e-17,
624 1.91294189103582677e-20,
625 1.35344656764205340e-23,
626 5.59535712428588720e-27,
627 1.35164257972401769e-30,
628 1.90784582843501167e-34,
629 1.57351920291442930e-38,
630 7.58312432328032845e-43,
631 2.13536275438697082e-47,
632 3.51352063787195769e-52,
633 3.37800830266396920e-57,
634 1.89769439468301000e-62,
635 6.22929926072668851e-68,
636 1.19481172006938722e-73,
637 1.33908181133005953e-79,
638 8.76924303483223939e-86,
639 3.35555576166254986e-92,
640 7.50264110688173024e-99,
641 9.80192200745410268e-106,
642 7.48265412822268959e-113,
643 3.33770122566809425e-120,
644 8.69934598159861140e-128,
645 1.32486951484088852e-135,
646 1.17898144201315253e-143,
647 6.13039120236180012e-152,
648 1.86258785950822098e-160,
649 3.30668408201432783e-169,
650 3.43017280887946235e-178,
651 2.07915397775808219e-187,
652 7.36384545323984966e-197,
653 1.52394760394085741e-206,
654 1.84281935046532100e-216,
655 1.30209553802992923e-226,
656 5.37588903521080531e-237,
657 1.29689584599763145e-247,
658 1.82813078022866562e-258,
659 1.50576355348684241e-269,
660 7.24692320799294194e-281,
661 2.03797051314726829e-292,
662 3.34880215927873807e-304,
670 return C(FADDEEVA_RE(
erfcx)(cimag(z)),
672 else if (cimag(z) == 0)
673 return C(exp(-
sqr(creal(z))), FADDEEVA(
w_im)(creal(z)));
676 if (relerr <= DBL_EPSILON) {
677 relerr = DBL_EPSILON;
678 a = 0.518321480430085929872;
679 c = 0.329973702884629072537;
680 a2 = 0.268657157075235951582;
682 const double pi = 3.14159265358979323846264338327950288419716939937510582;
685 a = pi / sqrt(-log(relerr * 0.5));
689 const double x = fabs(creal(z));
690 const double y = cimag(z), ya = fabs(y);
694 double sum1 = 0, sum2 = 0, sum3 = 0, sum4 = 0, sum5 = 0;
696 #define USE_CONTINUED_FRACTION 1 // 1 to use continued fraction for large |z| 698 #if USE_CONTINUED_FRACTION 704 && (ya > 0.1 || (x > 8 && ya > 1e-10) || x > 28))) {
716 const double ispi = 0.56418958354775628694807945156;
717 double xs = y < 0 ? -creal(z) : creal(z);
722 double yax = ya / xs;
723 double denom = ispi / (xs + yax * ya);
724 ret = C(denom * yax, denom);
725 }
else if (isinf(ya))
726 return ((isnan(x) || y < 0) ? C(
NaN,
NaN) : C(0, 0));
728 double xya = xs / ya;
729 double denom = ispi / (xya * xs + ya);
730 ret = C(denom, denom * xya);
733 double dr = xs * xs - ya * ya - 0.5, di = 2 * xs * ya;
734 double denom = ispi / (dr * dr + di * di);
735 ret = C(denom * (xs * di - ya * dr), denom * (xs * dr + ya * di));
738 const double c0 = 3.9, c1 = 11.398, c2 = 0.08254, c3 = 0.1421,
740 double nu = floor(c0 + c1 / (c2 * x + c3 * ya + c4));
741 double wr = xs, wi = ya;
742 for (nu = 0.5 * (nu - 1); nu > 0.4; nu -= 0.5) {
744 double denom = nu / (wr * wr + wi * wi);
745 wr = xs - wr * denom;
746 wi = ya + wi * denom;
749 double denom = ispi / (wr * wr + wi * wi);
750 ret = C(denom * wi, denom * wr);
757 return 2.0 * cexp(C((ya - xs) * (xs + ya), 2 * xs * y)) - ret;
761 #else // !USE_CONTINUED_FRACTION 763 const double ispi = 0.56418958354775628694807945156;
764 double xs = y < 0 ? -creal(z) : creal(z);
767 double yax = ya / xs;
768 double denom = ispi / (xs + yax * ya);
769 ret = C(denom * yax, denom);
771 double xya = xs / ya;
772 double denom = ispi / (xya * xs + ya);
773 ret = C(denom, denom * xya);
779 return 2.0 * cexp(C((ya - xs) * (xs + ya), 2 * xs * y)) - ret;
783 #endif // !USE_CONTINUED_FRACTION 798 double prod2ax = 1, prodm2ax = 1;
808 if (relerr == DBL_EPSILON) {
810 const double x2 = x * x;
811 expx2 = 1 - x2 * (1 - 0.5 * x2);
813 const double ax2 = 1.036642960860171859744 * x;
814 const double exp2ax =
815 1 + ax2 * (1 + ax2 * (0.5 + 0.166666666666666666667 * ax2));
816 const double expm2ax =
817 1 - ax2 * (1 - ax2 * (0.5 - 0.166666666666666666667 * ax2));
818 for (
int n = 1; 1; ++n) {
819 const double coef =
expa2n2[n - 1] * expx2 / (a2 * (n * n) + y * y);
823 sum2 += coef * prodm2ax;
824 sum3 += coef * prod2ax;
827 sum5 += coef * (2 * a) * n *
sinh_taylor((2 * a) * n * x);
830 if (coef * prod2ax < relerr * sum3)
835 const double exp2ax = exp((2 * a) * x), expm2ax = 1 / exp2ax;
836 for (
int n = 1; 1; ++n) {
837 const double coef =
expa2n2[n - 1] * expx2 / (a2 * (n * n) + y * y);
841 sum2 += coef * prodm2ax;
842 sum4 += (coef * prodm2ax) * (a * n);
843 sum3 += coef * prod2ax;
844 sum5 += (coef * prod2ax) * (a * n);
846 if ((coef * prod2ax) * (a * n) < relerr * sum5)
851 const double exp2ax = exp((2 * a) * x), expm2ax = 1 / exp2ax;
853 const double x2 = x * x;
854 expx2 = 1 - x2 * (1 - 0.5 * x2);
855 for (
int n = 1; 1; ++n) {
857 exp(-a2 * (n * n)) * expx2 / (a2 * (n * n) + y * y);
861 sum2 += coef * prodm2ax;
862 sum3 += coef * prod2ax;
865 sum5 += coef * (2 * a) * n *
sinh_taylor((2 * a) * n * x);
868 if (coef * prod2ax < relerr * sum3)
873 for (
int n = 1; 1; ++n) {
875 exp(-a2 * (n * n)) * expx2 / (a2 * (n * n) + y * y);
879 sum2 += coef * prodm2ax;
880 sum4 += (coef * prodm2ax) * (a * n);
881 sum3 += coef * prod2ax;
882 sum5 += (coef * prod2ax) * (a * n);
884 if ((coef * prod2ax) * (a * n) < relerr * sum5)
889 const double expx2erfcxy =
891 ? expx2 * FADDEEVA_RE(
erfcx)(y)
892 : 2 * exp(y * y - x * x);
894 const double sinxy = sin(x * y);
895 ret = (expx2erfcxy - c * y * sum1) * cos(2 * x * y) +
896 (c * x * expx2) * sinxy *
sinc(x * y, sinxy);
898 double xs = creal(z);
899 const double sinxy = sin(xs * y);
900 const double sin2xy = sin(2 * xs * y), cos2xy = cos(2 * xs * y);
901 const double coef1 = expx2erfcxy - c * y * sum1;
902 const double coef2 = c * xs * expx2;
903 ret = C(coef1 * cos2xy + coef2 * sinxy *
sinc(xs * y, sinxy),
904 coef2 *
sinc(2 * xs * y, sin2xy) - coef1 * sin2xy);
912 #if USE_CONTINUED_FRACTION 921 ret =
cpolar(2 * exp(y * y - x * x) - exp(-x * x), -2 * creal(z) * y);
926 double n0 = floor(x / a + 0.5);
927 double dx = a * n0 - x;
928 sum3 = exp(-dx * dx) / (a2 * (n0 * n0) + y * y);
929 sum5 = a * n0 * sum3;
930 double exp1 = exp(4 * a * dx), exp1dn = 1;
932 for (dn = 1; n0 - dn > 0; ++dn) {
933 double np = n0 + dn, nm = n0 - dn;
934 double tp = exp(-
sqr(a * dn + dx));
935 double tm = tp * (exp1dn *= exp1);
936 tp /= (a2 * (np * np) + y * y);
937 tm /= (a2 * (nm * nm) + y * y);
939 sum5 += a * (np * tp + nm * tm);
940 if (a * (np * tp + nm * tm) < relerr * sum5)
944 double np = n0 + dn++;
945 double tp = exp(-
sqr(a * dn + dx)) / (a2 * (np * np) + y * y);
948 if (a * np * tp < relerr * sum5)
953 return ret + C((0.5 * c) * y * (sum2 + sum3),
954 (0.5 * c) * copysign(sum5 - sum4, creal(z)));
999 double t = 2 * y100 - 1;
1000 return 0.70878032454106438663e-3 +
1001 (0.71234091047026302958e-3 +
1002 (0.35779077297597742384e-5 +
1003 (0.17403143962587937815e-7 +
1004 (0.81710660047307788845e-10 +
1005 (0.36885022360434957634e-12 + 0.15917038551111111111e-14 * t) *
1013 double t = 2 * y100 - 3;
1014 return 0.21479143208285144230e-2 +
1015 (0.72686402367379996033e-3 +
1016 (0.36843175430938995552e-5 +
1017 (0.18071841272149201685e-7 +
1018 (0.85496449296040325555e-10 +
1019 (0.38852037518534291510e-12 + 0.16868473576888888889e-14 * t) *
1027 double t = 2 * y100 - 5;
1028 return 0.36165255935630175090e-2 +
1029 (0.74182092323555510862e-3 +
1030 (0.37948319957528242260e-5 +
1031 (0.18771627021793087350e-7 +
1032 (0.89484715122415089123e-10 +
1033 (0.40935858517772440862e-12 + 0.17872061464888888889e-14 * t) *
1041 double t = 2 * y100 - 7;
1042 return 0.51154983860031979264e-2 +
1043 (0.75722840734791660540e-3 +
1044 (0.39096425726735703941e-5 +
1045 (0.19504168704300468210e-7 +
1046 (0.93687503063178993915e-10 +
1047 (0.43143925959079664747e-12 + 0.18939926435555555556e-14 * t) *
1055 double t = 2 * y100 - 9;
1056 return 0.66457513172673049824e-2 +
1057 (0.77310406054447454920e-3 +
1058 (0.40289510589399439385e-5 +
1059 (0.20271233238288381092e-7 +
1060 (0.98117631321709100264e-10 +
1061 (0.45484207406017752971e-12 + 0.20076352213333333333e-14 * t) *
1069 double t = 2 * y100 - 11;
1070 return 0.82082389970241207883e-2 +
1071 (0.78946629611881710721e-3 +
1072 (0.41529701552622656574e-5 +
1073 (0.21074693344544655714e-7 +
1074 (0.10278874108587317989e-9 +
1075 (0.47965201390613339638e-12 + 0.21285907413333333333e-14 * t) *
1083 double t = 2 * y100 - 13;
1084 return 0.98039537275352193165e-2 +
1085 (0.80633440108342840956e-3 +
1086 (0.42819241329736982942e-5 +
1087 (0.21916534346907168612e-7 +
1088 (0.10771535136565470914e-9 +
1089 (0.50595972623692822410e-12 + 0.22573462684444444444e-14 * t) *
1097 double t = 2 * y100 - 15;
1098 return 0.11433927298290302370e-1 +
1099 (0.82372858383196561209e-3 +
1100 (0.44160495311765438816e-5 +
1101 (0.22798861426211986056e-7 +
1102 (0.11291291745879239736e-9 +
1103 (0.53386189365816880454e-12 + 0.23944209546666666667e-14 * t) *
1111 double t = 2 * y100 - 17;
1112 return 0.13099232878814653979e-1 +
1113 (0.84167002467906968214e-3 +
1114 (0.45555958988457506002e-5 +
1115 (0.23723907357214175198e-7 +
1116 (0.11839789326602695603e-9 +
1117 (0.56346163067550237877e-12 + 0.25403679644444444444e-14 * t) *
1125 double t = 2 * y100 - 19;
1126 return 0.14800987015587535621e-1 +
1127 (0.86018092946345943214e-3 +
1128 (0.47008265848816866105e-5 +
1129 (0.24694040760197315333e-7 +
1130 (0.12418779768752299093e-9 +
1131 (0.59486890370320261949e-12 + 0.26957764568888888889e-14 * t) *
1139 double t = 2 * y100 - 21;
1140 return 0.16540351739394069380e-1 +
1141 (0.87928458641241463952e-3 +
1142 (0.48520195793001753903e-5 +
1143 (0.25711774900881709176e-7 +
1144 (0.13030128534230822419e-9 +
1145 (0.62820097586874779402e-12 + 0.28612737351111111111e-14 * t) *
1153 double t = 2 * y100 - 23;
1154 return 0.18318536789842392647e-1 +
1155 (0.89900542647891721692e-3 +
1156 (0.50094684089553365810e-5 +
1157 (0.26779777074218070482e-7 +
1158 (0.13675822186304615566e-9 +
1159 (0.66358287745352705725e-12 + 0.30375273884444444444e-14 * t) *
1167 double t = 2 * y100 - 25;
1168 return 0.20136801964214276775e-1 +
1169 (0.91936908737673676012e-3 +
1170 (0.51734830914104276820e-5 +
1171 (0.27900878609710432673e-7 +
1172 (0.14357976402809042257e-9 +
1173 (0.70114790311043728387e-12 + 0.32252476000000000000e-14 * t) *
1181 double t = 2 * y100 - 27;
1182 return 0.21996459598282740954e-1 +
1183 (0.94040248155366777784e-3 +
1184 (0.53443911508041164739e-5 +
1185 (0.29078085538049374673e-7 +
1186 (0.15078844500329731137e-9 +
1187 (0.74103813647499204269e-12 + 0.34251892320000000000e-14 * t) *
1195 double t = 2 * y100 - 29;
1196 return 0.23898877187226319502e-1 +
1197 (0.96213386835900177540e-3 +
1198 (0.55225386998049012752e-5 +
1199 (0.30314589961047687059e-7 +
1200 (0.15840826497296335264e-9 +
1201 (0.78340500472414454395e-12 + 0.36381553564444444445e-14 * t) *
1209 double t = 2 * y100 - 31;
1210 return 0.25845480155298518485e-1 +
1211 (0.98459293067820123389e-3 +
1212 (0.57082915920051843672e-5 +
1213 (0.31613782169164830118e-7 +
1214 (0.16646478745529630813e-9 +
1215 (0.82840985928785407942e-12 + 0.38649975768888888890e-14 * t) *
1223 double t = 2 * y100 - 33;
1224 return 0.27837754783474696598e-1 +
1225 (0.10078108563256892757e-2 +
1226 (0.59020366493792212221e-5 +
1227 (0.32979263553246520417e-7 +
1228 (0.17498524159268458073e-9 +
1229 (0.87622459124842525110e-12 + 0.41066206488888888890e-14 * t) *
1237 double t = 2 * y100 - 35;
1238 return 0.29877251304899307550e-1 +
1239 (0.10318204245057349310e-2 +
1240 (0.61041829697162055093e-5 +
1241 (0.34414860359542720579e-7 +
1242 (0.18399863072934089607e-9 +
1243 (0.92703227366365046533e-12 + 0.43639844053333333334e-14 * t) *
1251 double t = 2 * y100 - 37;
1252 return 0.31965587178596443475e-1 +
1253 (0.10566560976716574401e-2 +
1254 (0.63151633192414586770e-5 +
1255 (0.35924638339521924242e-7 +
1256 (0.19353584758781174038e-9 +
1257 (0.98102783859889264382e-12 + 0.46381060817777777779e-14 * t) *
1265 double t = 2 * y100 - 39;
1266 return 0.34104450552588334840e-1 +
1267 (0.10823541191350532574e-2 +
1268 (0.65354356159553934436e-5 +
1269 (0.37512918348533521149e-7 +
1270 (0.20362979635817883229e-9 +
1271 (0.10384187833037282363e-11 + 0.49300625262222222221e-14 * t) *
1279 double t = 2 * y100 - 41;
1280 return 0.36295603928292425716e-1 +
1281 (0.11089526167995268200e-2 +
1282 (0.67654845095518363577e-5 +
1283 (0.39184292949913591646e-7 +
1284 (0.21431552202133775150e-9 +
1285 (0.10994259106646731797e-11 + 0.52409949102222222221e-14 * t) *
1293 double t = 2 * y100 - 43;
1294 return 0.38540888038840509795e-1 +
1295 (0.11364917134175420009e-2 +
1296 (0.70058230641246312003e-5 +
1297 (0.40943644083718586939e-7 +
1298 (0.22563034723692881631e-9 +
1299 (0.11642841011361992885e-11 + 0.55721092871111111110e-14 * t) *
1307 double t = 2 * y100 - 45;
1308 return 0.40842225954785960651e-1 +
1309 (0.11650136437945673891e-2 +
1310 (0.72569945502343006619e-5 +
1311 (0.42796161861855042273e-7 +
1312 (0.23761401711005024162e-9 +
1313 (0.12332431172381557035e-11 + 0.59246802364444444445e-14 * t) *
1321 double t = 2 * y100 - 47;
1322 return 0.43201627431540222422e-1 +
1323 (0.11945628793917272199e-2 +
1324 (0.75195743532849206263e-5 +
1325 (0.44747364553960993492e-7 +
1326 (0.25030885216472953674e-9 +
1327 (0.13065684400300476484e-11 + 0.63000532853333333334e-14 * t) *
1335 double t = 2 * y100 - 49;
1336 return 0.45621193513810471438e-1 +
1337 (0.12251862608067529503e-2 +
1338 (0.77941720055551920319e-5 +
1339 (0.46803119830954460212e-7 +
1340 (0.26375990983978426273e-9 +
1341 (0.13845421370977119765e-11 + 0.66996477404444444445e-14 * t) *
1349 double t = 2 * y100 - 51;
1350 return 0.48103121413299865517e-1 +
1351 (0.12569331386432195113e-2 +
1352 (0.80814333496367673980e-5 +
1353 (0.48969667335682018324e-7 +
1354 (0.27801515481905748484e-9 +
1355 (0.14674637611609884208e-11 + 0.71249589351111111110e-14 * t) *
1363 double t = 2 * y100 - 53;
1364 return 0.50649709676983338501e-1 +
1365 (0.12898555233099055810e-2 +
1366 (0.83820428414568799654e-5 +
1367 (0.51253642652551838659e-7 +
1368 (0.29312563849675507232e-9 +
1369 (0.15556512782814827846e-11 + 0.75775607822222222221e-14 * t) *
1377 double t = 2 * y100 - 55;
1378 return 0.53263363664388864181e-1 +
1379 (0.13240082443256975769e-2 +
1380 (0.86967260015007658418e-5 +
1381 (0.53662102750396795566e-7 +
1382 (0.30914568786634796807e-9 +
1383 (0.16494420240828493176e-11 + 0.80591079644444444445e-14 * t) *
1391 double t = 2 * y100 - 57;
1392 return 0.55946601353500013794e-1 +
1393 (0.13594491197408190706e-2 +
1394 (0.90262520233016380987e-5 +
1395 (0.56202552975056695376e-7 +
1396 (0.32613310410503135996e-9 +
1397 (0.17491936862246367398e-11 + 0.85713381688888888890e-14 * t) *
1405 double t = 2 * y100 - 59;
1406 return 0.58702059496154081813e-1 +
1407 (0.13962391363223647892e-2 +
1408 (0.93714365487312784270e-5 +
1409 (0.58882975670265286526e-7 +
1410 (0.34414937110591753387e-9 +
1411 (0.18552853109751857859e-11 + 0.91160736711111111110e-14 * t) *
1419 double t = 2 * y100 - 61;
1420 return 0.61532500145144778048e-1 +
1421 (0.14344426411912015247e-2 +
1422 (0.97331446201016809696e-5 +
1423 (0.61711860507347175097e-7 +
1424 (0.36325987418295300221e-9 +
1425 (0.19681183310134518232e-11 + 0.96952238400000000000e-14 * t) *
1433 double t = 2 * y100 - 63;
1434 return 0.64440817576653297993e-1 +
1435 (0.14741275456383131151e-2 +
1436 (0.10112293819576437838e-4 +
1437 (0.64698236605933246196e-7 +
1438 (0.38353412915303665586e-9 +
1439 (0.20881176114385120186e-11 + 0.10310784480000000000e-13 * t) *
1447 double t = 2 * y100 - 65;
1448 return 0.67430045633130393282e-1 +
1449 (0.15153655418916540370e-2 +
1450 (0.10509857606888328667e-4 +
1451 (0.67851706529363332855e-7 +
1452 (0.40504602194811140006e-9 +
1453 (0.22157325110542534469e-11 + 0.10964842115555555556e-13 * t) *
1461 double t = 2 * y100 - 67;
1462 return 0.70503365513338850709e-1 +
1463 (0.15582323336495709827e-2 +
1464 (0.10926868866865231089e-4 +
1465 (0.71182482239613507542e-7 +
1466 (0.42787405890153386710e-9 +
1467 (0.23514379522274416437e-11 + 0.11659571751111111111e-13 * t) *
1475 double t = 2 * y100 - 69;
1476 return 0.73664114037944596353e-1 +
1477 (0.16028078812438820413e-2 +
1478 (0.11364423678778207991e-4 +
1479 (0.74701423097423182009e-7 +
1480 (0.45210162777476488324e-9 +
1481 (0.24957355004088569134e-11 + 0.12397238257777777778e-13 * t) *
1489 double t = 2 * y100 - 71;
1490 return 0.76915792420819562379e-1 +
1491 (0.16491766623447889354e-2 +
1492 (0.11823685320041302169e-4 +
1493 (0.78420075993781544386e-7 +
1494 (0.47781726956916478925e-9 +
1495 (0.26491544403815724749e-11 + 0.13180196462222222222e-13 * t) *
1503 double t = 2 * y100 - 73;
1504 return 0.80262075578094612819e-1 +
1505 (0.16974279491709504117e-2 +
1506 (0.12305888517309891674e-4 +
1507 (0.82350717698979042290e-7 +
1508 (0.50511496109857113929e-9 +
1509 (0.28122528497626897696e-11 + 0.14010889635555555556e-13 * t) *
1517 double t = 2 * y100 - 75;
1518 return 0.83706822008980357446e-1 +
1519 (0.17476561032212656962e-2 +
1520 (0.12812343958540763368e-4 +
1521 (0.86506399515036435592e-7 +
1522 (0.53409440823869467453e-9 +
1523 (0.29856186620887555043e-11 + 0.14891851591111111111e-13 * t) *
1531 double t = 2 * y100 - 77;
1532 return 0.87254084284461718231e-1 +
1533 (0.17999608886001962327e-2 +
1534 (0.13344443080089492218e-4 +
1535 (0.90900994316429008631e-7 +
1536 (0.56486134972616465316e-9 +
1537 (0.31698707080033956934e-11 + 0.15825697795555555556e-13 * t) *
1545 double t = 2 * y100 - 79;
1546 return 0.90908120182172748487e-1 +
1547 (0.18544478050657699758e-2 +
1548 (0.13903663143426120077e-4 +
1549 (0.95549246062549906177e-7 +
1550 (0.59752787125242054315e-9 +
1551 (0.33656597366099099413e-11 + 0.16815130613333333333e-13 * t) *
1559 double t = 2 * y100 - 81;
1560 return 0.94673404508075481121e-1 +
1561 (0.19112284419887303347e-2 +
1562 (0.14491572616545004930e-4 +
1563 (0.10046682186333613697e-6 +
1564 (0.63221272959791000515e-9 +
1565 (0.35736693975589130818e-11 + 0.17862931591111111111e-13 * t) *
1573 double t = 2 * y100 - 83;
1574 return 0.98554641648004456555e-1 +
1575 (0.19704208544725622126e-2 +
1576 (0.15109836875625443935e-4 +
1577 (0.10567036667675984067e-6 +
1578 (0.66904168640019354565e-9 +
1579 (0.37946171850824333014e-11 + 0.18971959040000000000e-13 * t) *
1587 double t = 2 * y100 - 85;
1588 return 0.10255677889470089531e0 +
1589 (0.20321499629472857418e-2 +
1590 (0.15760224242962179564e-4 +
1591 (0.11117756071353507391e-6 +
1592 (0.70814785110097658502e-9 +
1593 (0.40292553276632563925e-11 + 0.20145143075555555556e-13 * t) *
1601 double t = 2 * y100 - 87;
1602 return 0.10668502059865093318e0 +
1603 (0.20965479776148731610e-2 +
1604 (0.16444612377624983565e-4 +
1605 (0.11700717962026152749e-6 +
1606 (0.74967203250938418991e-9 +
1607 (0.42783716186085922176e-11 + 0.21385479360000000000e-13 * t) *
1615 double t = 2 * y100 - 89;
1616 return 0.11094484319386444474e0 +
1617 (0.21637548491908170841e-2 +
1618 (0.17164995035719657111e-4 +
1619 (0.12317915750735938089e-6 +
1620 (0.79376309831499633734e-9 +
1621 (0.45427901763106353914e-11 + 0.22696025653333333333e-13 * t) *
1629 double t = 2 * y100 - 91;
1630 return 0.11534201115268804714e0 +
1631 (0.22339187474546420375e-2 +
1632 (0.17923489217504226813e-4 +
1633 (0.12971465288245997681e-6 +
1634 (0.84057834180389073587e-9 +
1635 (0.48233721206418027227e-11 + 0.24079890062222222222e-13 * t) *
1643 double t = 2 * y100 - 93;
1644 return 0.11988259392684094740e0 +
1645 (0.23071965691918689601e-2 +
1646 (0.18722342718958935446e-4 +
1647 (0.13663611754337957520e-6 +
1648 (0.89028385488493287005e-9 +
1649 (0.51210161569225846701e-11 + 0.25540227111111111111e-13 * t) *
1657 double t = 2 * y100 - 95;
1658 return 0.12457298393509812907e0 +
1659 (0.23837544771809575380e-2 +
1660 (0.19563942105711612475e-4 +
1661 (0.14396736847739470782e-6 +
1662 (0.94305490646459247016e-9 +
1663 (0.54366590583134218096e-11 + 0.27080225920000000000e-13 * t) *
1671 double t = 2 * y100 - 97;
1672 return 0.12941991566142438816e0 +
1673 (0.24637684719508859484e-2 +
1674 (0.20450821127475879816e-4 +
1675 (0.15173366280523906622e-6 +
1676 (0.99907632506389027739e-9 +
1677 (0.57712760311351625221e-11 + 0.28703099555555555556e-13 * t) *
1685 double t = 2 * y100 - 99;
1686 return 0.13443048593088696613e0 +
1687 (0.25474249981080823877e-2 +
1688 (0.21385669591362915223e-4 +
1689 (0.15996177579900443030e-6 +
1690 (0.10585428844575134013e-8 +
1691 (0.61258809536787882989e-11 + 0.30412080142222222222e-13 * t) *
1699 double t = 2 * y100 - 101;
1700 return 0.13961217543434561353e0 +
1701 (0.26349215871051761416e-2 +
1702 (0.22371342712572567744e-4 +
1703 (0.16868008199296822247e-6 +
1704 (0.11216596910444996246e-8 +
1705 (0.65015264753090890662e-11 + 0.32210394506666666666e-13 * t) *
1713 double t = 2 * y100 - 103;
1714 return 0.14497287157673800690e0 +
1715 (0.27264675383982439814e-2 +
1716 (0.23410870961050950197e-4 +
1717 (0.17791863939526376477e-6 +
1718 (0.11886425714330958106e-8 +
1719 (0.68993039665054288034e-11 + 0.34101266222222222221e-13 * t) *
1727 double t = 2 * y100 - 105;
1728 return 0.15052089272774618151e0 +
1729 (0.28222846410136238008e-2 +
1730 (0.24507470422713397006e-4 +
1731 (0.18770927679626136909e-6 +
1732 (0.12597184587583370712e-8 +
1733 (0.73203433049229821618e-11 + 0.36087889048888888890e-13 * t) *
1741 double t = 2 * y100 - 107;
1742 return 0.15626501395774612325e0 +
1743 (0.29226079376196624949e-2 +
1744 (0.25664553693768450545e-4 +
1745 (0.19808568415654461964e-6 +
1746 (0.13351257759815557897e-8 +
1747 (0.77658124891046760667e-11 + 0.38173420035555555555e-13 * t) *
1755 double t = 2 * y100 - 109;
1756 return 0.16221449434620737567e0 +
1757 (0.30276865332726475672e-2 +
1758 (0.26885741326534564336e-4 +
1759 (0.20908350604346384143e-6 +
1760 (0.14151148144240728728e-8 +
1761 (0.82369170665974313027e-11 + 0.40360957457777777779e-13 * t) *
1769 double t = 2 * y100 - 111;
1770 return 0.16837910595412130659e0 +
1771 (0.31377844510793082301e-2 +
1772 (0.28174873844911175026e-4 +
1773 (0.22074043807045782387e-6 +
1774 (0.14999481055996090039e-8 +
1775 (0.87348993661930809254e-11 + 0.42653528977777777779e-13 * t) *
1783 double t = 2 * y100 - 113;
1784 return 0.17476916455659369953e0 +
1785 (0.32531815370903068316e-2 +
1786 (0.29536024347344364074e-4 +
1787 (0.23309632627767074202e-6 +
1788 (0.15899007843582444846e-8 +
1789 (0.92610375235427359475e-11 + 0.45054073102222222221e-13 * t) *
1797 double t = 2 * y100 - 115;
1798 return 0.18139556223643701364e0 +
1799 (0.33741744168096996041e-2 +
1800 (0.30973511714709500836e-4 +
1801 (0.24619326937592290996e-6 +
1802 (0.16852609412267750744e-8 +
1803 (0.98166442942854895573e-11 + 0.47565418097777777779e-13 * t) *
1811 double t = 2 * y100 - 117;
1812 return 0.18826980194443664549e0 +
1813 (0.35010775057740317997e-2 +
1814 (0.32491914440014267480e-4 +
1815 (0.26007572375886319028e-6 +
1816 (0.17863299617388376116e-8 +
1817 (0.10403065638343878679e-10 + 0.50190265831111111110e-13 * t) *
1825 double t = 2 * y100 - 119;
1826 return 0.19540403413693967350e0 +
1827 (0.36342240767211326315e-2 +
1828 (0.34096085096200907289e-4 +
1829 (0.27479061117017637474e-6 +
1830 (0.18934228504790032826e-8 +
1831 (0.11021679075323598664e-10 + 0.52931171733333333334e-13 * t) *
1839 double t = 2 * y100 - 121;
1840 return 0.20281109560651886959e0 +
1841 (0.37739673859323597060e-2 +
1842 (0.35791165457592409054e-4 +
1843 (0.29038742889416172404e-6 +
1844 (0.20068685374849001770e-8 +
1845 (0.11673891799578381999e-10 + 0.55790523093333333334e-13 * t) *
1853 double t = 2 * y100 - 123;
1854 return 0.21050455062669334978e0 +
1855 (0.39206818613925652425e-2 +
1856 (0.37582602289680101704e-4 +
1857 (0.30691836231886877385e-6 +
1858 (0.21270101645763677824e-8 +
1859 (0.12361138551062899455e-10 + 0.58770520160000000000e-13 * t) *
1867 double t = 2 * y100 - 125;
1868 return 0.21849873453703332479e0 +
1869 (0.40747643554689586041e-2 +
1870 (0.39476163820986711501e-4 +
1871 (0.32443839970139918836e-6 +
1872 (0.22542053491518680200e-8 +
1873 (0.13084879235290858490e-10 + 0.61873153262222222221e-13 * t) *
1881 double t = 2 * y100 - 127;
1882 return 0.22680879990043229327e0 +
1883 (0.42366354648628516935e-2 +
1884 (0.41477956909656896779e-4 +
1885 (0.34300544894502810002e-6 +
1886 (0.23888264229264067658e-8 +
1887 (0.13846596292818514601e-10 + 0.65100183751111111110e-13 * t) *
1895 double t = 2 * y100 - 129;
1896 return 0.23545076536988703937e0 +
1897 (0.44067409206365170888e-2 +
1898 (0.43594444916224700881e-4 +
1899 (0.36268045617760415178e-6 +
1900 (0.25312606430853202748e-8 +
1901 (0.14647791812837903061e-10 + 0.68453122631111111110e-13 * t) *
1909 double t = 2 * y100 - 131;
1910 return 0.24444156740777432838e0 +
1911 (0.45855530511605787178e-2 +
1912 (0.45832466292683085475e-4 +
1913 (0.38352752590033030472e-6 +
1914 (0.26819103733055603460e-8 +
1915 (0.15489984390884756993e-10 + 0.71933206364444444445e-13 * t) *
1923 double t = 2 * y100 - 133;
1924 return 0.25379911500634264643e0 +
1925 (0.47735723208650032167e-2 +
1926 (0.48199253896534185372e-4 +
1927 (0.40561404245564732314e-6 +
1928 (0.28411932320871165585e-8 +
1929 (0.16374705736458320149e-10 + 0.75541379822222222221e-13 * t) *
1937 double t = 2 * y100 - 135;
1938 return 0.26354234756393613032e0 +
1939 (0.49713289477083781266e-2 +
1940 (0.50702455036930367504e-4 +
1941 (0.42901079254268185722e-6 +
1942 (0.30095422058900481753e-8 +
1943 (0.17303497025347342498e-10 + 0.79278273368888888890e-13 * t) *
1951 double t = 2 * y100 - 137;
1952 return 0.27369129607732343398e0 +
1953 (0.51793846023052643767e-2 +
1954 (0.53350152258326602629e-4 +
1955 (0.45379208848865015485e-6 +
1956 (0.31874057245814381257e-8 +
1957 (0.18277905010245111046e-10 + 0.83144182364444444445e-13 * t) *
1965 double t = 2 * y100 - 139;
1966 return 0.28426714781640316172e0 +
1967 (0.53983341916695141966e-2 +
1968 (0.56150884865255810638e-4 +
1969 (0.48003589196494734238e-6 +
1970 (0.33752476967570796349e-8 +
1971 (0.19299477888083469086e-10 + 0.87139049137777777779e-13 * t) *
1979 double t = 2 * y100 - 141;
1980 return 0.29529231465348519920e0 +
1981 (0.56288077305420795663e-2 +
1982 (0.59113671189913307427e-4 +
1983 (0.50782393781744840482e-6 +
1984 (0.35735475025851713168e-8 +
1985 (0.20369760937017070382e-10 + 0.91262442613333333334e-13 * t) *
1993 double t = 2 * y100 - 143;
1994 return 0.30679050522528838613e0 +
1995 (0.58714723032745403331e-2 +
1996 (0.62248031602197686791e-4 +
1997 (0.53724185766200945789e-6 +
1998 (0.37827999418960232678e-8 +
1999 (0.21490291930444538307e-10 + 0.95513539182222222221e-13 * t) *
2007 double t = 2 * y100 - 145;
2008 return 0.31878680111173319425e0 +
2009 (0.61270341192339103514e-2 +
2010 (0.65564012259707640976e-4 +
2011 (0.56837930287837738996e-6 +
2012 (0.40035151353392378882e-8 +
2013 (0.22662596341239294792e-10 + 0.99891109760000000000e-13 * t) *
2021 double t = 2 * y100 - 147;
2022 return 0.33130773722152622027e0 +
2023 (0.63962406646798080903e-2 +
2024 (0.69072209592942396666e-4 +
2025 (0.60133006661885941812e-6 +
2026 (0.42362183765883466691e-8 +
2027 (0.23888182347073698382e-10 + 0.10439349811555555556e-12 * t) *
2035 double t = 2 * y100 - 149;
2036 return 0.34438138658041336523e0 +
2037 (0.66798829540414007258e-2 +
2038 (0.72783795518603561144e-4 +
2039 (0.63619220443228800680e-6 +
2040 (0.44814499336514453364e-8 +
2041 (0.25168535651285475274e-10 + 0.10901861383111111111e-12 * t) *
2049 double t = 2 * y100 - 151;
2050 return 0.35803744972380175583e0 +
2051 (0.69787978834882685031e-2 +
2052 (0.76710543371454822497e-4 +
2053 (0.67306815308917386747e-6 +
2054 (0.47397647975845228205e-8 +
2055 (0.26505114141143050509e-10 + 0.11376390933333333333e-12 * t) *
2063 double t = 2 * y100 - 153;
2064 return 0.37230734890119724188e0 +
2065 (0.72938706896461381003e-2 +
2066 (0.80864854542670714092e-4 +
2067 (0.71206484718062688779e-6 +
2068 (0.50117323769745883805e-8 +
2069 (0.27899342394100074165e-10 + 0.11862637614222222222e-12 * t) *
2077 double t = 2 * y100 - 155;
2078 return 0.38722432730555448223e0 +
2079 (0.76260375162549802745e-2 +
2080 (0.85259785810004603848e-4 +
2081 (0.75329383305171327677e-6 +
2082 (0.52979361368388119355e-8 +
2083 (0.29352606054164086709e-10 + 0.12360253370666666667e-12 * t) *
2091 double t = 2 * y100 - 157;
2092 return 0.40282355354616940667e0 +
2093 (0.79762880915029728079e-2 +
2094 (0.89909077342438246452e-4 +
2095 (0.79687137961956194579e-6 +
2096 (0.55989731807360403195e-8 +
2097 (0.30866246101464869050e-10 + 0.12868841946666666667e-12 * t) *
2105 double t = 2 * y100 - 159;
2106 return 0.41914223158913787649e0 +
2107 (0.83456685186950463538e-2 +
2108 (0.94827181359250161335e-4 +
2109 (0.84291858561783141014e-6 +
2110 (0.59154537751083485684e-8 +
2111 (0.32441553034347469291e-10 + 0.13387957943111111111e-12 * t) *
2119 double t = 2 * y100 - 161;
2120 return 0.43621971639463786896e0 +
2121 (0.87352841828289495773e-2 +
2122 (0.10002929142066799966e-3 +
2123 (0.89156148280219880024e-6 +
2124 (0.62480008150788597147e-8 +
2125 (0.34079760983458878910e-10 + 0.13917107176888888889e-12 * t) *
2133 double t = 2 * y100 - 163;
2134 return 0.45409763548534330981e0 +
2135 (0.91463027755548240654e-2 +
2136 (0.10553137232446167258e-3 +
2137 (0.94293113464638623798e-6 +
2138 (0.65972492312219959885e-8 +
2139 (0.35782041795476563662e-10 + 0.14455745872000000000e-12 * t) *
2147 double t = 2 * y100 - 165;
2148 return 0.47282001668512331468e0 +
2149 (0.95799574408860463394e-2 +
2150 (0.11135019058000067469e-3 +
2151 (0.99716373005509038080e-6 +
2152 (0.69638453369956970347e-8 +
2153 (0.37549499088161345850e-10 + 0.15003280712888888889e-12 * t) *
2161 double t = 2 * y100 - 167;
2162 return 0.49243342227179841649e0 +
2163 (0.10037550043909497071e-1 +
2164 (0.11750334542845234952e-3 +
2165 (0.10544006716188967172e-5 +
2166 (0.73484461168242224872e-8 +
2167 (0.39383162326435752965e-10 + 0.15559069118222222222e-12 * t) *
2175 double t = 2 * y100 - 169;
2176 return 0.51298708979209258326e0 +
2177 (0.10520454564612427224e-1 +
2178 (0.12400930037494996655e-3 +
2179 (0.11147886579371265246e-5 +
2180 (0.77517184550568711454e-8 +
2181 (0.41283980931872622611e-10 + 0.16122419680000000000e-12 * t) *
2189 double t = 2 * y100 - 171;
2190 return 0.53453307979101369843e0 +
2191 (0.11030120618800726938e-1 +
2192 (0.13088741519572269581e-3 +
2193 (0.11784797595374515432e-5 +
2194 (0.81743383063044825400e-8 +
2195 (0.43252818449517081051e-10 + 0.16692592640000000000e-12 * t) *
2203 double t = 2 * y100 - 173;
2204 return 0.55712643071169299478e0 +
2205 (0.11568077107929735233e-1 +
2206 (0.13815797838036651289e-3 +
2207 (0.12456314879260904558e-5 +
2208 (0.86169898078969313597e-8 +
2209 (0.45290446811539652525e-10 + 0.17268801084444444444e-12 * t) *
2217 double t = 2 * y100 - 175;
2218 return 0.58082532122519320968e0 +
2219 (0.12135935999503877077e-1 +
2220 (0.14584223996665838559e-3 +
2221 (0.13164068573095710742e-5 +
2222 (0.90803643355106020163e-8 +
2223 (0.47397540713124619155e-10 + 0.17850211608888888889e-12 * t) *
2231 double t = 2 * y100 - 177;
2232 return 0.60569124025293375554e0 +
2233 (0.12735396239525550361e-1 +
2234 (0.15396244472258863344e-3 +
2235 (0.13909744385382818253e-5 +
2236 (0.95651595032306228245e-8 +
2237 (0.49574672127669041550e-10 + 0.18435945564444444444e-12 * t) *
2245 double t = 2 * y100 - 179;
2246 return 0.63178916494715716894e0 +
2247 (0.13368247798287030927e-1 +
2248 (0.16254186562762076141e-3 +
2249 (0.14695084048334056083e-5 +
2250 (0.10072078109604152350e-7 +
2251 (0.51822304995680707483e-10 + 0.19025081422222222222e-12 * t) *
2259 double t = 2 * y100 - 181;
2260 return 0.65918774689725319200e0 +
2261 (0.14036375850601992063e-1 +
2262 (0.17160483760259706354e-3 +
2263 (0.15521885688723188371e-5 +
2264 (0.10601827031535280590e-7 +
2265 (0.54140790105837520499e-10 + 0.19616655146666666667e-12 * t) *
2273 double t = 2 * y100 - 183;
2274 return 0.68795950683174433822e0 +
2275 (0.14741765091365869084e-1 +
2276 (0.18117679143520433835e-3 +
2277 (0.16392004108230585213e-5 +
2278 (0.11155116068018043001e-7 +
2279 (0.56530360194925690374e-10 + 0.20209663662222222222e-12 * t) *
2287 double t = 2 * y100 - 185;
2288 return 0.71818103808729967036e0 +
2289 (0.15486504187117112279e-1 +
2290 (0.19128428784550923217e-3 +
2291 (0.17307350969359975848e-5 +
2292 (0.11732656736113607751e-7 +
2293 (0.58991125287563833603e-10 + 0.20803065333333333333e-12 * t) *
2301 double t = 2 * y100 - 187;
2302 return 0.74993321911726254661e0 +
2303 (0.16272790364044783382e-1 +
2304 (0.20195505163377912645e-3 +
2305 (0.18269894883203346953e-5 +
2306 (0.12335161021630225535e-7 +
2307 (0.61523068312169087227e-10 + 0.21395783431111111111e-12 * t) *
2315 double t = 2 * y100 - 189;
2316 return 0.78330143531283492729e0 +
2317 (0.17102934132652429240e-1 +
2318 (0.21321800585063327041e-3 +
2319 (0.19281661395543913713e-5 +
2320 (0.12963340087354341574e-7 +
2321 (0.64126040998066348872e-10 + 0.21986708942222222222e-12 * t) *
2329 double t = 2 * y100 - 191;
2330 return 0.81837581041023811832e0 +
2331 (0.17979364149044223802e-1 +
2332 (0.22510330592753129006e-3 +
2333 (0.20344732868018175389e-5 +
2334 (0.13617902941839949718e-7 +
2335 (0.66799760083972474642e-10 + 0.22574701262222222222e-12 * t) *
2343 double t = 2 * y100 - 193;
2344 return 0.85525144775685126237e0 +
2345 (0.18904632212547561026e-1 +
2346 (0.23764237370371255638e-3 +
2347 (0.21461248251306387979e-5 +
2348 (0.14299555071870523786e-7 +
2349 (0.69543803864694171934e-10 + 0.23158593688888888889e-12 * t) *
2357 double t = 2 * y100 - 195;
2358 return 0.89402868170849933734e0 +
2359 (0.19881418399127202569e-1 +
2360 (0.25086793128395995798e-3 +
2361 (0.22633402747585233180e-5 +
2362 (0.15008997042116532283e-7 +
2363 (0.72357609075043941261e-10 + 0.23737194737777777778e-12 * t) *
2371 double t = 2 * y100 - 197;
2372 return 0.93481333942870796363e0 +
2373 (0.20912536329780368893e-1 +
2374 (0.26481403465998477969e-3 +
2375 (0.23863447359754921676e-5 +
2376 (0.15746923065472184451e-7 +
2377 (0.75240468141720143653e-10 + 0.24309291271111111111e-12 * t) *
2385 double t = 2 * y100 - 199;
2386 return 0.97771701335885035464e0 +
2387 (0.22000938572830479551e-1 +
2388 (0.27951610702682383001e-3 +
2389 (0.25153688325245314530e-5 +
2390 (0.16514019547822821453e-7 +
2391 (0.78191526829368231251e-10 + 0.24873652355555555556e-12 * t) *
2407 const double ispi = 0.56418958354775628694807945156;
2412 return ispi * ((x * x) * (x * x + 4.5) + 2) /
2413 (x * ((x * x) * (x * x + 5) + 3.75));
2417 return x < -26.7 ? HUGE_VAL
2418 : (x < -6.1 ? 2 * exp(x * x)
2419 : 2 * exp(x * x) -
erfcx_y100(400 / (4 - x)));
2441 switch ((
int)y100) {
2443 double t = 2 * y100 - 1;
2444 return 0.28351593328822191546e-2 +
2445 (0.28494783221378400759e-2 +
2446 (0.14427470563276734183e-4 +
2447 (0.10939723080231588129e-6 +
2448 (0.92474307943275042045e-9 +
2449 (0.89128907666450075245e-11 + 0.92974121935111111110e-13 * t) *
2457 double t = 2 * y100 - 3;
2458 return 0.85927161243940350562e-2 +
2459 (0.29085312941641339862e-2 +
2460 (0.15106783707725582090e-4 +
2461 (0.11716709978531327367e-6 +
2462 (0.10197387816021040024e-8 +
2463 (0.10122678863073360769e-10 + 0.10917479678400000000e-12 * t) *
2471 double t = 2 * y100 - 5;
2472 return 0.14471159831187703054e-1 +
2473 (0.29703978970263836210e-2 +
2474 (0.15835096760173030976e-4 +
2475 (0.12574803383199211596e-6 +
2476 (0.11278672159518415848e-8 +
2477 (0.11547462300333495797e-10 + 0.12894535335111111111e-12 * t) *
2485 double t = 2 * y100 - 7;
2486 return 0.20476320420324610618e-1 +
2487 (0.30352843012898665856e-2 +
2488 (0.16617609387003727409e-4 +
2489 (0.13525429711163116103e-6 +
2490 (0.12515095552507169013e-8 +
2491 (0.13235687543603382345e-10 + 0.15326595042666666667e-12 * t) *
2499 double t = 2 * y100 - 9;
2500 return 0.26614461952489004566e-1 +
2501 (0.31034189276234947088e-2 +
2502 (0.17460268109986214274e-4 +
2503 (0.14582130824485709573e-6 +
2504 (0.13935959083809746345e-8 +
2505 (0.15249438072998932900e-10 + 0.18344741882133333333e-12 * t) *
2513 double t = 2 * y100 - 11;
2514 return 0.32892330248093586215e-1 +
2515 (0.31750557067975068584e-2 +
2516 (0.18369907582308672632e-4 +
2517 (0.15761063702089457882e-6 +
2518 (0.15577638230480894382e-8 +
2519 (0.17663868462699097951e-10 +
2520 (0.22126732680711111111e-12 + 0.30273474177737853668e-14 * t) *
2529 double t = 2 * y100 - 13;
2530 return 0.39317207681134336024e-1 +
2531 (0.32504779701937539333e-2 +
2532 (0.19354426046513400534e-4 +
2533 (0.17081646971321290539e-6 +
2534 (0.17485733959327106250e-8 +
2535 (0.20593687304921961410e-10 +
2536 (0.26917401949155555556e-12 + 0.38562123837725712270e-14 * t) *
2545 double t = 2 * y100 - 15;
2546 return 0.45896976511367738235e-1 +
2547 (0.33300031273110976165e-2 +
2548 (0.20423005398039037313e-4 +
2549 (0.18567412470376467303e-6 +
2550 (0.19718038363586588213e-8 +
2551 (0.24175006536781219807e-10 +
2552 (0.33059982791466666666e-12 + 0.49756574284439426165e-14 * t) *
2561 double t = 2 * y100 - 17;
2562 return 0.52640192524848962855e-1 +
2563 (0.34139883358846720806e-2 +
2564 (0.21586390240603337337e-4 +
2565 (0.20247136501568904646e-6 +
2566 (0.22348696948197102935e-8 +
2567 (0.28597516301950162548e-10 +
2568 (0.41045502119111111110e-12 + 0.65151614515238361946e-14 * t) *
2577 double t = 2 * y100 - 19;
2578 return 0.59556171228656770456e-1 +
2579 (0.35028374386648914444e-2 +
2580 (0.22857246150998562824e-4 +
2581 (0.22156372146525190679e-6 +
2582 (0.25474171590893813583e-8 +
2583 (0.34122390890697400584e-10 +
2584 (0.51593189879111111110e-12 + 0.86775076853908006938e-14 * t) *
2593 double t = 2 * y100 - 21;
2594 return 0.66655089485108212551e-1 +
2595 (0.35970095381271285568e-2 +
2596 (0.24250626164318672928e-4 +
2597 (0.24339561521785040536e-6 +
2598 (0.29221990406518411415e-8 +
2599 (0.41117013527967776467e-10 +
2600 (0.65786450716444444445e-12 + 0.11791885745450623331e-13 * t) *
2609 double t = 2 * y100 - 23;
2610 return 0.73948106345519174661e-1 +
2611 (0.36970297216569341748e-2 +
2612 (0.25784588137312868792e-4 +
2613 (0.26853012002366752770e-6 +
2614 (0.33763958861206729592e-8 +
2615 (0.50111549981376976397e-10 +
2616 (0.85313857496888888890e-12 + 0.16417079927706899860e-13 * t) *
2625 double t = 2 * y100 - 25;
2626 return 0.81447508065002963203e-1 +
2627 (0.38035026606492705117e-2 +
2628 (0.27481027572231851896e-4 +
2629 (0.29769200731832331364e-6 +
2630 (0.39336816287457655076e-8 +
2631 (0.61895471132038157624e-10 +
2632 (0.11292303213511111111e-11 + 0.23558532213703884304e-13 * t) *
2641 double t = 2 * y100 - 27;
2642 return 0.89166884027582716628e-1 +
2643 (0.39171301322438946014e-2 +
2644 (0.29366827260422311668e-4 +
2645 (0.33183204390350724895e-6 +
2646 (0.46276006281647330524e-8 +
2647 (0.77692631378169813324e-10 +
2648 (0.15335153258844444444e-11 + 0.35183103415916026911e-13 * t) *
2657 double t = 2 * y100 - 29;
2658 return 0.97121342888032322019e-1 +
2659 (0.40387340353207909514e-2 +
2660 (0.31475490395950776930e-4 +
2661 (0.37222714227125135042e-6 +
2662 (0.55074373178613809996e-8 +
2663 (0.99509175283990337944e-10 +
2664 (0.21552645758222222222e-11 + 0.55728651431872687605e-13 * t) *
2673 double t = 2 * y100 - 31;
2674 return 0.10532778218603311137e0 +
2675 (0.41692873614065380607e-2 +
2676 (0.33849549774889456984e-4 +
2677 (0.42064596193692630143e-6 +
2678 (0.66494579697622432987e-8 +
2679 (0.13094103581931802337e-9 +
2680 (0.31896187409777777778e-11 + 0.97271974184476560742e-13 * t) *
2689 double t = 2 * y100 - 33;
2690 return 0.11380523107427108222e0 +
2691 (0.43099572287871821013e-2 +
2692 (0.36544324341565929930e-4 +
2693 (0.47965044028581857764e-6 +
2694 (0.81819034238463698796e-8 +
2695 (0.17934133239549647357e-9 +
2696 (0.50956666166186293627e-11 +
2697 (0.18850487318190638010e-12 + 0.79697813173519853340e-14 * t) *
2707 double t = 2 * y100 - 35;
2708 return 0.12257529703447467345e0 +
2709 (0.44621675710026986366e-2 +
2710 (0.39634304721292440285e-4 +
2711 (0.55321553769873381819e-6 +
2712 (0.10343619428848520870e-7 +
2713 (0.26033830170470368088e-9 +
2714 (0.87743837749108025357e-11 +
2715 (0.34427092430230063401e-12 + 0.10205506615709843189e-13 * t) *
2725 double t = 2 * y100 - 37;
2726 return 0.13166276955656699478e0 +
2727 (0.46276970481783001803e-2 +
2728 (0.43225026380496399310e-4 +
2729 (0.64799164020016902656e-6 +
2730 (0.13580082794704641782e-7 +
2731 (0.39839800853954313927e-9 +
2732 (0.14431142411840000000e-10 + 0.42193457308830027541e-12 * t) *
2741 double t = 2 * y100 - 39;
2742 return 0.14109647869803356475e0 +
2743 (0.48088424418545347758e-2 +
2744 (0.47474504753352150205e-4 +
2745 (0.77509866468724360352e-6 +
2746 (0.18536851570794291724e-7 +
2747 (0.60146623257887570439e-9 +
2748 (0.18533978397305276318e-10 +
2749 (0.41033845938901048380e-13 - 0.46160680279304825485e-13 * t) *
2759 double t = 2 * y100 - 41;
2760 return 0.15091057940548936603e0 +
2761 (0.50086864672004685703e-2 +
2762 (0.52622482832192230762e-4 +
2763 (0.95034664722040355212e-6 +
2764 (0.25614261331144718769e-7 +
2765 (0.80183196716888606252e-9 +
2766 (0.12282524750534352272e-10 + (-0.10531774117332273617e-11 -
2767 0.86157181395039646412e-13 * t) *
2777 double t = 2 * y100 - 43;
2778 return 0.16114648116017010770e0 +
2779 (0.52314661581655369795e-2 +
2780 (0.59005534545908331315e-4 +
2781 (0.11885518333915387760e-5 +
2782 (0.33975801443239949256e-7 +
2783 (0.82111547144080388610e-9 + (-0.12357674017312854138e-10 +
2784 (-0.24355112256914479176e-11 -
2785 0.75155506863572930844e-13 * t) *
2795 double t = 2 * y100 - 45;
2796 return 0.17185551279680451144e0 +
2797 (0.54829002967599420860e-2 +
2798 (0.67013226658738082118e-4 +
2799 (0.14897400671425088807e-5 +
2800 (0.40690283917126153701e-7 +
2801 (0.44060872913473778318e-9 +
2802 (-0.52641873433280000000e-10 - 0.30940587864543343124e-11 * t) *
2811 double t = 2 * y100 - 47;
2812 return 0.18310194559815257381e0 +
2813 (0.57701559375966953174e-2 +
2814 (0.76948789401735193483e-4 +
2815 (0.18227569842290822512e-5 +
2816 (0.41092208344387212276e-7 +
2817 (-0.44009499965694442143e-9 + (-0.92195414685628803451e-10 +
2818 (-0.22657389705721753299e-11 +
2819 0.10004784908106839254e-12 * t) *
2829 double t = 2 * y100 - 49;
2830 return 0.19496527191546630345e0 +
2831 (0.61010853144364724856e-2 +
2832 (0.88812881056342004864e-4 +
2833 (0.21180686746360261031e-5 +
2834 (0.30652145555130049203e-7 +
2835 (-0.16841328574105890409e-8 +
2836 (-0.11008129460612823934e-9 + (-0.12180794204544515779e-12 +
2837 0.15703325634590334097e-12 * t) *
2847 double t = 2 * y100 - 51;
2848 return 0.20754006813966575720e0 +
2849 (0.64825787724922073908e-2 +
2850 (0.10209599627522311893e-3 +
2851 (0.22785233392557600468e-5 +
2852 (0.73495224449907568402e-8 +
2853 (-0.29442705974150112783e-8 +
2854 (-0.94082603434315016546e-10 +
2855 (0.23609990400179321267e-11 + 0.14141908654269023788e-12 * t) *
2865 double t = 2 * y100 - 53;
2866 return 0.22093185554845172146e0 +
2867 (0.69182878150187964499e-2 +
2868 (0.11568723331156335712e-3 +
2869 (0.22060577946323627739e-5 +
2870 (-0.26929730679360840096e-7 +
2871 (-0.38176506152362058013e-8 +
2872 (-0.47399503861054459243e-10 +
2873 (0.40953700187172127264e-11 + 0.69157730376118511127e-13 * t) *
2883 double t = 2 * y100 - 55;
2884 return 0.23524827304057813918e0 +
2885 (0.74063350762008734520e-2 +
2886 (0.12796333874615790348e-3 +
2887 (0.18327267316171054273e-5 +
2888 (-0.66742910737957100098e-7 +
2889 (-0.40204740975496797870e-8 +
2890 (0.14515984139495745330e-10 +
2891 (0.44921608954536047975e-11 - 0.18583341338983776219e-13 * t) *
2901 double t = 2 * y100 - 57;
2902 return 0.25058626331812744775e0 +
2903 (0.79377285151602061328e-2 +
2904 (0.13704268650417478346e-3 +
2905 (0.11427511739544695861e-5 +
2906 (-0.10485442447768377485e-6 +
2907 (-0.34850364756499369763e-8 +
2908 (0.72656453829502179208e-10 +
2909 (0.36195460197779299406e-11 - 0.84882136022200714710e-13 * t) *
2919 double t = 2 * y100 - 59;
2920 return 0.26701724900280689785e0 +
2921 (0.84959936119625864274e-2 +
2922 (0.14112359443938883232e-3 +
2923 (0.17800427288596909634e-6 +
2924 (-0.13443492107643109071e-6 +
2925 (-0.23512456315677680293e-8 +
2926 (0.11245846264695936769e-9 +
2927 (0.19850501334649565404e-11 - 0.11284666134635050832e-12 * t) *
2937 double t = 2 * y100 - 61;
2938 return 0.28457293586253654144e0 +
2939 (0.90581563892650431899e-2 +
2940 (0.13880520331140646738e-3 +
2941 (-0.97262302362522896157e-6 +
2942 (-0.15077100040254187366e-6 +
2943 (-0.88574317464577116689e-9 +
2944 (0.12760311125637474581e-9 +
2945 (0.20155151018282695055e-12 - 0.10514169375181734921e-12 * t) *
2955 double t = 2 * y100 - 63;
2956 return 0.30323425595617385705e0 +
2957 (0.95968346790597422934e-2 +
2958 (0.12931067776725883939e-3 +
2959 (-0.21938741702795543986e-5 +
2960 (-0.15202888584907373963e-6 +
2961 (0.61788350541116331411e-9 +
2962 (0.11957835742791248256e-9 + (-0.12598179834007710908e-11 -
2963 0.75151817129574614194e-13 * t) *
2973 double t = 2 * y100 - 65;
2974 return 0.32292521181517384379e0 +
2975 (0.10082957727001199408e-1 +
2976 (0.11257589426154962226e-3 +
2977 (-0.33670890319327881129e-5 +
2978 (-0.13910529040004008158e-6 +
2979 (0.19170714373047512945e-8 +
2980 (0.94840222377720494290e-10 + (-0.21650018351795353201e-11 -
2981 0.37875211678024922689e-13 * t) *
2991 double t = 2 * y100 - 67;
2992 return 0.34351233557911753862e0 +
2993 (0.10488575435572745309e-1 +
2994 (0.89209444197248726614e-4 +
2995 (-0.43893459576483345364e-5 +
2996 (-0.11488595830450424419e-6 +
2997 (0.28599494117122464806e-8 +
2998 (0.61537542799857777779e-10 - 0.24935749227658002212e-11 * t) *
3007 double t = 2 * y100 - 69;
3008 return 0.36480946642143669093e0 +
3009 (0.10789304203431861366e-1 +
3010 (0.60357993745283076834e-4 +
3011 (-0.51855862174130669389e-5 +
3012 (-0.83291664087289801313e-7 +
3013 (0.33898011178582671546e-8 +
3014 (0.27082948188277716482e-10 + (-0.23603379397408694974e-11 +
3015 0.19328087692252869842e-13 * t) *
3025 double t = 2 * y100 - 71;
3026 return 0.38658679935694939199e0 +
3027 (0.10966119158288804999e-1 +
3028 (0.27521612041849561426e-4 +
3029 (-0.57132774537670953638e-5 +
3030 (-0.48404772799207914899e-7 +
3031 (0.35268354132474570493e-8 + (-0.32383477652514618094e-11 +
3032 (-0.19334202915190442501e-11 +
3033 0.32333189861286460270e-13 * t) *
3043 double t = 2 * y100 - 73;
3044 return 0.40858275583808707870e0 +
3045 (0.11006378016848466550e-1 +
3046 (-0.76396376685213286033e-5 +
3047 (-0.59609835484245791439e-5 +
3048 (-0.13834610033859313213e-7 +
3049 (0.33406952974861448790e-8 + (-0.26474915974296612559e-10 +
3050 (-0.13750229270354351983e-11 +
3051 0.36169366979417390637e-13 * t) *
3061 double t = 2 * y100 - 75;
3062 return 0.43051714914006682977e0 +
3063 (0.10904106549500816155e-1 +
3064 (-0.43477527256787216909e-4 +
3065 (-0.59429739547798343948e-5 +
3066 (0.17639200194091885949e-7 +
3067 (0.29235991689639918688e-8 + (-0.41718791216277812879e-10 +
3068 (-0.81023337739508049606e-12 +
3069 0.33618915934461994428e-13 * t) *
3079 double t = 2 * y100 - 77;
3080 return 0.45210428135559607406e0 +
3081 (0.10659670756384400554e-1 +
3082 (-0.78488639913256978087e-4 +
3083 (-0.56919860886214735936e-5 +
3084 (0.44181850467477733407e-7 +
3085 (0.23694306174312688151e-8 + (-0.49492621596685443247e-10 +
3086 (-0.31827275712126287222e-12 +
3087 0.27494438742721623654e-13 * t) *
3097 double t = 2 * y100 - 79;
3098 return 0.47306491195005224077e0 +
3099 (0.10279006119745977570e-1 +
3100 (-0.11140268171830478306e-3 +
3101 (-0.52518035247451432069e-5 +
3102 (0.64846898158889479518e-7 +
3103 (0.17603624837787337662e-8 +
3104 (-0.51129481592926104316e-10 +
3105 (0.62674584974141049511e-13 + 0.20055478560829935356e-13 * t) *
3115 double t = 2 * y100 - 81;
3116 return 0.49313638965719857647e0 +
3117 (0.97725799114772017662e-2 +
3118 (-0.14122854267291533334e-3 +
3119 (-0.46707252568834951907e-5 +
3120 (0.79421347979319449524e-7 +
3121 (0.11603027184324708643e-8 +
3122 (-0.48269605844397175946e-10 +
3123 (0.32477251431748571219e-12 + 0.12831052634143527985e-13 * t) *
3133 double t = 2 * y100 - 83;
3134 return 0.51208057433416004042e0 +
3135 (0.91542422354009224951e-2 +
3136 (-0.16726530230228647275e-3 +
3137 (-0.39964621752527649409e-5 +
3138 (0.88232252903213171454e-7 +
3139 (0.61343113364949928501e-9 +
3140 (-0.42516755603130443051e-10 +
3141 (0.47910437172240209262e-12 + 0.66784341874437478953e-14 * t) *
3151 double t = 2 * y100 - 85;
3152 return 0.52968945458607484524e0 +
3153 (0.84400880445116786088e-2 +
3154 (-0.18908729783854258774e-3 +
3155 (-0.32725905467782951931e-5 +
3156 (0.91956190588652090659e-7 +
3157 (0.14593989152420122909e-9 +
3158 (-0.35239490687644444445e-10 + 0.54613829888448694898e-12 * t) *
3167 double t = 2 * y100 - 87;
3168 return 0.54578857454330070965e0 +
3169 (0.76474155195880295311e-2 +
3170 (-0.20651230590808213884e-3 +
3171 (-0.25364339140543131706e-5 +
3172 (0.91455367999510681979e-7 +
3173 (-0.23061359005297528898e-9 +
3174 (-0.27512928625244444444e-10 + 0.54895806008493285579e-12 * t) *
3183 double t = 2 * y100 - 89;
3184 return 0.56023851910298493910e0 +
3185 (0.67938321739997196804e-2 +
3186 (-0.21956066613331411760e-3 +
3187 (-0.18181127670443266395e-5 +
3188 (0.87650335075416845987e-7 +
3189 (-0.51548062050366615977e-9 +
3190 (-0.20068462174044444444e-10 + 0.50912654909758187264e-12 * t) *
3199 double t = 2 * y100 - 91;
3200 return 0.57293478057455721150e0 +
3201 (0.58965321010394044087e-2 +
3202 (-0.22841145229276575597e-3 +
3203 (-0.11404605562013443659e-5 +
3204 (0.81430290992322326296e-7 +
3205 (-0.71512447242755357629e-9 +
3206 (-0.13372664928000000000e-10 + 0.44461498336689298148e-12 * t) *
3215 double t = 2 * y100 - 93;
3216 return 0.58380635448407827360e0 +
3217 (0.49717469530842831182e-2 +
3218 (-0.23336001540009645365e-3 +
3219 (-0.51952064448608850822e-6 +
3220 (0.73596577815411080511e-7 +
3221 (-0.84020916763091566035e-9 +
3222 (-0.76700972702222222221e-11 + 0.36914462807972467044e-12 * t) *
3231 double t = 2 * y100 - 95;
3232 return 0.59281340237769489597e0 +
3233 (0.40343592069379730568e-2 +
3234 (-0.23477963738658326185e-3 +
3235 (0.34615944987790224234e-7 +
3236 (0.64832803248395814574e-7 +
3237 (-0.90329163587627007971e-9 +
3238 (-0.30421940400000000000e-11 + 0.29237386653743536669e-12 * t) *
3247 double t = 2 * y100 - 97;
3248 return 0.59994428743114271918e0 +
3249 (0.30976579788271744329e-2 +
3250 (-0.23308875765700082835e-3 +
3251 (0.51681681023846925160e-6 +
3252 (0.55694594264948268169e-7 +
3253 (-0.91719117313243464652e-9 +
3254 (0.53982743680000000000e-12 + 0.22050829296187771142e-12 * t) *
3263 double t = 2 * y100 - 99;
3264 return 0.60521224471819875444e0 +
3265 (0.21732138012345456060e-2 +
3266 (-0.22872428969625997456e-3 +
3267 (0.92588959922653404233e-6 +
3268 (0.46612665806531930684e-7 +
3269 (-0.89393722514414153351e-9 +
3270 (0.31718550353777777778e-11 + 0.15705458816080549117e-12 * t) *
3279 double t = 2 * y100 - 101;
3280 return 0.60865189969791123620e0 +
3281 (0.12708480848877451719e-2 +
3282 (-0.22212090111534847166e-3 +
3283 (0.12636236031532793467e-5 +
3284 (0.37904037100232937574e-7 +
3285 (-0.84417089968101223519e-9 +
3286 (0.49843180828444444445e-11 + 0.10355439441049048273e-12 * t) *
3295 double t = 2 * y100 - 103;
3296 return 0.61031580103499200191e0 +
3297 (0.39867436055861038223e-3 +
3298 (-0.21369573439579869291e-3 +
3299 (0.15339402129026183670e-5 +
3300 (0.29787479206646594442e-7 +
3301 (-0.77687792914228632974e-9 +
3302 (0.61192452741333333334e-11 + 0.60216691829459295780e-13 * t) *
3311 double t = 2 * y100 - 105;
3312 return 0.61027109047879835868e0 +
3313 (-0.43680904508059878254e-3 +
3314 (-0.20383783788303894442e-3 +
3315 (0.17421743090883439959e-5 +
3316 (0.22400425572175715576e-7 +
3317 (-0.69934719320045128997e-9 +
3318 (0.67152759655111111110e-11 + 0.26419960042578359995e-13 * t) *
3327 double t = 2 * y100 - 107;
3328 return 0.60859639489217430521e0 +
3329 (-0.12305921390962936873e-2 +
3330 (-0.19290150253894682629e-3 +
3331 (0.18944904654478310128e-5 +
3332 (0.15815530398618149110e-7 +
3333 (-0.61726850580964876070e-9 + 0.68987888999111111110e-11 * t) *
3341 double t = 2 * y100 - 109;
3342 return 0.60537899426486075181e0 +
3343 (-0.19790062241395705751e-2 +
3344 (-0.18120271393047062253e-3 +
3345 (0.19974264162313241405e-5 +
3346 (0.10055795094298172492e-7 +
3347 (-0.53491997919318263593e-9 +
3348 (0.67794550295111111110e-11 - 0.17059208095741511603e-13 * t) *
3357 double t = 2 * y100 - 111;
3358 return 0.60071229457904110537e0 +
3359 (-0.26795676776166354354e-2 +
3360 (-0.16901799553627508781e-3 +
3361 (0.20575498324332621581e-5 +
3362 (0.51077165074461745053e-8 +
3363 (-0.45536079828057221858e-9 +
3364 (0.64488005516444444445e-11 - 0.29311677573152766338e-13 * t) *
3373 double t = 2 * y100 - 113;
3374 return 0.59469361520112714738e0 +
3375 (-0.33308208190600993470e-2 +
3376 (-0.15658501295912405679e-3 +
3377 (0.20812116912895417272e-5 +
3378 (0.93227468760614182021e-9 +
3379 (-0.38066673740116080415e-9 +
3380 (0.59806790359111111110e-11 - 0.36887077278950440597e-13 * t) *
3389 double t = 2 * y100 - 115;
3390 return 0.58742228631775388268e0 +
3391 (-0.39321858196059227251e-2 +
3392 (-0.14410441141450122535e-3 +
3393 (0.20743790018404020716e-5 +
3394 (-0.25261903811221913762e-8 +
3395 (-0.31212416519526924318e-9 +
3396 (0.54328422462222222221e-11 - 0.40864152484979815972e-13 * t) *
3405 double t = 2 * y100 - 117;
3406 return 0.57899804200033018447e0 +
3407 (-0.44838157005618913447e-2 +
3408 (-0.13174245966501437965e-3 +
3409 (0.20425306888294362674e-5 +
3410 (-0.53330296023875447782e-8 +
3411 (-0.25041289435539821014e-9 +
3412 (0.48490437205333333334e-11 - 0.42162206939169045177e-13 * t) *
3421 double t = 2 * y100 - 119;
3422 return 0.56951968796931245974e0 +
3423 (-0.49864649488074868952e-2 +
3424 (-0.11963416583477567125e-3 +
3425 (0.19906021780991036425e-5 +
3426 (-0.75580140299436494248e-8 +
3427 (-0.19576060961919820491e-9 +
3428 (0.42613011928888888890e-11 - 0.41539443304115604377e-13 * t) *
3437 double t = 2 * y100 - 121;
3438 return 0.55908401930063918964e0 +
3439 (-0.54413711036826877753e-2 +
3440 (-0.10788661102511914628e-3 +
3441 (0.19229663322982839331e-5 +
3442 (-0.92714731195118129616e-8 +
3443 (-0.14807038677197394186e-9 +
3444 (0.36920870298666666666e-11 - 0.39603726688419162617e-13 * t) *
3453 double t = 2 * y100 - 123;
3454 return 0.54778496152925675315e0 +
3455 (-0.58501497933213396670e-2 +
3456 (-0.96582314317855227421e-4 +
3457 (0.18434405235069270228e-5 +
3458 (-0.10541580254317078711e-7 +
3459 (-0.10702303407788943498e-9 +
3460 (0.31563175582222222222e-11 - 0.36829748079110481422e-13 * t) *
3469 double t = 2 * y100 - 125;
3470 return 0.53571290831682823999e0 +
3471 (-0.62147030670760791791e-2 +
3472 (-0.85782497917111760790e-4 +
3473 (0.17553116363443470478e-5 +
3474 (-0.11432547349815541084e-7 +
3475 (-0.72157091369041330520e-10 +
3476 (0.26630811607111111111e-11 - 0.33578660425893164084e-13 * t) *
3485 double t = 2 * y100 - 127;
3486 return 0.52295422962048434978e0 +
3487 (-0.65371404367776320720e-2 +
3488 (-0.75530164941473343780e-4 +
3489 (0.16613725797181276790e-5 +
3490 (-0.12003521296598910761e-7 +
3491 (-0.42929753689181106171e-10 +
3492 (0.22170894940444444444e-11 - 0.30117697501065110505e-13 * t) *
3501 double t = 2 * y100 - 129;
3502 return 0.50959092577577886140e0 +
3503 (-0.68197117603118591766e-2 +
3504 (-0.65852936198953623307e-4 +
3505 (0.15639654113906716939e-5 +
3506 (-0.12308007991056524902e-7 +
3507 (-0.18761997536910939570e-10 +
3508 (0.18198628922666666667e-11 - 0.26638355362285200932e-13 * t) *
3517 double t = 2 * y100 - 131;
3518 return 0.49570040481823167970e0 +
3519 (-0.70647509397614398066e-2 +
3520 (-0.56765617728962588218e-4 +
3521 (0.14650274449141448497e-5 +
3522 (-0.12393681471984051132e-7 +
3523 (0.92904351801168955424e-12 +
3524 (0.14706755960177777778e-11 - 0.23272455351266325318e-13 * t) *
3533 double t = 2 * y100 - 133;
3534 return 0.48135536250935238066e0 +
3535 (-0.72746293327402359783e-2 +
3536 (-0.48272489495730030780e-4 +
3537 (0.13661377309113939689e-5 +
3538 (-0.12302464447599382189e-7 +
3539 (0.16707760028737074907e-10 +
3540 (0.11672928324444444444e-11 - 0.20105801424709924499e-13 * t) *
3549 double t = 2 * y100 - 135;
3550 return 0.46662374675511439448e0 +
3551 (-0.74517177649528487002e-2 +
3552 (-0.40369318744279128718e-4 +
3553 (0.12685621118898535407e-5 +
3554 (-0.12070791463315156250e-7 +
3555 (0.29105507892605823871e-10 +
3556 (0.90653314645333333334e-12 - 0.17189503312102982646e-13 * t) *
3565 double t = 2 * y100 - 137;
3566 return 0.45156879030168268778e0 +
3567 (-0.75983560650033817497e-2 +
3568 (-0.33045110380705139759e-4 +
3569 (0.11732956732035040896e-5 +
3570 (-0.11729986947158201869e-7 +
3571 (0.38611905704166441308e-10 +
3572 (0.68468768305777777779e-12 - 0.14549134330396754575e-13 * t) *
3581 double t = 2 * y100 - 139;
3582 return 0.43624909769330896904e0 +
3583 (-0.77168291040309554679e-2 +
3584 (-0.26283612321339907756e-4 +
3585 (0.10811018836893550820e-5 +
3586 (-0.11306707563739851552e-7 +
3587 (0.45670446788529607380e-10 +
3588 (0.49782492549333333334e-12 - 0.12191983967561779442e-13 * t) *
3597 double t = 2 * y100 - 141;
3598 return 0.42071877443548481181e0 +
3599 (-0.78093484015052730097e-2 +
3600 (-0.20064596897224934705e-4 +
3601 (0.99254806680671890766e-6 +
3602 (-0.10823412088884741451e-7 +
3603 (0.50677203326904716247e-10 +
3604 (0.34200547594666666666e-12 - 0.10112698698356194618e-13 * t) *
3613 double t = 2 * y100 - 143;
3614 return 0.40502758809710844280e0 +
3615 (-0.78780384460872937555e-2 +
3616 (-0.14364940764532853112e-4 +
3617 (0.90803709228265217384e-6 +
3618 (-0.10298832847014466907e-7 +
3619 (0.53981671221969478551e-10 +
3620 (0.21342751381333333333e-12 - 0.82975901848387729274e-14 * t) *
3629 double t = 2 * y100 - 145;
3630 return 0.38922115269731446690e0 +
3631 (-0.79249269708242064120e-2 +
3632 (-0.91595258799106970453e-5 +
3633 (0.82783535102217576495e-6 +
3634 (-0.97484311059617744437e-8 +
3635 (0.55889029041660225629e-10 +
3636 (0.10851981336888888889e-12 - 0.67278553237853459757e-14 * t) *
3645 double t = 2 * y100 - 147;
3646 return 0.37334112915460307335e0 +
3647 (-0.79519385109223148791e-2 +
3648 (-0.44219833548840469752e-5 +
3649 (0.75209719038240314732e-6 +
3650 (-0.91848251458553190451e-8 +
3651 (0.56663266668051433844e-10 +
3652 (0.23995894257777777778e-13 - 0.53819475285389344313e-14 * t) *
3661 double t = 2 * y100 - 149;
3662 return 0.35742543583374223085e0 +
3663 (-0.79608906571527956177e-2 +
3664 (-0.12530071050975781198e-6 +
3665 (0.68088605744900552505e-6 +
3666 (-0.86181844090844164075e-8 +
3667 (0.56530784203816176153e-10 +
3668 (-0.43120012248888888890e-13 - 0.42372603392496813810e-14 * t) *
3677 double t = 2 * y100 - 151;
3678 return 0.34150846431979618536e0 +
3679 (-0.79534924968773806029e-2 +
3680 (0.37576885610891515813e-5 +
3681 (0.61419263633090524326e-6 +
3682 (-0.80565865409945960125e-8 +
3683 (0.55684175248749269411e-10 +
3684 (-0.95486860764444444445e-13 - 0.32712946432984510595e-14 * t) *
3693 double t = 2 * y100 - 153;
3694 return 0.32562129649136346824e0 +
3695 (-0.79313448067948884309e-2 +
3696 (0.72539159933545300034e-5 +
3697 (0.55195028297415503083e-6 +
3698 (-0.75063365335570475258e-8 +
3699 (0.54281686749699595941e-10 - 0.13545424295111111111e-12 * t) *
3707 double t = 2 * y100 - 155;
3708 return 0.30979191977078391864e0 +
3709 (-0.78959416264207333695e-2 +
3710 (0.10389774377677210794e-4 +
3711 (0.49404804463196316464e-6 +
3712 (-0.69722488229411164685e-8 +
3713 (0.52469254655951393842e-10 - 0.16507860650666666667e-12 * t) *
3721 double t = 2 * y100 - 157;
3722 return 0.29404543811214459904e0 +
3723 (-0.78486728990364155356e-2 +
3724 (0.13190885683106990459e-4 +
3725 (0.44034158861387909694e-6 +
3726 (-0.64578942561562616481e-8 +
3727 (0.50354306498006928984e-10 - 0.18614473550222222222e-12 * t) *
3735 double t = 2 * y100 - 159;
3736 return 0.27840427686253660515e0 +
3737 (-0.77908279176252742013e-2 +
3738 (0.15681928798708548349e-4 +
3739 (0.39066226205099807573e-6 +
3740 (-0.59658144820660420814e-8 +
3741 (0.48030086420373141763e-10 - 0.20018995173333333333e-12 * t) *
3749 double t = 2 * y100 - 161;
3750 return 0.26288838011163800908e0 +
3751 (-0.77235993576119469018e-2 +
3752 (0.17886516796198660969e-4 +
3753 (0.34482457073472497720e-6 +
3754 (-0.54977066551955420066e-8 +
3755 (0.45572749379147269213e-10 - 0.20852924954666666667e-12 * t) *
3763 double t = 2 * y100 - 163;
3764 return 0.24751539954181029717e0 +
3765 (-0.76480877165290370975e-2 +
3766 (0.19827114835033977049e-4 +
3767 (0.30263228619976332110e-6 +
3768 (-0.50545814570120129947e-8 +
3769 (0.43043879374212005966e-10 - 0.21228012028444444444e-12 * t) *
3777 double t = 2 * y100 - 165;
3778 return 0.23230087411688914593e0 +
3779 (-0.75653060136384041587e-2 +
3780 (0.21524991113020016415e-4 +
3781 (0.26388338542539382413e-6 +
3782 (-0.46368974069671446622e-8 +
3783 (0.40492715758206515307e-10 - 0.21238627815111111111e-12 * t) *
3791 double t = 2 * y100 - 167;
3792 return 0.21725840021297341931e0 +
3793 (-0.74761846305979730439e-2 +
3794 (0.23000194404129495243e-4 +
3795 (0.22837400135642906796e-6 +
3796 (-0.42446743058417541277e-8 +
3797 (0.37958104071765923728e-10 - 0.20963978568888888889e-12 * t) *
3805 double t = 2 * y100 - 169;
3806 return 0.20239979200788191491e0 +
3807 (-0.73815761980493466516e-2 +
3808 (0.24271552727631854013e-4 +
3809 (0.19590154043390012843e-6 +
3810 (-0.38775884642456551753e-8 +
3811 (0.35470192372162901168e-10 - 0.20470131678222222222e-12 * t) *
3819 double t = 2 * y100 - 171;
3820 return 0.18773523211558098962e0 +
3821 (-0.72822604530339834448e-2 +
3822 (0.25356688567841293697e-4 +
3823 (0.16626710297744290016e-6 +
3824 (-0.35350521468015310830e-8 +
3825 (0.33051896213898864306e-10 - 0.19811844544000000000e-12 * t) *
3833 double t = 2 * y100 - 173;
3834 return 0.17327341258479649442e0 +
3835 (-0.71789490089142761950e-2 +
3836 (0.26272046822383820476e-4 +
3837 (0.13927732375657362345e-6 +
3838 (-0.32162794266956859603e-8 +
3839 (0.30720156036105652035e-10 - 0.19034196304000000000e-12 * t) *
3847 double t = 2 * y100 - 175;
3848 return 0.15902166648328672043e0 +
3849 (-0.70722899934245504034e-2 +
3850 (0.27032932310132226025e-4 +
3851 (0.11474573347816568279e-6 +
3852 (-0.29203404091754665063e-8 +
3853 (0.28487010262547971859e-10 - 0.18174029063111111111e-12 * t) *
3861 double t = 2 * y100 - 177;
3862 return 0.14498609036610283865e0 +
3863 (-0.69628725220045029273e-2 +
3864 (0.27653554229160596221e-4 +
3865 (0.92493727167393036470e-7 +
3866 (-0.26462055548683583849e-8 +
3867 (0.26360506250989943739e-10 - 0.17261211260444444444e-12 * t) *
3875 double t = 2 * y100 - 179;
3876 return 0.13117165798208050667e0 +
3877 (-0.68512309830281084723e-2 +
3878 (0.28147075431133863774e-4 +
3879 (0.72351212437979583441e-7 +
3880 (-0.23927816200314358570e-8 +
3881 (0.24345469651209833155e-10 - 0.16319736960000000000e-12 * t) *
3889 double t = 2 * y100 - 181;
3890 return 0.11758232561160626306e0 +
3891 (-0.67378491192463392927e-2 +
3892 (0.28525664781722907847e-4 +
3893 (0.54156999310046790024e-7 +
3894 (-0.21589405340123827823e-8 +
3895 (0.22444150951727334619e-10 - 0.15368675584000000000e-12 * t) *
3903 double t = 2 * y100 - 183;
3904 return 0.10422112945361673560e0 +
3905 (-0.66231638959845581564e-2 +
3906 (0.28800551216363918088e-4 +
3907 (0.37758983397952149613e-7 +
3908 (-0.19435423557038933431e-8 +
3909 (0.20656766125421362458e-10 - 0.14422990012444444444e-12 * t) *
3917 double t = 2 * y100 - 185;
3918 return 0.91090275493541084785e-1 +
3919 (-0.65075691516115160062e-2 +
3920 (0.28982078385527224867e-4 +
3921 (0.23014165807643012781e-7 +
3922 (-0.17454532910249875958e-8 +
3923 (0.18981946442680092373e-10 - 0.13494234691555555556e-12 * t) *
3931 double t = 2 * y100 - 187;
3932 return 0.78191222288771379358e-1 +
3933 (-0.63914190297303976434e-2 +
3934 (0.29079759021299682675e-4 +
3935 (0.97885458059415717014e-8 +
3936 (-0.15635596116134296819e-8 +
3937 (0.17417110744051331974e-10 - 0.12591151763555555556e-12 * t) *
3945 double t = 2 * y100 - 189;
3946 return 0.65524757106147402224e-1 +
3947 (-0.62750311956082444159e-2 +
3948 (0.29102328354323449795e-4 +
3949 (-0.20430838882727954582e-8 +
3950 (-0.13967781903855367270e-8 +
3951 (0.15958771833747057569e-10 - 0.11720175765333333333e-12 * t) *
3959 double t = 2 * y100 - 191;
3960 return 0.53091065838453612773e-1 +
3961 (-0.61586898417077043662e-2 +
3962 (0.29057796072960100710e-4 +
3963 (-0.12597414620517987536e-7 +
3964 (-0.12440642607426861943e-8 +
3965 (0.14602787128447932137e-10 - 0.10885859114666666667e-12 * t) *
3973 double t = 2 * y100 - 193;
3974 return 0.40889797115352738582e-1 +
3975 (-0.60426484889413678200e-2 +
3976 (0.28953496450191694606e-4 +
3977 (-0.21982952021823718400e-7 +
3978 (-0.11044169117553026211e-8 +
3979 (0.13344562332430552171e-10 - 0.10091231402844444444e-12 * t) *
3992 return x * (1.1283791670955125739 -
3993 x2 * (0.75225277806367504925 -
3994 x2 * (0.30090111122547001970 -
3995 x2 * (0.085971746064420005629 -
3996 x2 * 0.016931216931216931217))));
4007 const double ispi = 0.56418958354775628694807945156;
4012 return ispi * ((x * x) * (x * x - 4.5) + 2) /
4013 (x * ((x * x) * (x * x - 5) + 3.75));
4018 const double ispi = 0.56418958354775628694807945156;
4023 return ispi * ((x * x) * (x * x - 4.5) + 2) /
4024 (x * ((x * x) * (x * x - 5) + 3.75));
4033 #ifdef TEST_FADDEEVA 4042 static double relerr(
double a,
double b) {
4043 if (isnan(a) || isnan(b) || isinf(a) || isinf(b)) {
4044 if ((isnan(a) && !isnan(b)) || (!isnan(a) && isnan(b)) ||
4045 (isinf(a) && !isinf(b)) || (!isinf(a) && isinf(b)) ||
4046 (isinf(a) && isinf(b) && a * b < 0))
4051 return b == 0 ? 0 :
Inf;
4053 return fabs((b - a) / a);
4057 double errmax_all = 0;
4059 printf(
"############# w(z) tests #############\n");
4060 #define NTST 57 // define instead of const for C compatibility 4061 cmplx z[NTST] = {C(624.2, -0.26123),
4067 C(-0.0000000234545, 1.1234),
4081 C(2.611780000000000e+01, 4.540909610972489e+03),
4097 C(10.01, -0.99e-10),
4125 C(-3.78270245518980507452677445620103199303131110e-7,
4126 0.000903861276433172057331093754199933411710053155),
4127 C(0.1764906227004816847297495349730234591778719532788,
4128 -0.02146550539468457616788719893991501311573031095617),
4129 C(0.2410250715772692146133539023007113781272362309451,
4130 0.06087579663428089745895459735240964093522265589350),
4131 C(0.30474420525691259245713884106959496013413834051768,
4132 -0.20821893820283162728743734725471561394145872072738),
4133 C(7.317131068972378096865595229600561710140617977e34,
4134 8.321873499714402777186848353320412813066170427e34),
4135 C(0.0615698507236323685519612934241429530190806818395,
4136 -0.00676005783716575013073036218018565206070072304635),
4137 C(0.3960793007699874918961319170187598400134746631,
4138 -5.593152259116644920546186222529802777409274656e-9),
4139 C(0.08217199226739447943295069917990417630675021771804,
4140 -0.04701291087643609891018366143118110965272615832184),
4141 C(0.00457246000350281640952328010227885008541748668738,
4142 -0.00804900791411691821818731763401840373998654987934),
4143 C(0.8746342859608052666092782112565360755791467973338452, 0.),
4144 C(0.00468190164965444174367477874864366058339647648741,
4145 0.0510735563901306197993676329845149741675029197050),
4146 C(-0.0023193175200187620902125853834909543869428763219,
4147 -0.025460054739731556004902057663500272721780776336),
4148 C(9.11463368405637174660562096516414499772662584e304,
4149 3.97101807145263333769664875189354358563218932e305),
4150 C(-4.4927207857715598976165541011143706155432296e281,
4151 -2.8019591213423077494444700357168707775769028e281),
4152 C(2.820947917809305132678577516325951485807107151e-6,
4153 2.820947917668257736791638444590253942253354058e-6),
4154 C(2.82094791773878143474039725787438662716372268e-15,
4155 2.82094791773878143474039725773333923127678361e-15),
4156 C(-0.0000563851289696244350147899376081488003110150498,
4157 -0.000169211755126812174631861529808288295454992688),
4158 C(-5.586035480670854326218608431294778077663867e-162,
4159 5.586035480670854326218608431294778077663867e-161),
4160 C(0.00016318325137140451888255634399123461580248456,
4161 -0.095232456573009287370728788146686162555021209999),
4162 C(0.69504753678406939989115375989939096800793577783885,
4163 -1.8916411171103639136680830887017670616339912024317),
4164 C(0.0001242418269653279656612334210746733213167234822,
4165 7.145975826320186888508563111992099992116786763e-7),
4166 C(2.318587329648353318615800865959225429377529825e-8,
4167 6.182899545728857485721417893323317843200933380e-8),
4168 C(-0.0133426877243506022053521927604277115767311800303,
4169 -0.0148087097143220769493341484176979826888871576145),
4170 C(1.00000000000000012412170838050638522857747934,
4171 1.12837916709551279389615890312156495593616433e-16),
4172 C(0.9999999853310704677583504063775310832036830015,
4173 2.595272024519678881897196435157270184030360773e-8),
4174 C(-1.4731421795638279504242963027196663601154624e-15,
4175 0.090727659684127365236479098488823462473074709),
4176 C(5.79246077884410284575834156425396800754409308e-18,
4177 0.0907276596841273652364790985059772809093822374),
4178 C(0.0884658993528521953466533278764830881245144368,
4179 1.37088352495749125283269718778582613192166760e-22),
4180 C(0.0345480845419190424370085249304184266813447878,
4181 2.11161102895179044968099038990446187626075258e-23),
4182 C(6.63967719958073440070225527042829242391918213e-36,
4183 0.0630820900592582863713653132559743161572639353),
4184 C(0.00179435233208702644891092397579091030658500743634,
4185 0.0951983814805270647939647438459699953990788064762),
4186 C(9.09760377102097999924241322094863528771095448e-13,
4187 0.0709979210725138550986782242355007611074966717),
4188 C(7.2049510279742166460047102593255688682910274423e-304,
4189 0.0201552956479526953866611812593266285000876784321),
4190 C(3.04543604652250734193622967873276113872279682e-44,
4191 0.0566481651760675042930042117726713294607499165),
4192 C(3.04543604652250734193622967873276113872279682e-44,
4193 0.0566481651760675042930042117726713294607499165),
4194 C(0.5659928732065273429286988428080855057102069081e-12,
4195 0.056648165176067504292998527162143030538756683302),
4196 C(-0.56599287320652734292869884280802459698927645e-12,
4197 0.0566481651760675042929985271621430305387566833029),
4198 C(0.0796884251721652215687859778119964009569455462,
4199 1.11474461817561675017794941973556302717225126e-22),
4200 C(0.07817195821247357458545539935996687005781943386550,
4201 -0.01093913670103576690766705513142246633056714279654),
4202 C(0.04670032980990449912809326141164730850466208439937,
4203 0.03944038961933534137558064191650437353429669886545),
4204 C(0.36787944117144232159552377016146086744581113103176,
4205 0.60715770584139372911503823580074492116122092866515),
4206 C(0, 0.010259688805536830986089913987516716056946786526145),
4207 C(0.99004983374916805357390597718003655777207908125383,
4208 -0.11208866436449538036721343053869621153527769495574),
4209 C(0.99999999999999999999999999999999999999990000,
4210 1.12837916709551257389615890312154517168802603e-20),
4211 C(0.999999999999943581041645226871305192054749891144158, 0),
4212 C(0.0110604154853277201542582159216317923453996211744250, 0),
4225 for (
int i = 0; i < NTST; ++i) {
4226 cmplx fw = FADDEEVA(w)(z[i], 0.);
4227 double re_err = relerr(creal(w[i]), creal(fw));
4228 double im_err = relerr(cimag(w[i]), cimag(fw));
4230 "w(%g%+gi) = %g%+gi (vs. %g%+gi), re/im rel. err. = %0.2g/%0.2g)\n",
4231 creal(z[i]), cimag(z[i]), creal(fw), cimag(fw), creal(w[i]),
4232 cimag(w[i]), re_err, im_err);
4233 if (re_err > errmax)
4235 if (im_err > errmax)
4238 if (errmax > 1e-13) {
4239 printf(
"FAILURE -- relative error %g too large!\n", errmax);
4242 printf(
"SUCCESS (max relative error = %g)\n", errmax);
4243 if (errmax > errmax_all)
4244 errmax_all = errmax;
4248 #define NTST 41 // define instead of const for C compatibility 4249 cmplx z[NTST] = {C(1, 2),
4259 C(-4.9e-3, 4.95e-3),
4291 C(-0.5366435657785650339917955593141927494421,
4292 -5.049143703447034669543036958614140565553),
4293 C(0.5366435657785650339917955593141927494421,
4294 -5.049143703447034669543036958614140565553),
4295 C(-0.5366435657785650339917955593141927494421,
4296 5.049143703447034669543036958614140565553),
4297 C(0.5366435657785650339917955593141927494421,
4298 5.049143703447034669543036958614140565553),
4299 C(0.3359473673830576996788000505817956637777e304,
4300 -0.1999896139679880888755589794455069208455e304),
4301 C(0.3584459971462946066523939204836760283645e278,
4302 0.3818954885257184373734213077678011282505e280),
4303 C(0.9996020422657148639102150147542224526887,
4304 0.00002801044116908227889681753993542916894856),
4307 C(0.005754683859034800134412990541076554934877,
4308 0.1128349818335058741511924929801267822634e-7),
4309 C(-0.005529149142341821193633460286828381876955,
4310 0.005585388387864706679609092447916333443570),
4311 C(0.007099365669981359632319829148438283865814,
4312 0.6149347012854211635026981277569074001219),
4313 C(0.3981176338702323417718189922039863062440e8,
4314 -0.8298176341665249121085423917575122140650e10),
4316 C(0.007389128308257135427153919483147229573895,
4317 0.6149332524601658796226417164791221815139),
4318 C(0.4143671923267934479245651547534414976991e8,
4319 -0.8298168216818314211557046346850921446950e10),
4321 C(0.1128379167099649964175513742247082845155e-5,
4322 0.2256758334191777400570377193451519478895e-5),
4323 C(0, 0.2256758334194034158904576117253481476197e-5),
4324 C(0, 18.56480241457555259870429191324101719886),
4325 C(0, 0.1474797539628786202447733153131835124599e173),
4339 C(0.07924380404615782687930591956705225541145,
4340 0.07872776218046681145537914954027729115247),
4341 C(0.07885775828512276968931773651224684454495,
4342 -0.0007860046704118224342390725280161272277506),
4343 C(-0.1012806432747198859687963080684978759881,
4344 0.0007834934747022035607566216654982820299469),
4345 C(-0.1020998418798097910247132140051062512527,
4346 0.1010030778892310851309082083238896270340),
4347 C(-0.0007962891763147907785684591823889484764272,
4348 0.1018289385936278171741809237435404896152),
4349 C(0.07886408666470478681566329888615410479530,
4350 0.01010604288780868961492224347707949372245),
4351 C(0.07886723099940260286824654364807981336591,
4352 0.01235199327873258197931147306290916629654)};
4353 #define TST(f, isc) \ 4354 printf("############# " #f "(z) tests #############\n"); \ 4355 double errmax = 0; \ 4356 for (int i = 0; i < NTST; ++i) { \ 4357 cmplx fw = FADDEEVA(f)(z[i], 0.); \ 4358 double re_err = relerr(creal(w[i]), creal(fw)); \ 4359 double im_err = relerr(cimag(w[i]), cimag(fw)); \ 4361 "(%g%+gi) = %g%+gi (vs. %g%+gi), re/im rel. err. = %0.2g/%0.2g)\n", \ 4362 creal(z[i]), cimag(z[i]), creal(fw), cimag(fw), creal(w[i]), \ 4363 cimag(w[i]), re_err, im_err); \ 4364 if (re_err > errmax) \ 4366 if (im_err > errmax) \ 4369 if (errmax > 1e-13) { \ 4370 printf("FAILURE -- relative error %g too large!\n", errmax); \ 4373 printf("Checking " #f "(x) special case...\n"); \ 4374 for (int i = 0; i < 10000; ++i) { \ 4375 double x = pow(10., -300. + i * 600. / (10000 - 1)); \ 4377 relerr(FADDEEVA_RE(f)(x), creal(FADDEEVA(f)(C(x, x * isc), 0.))); \ 4378 if (re_err > errmax) \ 4381 relerr(FADDEEVA_RE(f)(-x), creal(FADDEEVA(f)(C(-x, x * isc), 0.))); \ 4382 if (re_err > errmax) \ 4387 relerr(FADDEEVA_RE(f)(Inf), creal(FADDEEVA(f)(C(Inf, 0.), 0.))); \ 4388 if (re_err > errmax) \ 4391 relerr(FADDEEVA_RE(f)(-Inf), creal(FADDEEVA(f)(C(-Inf, 0.), 0.))); \ 4392 if (re_err > errmax) \ 4394 re_err = relerr(FADDEEVA_RE(f)(NaN), creal(FADDEEVA(f)(C(NaN, 0.), 0.))); \ 4395 if (re_err > errmax) \ 4398 if (errmax > 1e-13) { \ 4399 printf("FAILURE -- relative error %g too large!\n", errmax); \ 4402 printf("SUCCESS (max relative error = %g)\n", errmax); \ 4403 if (errmax > errmax_all) \ 4412 #define NTST 1 // define instead of const for C compatibility 4413 cmplx z[NTST] = {C(1.234, 0.5678)};
4415 C(1.081032284405373149432716643834106923212,
4416 1.926775520840916645838949402886591180834)};
4423 #define NTST 1 // define instead of const for C compatibility 4424 cmplx z[NTST] = {C(1.234, 0.5678)};
4426 C(0.3382187479799972294747793561190487832579,
4427 -0.1116077470811648467464927471872945833154)};
4432 #define NTST 30 // define instead of const for C compatibility 4434 C(1, 2), C(-1, 2), C(1, -2), C(-1, -2),
4435 C(9, -28), C(21, -33), C(1e3, 1e3), C(-3001, -1000),
4436 C(1e160, -1e159), C(5.1e-3, 1e-8), C(0, 2e-6), C(0, 2),
4437 C(0, 20), C(0, 200), C(2e-6, 0), C(2, 0),
4438 C(20, 0), C(200, 0), C(
Inf, 0), C(-
Inf, 0),
4443 C(1.536643565778565033991795559314192749442,
4444 5.049143703447034669543036958614140565553),
4445 C(0.4633564342214349660082044406858072505579,
4446 5.049143703447034669543036958614140565553),
4447 C(1.536643565778565033991795559314192749442,
4448 -5.049143703447034669543036958614140565553),
4449 C(0.4633564342214349660082044406858072505579,
4450 -5.049143703447034669543036958614140565553),
4451 C(-0.3359473673830576996788000505817956637777e304,
4452 0.1999896139679880888755589794455069208455e304),
4453 C(-0.3584459971462946066523939204836760283645e278,
4454 -0.3818954885257184373734213077678011282505e280),
4455 C(0.0003979577342851360897849852457775473112748,
4456 -0.00002801044116908227889681753993542916894856),
4459 C(0.9942453161409651998655870094589234450651,
4460 -0.1128349818335058741511924929801267822634e-7),
4461 C(1, -0.2256758334194034158904576117253481476197e-5),
4462 C(1, -18.56480241457555259870429191324101719886),
4463 C(1, -0.1474797539628786202447733153131835124599e173),
4465 C(0.9999977432416658119838633199332831406314, 0),
4466 C(0.004677734981047265837930743632747071389108, 0),
4467 C(0.5395865611607900928934999167905345604088e-175, 0),
4485 #define NTST 48 // define instead of const for C compatibility 4486 cmplx z[NTST] = {C(2, 1), C(-2, 1),
4487 C(2, -1), C(-2, -1),
4488 C(-28, 9), C(33, -21),
4489 C(1e3, 1e3), C(-1000, -3001),
4490 C(1e-8, 5.1e-3), C(4.95e-3, -4.9e-3),
4491 C(5.1e-3, 5.1e-3), C(0.5, 4.9e-3),
4492 C(-0.5e1, 4.9e-4), C(-0.5e2, -4.9e-5),
4493 C(0.5e3, 4.9e-6), C(0.5, 5.1e-3),
4494 C(-0.5e1, 5.1e-4), C(-0.5e2, -5.1e-5),
4495 C(1e-6, 2e-6), C(2e-6, 0),
4497 C(200, 0), C(0, 4.9e-3),
4498 C(0, -5.1e-3), C(0, 2e-6),
4500 C(0, -200), C(
Inf, 0),
4506 C(39, 6.4e-5), C(41, 6.09e-5),
4507 C(4.9e7, 5e-11), C(5.1e7, 4.8e-11),
4508 C(1e9, 2.4e-12), C(1e11, 2.4e-14),
4509 C(1e13, 2.4e-16), C(1e300, 2.4e-303)};
4511 C(0.1635394094345355614904345232875688576839,
4512 -0.1531245755371229803585918112683241066853),
4513 C(-0.1635394094345355614904345232875688576839,
4514 -0.1531245755371229803585918112683241066853),
4515 C(0.1635394094345355614904345232875688576839,
4516 0.1531245755371229803585918112683241066853),
4517 C(-0.1635394094345355614904345232875688576839,
4518 0.1531245755371229803585918112683241066853),
4519 C(-0.01619082256681596362895875232699626384420,
4520 -0.005210224203359059109181555401330902819419),
4521 C(0.01078377080978103125464543240346760257008,
4522 0.006866888783433775382193630944275682670599),
4523 C(-0.5808616819196736225612296471081337245459,
4524 0.6688593905505562263387760667171706325749),
4526 C(0.1000052020902036118082966385855563526705e-7,
4527 0.005100088434920073153418834680320146441685),
4528 C(0.004950156837581592745389973960217444687524,
4529 -0.004899838305155226382584756154100963570500),
4530 C(0.005100176864319675957314822982399286703798,
4531 0.005099823128319785355949825238269336481254),
4532 C(0.4244534840871830045021143490355372016428,
4533 0.002820278933186814021399602648373095266538),
4534 C(-0.1021340733271046543881236523269967674156,
4535 -0.00001045696456072005761498961861088944159916),
4536 C(-0.01000200120119206748855061636187197886859,
4537 0.9805885888237419500266621041508714123763e-8),
4538 C(0.001000002000012000023960527532953151819595,
4539 -0.9800058800588007290937355024646722133204e-11),
4540 C(0.4244549085628511778373438768121222815752,
4541 0.002935393851311701428647152230552122898291),
4542 C(-0.1021340732357117208743299813648493928105,
4543 -0.00001088377943049851799938998805451564893540),
4544 C(-0.01000200120119126652710792390331206563616,
4545 0.1020612612857282306892368985525393707486e-7),
4546 C(0.1000000000007333333333344266666666664457e-5,
4547 0.2000000000001333333333323199999999978819e-5),
4548 C(0.1999999999994666666666675199999999990248e-5, 0),
4549 C(0.3013403889237919660346644392864226952119, 0),
4550 C(0.02503136792640367194699495234782353186858, 0),
4551 C(0.002500031251171948248596912483183760683918, 0),
4552 C(0, 0.004900078433419939164774792850907128053308),
4553 C(0, -0.005100088434920074173454208832365950009419),
4554 C(0, 0.2000000000005333333333341866666666676419e-5),
4555 C(0, -48.16001211429122974789822893525016528191),
4556 C(0, 0.4627407029504443513654142715903005954668e174),
4569 C(0.01282473148489433743567240624939698290584,
4570 -0.2105957276516618621447832572909153498104e-7),
4571 C(0.01219875253423634378984109995893708152885,
4572 -0.1813040560401824664088425926165834355953e-7),
4573 C(0.1020408163265306334945473399689037886997e-7,
4574 -0.1041232819658476285651490827866174985330e-25),
4575 C(0.9803921568627452865036825956835185367356e-8,
4576 -0.9227220299884665067601095648451913375754e-26),
4577 C(0.5000000000000000002500000000000000003750e-9,
4578 -0.1200000000000000001800000188712838420241e-29),
4579 C(5.00000000000000000000025000000000000000000003e-12,
4580 -1.20000000000000000000018000000000000000000004e-36),
4581 C(5.00000000000000000000000002500000000000000000e-14,
4582 -1.20000000000000000000000001800000000000000000e-42),
4586 printf(
"#####################################\n");
4587 printf(
"SUCCESS (max relative error = %g)\n", errmax_all);
double FADDEEVA_RE() Dawson(double x)
double FADDEEVA_RE() erfc(double x)
cmplx FADDEEVA() w(cmplx z, double relerr)
static double sinc(double x, double sinx)
cmplx FADDEEVA() erfi(cmplx z, double relerr)
int main(int argc, char **argv)
Simple Dalitz plot fit of the channel J/psi -> gamma pi0 pi0.
static double w_im_y100(double y100, double x)
static double sinh_taylor(double x)
static double sqr(double x)
double FADDEEVA() w_im(double x)
static cmplx cpolar(double r, double t)
static const double expa2n2[]
static double erfcx_y100(double y100)
cmplx FADDEEVA() erfcx(cmplx z, double relerr)
double FADDEEVA_RE() erf(double x)