diff options
| -rw-r--r-- | libgcompat/math.c | 1092 |
1 files changed, 1090 insertions, 2 deletions
diff --git a/libgcompat/math.c b/libgcompat/math.c index db123f8..83abbf4 100644 --- a/libgcompat/math.c +++ b/libgcompat/math.c @@ -1,6 +1,77 @@ -#include <math.h> /* isfinite, isinf, isnan */ +#define _GNU_SOURCE /* Extra maths functions */ +#include <math.h> /* Literally everything */ -#include "alias.h" /* weak_alias */ +#include "alias.h" /* weak_alias */ +#include "internal.h" /* GCOMPAT__assert_with_reason */ + +/** + * Multiplies the first argument x by FLT_RADIX (probably 2) to the power of y. + */ +long double scalbl(long double x, long double y) +{ + /* + * XXX strictly not correct but: + * 1) Good Enough(TM) + * 2) scalbl is deprecated anyway + * */ + return scalblnl(x, (long int)y); +} + +/* + * The below require support for ynl/jnl which doesn't exist in musl and isn't + * implemented in gcompat yet + */ +#if 0 +/** + * Return Bessel functions of x of the first kind of order n. + */ +long double jnl(int n, long double x) +{ + /* TODO implement */ + return 0; +} + +/** + * Return Bessel functions of x of the first kind of order 0. + */ +long double j0l(long double n) +{ + return jnl(0, n); +} + +/** + * Return Bessel functions of x of the first kind of order 1. + */ +long double j1l(long double n) +{ + return jnl(1, n); +} + +/** + * Return Bessel functions of x of the second kind of order n. + */ +long double ynl(int n, long double x) +{ + /* TODO implement */ + return 0; +} + +/** + * Return Bessel functions of x of the second kind of order 0. + */ +long double y0l(long double n) +{ + return ynl(0, n); +} + +/** + * Return Bessel functions of x of the second kind of order 1. + */ +long double y1l(long double n) +{ + return ynl(1, n); +} +#endif /** * Test for finite value. @@ -73,6 +144,7 @@ int __isnan(double arg) } weak_alias(__isnan, isnan); + /** * Test for a NaN. * @@ -94,3 +166,1019 @@ int __isnanl(long double arg) return isnan(arg); } weak_alias(__isnanl, isnanl); + + +/* + * Finite specialisations of functions used by glibc, that aren't supposed to + * return infinity. + */ + +#define _ASSERT_FINITE(finite_fn, res) \ + GCOMPAT__assert_with_reason(finite_fn(res), \ + "infinite value returned in a function that returns a " \ + "finite result"); + +#define ASSERT_FINITEF(res) _ASSERT_FINITE(isinff, res) +#define ASSERT_FINITE(res) _ASSERT_FINITE(isinf, res) +#define ASSERT_FINITEL(res) _ASSERT_FINITE(isinfl, res) + +/** + * Returns the principal value of the arc cosine of x, expressed in radians. + */ +float __acosf_finite(float x) +{ + float res = acosf(x); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Returns the principal value of the arc cosine of x, expressed in radians. + */ +double __acos_finite(double x) +{ + double res = acos(x); + + ASSERT_FINITE(res); + + return res; +} + +/** + * Returns the principal value of the arc cosine of x, expressed in radians. + */ +long double __acosl_finite(long double x) +{ + long double res = acosl(x); + + ASSERT_FINITEL(res); + + return res; +} + +/** + * Returns the nonnegative area hyperbolic cosine of x. + */ +double __acosh_finite(double x) +{ + double res = acosh(x); + + ASSERT_FINITE(res); + + return res; +} + +/** + * Returns the nonnegative area hyperbolic cosine of x. + */ +float __acoshf_finite(float x) +{ + float res = acoshf(x); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Returns the nonnegative area hyperbolic cosine of x. + */ +long double __acoshl_finite(long double x) +{ + long double res = acoshl(x); + + ASSERT_FINITEL(res); + + return res; +} + +/** + * Returns the principal value of the arc sine of x, expressed in radians. + */ +float __asinf_finite(float x) +{ + float res = asinf(x); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Returns the principal value of the arc sine of x, expressed in radians. + */ +double __asin_finite(double x) +{ + double res = asin(x); + + ASSERT_FINITE(res); + + return res; +} + +/** + * Returns the principal value of the arc sine of x, expressed in radians. + */ +long double __asinl_finite(long double x) +{ + long double res = asinl(x); + + ASSERT_FINITEL(res); + + return res; +} + + +/** + * Returns the principal value of the arc tangent of x/y, expressed in radians. + */ +float __atan2f_finite(float x, float y) +{ + float res = atan2f(x, y); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Returns the principal value of the arc tangent of x/y, expressed in radians. + */ +double __atan2_finite(double x, double y) +{ + double res = atan2(x, y); + + ASSERT_FINITE(res); + + return res; +} + +/** + * Returns the principal value of the arc tangent of x/y, expressed in radians. + */ +long double __atan2l_finite(long double x, long double y) +{ + long double res = atan2l(x, y); + + ASSERT_FINITEL(res); + + return res; +} + +/** + * Returns the area hyperbolic tangent of x. + */ +float __atanhf_finite(float x) +{ + float res = atanhf(x); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Returns the area hyperbolic tangent of x. + */ +double __atanh_finite(double x) +{ + double res = atanh(x); + + ASSERT_FINITE(res); + + return res; +} + +/** + * Returns the area hyperbolic tangent of x. + */ +long double __atanhl_finite(long double x) +{ + long double res = atanhl(x); + + ASSERT_FINITEL(res); + + return res; +} + + +/** + * Returns the hyperbolic cosine of x. + */ +float __coshf_finite(float x) +{ + float res = coshf(x); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Returns the hyperbolic cosine of x. + */ +double __cosh_finite(double x) +{ + double res = cosh(x); + + ASSERT_FINITE(res); + + return res; +} + +/** + * Returns the hyperbolic cosine of x. + */ +long double __coshl_finite(long double x) +{ + long double res = coshl(x); + + ASSERT_FINITEL(res); + + return res; +} + +/** + * Return the value of 10 raised to the power of x. + */ +float __exp10f_finite(float x) +{ + float res = exp10f(x); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Return the value of 10 raised to the power of x. + */ +double __exp10_finite(double x) +{ + double res = exp10(x); + + ASSERT_FINITE(res); + + return res; +} + +/** + * Return the value of 10 raised to the power of x. + */ +long double __exp10l_finite(long double x) +{ + long double res = exp10l(x); + + ASSERT_FINITEL(res); + + return res; +} + +/** + * Returns the base-2 exponential function of x, which is 2 raised to the power x + */ +float __exp2f_finite(float x) +{ + float res = exp2f(x); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Returns the base-2 exponential function of x, which is 2 raised to the power x + */ +double __exp2_finite(double x) +{ + double res = exp2(x); + + ASSERT_FINITE(res); + + return res; +} + +/** + * Returns the base-2 exponential function of x, which is 2 raised to the power x + */ +long double __exp2l_finite(long double x) +{ + long double res = exp2l(x); + + ASSERT_FINITEL(res); + + return res; +} + +/** + * Returns the base-e exponential function of x, which is e raised to the power x + */ +float __expf_finite(float x) +{ + float res = expf(x); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Returns the base-e exponential function of x, which is e raised to the power x + */ +double __exp_finite(double x) +{ + double res = exp(x); + + ASSERT_FINITE(res); + + return res; +} + +/** + * Returns the base-e exponential function of x, which is e raised to the power x + */ +long double __expl_finite(long double x) +{ + long double res = expl(x); + + ASSERT_FINITEL(res); + + return res; +} + +/** + * Returns the floating-point remainder of x/y (rounded towards zero) + */ +float __fmodf_finite(float x, float y) +{ + float res = fmodf(x, y); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Returns the floating-point remainder of x/y (rounded towards zero) + */ +double __fmod_finite(double x, double y) +{ + double res = fmod(x, y); + + ASSERT_FINITE(res); + + return res; +} + +/** + * Returns the floating-point remainder of x/y (rounded towards zero) + */ +long double __fmodl_finite(long double x, long double y) +{ + long double res = fmodl(x, y); + + ASSERT_FINITEL(res); + + return res; +} + +/** + * Computes the square root of the sum of the squares of x and y, without undue + * overflow or underflow at intermediate stages of the computation. + */ +float __hypotf_finite(float x, float y) +{ + float res = hypotf(x, y); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Computes the square root of the sum of the squares of x and y, without undue + * overflow or underflow at intermediate stages of the computation. + */ +double __hypot_finite(double x, double y) +{ + double res = hypot(x, y); + + ASSERT_FINITE(res); + + return res; +} + +/** + * Computes the square root of the sum of the squares of x and y, without undue + * overflow or underflow at intermediate stages of the computation. + */ +long double __hypotl_finite(long double x, long double y) +{ + long double res = hypotl(x, y); + + ASSERT_FINITEL(res); + + return res; +} + +/** + * Return Bessel functions of x of the first kind of orders 0. + */ +float __j0f_finite(float x) +{ + float res = j0f(x); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Return Bessel functions of x of the first kind of orders 0. + */ +double __j0_finite(double x) +{ + double res = j0(x); + + ASSERT_FINITE(res); + + return res; +} + +/* The below requires support for j0l, see above */ +#if 0 +/** + * Return Bessel functions of x of the first kind of orders 0. + */ +long double __j0l_finite(long double x) +{ + long double res = j0l(x); + + ASSERT_FINITEL(res); + + return res; +} +#endif + +/** + * Return Bessel functions of x of the first kind of orders 1. + */ +float __j1f_finite(float x) +{ + float res = j1f(x); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Return Bessel functions of x of the first kind of orders 1. + */ +double __j1_finite(double x) +{ + double res = j1(x); + + ASSERT_FINITE(res); + + return res; +} + +/* The below requires support for j1l, see above */ +#if 0 +/** + * Return Bessel functions of x of the first kind of orders 1. + */ +long double __j1l_finite(long double x) +{ + long double res = j1l(x); + + ASSERT_FINITEL(res); + + return res; +} +#endif + +/** + * Return the Bessel function of x of the first kind of order n. + */ +float __jnf_finite(int n, float x) +{ + float res = jnf(n, x); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Return the Bessel function of x of the first kind of order n. + */ +double __jn_finite(int n, double x) +{ + double res = jn(n, x); + + ASSERT_FINITE(res); + + return res; +} + +/* The below requires support for jnl, see above */ +#if 0 +/** + * Return the Bessel function of x of the first kind of order n. + */ +long double __jnl_finite(int n, long double x) +{ + long double res = jnl(n, x); + + ASSERT_FINITEL(res); + + return res; +} +#endif + +/** + * Returns the natural logarithm of the absolute value of the Gamma function. + */ +float __lgammaf_finite(float x) +{ + float res = lgammaf(x); + + ASSERT_FINITEF(res); + + return res; +} +alias(__lgammaf_finite, __gammaf_finite); + +/** + * Returns the natural logarithm of the absolute value of the Gamma function. + */ +double __lgamma_finite(double x) +{ + double res = lgamma(x); + + ASSERT_FINITE(res); + + return res; +} +alias(__lgamma_finite, __gamma_finite); + +/** + * Returns the natural logarithm of the absolute value of the Gamma function. + */ +long double __lgammal_finite(long double x) +{ + long double res = lgammal(x); + + ASSERT_FINITEL(res); + + return res; +} +alias(__lgammal_finite, __gammal_finite); + +/** + * Returns the natural logarithm of the absolute value of the Gamma function. + */ +float __lgammaf_r_finite(float x, int *p) +{ + float res = lgammaf_r(x, p); + + ASSERT_FINITEF(res); + + return res; +} +alias(__lgammaf_r_finite, __gammaf_r_finite); + +/** + * Returns the natural logarithm of the absolute value of the Gamma function. + */ +double __lgamma_r_finite(double x, int *p) +{ + double res = lgamma_r(x, p); + + ASSERT_FINITE(res); + + return res; +} +alias(__lgamma_r_finite, __gamma_r_finite); + +/** + * Returns the natural logarithm of the absolute value of the Gamma function. + */ +long double __lgammal_r_finite(long double x, int *p) +{ + long double res = lgammal_r(x, p); + + ASSERT_FINITEL(res); + + return res; +} +alias(__lgammal_r_finite, __gammal_r_finite); + +/** + * Returns the common (base-10) logarithm of x. + */ +float __log10f_finite(float x) +{ + float res = log10f(x); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Returns the common (base-10) logarithm of x. + */ +double __log10_finite(double x) +{ + double res = log10(x); + + ASSERT_FINITE(res); + + return res; +} + +/** + * Returns the common (base-10) logarithm of x. + */ +long double __log10l_finite(long double x) +{ + long double res = log10l(x); + + ASSERT_FINITEL(res); + + return res; +} + +/** + * Returns the binary (base-2) logarithm of x. + */ +float __log2f_finite(float x) +{ + float res = log2f(x); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Returns the binary (base-2) logarithm of x. + */ +double __log2_finite(double x) +{ + double res = log2(x); + + ASSERT_FINITE(res); + + return res; +} + +/** + * Returns the binary (base-2) logarithm of x. + */ +long double __log2l_finite(long double x) +{ + long double res = log2l(x); + + ASSERT_FINITEL(res); + + return res; +} + +/** + * Returns the natural logarithm of x. + */ +float __logf_finite(float x) +{ + float res = logf(x); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Returns the natural logarithm of x. + */ +double __log_finite(double x) +{ + double res = log(x); + + ASSERT_FINITE(res); + + return res; +} + +/** + * Returns the natural logarithm of x. + */ +long double __logl_finite(long double x) +{ + long double res = logl(x); + + ASSERT_FINITEL(res); + + return res; +} + +/** + * Returns x raised to the y exponent. + */ +float __powf_finite(float x, float y) +{ + float res = powf(x, y); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Returns x raised to the y exponent. + */ +double __pow_finite(double x, double y) +{ + double res = pow(x, y); + + ASSERT_FINITE(res); + + return res; +} + +/** + * Returns x raised to the y exponent. + */ +long double __powl_finite(long double x, long double y) +{ + long double res = powl(x, y); + + ASSERT_FINITEL(res); + + return res; +} + +/** + * Returns the floating-point remainder of x/y (rounded to nearest). + */ +float __remainderf_finite(float x, float y) +{ + float res = remainderf(x, y); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Returns the floating-point remainder of x/y (rounded to nearest). + */ +double __remainder_finite(double x, double y) +{ + double res = remainder(x, y); + + ASSERT_FINITE(res); + + return res; +} + +/** + * Returns the floating-point remainder of x/y (rounded to nearest). + */ +long double __remainderl_finite(long double x, long double y) +{ + long double res = remainderl(x, y); + + ASSERT_FINITEL(res); + + return res; +} + +/** + * Multiplies the first argument x by FLT_RADIX (probably 2) to the power of y. + */ +float __scalbf_finite(float x, float y) +{ + float res = scalbf(x, y); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Multiplies the first argument x by FLT_RADIX (probably 2) to the power of y. + */ +double __scalb_finite(double x, double y) +{ + double res = scalb(x, y); + + ASSERT_FINITE(res); + + return res; +} + +/** + * Multiplies the first argument x by FLT_RADIX (probably 2) to the power of y. + */ +long double __scalbl_finite(long double x, long double y) +{ + long double res = scalbl(x, y); + + ASSERT_FINITEL(res); + + return res; +} + +/** + * Returns the hyperbolic sine of x. + */ +float __sinhf_finite(float x) +{ + float res = sinhf(x); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Returns the hyperbolic sine of x. + */ +double __sinh_finite(double x) +{ + double res = sinh(x); + + ASSERT_FINITE(res); + + return res; +} + +/** + * Returns the hyperbolic sine of x. + */ +long double __sinhl_finite(long double x) +{ + long double res = sinhl(x); + + ASSERT_FINITEL(res); + + return res; +} + +/** + * Returns the square root of x. + */ +float __sqrtf_finite(float x) +{ + float res = sqrtf(x); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Returns the square root of x. + */ +double __sqrt_finite(double x) +{ + double res = sqrt(x); + + ASSERT_FINITE(res); + + return res; +} + +/** + * Returns the square root of x. + */ +long double __sqrtl_finite(long double x) +{ + long double res = sqrtl(x); + + ASSERT_FINITEL(res); + + return res; +} + +/** + * Return Bessel functions of x of the second kind of order 0. + */ +float __y0f_finite(float x) +{ + float res = y0f(x); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Return Bessel functions of x of the second kind of order 0. + */ +double __y0_finite(double x) +{ + double res = y0(x); + + ASSERT_FINITE(res); + + return res; +} + +/* The below requires support for y0l, see above */ +#if 0 +/** + * Return Bessel functions of x of the second kind of order 0. + */ +long double __y0l_finite(long double x) +{ + long double res = y0l(x); + + ASSERT_FINITEL(res); + + return res; +} +#endif + +/** + * Return Bessel functions of x of the second kind of order 1. + */ +float __y1f_finite(float x) +{ + float res = y1f(x); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Return Bessel functions of x of the second kind of order 1. + */ +double __y1_finite(double x) +{ + double res = y1(x); + + ASSERT_FINITE(res); + + return res; +} + +/* The below requires support for y1l, see above */ +#if 0 +/** + * Return Bessel functions of x of the second kind of order 1. + */ +long double __y1l_finite(long double x) +{ + long double res = y1l(x); + + ASSERT_FINITEL(res); + + return res; +} +#endif + +/** + * Return Bessel functions of x of the second kind of order n. + */ +float __ynf_finite(int n, float x) +{ + float res = ynf(n, x); + + ASSERT_FINITEF(res); + + return res; +} + +/** + * Return Bessel functions of x of the second kind of order n. + */ +double __yn_finite(int n, double x) +{ + double res = yn(n, x); + + ASSERT_FINITE(res); + + return res; +} + +/* The below requires support for ynl, see above */ +#if 0 +/** + * Return Bessel functions of x of the second kind of order n. + */ +long double __ynl_finite(int n, long double x) +{ + long double res = ynl(n, x); + + ASSERT_FINITEL(res); + + return res; +} +#endif |