Piotr Caban : msvcrt: Import tgamma implementation from musl.
Alexandre Julliard
julliard at winehq.org
Tue Jun 8 16:34:27 CDT 2021
Module: wine
Branch: master
Commit: 16fa6b78460c83c08d28aa0a77d6077ea5010650
URL: https://source.winehq.org/git/wine.git/?a=commit;h=16fa6b78460c83c08d28aa0a77d6077ea5010650
Author: Piotr Caban <piotr at codeweavers.com>
Date: Tue Jun 8 21:19:27 2021 +0200
msvcrt: Import tgamma implementation from musl.
Signed-off-by: Piotr Caban <piotr at codeweavers.com>
Signed-off-by: Alexandre Julliard <julliard at winehq.org>
---
configure | 1 -
configure.ac | 1 -
dlls/msvcrt/math.c | 122 ++++++++++++++++++++++++++++++++++++++++++++++++--
dlls/msvcrt/unixlib.c | 14 ------
dlls/msvcrt/unixlib.h | 1 -
include/config.h.in | 3 --
6 files changed, 118 insertions(+), 24 deletions(-)
diff --git a/configure b/configure
index 9acc579ca11..c0d24695456 100755
--- a/configure
+++ b/configure
@@ -19623,7 +19623,6 @@ for ac_func in \
exp2 \
exp2f \
fmaf \
- tgamma \
tgammaf
do :
diff --git a/configure.ac b/configure.ac
index bf7efa0690c..748913f7216 100644
--- a/configure.ac
+++ b/configure.ac
@@ -2662,7 +2662,6 @@ AC_CHECK_FUNCS(\
exp2 \
exp2f \
fmaf \
- tgamma \
tgammaf
)
LIBS="$ac_save_LIBS"
diff --git a/dlls/msvcrt/math.c b/dlls/msvcrt/math.c
index a471f07f046..d419c34a947 100644
--- a/dlls/msvcrt/math.c
+++ b/dlls/msvcrt/math.c
@@ -8647,16 +8647,15 @@ end:
/* sin(pi*x) assuming x > 2^-100, if sin(pi*x)==0 the sign is arbitrary */
static double sin_pi(double x)
{
- static const double pi = 3.14159265358979311600e+00;
int n;
/* spurious inexact if odd int */
x = 2.0 * (x * 0.5 - floor(x * 0.5)); /* x mod 2.0 */
- n = (int)(x * 4.0);
+ n = x * 4.0;
n = (n + 1) / 2;
x -= n * 0.5f;
- x *= pi;
+ x *= M_PI;
switch (n) {
default: /* case 4: */
@@ -9047,12 +9046,127 @@ float CDECL lgammaf(float x)
return r;
}
+static double tgamma_S(double x)
+{
+ static const double Snum[] = {
+ 23531376880.410759688572007674451636754734846804940,
+ 42919803642.649098768957899047001988850926355848959,
+ 35711959237.355668049440185451547166705960488635843,
+ 17921034426.037209699919755754458931112671403265390,
+ 6039542586.3520280050642916443072979210699388420708,
+ 1439720407.3117216736632230727949123939715485786772,
+ 248874557.86205415651146038641322942321632125127801,
+ 31426415.585400194380614231628318205362874684987640,
+ 2876370.6289353724412254090516208496135991145378768,
+ 186056.26539522349504029498971604569928220784236328,
+ 8071.6720023658162106380029022722506138218516325024,
+ 210.82427775157934587250973392071336271166969580291,
+ 2.5066282746310002701649081771338373386264310793408,
+ };
+ static const double Sden[] = {
+ 0, 39916800, 120543840, 150917976, 105258076, 45995730, 13339535,
+ 2637558, 357423, 32670, 1925, 66, 1,
+ };
+
+ double num = 0, den = 0;
+ int i;
+
+ /* to avoid overflow handle large x differently */
+ if (x < 8)
+ for (i = ARRAY_SIZE(Snum) - 1; i >= 0; i--) {
+ num = num * x + Snum[i];
+ den = den * x + Sden[i];
+ }
+ else
+ for (i = 0; i < ARRAY_SIZE(Snum); i++) {
+ num = num / x + Snum[i];
+ den = den / x + Sden[i];
+ }
+ return num / den;
+}
+
/*********************************************************************
* tgamma (MSVCR120.@)
+ *
+ * Copied from musl: src/math/tgamma.c
*/
double CDECL tgamma(double x)
{
- return unix_funcs->tgamma( x );
+ static const double gmhalf = 5.524680040776729583740234375;
+ static const double fact[] = {
+ 1, 1, 2, 6, 24, 120, 720, 5040.0, 40320.0, 362880.0, 3628800.0, 39916800.0,
+ 479001600.0, 6227020800.0, 87178291200.0, 1307674368000.0, 20922789888000.0,
+ 355687428096000.0, 6402373705728000.0, 121645100408832000.0,
+ 2432902008176640000.0, 51090942171709440000.0, 1124000727777607680000.0,
+ };
+
+ union {double f; UINT64 i;} u = {x};
+ double absx, y, dy, z, r;
+ UINT32 ix = u.i >> 32 & 0x7fffffff;
+ int sign = u.i >> 63;
+
+ /* special cases */
+ if (ix >= 0x7ff00000) {
+ /* tgamma(nan)=nan, tgamma(inf)=inf, tgamma(-inf)=nan with invalid */
+ if (u.i == 0xfff0000000000000ULL)
+ *_errno() = EDOM;
+ return x + INFINITY;
+ }
+ if (ix < (0x3ff - 54) << 20) {
+ /* |x| < 2^-54: tgamma(x) ~ 1/x, +-0 raises div-by-zero */
+ if (x == 0.0)
+ *_errno() = ERANGE;
+ return 1 / x;
+ }
+
+ /* integer arguments */
+ /* raise inexact when non-integer */
+ if (x == floor(x)) {
+ if (sign) {
+ *_errno() = EDOM;
+ return 0 / (x - x);
+ }
+ if (x <= ARRAY_SIZE(fact))
+ return fact[(int)x - 1];
+ }
+
+ /* x >= 172: tgamma(x)=inf with overflow */
+ /* x =< -184: tgamma(x)=+-0 with underflow */
+ if (ix >= 0x40670000) { /* |x| >= 184 */
+ *_errno() = ERANGE;
+ if (sign) {
+ fp_barrierf(0x1p-126 / x);
+ return 0;
+ }
+ x *= 0x1p1023;
+ return x;
+ }
+
+ absx = sign ? -x : x;
+
+ /* handle the error of x + g - 0.5 */
+ y = absx + gmhalf;
+ if (absx > gmhalf) {
+ dy = y - absx;
+ dy -= gmhalf;
+ } else {
+ dy = y - gmhalf;
+ dy -= absx;
+ }
+
+ z = absx - 0.5;
+ r = tgamma_S(absx) * exp(-y);
+ if (x < 0) {
+ /* reflection formula for negative x */
+ /* sinpi(absx) is not 0, integers are already handled */
+ r = -M_PI / (sin_pi(absx) * absx * r);
+ dy = -dy;
+ z = -z;
+ }
+ r += dy * (gmhalf + 0.5) * r / y;
+ z = pow(y, 0.5 * z);
+ y = r * z * z;
+ return y;
}
/*********************************************************************
diff --git a/dlls/msvcrt/unixlib.c b/dlls/msvcrt/unixlib.c
index 012ec6f6b85..26af9812f88 100644
--- a/dlls/msvcrt/unixlib.c
+++ b/dlls/msvcrt/unixlib.c
@@ -110,19 +110,6 @@ static float CDECL unix_powf( float x, float y )
return powf( x, y );
}
-/*********************************************************************
- * tgamma
- */
-static double CDECL unix_tgamma(double x)
-{
-#ifdef HAVE_TGAMMA
- return tgamma(x);
-#else
- FIXME( "not implemented\n" );
- return 0;
-#endif
-}
-
/*********************************************************************
* tgammaf
*/
@@ -145,7 +132,6 @@ static const struct unix_funcs funcs =
unix_fmaf,
unix_pow,
unix_powf,
- unix_tgamma,
unix_tgammaf,
};
diff --git a/dlls/msvcrt/unixlib.h b/dlls/msvcrt/unixlib.h
index a9c569d6f63..f119a0a35d0 100644
--- a/dlls/msvcrt/unixlib.h
+++ b/dlls/msvcrt/unixlib.h
@@ -30,7 +30,6 @@ struct unix_funcs
float (CDECL *fmaf)(float x, float y, float z);
double (CDECL *pow)(double x, double y);
float (CDECL *powf)(float x, float y);
- double (CDECL *tgamma)(double x);
float (CDECL *tgammaf)(float x);
};
diff --git a/include/config.h.in b/include/config.h.in
index 537dda5f68d..aa74e710e06 100644
--- a/include/config.h.in
+++ b/include/config.h.in
@@ -929,9 +929,6 @@
/* Define to 1 if you have the <termios.h> header file. */
#undef HAVE_TERMIOS_H
-/* Define to 1 if you have the `tgamma' function. */
-#undef HAVE_TGAMMA
-
/* Define to 1 if you have the `tgammaf' function. */
#undef HAVE_TGAMMAF
More information about the wine-cvs
mailing list