Lines Matching refs:real
33 enum real SQRT2PI = 2.50662827463100050242E0L; // sqrt(2pi)
34 immutable real EULERGAMMA = 0.57721_56649_01532_86060_65120_90082_40243_10421_59335_93992L; /** Eul…
38 immutable real[8] GammaNumeratorCoeffs = [ 1.0L, variable
44 immutable real[9] GammaDenominatorCoeffs = [ 1.0L, variable
50 immutable real[9] GammaSmallCoeffs = [ 1.0L, variable
56 immutable real[9] GammaSmallNegCoeffs = [ -1.0L, variable
62 immutable real[7] logGammaStirlingCoeffs = [ variable
68 immutable real[7] logGammaNumerator = [ variable
74 immutable real[8] logGammaDenominator = [ variable
88 real gammaStirling(real x) in gammaStirling()
92 static immutable real[9] SmallStirlingCoeffs = [ in gammaStirling() local
98 static immutable real[7] LargeStirlingCoeffs = [ 1.0L, in gammaStirling() local
104 real w = 1.0L/x; in gammaStirling()
105 real y = exp(x); in gammaStirling()
111 real v = pow( x, 0.5L * x - 0.25L ); in gammaStirling()
117 static if (floatTraits!(real).realFormat == RealFormat.ieeeDouble) in gammaStirling()
122 real v = pow( x, 0.5 * x - 0.25 ); in gammaStirling()
145 real igammaTemmeLarge(real a, real x) in igammaTemmeLarge()
147 static immutable real[][13] coef = [ in igammaTemmeLarge() local
221 if (x == real.infinity && a != real.infinity) in igammaTemmeLarge()
224 if (x != real.infinity && a == real.infinity) in igammaTemmeLarge()
227 real sigma = (x - a) / a; in igammaTemmeLarge()
228 real phi = sigma - log(sigma + 1); in igammaTemmeLarge()
230 real y = a * phi; in igammaTemmeLarge()
231 real z = sqrt(2 * phi); in igammaTemmeLarge()
235 real[13] workspace; in igammaTemmeLarge()
239 real result = poly(1 / a, workspace); in igammaTemmeLarge()
253 static if (floatTraits!(real).realFormat == RealFormat.ieeeQuadruple)
254 enum real MAXGAMMA = 1755.5483429L;
255 else static if (floatTraits!(real).realFormat == RealFormat.ieeeExtended)
256 enum real MAXGAMMA = 1755.5483429L;
257 else static if (floatTraits!(real).realFormat == RealFormat.ieeeExtended53)
258 enum real MAXGAMMA = 1755.5483429L;
259 else static if (floatTraits!(real).realFormat == RealFormat.ieeeDouble)
260 enum real MAXGAMMA = 171.6243769L;
285 real gamma(real x) in gamma()
297 real q, z; in gamma()
299 if (x == -x.infinity) return real.nan; in gamma()
300 if ( fabs(x) > MAXGAMMA ) return real.infinity; in gamma()
313 if (x < -1/real.epsilon) in gamma()
316 return real.nan; in gamma()
321 return real.nan; // poles for all integers <0. in gamma()
322 real p = intpart; in gamma()
333 if ( z <= PI/real.max ) return sgngam * real.infinity; in gamma()
362 return real.nan; in gamma()
391 real fact = 1.0L;
392 for (int i=1; fact<real.max; ++i)
396 assert(feqrel(gamma(i*1.0L), fact) >= real.mant_dig-15);
399 assert(gamma(0.0) == real.infinity);
400 assert(gamma(-0.0) == -real.infinity);
404 assert(gamma(real.infinity) == real.infinity);
405 assert(gamma(real.max) == real.infinity);
406 assert(isNaN(gamma(-real.infinity)));
407 assert(gamma(real.min_normal*real.epsilon) == real.infinity);
408 assert(gamma(MAXGAMMA)< real.infinity);
409 assert(gamma(MAXGAMMA*2) == real.infinity);
412 real SQRT_PI = 1.77245385090551602729816748334114518279754945612238L;
415 assert(feqrel(gamma(0.5L), SQRT_PI) >= real.mant_dig-1);
416 assert(feqrel(gamma(17.25L), 4.224986665692703551570937158682064589938e13L) >= real.mant_dig-4);
418 …rel(gamma(1.0 / 3.0L), 2.67893853470774763365569294097467764412868937795730L) >= real.mant_dig-2);
420 3.62560990822190831193068515586767200299516768288006L) >= real.mant_dig-1);
422 4.59084371199880305320475827592915200343410999829340L) >= real.mant_dig-1);
440 real logGamma(real x) in logGamma()
453 real q, w, z, f, nx; in logGamma()
462 real p = floor(q); in logGamma()
464 return real.infinity; in logGamma()
466 real sgngam = 1; in logGamma()
477 return sgngam * real.infinity; in logGamma()
499 return real.infinity; in logGamma()
516 real p = x * rationalPoly( x, logGammaNumerator, logGammaDenominator); in logGamma()
523 const real LOGSQRT2PI = 0.91893853320467274178L; // log( sqrt( 2*pi ) ) in logGamma()
527 real p = 1.0L / (x*x); in logGamma()
535 assert(logGamma(real.infinity) == real.infinity);
536 assert(logGamma(-1.0) == real.infinity);
537 assert(logGamma(0.0) == real.infinity);
538 assert(logGamma(-50.0) == real.infinity);
541 assert(logGamma(real.min_normal*real.epsilon) == real.infinity);
542 assert(logGamma(-real.min_normal*real.epsilon) == real.infinity);
545 immutable static real[] testpoints = [ variable
566 assert( feqrel(logGamma(testpoints[i]), testpoints[i+1]) > real.mant_dig-5);
569 assert( feqrel(log(fabs(gamma(testpoints[i]))), testpoints[i+1]) > real.mant_dig-5);
572 assert(feqrel(logGamma(-50.2L),log(fabs(gamma(-50.2L)))) > real.mant_dig-2);
573 assert(feqrel(logGamma(-0.008L),log(fabs(gamma(-0.008L)))) > real.mant_dig-2);
574 assert(feqrel(logGamma(-38.8L),log(fabs(gamma(-38.8L)))) > real.mant_dig-4);
575 static if (real.mant_dig >= 64) // incl. 80-bit reals
576 assert(feqrel(logGamma(1500.0L),log(gamma(1500.0L))) > real.mant_dig-2);
577 else static if (real.mant_dig >= 53) // incl. 64-bit reals
578 assert(feqrel(logGamma(150.0L),log(gamma(150.0L))) > real.mant_dig-2);
599 static if (floatTraits!(real).realFormat == RealFormat.ieeeQuadruple)
601 enum real MAXLOG = 0x1.62e42fefa39ef35793c7673007e6p+13L; // log(real.max)
602 …enum real MINLOG = -0x1.6546282207802c89d24d65e96274p+13L; // log(real.min_normal*real.epsilon) = …
604 else static if (floatTraits!(real).realFormat == RealFormat.ieeeExtended)
606 enum real MAXLOG = 0x1.62e42fefa39ef358p+13L; // log(real.max)
607 …enum real MINLOG = -0x1.6436716d5406e6d8p+13L; // log(real.min_normal*real.epsilon) = log(smallest…
609 else static if (floatTraits!(real).realFormat == RealFormat.ieeeExtended53)
611 enum real MAXLOG = 0x1.62e42fefa39ef358p+13L; // log(real.max)
612 …enum real MINLOG = -0x1.6436716d5406e6d8p+13L; // log(real.min_normal*real.epsilon) = log(smallest…
614 else static if (floatTraits!(real).realFormat == RealFormat.ieeeDouble)
616 enum real MAXLOG = 0x1.62e42fefa39efp+9L; // log(real.max)
617 …enum real MINLOG = -0x1.74385446d71c3p+9L; // log(real.min_normal*real.epsilon) = log(smallest den…
622 enum real BETA_BIG = 9.223372036854775808e18L;
623 enum real BETA_BIGINV = 1.084202172485504434007e-19L;
646 real betaIncomplete(real aa, real bb, real xx ) in betaIncomplete()
652 return real.nan; // domain error in betaIncomplete()
659 return real.nan; // domain error in betaIncomplete()
665 real x; in betaIncomplete()
666 real xc; // = 1 - x in betaIncomplete()
668 real a, b; in betaIncomplete()
695 real w; in betaIncomplete()
699 real y = x * (a+b-2.0L) - (a-1.0L); in betaIncomplete()
714 real t = b * log(xc); in betaIncomplete()
746 if ( t <= real.epsilon ) in betaIncomplete()
748 t = 1.0L - real.epsilon; in betaIncomplete()
764 real betaIncompleteInv(real aa, real bb, real yy0 ) in betaIncompleteInv()
766 real a, b, y0, d, y, x, x0, x1, lgm, yp, di, dithresh, yl, yh, xt; in betaIncompleteInv()
847 x = 1.0L - real.epsilon; in betaIncompleteInv()
907 if ( rflg == 1 && x1 < real.epsilon ) in betaIncompleteInv()
930 x = 1.0L - real.epsilon; in betaIncompleteInv()
937 x = real.min_normal * real.min_normal; in betaIncompleteInv()
1006 if ( fabs(d/x) < (128.0L * real.epsilon) ) in betaIncompleteInv()
1010 dithresh = 256.0L * real.epsilon; in betaIncompleteInv()
1016 if ( x <= real.epsilon ) in betaIncompleteInv()
1017 x = 1.0L - real.epsilon; in betaIncompleteInv()
1043 assert(feqrel(betaIncomplete(8, 10, 0.2L), 0.010_934_315_234_099_2L) >= real.mant_dig - 5);
1044 …assert(feqrel(betaIncomplete(2, 2.5L, 0.9L), 0.989_722_597_604_452_767_171_003_59L) >= real.mant_d…
1045 static if (real.mant_dig >= 64) // incl. 80-bit reals
1046 …(betaIncomplete(1000, 800, 0.5L), 1.179140859734704555102808541457164E-06L) >= real.mant_dig - 13);
1048 …(betaIncomplete(1000, 800, 0.5L), 1.179140859734704555102808541457164E-06L) >= real.mant_dig - 14);
1049 …assert(feqrel(betaIncomplete(0.0001, 10000, 0.0001L), 0.999978059362107134278786L) >= real.mant_di…
1051 assert(feqrel(betaIncompleteInv(8, 10, 0.010_934_315_234_099_2L), 0.2L) >= real.mant_dig - 2);
1052 …assert(feqrel(betaIncomplete(0.01L, 498.437L, 0.0121433L), 0.99999664562033077636065L) >= real.man…
1053 …assert(feqrel(betaIncompleteInv(5, 10, 0.2000002972865658842L), 0.229121208190918L) >= real.mant_d…
1054 …assert(feqrel(betaIncompleteInv(4, 7, 0.8000002209179505L), 0.483657360076904L) >= real.mant_dig -…
1062 static if (real.mant_dig == 64) // 80-bit reals
1064 assert(betaIncompleteInv(0.01L, 8e-48L, 5.45464e-20L) == 1-real.epsilon);
1065 assert(betaIncompleteInv(0.01L, 8e-48L, 9e-26L) == 1-real.epsilon);
1076 0x1.f97749d90c7adba8p-63L) >= real.mant_dig - 39);
1077 real a1 = 3.40483L;
1079 real b1 = 2.82847e-25L;
1080 …assert(feqrel(betaIncompleteInv(0.01L, b1, 9e-26L), 0x1.549696104490aa9p-830L) >= real.mant_dig-10…
1087 assert(betaIncomplete(1.16251e20L, 2.18e39L, 5.45e-20L) == -real.infinity);
1099 real betaDistExpansion1(real a, real b, real x ) in betaDistExpansion1()
1101 real xk, pk, pkm1, pkm2, qk, qkm1, qkm2; in betaDistExpansion1()
1102 real k1, k2, k3, k4, k5, k6, k7, k8; in betaDistExpansion1()
1103 real r, t, ans; in betaDistExpansion1()
1122 const real thresh = 3.0L * real.epsilon; in betaDistExpansion1()
1188 real betaDistExpansion2(real a, real b, real x ) in betaDistExpansion2()
1190 real xk, pk, pkm1, pkm2, qk, qkm1, qkm2; in betaDistExpansion2()
1191 real k1, k2, k3, k4, k5, k6, k7, k8; in betaDistExpansion2()
1192 real r, t, ans, z; in betaDistExpansion2()
1211 const real thresh = 3.0L * real.epsilon; in betaDistExpansion2()
1272 real betaDistPowerSeries(real a, real b, real x ) in betaDistPowerSeries()
1274 real ai = 1.0L / a; in betaDistPowerSeries()
1275 real u = (1.0L - b) * x; in betaDistPowerSeries()
1276 real v = u / (a + 1.0L); in betaDistPowerSeries()
1277 real t1 = v; in betaDistPowerSeries()
1278 real t = u; in betaDistPowerSeries()
1279 real n = 2.0L; in betaDistPowerSeries()
1280 real s = 0.0L; in betaDistPowerSeries()
1281 real z = real.epsilon * ai; in betaDistPowerSeries()
1329 real gammaIncomplete(real a, real x ) in gammaIncomplete()
1352 real ax = a * log(x) - x - logGamma(a);
1360 real r = a;
1361 real c = 1.0L;
1362 real ans = 1.0L;
1369 } while ( c/ans > real.epsilon );
1375 real gammaIncompleteCompl(real a, real x ) in gammaIncompleteCompl()
1391 const real MAXLOGL = 1.1356523406294143949492E4L;
1395 real ax = a * log(x) - x - logGamma(a);
1402 real y = 1.0L - a;
1403 real z = x + y + 1.0L;
1404 real c = 0.0L;
1406 real pk, qk, t;
1408 real pkm2 = 1.0L;
1409 real qkm2 = x;
1410 real pkm1 = x + 1.0L;
1411 real qkm1 = z * x;
1412 real ans = pkm1/qkm1;
1419 real yc = y * c;
1424 real r = pk/qk;
1437 const real BIG = 9.223372036854775808e18L;
1446 } while ( t > real.epsilon );
1463 real gammaIncompleteComplInv(real a, real p) in gammaIncompleteComplInv()
1471 if (p == 0) return real.infinity;
1473 real y0 = p;
1474 const real MAXLOGL = 1.1356523406294143949492E4L;
1475 real x0, x1, x, yl, yh, y, d, lgm, dithresh;
1479 x0 = real.max;
1483 dithresh = 4.0 * real.epsilon;
1524 if ( x0 == real.max )
1528 while ( x0 == real.max )
1600 assert(gammaIncomplete(4545, real.infinity)==1);
1607 assert(gammaIncompleteCompl(384, real.infinity)==0);
1608 assert(gammaIncompleteComplInv(3, 0)==real.infinity);
1612 static if (real.mant_dig >= 64) // incl. 80-bit reals
1620 immutable real [7] Bn_n = [ variable
1635 real digamma(real x) in digamma()
1639 real p, q, nz, s, w, y, z; in digamma()
1653 return real.nan; // singularity. in digamma()
1719 assert(feqrel(digamma(0.25), -PI/2 - 3* LN2 - EULERGAMMA) >= real.mant_dig-7);
1720 …assert(feqrel(digamma(1.0L/6), -PI/2 *sqrt(3.0L) - 2* LN2 -1.5*log(3.0L) - EULERGAMMA) >= real.man…
1722 assert(feqrel(digamma(2.5), -EULERGAMMA - 2*LN2 + 2.0 + 2.0L/3) >= real.mant_dig-9);
1727 real y=0;
1732 assert(feqrel(digamma(k+1.0), -EULERGAMMA + y) >= real.mant_dig-2);
1745 real logmdigamma(real x) in logmdigamma()
1751 return real.infinity; in logmdigamma()
1753 return real.nan; in logmdigamma()
1756 real s = x; in logmdigamma()
1757 real w = 0.0; in logmdigamma()
1764 real y; in logmdigamma()
1767 immutable real z = 1.0/(s * s); in logmdigamma()
1779 assert(logmdigamma(0.0) == real.infinity);
1797 real logmdigammaInverse(real y) in logmdigammaInverse()
1802 immutable maxY = logmdigamma(real.min_normal); in logmdigammaInverse()
1803 assert(maxY > 0 && maxY <= real.max); in logmdigammaInverse()
1812 return real.nan; in logmdigammaInverse()
1814 if (y < real.min_normal) in logmdigammaInverse()
1823 return 1 / findRoot((real x) => logmdigamma(1 / x) - y, y, 2*y); in logmdigammaInverse()
1832 immutable Tuple!(real, real)[5] testData = [
1843 assert(isClose(logmdigamma(logmdigammaInverse(real.min_normal)), real.min_normal, 1e-15L));
1844 assert(isClose(logmdigamma(logmdigammaInverse(real.max/2)), real.max/2, 1e-15L));
1846 assert(isClose(logmdigammaInverse(logmdigamma(real.min_normal)), real.min_normal, 1e-15L));
1847 assert(isClose(logmdigammaInverse(logmdigamma(real.max/2)), real.max/2, 1e-15L));