summaryrefslogtreecommitdiffstats
path: root/modules/fdlibm/patches
diff options
context:
space:
mode:
Diffstat (limited to 'modules/fdlibm/patches')
-rw-r--r--modules/fdlibm/patches/01_remove_unused_declarations_from_fdlibm_h.patch511
-rw-r--r--modules/fdlibm/patches/02_change_include_guard_in_fdlibm_h.patch37
-rw-r--r--modules/fdlibm/patches/03_put_fdlibm_functions_into_fdlibm_namespace.patch36
-rw-r--r--modules/fdlibm/patches/04_include_fdlibm_h_from_math_private_h.patch1214
-rw-r--r--modules/fdlibm/patches/05_include_stdint_h_in_math_private_h.patch21
-rw-r--r--modules/fdlibm/patches/06_use_mfbt_endian_h_in_math_private_h.patch116
-rw-r--r--modules/fdlibm/patches/07_add_fdlibm_namespace_to_functions_defined_and_used_in_fdlibm.patch71
-rw-r--r--modules/fdlibm/patches/08_remove_weak_reference_macro.patch459
-rw-r--r--modules/fdlibm/patches/09_comment_out_rcsid_variable.patch1362
-rw-r--r--modules/fdlibm/patches/10_remove_unused_function_from_k_exp_cpp.patch107
-rw-r--r--modules/fdlibm/patches/11_include_cfloat_to_use_flt_eval_method.patch21
-rw-r--r--modules/fdlibm/patches/12_define_u_int32_t_and_u_int64_t_on_windows.patch27
-rw-r--r--modules/fdlibm/patches/13_define_strict_assign_even_if_flt_eval_method_is_not_defined.patch31
-rw-r--r--modules/fdlibm/patches/14_do_not_use_hexadecimal_floating_point_number.patch68
-rw-r--r--modules/fdlibm/patches/15_remove_unused_rintl_function_from_s_nearbyint_cpp.patch13
-rw-r--r--modules/fdlibm/patches/16_use_safer_strict_assign_on_visual_studio.patch22
-rw-r--r--modules/fdlibm/patches/17_exp_exact_result_for_positive_one.patch40
-rw-r--r--modules/fdlibm/patches/18_use_stdlib_sqrt.patch289
-rw-r--r--modules/fdlibm/patches/19_remove_unneeded_round_to_integer_helpers.patch50
-rw-r--r--modules/fdlibm/patches/20_emulate_freebsd_internal_double_types.patch28
-rw-r--r--modules/fdlibm/patches/21_rem_pio_fixups.patch226
-rw-r--r--modules/fdlibm/patches/22_add_float_types.patch66
22 files changed, 4815 insertions, 0 deletions
diff --git a/modules/fdlibm/patches/01_remove_unused_declarations_from_fdlibm_h.patch b/modules/fdlibm/patches/01_remove_unused_declarations_from_fdlibm_h.patch
new file mode 100644
index 0000000000..6365f71488
--- /dev/null
+++ b/modules/fdlibm/patches/01_remove_unused_declarations_from_fdlibm_h.patch
@@ -0,0 +1,511 @@
+diff --git a/fdlibm.h b/fdlibm.h
+--- a/fdlibm.h
++++ b/fdlibm.h
+@@ -12,504 +12,66 @@
+ /*
+ * from: @(#)fdlibm.h 5.1 93/09/24
+ * $FreeBSD$
+ */
+
+ #ifndef _MATH_H_
+ #define _MATH_H_
+
+-#include <sys/cdefs.h>
+-#include <sys/_types.h>
+-#include <machine/_limits.h>
+-
+-/*
+- * ANSI/POSIX
+- */
+-extern const union __infinity_un {
+- unsigned char __uc[8];
+- double __ud;
+-} __infinity;
+-
+-extern const union __nan_un {
+- unsigned char __uc[sizeof(float)];
+- float __uf;
+-} __nan;
+-
+-#if __GNUC_PREREQ__(3, 3)
+-#define __MATH_BUILTIN_CONSTANTS
+-#endif
+-
+-#if __GNUC_PREREQ__(3, 0)
+-#define __MATH_BUILTIN_RELOPS
+-#endif
+-
+-#ifdef __MATH_BUILTIN_CONSTANTS
+-#define HUGE_VAL __builtin_huge_val()
+-#else
+-#define HUGE_VAL (__infinity.__ud)
+-#endif
+-
+-#if __ISO_C_VISIBLE >= 1999
+-#define FP_ILOGB0 (-__INT_MAX)
+-#define FP_ILOGBNAN __INT_MAX
+-
+-#ifdef __MATH_BUILTIN_CONSTANTS
+-#define HUGE_VALF __builtin_huge_valf()
+-#define HUGE_VALL __builtin_huge_vall()
+-#define INFINITY __builtin_inff()
+-#define NAN __builtin_nanf("")
+-#else
+-#define HUGE_VALF (float)HUGE_VAL
+-#define HUGE_VALL (long double)HUGE_VAL
+-#define INFINITY HUGE_VALF
+-#define NAN (__nan.__uf)
+-#endif /* __MATH_BUILTIN_CONSTANTS */
+-
+-#define MATH_ERRNO 1
+-#define MATH_ERREXCEPT 2
+-#define math_errhandling MATH_ERREXCEPT
+-
+-#define FP_FAST_FMAF 1
+-
+-/* Symbolic constants to classify floating point numbers. */
+-#define FP_INFINITE 0x01
+-#define FP_NAN 0x02
+-#define FP_NORMAL 0x04
+-#define FP_SUBNORMAL 0x08
+-#define FP_ZERO 0x10
+-
+-#if __STDC_VERSION__ >= 201112L || __has_extension(c_generic_selections)
+-#define __fp_type_select(x, f, d, ld) __extension__ _Generic((x), \
+- float: f(x), \
+- double: d(x), \
+- long double: ld(x))
+-#elif __GNUC_PREREQ__(3, 1) && !defined(__cplusplus)
+-#define __fp_type_select(x, f, d, ld) __builtin_choose_expr( \
+- __builtin_types_compatible_p(__typeof(x), long double), ld(x), \
+- __builtin_choose_expr( \
+- __builtin_types_compatible_p(__typeof(x), double), d(x), \
+- __builtin_choose_expr( \
+- __builtin_types_compatible_p(__typeof(x), float), f(x), (void)0)))
+-#else
+-#define __fp_type_select(x, f, d, ld) \
+- ((sizeof(x) == sizeof(float)) ? f(x) \
+- : (sizeof(x) == sizeof(double)) ? d(x) \
+- : ld(x))
+-#endif
+-
+-#define fpclassify(x) \
+- __fp_type_select(x, __fpclassifyf, __fpclassifyd, __fpclassifyl)
+-#define isfinite(x) __fp_type_select(x, __isfinitef, __isfinite, __isfinitel)
+-#define isinf(x) __fp_type_select(x, __isinff, __isinf, __isinfl)
+-#define isnan(x) \
+- __fp_type_select(x, __inline_isnanf, __inline_isnan, __inline_isnanl)
+-#define isnormal(x) __fp_type_select(x, __isnormalf, __isnormal, __isnormall)
+-
+-#ifdef __MATH_BUILTIN_RELOPS
+-#define isgreater(x, y) __builtin_isgreater((x), (y))
+-#define isgreaterequal(x, y) __builtin_isgreaterequal((x), (y))
+-#define isless(x, y) __builtin_isless((x), (y))
+-#define islessequal(x, y) __builtin_islessequal((x), (y))
+-#define islessgreater(x, y) __builtin_islessgreater((x), (y))
+-#define isunordered(x, y) __builtin_isunordered((x), (y))
+-#else
+-#define isgreater(x, y) (!isunordered((x), (y)) && (x) > (y))
+-#define isgreaterequal(x, y) (!isunordered((x), (y)) && (x) >= (y))
+-#define isless(x, y) (!isunordered((x), (y)) && (x) < (y))
+-#define islessequal(x, y) (!isunordered((x), (y)) && (x) <= (y))
+-#define islessgreater(x, y) (!isunordered((x), (y)) && \
+- ((x) > (y) || (y) > (x)))
+-#define isunordered(x, y) (isnan(x) || isnan(y))
+-#endif /* __MATH_BUILTIN_RELOPS */
+-
+-#define signbit(x) __fp_type_select(x, __signbitf, __signbit, __signbitl)
+-
+-typedef __double_t double_t;
+-typedef __float_t float_t;
+-#endif /* __ISO_C_VISIBLE >= 1999 */
+-
+-/*
+- * XOPEN/SVID
+- */
+-#if __BSD_VISIBLE || __XSI_VISIBLE
+-#define M_E 2.7182818284590452354 /* e */
+-#define M_LOG2E 1.4426950408889634074 /* log 2e */
+-#define M_LOG10E 0.43429448190325182765 /* log 10e */
+-#define M_LN2 0.69314718055994530942 /* log e2 */
+-#define M_LN10 2.30258509299404568402 /* log e10 */
+-#define M_PI 3.14159265358979323846 /* pi */
+ #define M_PI_2 1.57079632679489661923 /* pi/2 */
+-#define M_PI_4 0.78539816339744830962 /* pi/4 */
+-#define M_1_PI 0.31830988618379067154 /* 1/pi */
+-#define M_2_PI 0.63661977236758134308 /* 2/pi */
+-#define M_2_SQRTPI 1.12837916709551257390 /* 2/sqrt(pi) */
+-#define M_SQRT2 1.41421356237309504880 /* sqrt(2) */
+-#define M_SQRT1_2 0.70710678118654752440 /* 1/sqrt(2) */
+-
+-#define MAXFLOAT ((float)3.40282346638528860e+38)
+-extern int signgam;
+-#endif /* __BSD_VISIBLE || __XSI_VISIBLE */
+-
+-#if __BSD_VISIBLE
+-#if 0
+-/* Old value from 4.4BSD-Lite math.h; this is probably better. */
+-#define HUGE HUGE_VAL
+-#else
+-#define HUGE MAXFLOAT
+-#endif
+-#endif /* __BSD_VISIBLE */
+-
+-/*
+- * Most of these functions depend on the rounding mode and have the side
+- * effect of raising floating-point exceptions, so they are not declared
+- * as __pure2. In C99, FENV_ACCESS affects the purity of these functions.
+- */
+-__BEGIN_DECLS
+-/*
+- * ANSI/POSIX
+- */
+-int __fpclassifyd(double) __pure2;
+-int __fpclassifyf(float) __pure2;
+-int __fpclassifyl(long double) __pure2;
+-int __isfinitef(float) __pure2;
+-int __isfinite(double) __pure2;
+-int __isfinitel(long double) __pure2;
+-int __isinff(float) __pure2;
+-int __isinf(double) __pure2;
+-int __isinfl(long double) __pure2;
+-int __isnormalf(float) __pure2;
+-int __isnormal(double) __pure2;
+-int __isnormall(long double) __pure2;
+-int __signbit(double) __pure2;
+-int __signbitf(float) __pure2;
+-int __signbitl(long double) __pure2;
+-
+-static __inline int
+-__inline_isnan(__const double __x)
+-{
+-
+- return (__x != __x);
+-}
+-
+-static __inline int
+-__inline_isnanf(__const float __x)
+-{
+-
+- return (__x != __x);
+-}
+-
+-static __inline int
+-__inline_isnanl(__const long double __x)
+-{
+-
+- return (__x != __x);
+-}
+-
+-/*
+- * Define the following aliases, for compatibility with glibc and CUDA.
+- */
+-#define __isnan __inline_isnan
+-#define __isnanf __inline_isnanf
+-
+-/*
+- * Version 2 of the Single UNIX Specification (UNIX98) defined isnan() and
+- * isinf() as functions taking double. C99, and the subsequent POSIX revisions
+- * (SUSv3, POSIX.1-2001, define it as a macro that accepts any real floating
+- * point type. If we are targeting SUSv2 and C99 or C11 (or C++11) then we
+- * expose the newer definition, assuming that the language spec takes
+- * precedence over the operating system interface spec.
+- */
+-#if __XSI_VISIBLE > 0 && __XSI_VISIBLE < 600 && __ISO_C_VISIBLE < 1999
+-#undef isinf
+-#undef isnan
+-int isinf(double);
+-int isnan(double);
+-#endif
+-
+ double acos(double);
+ double asin(double);
+ double atan(double);
+ double atan2(double, double);
+ double cos(double);
+ double sin(double);
+ double tan(double);
+
+ double cosh(double);
+ double sinh(double);
+ double tanh(double);
+
+ double exp(double);
+-double frexp(double, int *); /* fundamentally !__pure2 */
+-double ldexp(double, int);
+ double log(double);
+ double log10(double);
+-double modf(double, double *); /* fundamentally !__pure2 */
+
+ double pow(double, double);
+-double sqrt(double);
+
+ double ceil(double);
+-double fabs(double) __pure2;
++double fabs(double);
+ double floor(double);
+-double fmod(double, double);
+
+-/*
+- * These functions are not in C90.
+- */
+-#if __BSD_VISIBLE || __ISO_C_VISIBLE >= 1999 || __XSI_VISIBLE
+ double acosh(double);
+ double asinh(double);
+ double atanh(double);
+ double cbrt(double);
+-double erf(double);
+-double erfc(double);
+ double exp2(double);
+ double expm1(double);
+-double fma(double, double, double);
+ double hypot(double, double);
+-int ilogb(double) __pure2;
+-double lgamma(double);
+-long long llrint(double);
+-long long llround(double);
+ double log1p(double);
+ double log2(double);
+-double logb(double);
+-long lrint(double);
+-long lround(double);
+-double nan(const char *) __pure2;
+-double nextafter(double, double);
+-double remainder(double, double);
+-double remquo(double, double, int *);
+ double rint(double);
+-#endif /* __BSD_VISIBLE || __ISO_C_VISIBLE >= 1999 || __XSI_VISIBLE */
+-
+-#if __BSD_VISIBLE || __XSI_VISIBLE
+-double j0(double);
+-double j1(double);
+-double jn(int, double);
+-double y0(double);
+-double y1(double);
+-double yn(int, double);
+-
+-#if __XSI_VISIBLE <= 500 || __BSD_VISIBLE
+-double gamma(double);
+-#endif
+-
+-#if __XSI_VISIBLE <= 600 || __BSD_VISIBLE
+-double scalb(double, double);
+-#endif
+-#endif /* __BSD_VISIBLE || __XSI_VISIBLE */
+-
+-#if __BSD_VISIBLE || __ISO_C_VISIBLE >= 1999
+-double copysign(double, double) __pure2;
+-double fdim(double, double);
+-double fmax(double, double) __pure2;
+-double fmin(double, double) __pure2;
++double copysign(double, double);
+ double nearbyint(double);
+-double round(double);
+-double scalbln(double, long);
+ double scalbn(double, int);
+-double tgamma(double);
+ double trunc(double);
+-#endif
+-
+-/*
+- * BSD math library entry points
+- */
+-#if __BSD_VISIBLE
+-double drem(double, double);
+-int finite(double) __pure2;
+-int isnanf(float) __pure2;
+-
+-/*
+- * Reentrant version of gamma & lgamma; passes signgam back by reference
+- * as the second argument; user must allocate space for signgam.
+- */
+-double gamma_r(double, int *);
+-double lgamma_r(double, int *);
+-
+-/*
+- * IEEE Test Vector
+- */
+-double significand(double);
+-#endif /* __BSD_VISIBLE */
+-
+-/* float versions of ANSI/POSIX functions */
+-#if __ISO_C_VISIBLE >= 1999
+ float acosf(float);
+ float asinf(float);
+ float atanf(float);
+-float atan2f(float, float);
+ float cosf(float);
+ float sinf(float);
+ float tanf(float);
+-
+-float coshf(float);
+-float sinhf(float);
+-float tanhf(float);
+-
+ float exp2f(float);
+ float expf(float);
+-float expm1f(float);
+-float frexpf(float, int *); /* fundamentally !__pure2 */
+-int ilogbf(float) __pure2;
+-float ldexpf(float, int);
+-float log10f(float);
+-float log1pf(float);
+-float log2f(float);
+ float logf(float);
+-float modff(float, float *); /* fundamentally !__pure2 */
+-
+ float powf(float, float);
+ float sqrtf(float);
+
+ float ceilf(float);
+-float fabsf(float) __pure2;
++float fabsf(float);
+ float floorf(float);
+-float fmodf(float, float);
+-float roundf(float);
+-
+-float erff(float);
+-float erfcf(float);
+-float hypotf(float, float);
+-float lgammaf(float);
+-float tgammaf(float);
+-
+-float acoshf(float);
+-float asinhf(float);
+-float atanhf(float);
+-float cbrtf(float);
+-float logbf(float);
+-float copysignf(float, float) __pure2;
+-long long llrintf(float);
+-long long llroundf(float);
+-long lrintf(float);
+-long lroundf(float);
+-float nanf(const char *) __pure2;
+ float nearbyintf(float);
+-float nextafterf(float, float);
+-float remainderf(float, float);
+-float remquof(float, float, int *);
+ float rintf(float);
+-float scalblnf(float, long);
+ float scalbnf(float, int);
+ float truncf(float);
+
+-float fdimf(float, float);
+-float fmaf(float, float, float);
+-float fmaxf(float, float) __pure2;
+-float fminf(float, float) __pure2;
+-#endif
+-
+-/*
+- * float versions of BSD math library entry points
+- */
+-#if __BSD_VISIBLE
+-float dremf(float, float);
+-int finitef(float) __pure2;
+-float gammaf(float);
+-float j0f(float);
+-float j1f(float);
+-float jnf(int, float);
+-float scalbf(float, float);
+-float y0f(float);
+-float y1f(float);
+-float ynf(int, float);
+-
+-/*
+- * Float versions of reentrant version of gamma & lgamma; passes
+- * signgam back by reference as the second argument; user must
+- * allocate space for signgam.
+- */
+-float gammaf_r(float, int *);
+-float lgammaf_r(float, int *);
+-
+-/*
+- * float version of IEEE Test Vector
+- */
+-float significandf(float);
+-#endif /* __BSD_VISIBLE */
+-
+-/*
+- * long double versions of ISO/POSIX math functions
+- */
+-#if __ISO_C_VISIBLE >= 1999
+-long double acoshl(long double);
+-long double acosl(long double);
+-long double asinhl(long double);
+-long double asinl(long double);
+-long double atan2l(long double, long double);
+-long double atanhl(long double);
+-long double atanl(long double);
+-long double cbrtl(long double);
+-long double ceill(long double);
+-long double copysignl(long double, long double) __pure2;
+-long double coshl(long double);
+-long double cosl(long double);
+-long double erfcl(long double);
+-long double erfl(long double);
+-long double exp2l(long double);
+-long double expl(long double);
+-long double expm1l(long double);
+-long double fabsl(long double) __pure2;
+-long double fdiml(long double, long double);
+-long double floorl(long double);
+-long double fmal(long double, long double, long double);
+-long double fmaxl(long double, long double) __pure2;
+-long double fminl(long double, long double) __pure2;
+-long double fmodl(long double, long double);
+-long double frexpl(long double, int *); /* fundamentally !__pure2 */
+-long double hypotl(long double, long double);
+-int ilogbl(long double) __pure2;
+-long double ldexpl(long double, int);
+-long double lgammal(long double);
+-long long llrintl(long double);
+-long long llroundl(long double);
+-long double log10l(long double);
+-long double log1pl(long double);
+-long double log2l(long double);
+-long double logbl(long double);
+-long double logl(long double);
+-long lrintl(long double);
+-long lroundl(long double);
+-long double modfl(long double, long double *); /* fundamentally !__pure2 */
+-long double nanl(const char *) __pure2;
+-long double nearbyintl(long double);
+-long double nextafterl(long double, long double);
+-double nexttoward(double, long double);
+-float nexttowardf(float, long double);
+-long double nexttowardl(long double, long double);
+-long double powl(long double, long double);
+-long double remainderl(long double, long double);
+-long double remquol(long double, long double, int *);
+-long double rintl(long double);
+-long double roundl(long double);
+-long double scalblnl(long double, long);
+-long double scalbnl(long double, int);
+-long double sinhl(long double);
+-long double sinl(long double);
+-long double sqrtl(long double);
+-long double tanhl(long double);
+-long double tanl(long double);
+-long double tgammal(long double);
+-long double truncl(long double);
+-#endif /* __ISO_C_VISIBLE >= 1999 */
+-
+-#if __BSD_VISIBLE
+-long double lgammal_r(long double, int *);
+-void sincos(double, double *, double *);
+-void sincosf(float, float *, float *);
+-void sincosl(long double, long double *, long double *);
+-double cospi(double);
+-float cospif(float);
+-long double cospil(long double);
+-double sinpi(double);
+-float sinpif(float);
+-long double sinpil(long double);
+-double tanpi(double);
+-float tanpif(float);
+-long double tanpil(long double);
+-#endif
+-
+-__END_DECLS
+-
+ #endif /* !_MATH_H_ */
diff --git a/modules/fdlibm/patches/02_change_include_guard_in_fdlibm_h.patch b/modules/fdlibm/patches/02_change_include_guard_in_fdlibm_h.patch
new file mode 100644
index 0000000000..e89eadbe32
--- /dev/null
+++ b/modules/fdlibm/patches/02_change_include_guard_in_fdlibm_h.patch
@@ -0,0 +1,37 @@
+--- a/fdlibm.h 2022-12-13 14:27:29.092805897 -0500
++++ b/fdlibm.h 2022-12-13 14:27:37.444836016 -0500
+@@ -9,18 +9,21 @@
+ * ====================================================
+ */
+
+ /*
+ * from: @(#)fdlibm.h 5.1 93/09/24
+ * $FreeBSD$
+ */
+
+-#ifndef _MATH_H_
+-#define _MATH_H_
++#ifndef mozilla_imported_fdlibm_h
++#define mozilla_imported_fdlibm_h
+
++#ifndef M_PI_2
+ #define M_PI_2 1.57079632679489661923 /* pi/2 */
++#endif
++
+ double acos(double);
+ double asin(double);
+ double atan(double);
+ double atan2(double, double);
+ double cos(double);
+ double sin(double);
+@@ -69,9 +72,9 @@
+ float ceilf(float);
+ float fabsf(float);
+ float floorf(float);
+ float nearbyintf(float);
+ float rintf(float);
+ float scalbnf(float, int);
+ float truncf(float);
+
+-#endif /* !_MATH_H_ */
++#endif /* !mozilla_imported_fdlibm_h */
diff --git a/modules/fdlibm/patches/03_put_fdlibm_functions_into_fdlibm_namespace.patch b/modules/fdlibm/patches/03_put_fdlibm_functions_into_fdlibm_namespace.patch
new file mode 100644
index 0000000000..2a2b147381
--- /dev/null
+++ b/modules/fdlibm/patches/03_put_fdlibm_functions_into_fdlibm_namespace.patch
@@ -0,0 +1,36 @@
+--- a/fdlibm.h 2022-12-13 14:28:36.393051585 -0500
++++ b/fdlibm.h 2022-12-13 14:28:47.533092878 -0500
+@@ -12,19 +12,21 @@
+ /*
+ * from: @(#)fdlibm.h 5.1 93/09/24
+ * $FreeBSD$
+ */
+
+ #ifndef mozilla_imported_fdlibm_h
+ #define mozilla_imported_fdlibm_h
+
++namespace fdlibm {
++
+ #ifndef M_PI_2
+ #define M_PI_2 1.57079632679489661923 /* pi/2 */
+ #endif
+
+ double acos(double);
+ double asin(double);
+ double atan(double);
+ double atan2(double, double);
+ double cos(double);
+ double sin(double);
+ double tan(double);
+@@ -69,9 +71,11 @@
+ float ceilf(float);
+ float fabsf(float);
+ float floorf(float);
+ float nearbyintf(float);
+ float rintf(float);
+ float scalbnf(float, int);
+ float truncf(float);
+
++} /* namespace fdlibm */
++
+ #endif /* !mozilla_imported_fdlibm_h */
diff --git a/modules/fdlibm/patches/04_include_fdlibm_h_from_math_private_h.patch b/modules/fdlibm/patches/04_include_fdlibm_h_from_math_private_h.patch
new file mode 100644
index 0000000000..8531976f89
--- /dev/null
+++ b/modules/fdlibm/patches/04_include_fdlibm_h_from_math_private_h.patch
@@ -0,0 +1,1214 @@
+diff --git a/e_acos.cpp b/e_acos.cpp
+--- a/e_acos.cpp
++++ b/e_acos.cpp
+@@ -35,17 +35,16 @@ __FBSDID("$FreeBSD$");
+ * if x is NaN, return x itself;
+ * if |x|>1, return NaN with invalid signal.
+ *
+ * Function needed: sqrt
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const double
+ one= 1.00000000000000000000e+00, /* 0x3FF00000, 0x00000000 */
+ pi = 3.14159265358979311600e+00, /* 0x400921FB, 0x54442D18 */
+ pio2_hi = 1.57079632679489655800e+00; /* 0x3FF921FB, 0x54442D18 */
+ static volatile double
+ pio2_lo = 6.12323399573676603587e-17; /* 0x3C91A626, 0x33145C07 */
+diff --git a/e_acosh.cpp b/e_acosh.cpp
+--- a/e_acosh.cpp
++++ b/e_acosh.cpp
+@@ -26,17 +26,16 @@ __FBSDID("$FreeBSD$");
+ *
+ * Special cases:
+ * acosh(x) is NaN with signal if x<1.
+ * acosh(NaN) is NaN without signal.
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const double
+ one = 1.0,
+ ln2 = 6.93147180559945286227e-01; /* 0x3FE62E42, 0xFEFA39EF */
+
+ double
+ __ieee754_acosh(double x)
+diff --git a/e_asin.cpp b/e_asin.cpp
+--- a/e_asin.cpp
++++ b/e_asin.cpp
+@@ -41,17 +41,16 @@ __FBSDID("$FreeBSD$");
+ * Special cases:
+ * if x is NaN, return x itself;
+ * if |x|>1, return NaN with invalid signal.
+ *
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const double
+ one = 1.00000000000000000000e+00, /* 0x3FF00000, 0x00000000 */
+ huge = 1.000e+300,
+ pio2_hi = 1.57079632679489655800e+00, /* 0x3FF921FB, 0x54442D18 */
+ pio2_lo = 6.12323399573676603587e-17, /* 0x3C91A626, 0x33145C07 */
+ pio4_hi = 7.85398163397448278999e-01, /* 0x3FE921FB, 0x54442D18 */
+diff --git a/e_atan2.cpp b/e_atan2.cpp
+--- a/e_atan2.cpp
++++ b/e_atan2.cpp
+@@ -39,17 +39,16 @@ __FBSDID("$FreeBSD$");
+ * The hexadecimal values are the intended ones for the following
+ * constants. The decimal values may be used, provided that the
+ * compiler will convert from decimal to binary accurately enough
+ * to produce the hexadecimal values shown.
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static volatile double
+ tiny = 1.0e-300;
+ static const double
+ zero = 0.0,
+ pi_o_4 = 7.8539816339744827900E-01, /* 0x3FE921FB, 0x54442D18 */
+ pi_o_2 = 1.5707963267948965580E+00, /* 0x3FF921FB, 0x54442D18 */
+diff --git a/e_atanh.cpp b/e_atanh.cpp
+--- a/e_atanh.cpp
++++ b/e_atanh.cpp
+@@ -30,17 +30,16 @@ __FBSDID("$FreeBSD$");
+ * atanh(x) is NaN if |x| > 1 with signal;
+ * atanh(NaN) is that NaN with no signal;
+ * atanh(+-1) is +-INF with signal.
+ *
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const double one = 1.0, huge = 1e300;
+ static const double zero = 0.0;
+
+ double
+ __ieee754_atanh(double x)
+ {
+diff --git a/e_cosh.cpp b/e_cosh.cpp
+--- a/e_cosh.cpp
++++ b/e_cosh.cpp
+@@ -32,17 +32,16 @@ __FBSDID("$FreeBSD$");
+ *
+ * Special cases:
+ * cosh(x) is |x| if x is +INF, -INF, or NaN.
+ * only cosh(0)=1 is exact for finite x.
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const double one = 1.0, half=0.5, huge = 1.0e300;
+
+ double
+ __ieee754_cosh(double x)
+ {
+ double t,w;
+diff --git a/e_exp.cpp b/e_exp.cpp
+--- a/e_exp.cpp
++++ b/e_exp.cpp
+@@ -73,17 +73,16 @@ __FBSDID("$FreeBSD$");
+ * The hexadecimal values are the intended ones for the following
+ * constants. The decimal values may be used, provided that the
+ * compiler will convert from decimal to binary accurately enough
+ * to produce the hexadecimal values shown.
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const double
+ one = 1.0,
+ halF[2] = {0.5,-0.5,},
+ o_threshold= 7.09782712893383973096e+02, /* 0x40862E42, 0xFEFA39EF */
+ u_threshold= -7.45133219101941108420e+02, /* 0xc0874910, 0xD52D3051 */
+ ln2HI[2] ={ 6.93147180369123816490e-01, /* 0x3fe62e42, 0xfee00000 */
+diff --git a/e_hypot.cpp b/e_hypot.cpp
+--- a/e_hypot.cpp
++++ b/e_hypot.cpp
+@@ -43,17 +43,16 @@ __FBSDID("$FreeBSD$");
+ *
+ * Accuracy:
+ * hypot(x,y) returns sqrt(x^2+y^2) with error less
+ * than 1 ulps (units in the last place)
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ double
+ __ieee754_hypot(double x, double y)
+ {
+ double a,b,t1,t2,y1,y2,w;
+ int32_t j,k,ha,hb;
+
+diff --git a/e_log.cpp b/e_log.cpp
+--- a/e_log.cpp
++++ b/e_log.cpp
+@@ -62,17 +62,16 @@ __FBSDID("$FreeBSD$");
+ * The hexadecimal values are the intended ones for the following
+ * constants. The decimal values may be used, provided that the
+ * compiler will convert from decimal to binary accurately enough
+ * to produce the hexadecimal values shown.
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const double
+ ln2_hi = 6.93147180369123816490e-01, /* 3fe62e42 fee00000 */
+ ln2_lo = 1.90821492927058770002e-10, /* 3dea39ef 35793c76 */
+ two54 = 1.80143985094819840000e+16, /* 43500000 00000000 */
+ Lg1 = 6.666666666666735130e-01, /* 3FE55555 55555593 */
+ Lg2 = 3.999999999940941908e-01, /* 3FD99999 9997FA04 */
+diff --git a/e_log10.cpp b/e_log10.cpp
+--- a/e_log10.cpp
++++ b/e_log10.cpp
+@@ -19,17 +19,16 @@ __FBSDID("$FreeBSD$");
+ * comments.
+ *
+ * log10(x) = (f - 0.5*f*f + k_log1p(f)) / ln10 + k * log10(2)
+ * in not-quite-routine extra precision.
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+ #include "k_log.h"
+
+ static const double
+ two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */
+ ivln10hi = 4.34294481878168880939e-01, /* 0x3fdbcb7b, 0x15200000 */
+ ivln10lo = 2.50829467116452752298e-11, /* 0x3dbb9438, 0xca9aadd5 */
+ log10_2hi = 3.01029995663611771306e-01, /* 0x3FD34413, 0x509F6000 */
+diff --git a/e_log2.cpp b/e_log2.cpp
+--- a/e_log2.cpp
++++ b/e_log2.cpp
+@@ -21,17 +21,16 @@ __FBSDID("$FreeBSD$");
+ * This reduces x to {k, 1+f} exactly as in e_log.c, then calls the kernel,
+ * then does the combining and scaling steps
+ * log2(x) = (f - 0.5*f*f + k_log1p(f)) / ln2 + k
+ * in not-quite-routine extra precision.
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+ #include "k_log.h"
+
+ static const double
+ two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */
+ ivln2hi = 1.44269504072144627571e+00, /* 0x3ff71547, 0x65200000 */
+ ivln2lo = 1.67517131648865118353e-10; /* 0x3de705fc, 0x2eefa200 */
+
+diff --git a/e_pow.cpp b/e_pow.cpp
+--- a/e_pow.cpp
++++ b/e_pow.cpp
+@@ -53,17 +53,16 @@ __FBSDID("$FreeBSD$");
+ * Constants :
+ * The hexadecimal values are the intended ones for the following
+ * constants. The decimal values may be used, provided that the
+ * compiler will convert from decimal to binary accurately enough
+ * to produce the hexadecimal values shown.
+ */
+
+ #include <float.h>
+-#include "math.h"
+ #include "math_private.h"
+
+ static const double
+ bp[] = {1.0, 1.5,},
+ dp_h[] = { 0.0, 5.84962487220764160156e-01,}, /* 0x3FE2B803, 0x40000000 */
+ dp_l[] = { 0.0, 1.35003920212974897128e-08,}, /* 0x3E4CFDEB, 0x43CFD006 */
+ zero = 0.0,
+ half = 0.5,
+diff --git a/e_rem_pio2.cpp b/e_rem_pio2.cpp
+--- a/e_rem_pio2.cpp
++++ b/e_rem_pio2.cpp
+@@ -19,17 +19,16 @@ __FBSDID("$FreeBSD$");
+ /* __ieee754_rem_pio2(x,y)
+ *
+ * return the remainder of x rem pi/2 in y[0]+y[1]
+ * use __kernel_rem_pio2()
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ /*
+ * invpio2: 53 bits of 2/pi
+ * pio2_1: first 33 bit of pi/2
+ * pio2_1t: pi/2 - pio2_1
+ * pio2_2: second 33 bit of pi/2
+ * pio2_2t: pi/2 - (pio2_1+pio2_2)
+diff --git a/e_sinh.cpp b/e_sinh.cpp
+--- a/e_sinh.cpp
++++ b/e_sinh.cpp
+@@ -29,17 +29,16 @@ __FBSDID("$FreeBSD$");
+ *
+ * Special cases:
+ * sinh(x) is |x| if x is +INF, -INF, or NaN.
+ * only sinh(0)=0 is exact for finite x.
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const double one = 1.0, shuge = 1.0e307;
+
+ double
+ __ieee754_sinh(double x)
+ {
+ double t,h;
+diff --git a/k_cos.cpp b/k_cos.cpp
+--- a/k_cos.cpp
++++ b/k_cos.cpp
+@@ -48,17 +48,16 @@ __FBSDID("$FreeBSD$");
+ * and tmp having the same precision as x. If they have extra
+ * precision due to compiler bugs, then the extra precision is
+ * only good provided it is retained in all terms of the final
+ * expression for cos(). Retention happens in all cases tested
+ * under FreeBSD, so don't pessimize things by forcibly clipping
+ * any extra precision in w.
+ */
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const double
+ one = 1.00000000000000000000e+00, /* 0x3FF00000, 0x00000000 */
+ C1 = 4.16666666666666019037e-02, /* 0x3FA55555, 0x5555554C */
+ C2 = -1.38888888888741095749e-03, /* 0xBF56C16C, 0x16C15177 */
+ C3 = 2.48015872894767294178e-05, /* 0x3EFA01A0, 0x19CB1590 */
+ C4 = -2.75573143513906633035e-07, /* 0xBE927E4F, 0x809C52AD */
+diff --git a/k_exp.cpp b/k_exp.cpp
+--- a/k_exp.cpp
++++ b/k_exp.cpp
+@@ -26,17 +26,16 @@
+ * SUCH DAMAGE.
+ */
+
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+
+ #include <complex.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const uint32_t k = 1799; /* constant for reduction */
+ static const double kln2 = 1246.97177782734161156; /* k * ln2 */
+
+ /*
+ * Compute exp(x), scaled to avoid spurious overflow. An exponent is
+ * returned separately in 'expt'.
+diff --git a/k_rem_pio2.cpp b/k_rem_pio2.cpp
+--- a/k_rem_pio2.cpp
++++ b/k_rem_pio2.cpp
+@@ -126,17 +126,16 @@ __FBSDID("$FreeBSD$");
+ * The hexadecimal values are the intended ones for the following
+ * constants. The decimal values may be used, provided that the
+ * compiler will convert from decimal to binary accurately enough
+ * to produce the hexadecimal values shown.
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const int init_jk[] = {3,4,4,6}; /* initial value for jk */
+
+ /*
+ * Table of constants for 2/pi, 396 Hex digits (476 decimal) of 2/pi
+ *
+ * integer array, contains the (24*i)-th to (24*i+23)-th
+diff --git a/k_sin.cpp b/k_sin.cpp
+--- a/k_sin.cpp
++++ b/k_sin.cpp
+@@ -39,17 +39,16 @@ __FBSDID("$FreeBSD$");
+ * ~ sin(x) + (1-x*x/2)*y
+ * For better accuracy, let
+ * 3 2 2 2 2
+ * r = x *(S2+x *(S3+x *(S4+x *(S5+x *S6))))
+ * then 3 2
+ * sin(x) = x + (S1*x + (x *(r-y/2)+y))
+ */
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const double
+ half = 5.00000000000000000000e-01, /* 0x3FE00000, 0x00000000 */
+ S1 = -1.66666666666666324348e-01, /* 0xBFC55555, 0x55555549 */
+ S2 = 8.33333333332248946124e-03, /* 0x3F811111, 0x1110F8A6 */
+ S3 = -1.98412698298579493134e-04, /* 0xBF2A01A0, 0x19C161D5 */
+ S4 = 2.75573137070700676789e-06, /* 0x3EC71DE3, 0x57B1FE7D */
+diff --git a/k_tan.cpp b/k_tan.cpp
+--- a/k_tan.cpp
++++ b/k_tan.cpp
+@@ -44,17 +44,16 @@ __FBSDID("$FreeBSD$");
+ * 3 2
+ * tan(x+y) = x + (T1*x + (x *(r+y)+y))
+ *
+ * 4. For x in [0.67434,pi/4], let y = pi/4 - x, then
+ * tan(x) = tan(pi/4-y) = (1-tan(y))/(1+tan(y))
+ * = 1 - 2*(tan(y) - (tan(y)^2)/(1+tan(y)))
+ */
+
+-#include "math.h"
+ #include "math_private.h"
+ static const double xxx[] = {
+ 3.33333333333334091986e-01, /* 3FD55555, 55555563 */
+ 1.33333333333201242699e-01, /* 3FC11111, 1110FE7A */
+ 5.39682539762260521377e-02, /* 3FABA1BA, 1BB341FE */
+ 2.18694882948595424599e-02, /* 3F9664F4, 8406D637 */
+ 8.86323982359930005737e-03, /* 3F8226E3, E96E8493 */
+ 3.59207910759131235356e-03, /* 3F6D6D22, C9560328 */
+diff --git a/math_private.h b/math_private.h
+--- a/math_private.h
++++ b/math_private.h
+@@ -15,16 +15,18 @@
+ */
+
+ #ifndef _MATH_PRIVATE_H_
+ #define _MATH_PRIVATE_H_
+
+ #include <sys/types.h>
+ #include <machine/endian.h>
+
++#include "fdlibm.h"
++
+ /*
+ * The original fdlibm code used statements like:
+ * n0 = ((*(int*)&one)>>29)^1; * index of high word *
+ * ix0 = *(n0+(int*)&x); * high word of x *
+ * ix1 = *((1-n0)+(int*)&x); * low word of x *
+ * to dig two 32 bit words out of the 64 bit IEEE floating point
+ * value. That is non-ANSI, and, moreover, the gcc instruction
+ * scheduler gets it wrong. We instead use the following macros.
+diff --git a/s_asinh.cpp b/s_asinh.cpp
+--- a/s_asinh.cpp
++++ b/s_asinh.cpp
+@@ -21,17 +21,16 @@ __FBSDID("$FreeBSD$");
+ * asinh(x) := x if 1+x*x=1,
+ * := sign(x)*(log(x)+ln2)) for large |x|, else
+ * := sign(x)*log(2|x|+1/(|x|+sqrt(x*x+1))) if|x|>2, else
+ * := sign(x)*log1p(|x| + x^2/(1 + sqrt(1+x^2)))
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const double
+ one = 1.00000000000000000000e+00, /* 0x3FF00000, 0x00000000 */
+ ln2 = 6.93147180559945286227e-01, /* 0x3FE62E42, 0xFEFA39EF */
+ huge= 1.00000000000000000000e+300;
+
+ double
+diff --git a/s_atan.cpp b/s_atan.cpp
+--- a/s_atan.cpp
++++ b/s_atan.cpp
+@@ -30,17 +30,16 @@ __FBSDID("$FreeBSD$");
+ * The hexadecimal values are the intended ones for the following
+ * constants. The decimal values may be used, provided that the
+ * compiler will convert from decimal to binary accurately enough
+ * to produce the hexadecimal values shown.
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const double atanhi[] = {
+ 4.63647609000806093515e-01, /* atan(0.5)hi 0x3FDDAC67, 0x0561BB4F */
+ 7.85398163397448278999e-01, /* atan(1.0)hi 0x3FE921FB, 0x54442D18 */
+ 9.82793723247329054082e-01, /* atan(1.5)hi 0x3FEF730B, 0xD281F69B */
+ 1.57079632679489655800e+00, /* atan(inf)hi 0x3FF921FB, 0x54442D18 */
+ };
+diff --git a/s_cbrt.cpp b/s_cbrt.cpp
+--- a/s_cbrt.cpp
++++ b/s_cbrt.cpp
+@@ -11,17 +11,16 @@
+ *
+ * Optimized by Bruce D. Evans.
+ */
+
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+
+ #include <float.h>
+-#include "math.h"
+ #include "math_private.h"
+
+ /* cbrt(x)
+ * Return cube root of x
+ */
+ static const u_int32_t
+ B1 = 715094163, /* B1 = (1023-1023/3-0.03306235651)*2**20 */
+ B2 = 696219795; /* B2 = (1023-1023/3-54/3-0.03306235651)*2**20 */
+diff --git a/s_ceil.cpp b/s_ceil.cpp
+--- a/s_ceil.cpp
++++ b/s_ceil.cpp
+@@ -19,17 +19,16 @@ __FBSDID("$FreeBSD$");
+ * Method:
+ * Bit twiddling.
+ * Exception:
+ * Inexact flag raised if x not equal to ceil(x).
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const double huge = 1.0e300;
+
+ double
+ ceil(double x)
+ {
+ int32_t i0,i1,j0;
+diff --git a/s_ceilf.cpp b/s_ceilf.cpp
+--- a/s_ceilf.cpp
++++ b/s_ceilf.cpp
+@@ -11,17 +11,16 @@
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const float huge = 1.0e30;
+
+ float
+ ceilf(float x)
+ {
+ int32_t i0,j0;
+diff --git a/s_copysign.cpp b/s_copysign.cpp
+--- a/s_copysign.cpp
++++ b/s_copysign.cpp
+@@ -14,17 +14,16 @@
+ __FBSDID("$FreeBSD$");
+
+ /*
+ * copysign(double x, double y)
+ * copysign(x,y) returns a value with the magnitude of x and
+ * with the sign bit of y.
+ */
+
+-#include "math.h"
+ #include "math_private.h"
+
+ double
+ copysign(double x, double y)
+ {
+ u_int32_t hx,hy;
+ GET_HIGH_WORD(hx,x);
+ GET_HIGH_WORD(hy,y);
+diff --git a/s_cos.cpp b/s_cos.cpp
+--- a/s_cos.cpp
++++ b/s_cos.cpp
+@@ -41,17 +41,16 @@ __FBSDID("$FreeBSD$");
+ * trig(NaN) is that NaN;
+ *
+ * Accuracy:
+ * TRIG(x) returns trig(x) nearly rounded
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #define INLINE_REM_PIO2
+ #include "math_private.h"
+ #include "e_rem_pio2.c"
+
+ double
+ cos(double x)
+ {
+ double y[2],z=0.0;
+diff --git a/s_expm1.cpp b/s_expm1.cpp
+--- a/s_expm1.cpp
++++ b/s_expm1.cpp
+@@ -105,17 +105,16 @@ __FBSDID("$FreeBSD$");
+ * The hexadecimal values are the intended ones for the following
+ * constants. The decimal values may be used, provided that the
+ * compiler will convert from decimal to binary accurately enough
+ * to produce the hexadecimal values shown.
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const double
+ one = 1.0,
+ tiny = 1.0e-300,
+ o_threshold = 7.09782712893383973096e+02,/* 0x40862E42, 0xFEFA39EF */
+ ln2_hi = 6.93147180369123816490e-01,/* 0x3fe62e42, 0xfee00000 */
+ ln2_lo = 1.90821492927058770002e-10,/* 0x3dea39ef, 0x35793c76 */
+diff --git a/s_fabs.cpp b/s_fabs.cpp
+--- a/s_fabs.cpp
++++ b/s_fabs.cpp
+@@ -12,17 +12,16 @@
+
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+
+ /*
+ * fabs(x) returns the absolute value of x.
+ */
+
+-#include "math.h"
+ #include "math_private.h"
+
+ double
+ fabs(double x)
+ {
+ u_int32_t high;
+ GET_HIGH_WORD(high,x);
+ SET_HIGH_WORD(x,high&0x7fffffff);
+diff --git a/s_floor.cpp b/s_floor.cpp
+--- a/s_floor.cpp
++++ b/s_floor.cpp
+@@ -19,17 +19,16 @@ __FBSDID("$FreeBSD$");
+ * Method:
+ * Bit twiddling.
+ * Exception:
+ * Inexact flag raised if x not equal to floor(x).
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const double huge = 1.0e300;
+
+ double
+ floor(double x)
+ {
+ int32_t i0,i1,j0;
+diff --git a/s_floorf.cpp b/s_floorf.cpp
+--- a/s_floorf.cpp
++++ b/s_floorf.cpp
+@@ -20,17 +20,16 @@ __FBSDID("$FreeBSD$");
+ * floorf(x)
+ * Return x rounded toward -inf to integral value
+ * Method:
+ * Bit twiddling.
+ * Exception:
+ * Inexact flag raised if x not equal to floorf(x).
+ */
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const float huge = 1.0e30;
+
+ float
+ floorf(float x)
+ {
+ int32_t i0,j0;
+diff --git a/s_log1p.cpp b/s_log1p.cpp
+--- a/s_log1p.cpp
++++ b/s_log1p.cpp
+@@ -75,17 +75,16 @@ __FBSDID("$FreeBSD$");
+ * if(u==1.0) return x ; else
+ * return log(u)*(x/(u-1.0));
+ *
+ * See HP-15C Advanced Functions Handbook, p.193.
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const double
+ ln2_hi = 6.93147180369123816490e-01, /* 3fe62e42 fee00000 */
+ ln2_lo = 1.90821492927058770002e-10, /* 3dea39ef 35793c76 */
+ two54 = 1.80143985094819840000e+16, /* 43500000 00000000 */
+ Lp1 = 6.666666666666735130e-01, /* 3FE55555 55555593 */
+ Lp2 = 3.999999999940941908e-01, /* 3FD99999 9997FA04 */
+diff --git a/s_nearbyint.cpp b/s_nearbyint.cpp
+--- a/s_nearbyint.cpp
++++ b/s_nearbyint.cpp
+@@ -25,17 +25,18 @@
+ * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
+ * SUCH DAMAGE.
+ */
+
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+
+ #include <fenv.h>
+-#include <math.h>
++
++#include "math_private.h"
+
+ /*
+ * We save and restore the floating-point environment to avoid raising
+ * an inexact exception. We can get away with using fesetenv()
+ * instead of feclearexcept()/feupdateenv() to restore the environment
+ * because the only exception defined for rint() is overflow, and
+ * rounding can't overflow as long as emax >= p.
+ *
+diff --git a/s_rint.cpp b/s_rint.cpp
+--- a/s_rint.cpp
++++ b/s_rint.cpp
+@@ -20,17 +20,16 @@ __FBSDID("$FreeBSD$");
+ * Method:
+ * Using floating addition.
+ * Exception:
+ * Inexact flag raised if x not equal to rint(x).
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const double
+ TWO52[2]={
+ 4.50359962737049600000e+15, /* 0x43300000, 0x00000000 */
+ -4.50359962737049600000e+15, /* 0xC3300000, 0x00000000 */
+ };
+
+diff --git a/s_rintf.cpp b/s_rintf.cpp
+--- a/s_rintf.cpp
++++ b/s_rintf.cpp
+@@ -14,17 +14,16 @@
+ */
+
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+
+ #include <float.h>
+ #include <stdint.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const float
+ TWO23[2]={
+ 8.3886080000e+06, /* 0x4b000000 */
+ -8.3886080000e+06, /* 0xcb000000 */
+ };
+
+diff --git a/s_scalbn.cpp b/s_scalbn.cpp
+--- a/s_scalbn.cpp
++++ b/s_scalbn.cpp
+@@ -2,19 +2,20 @@
+ * Copyright (c) 2005-2020 Rich Felker, et al.
+ *
+ * SPDX-License-Identifier: MIT
+ *
+ * Please see https://git.musl-libc.org/cgit/musl/tree/COPYRIGHT
+ * for all contributors to musl.
+ */
+ #include <float.h>
+-#include <math.h>
+ #include <stdint.h>
+
++#include "math_private.h"
++
+ double scalbn(double x, int n)
+ {
+ union {double f; uint64_t i;} u;
+ double_t y = x;
+
+ if (n > 1023) {
+ y *= 0x1p1023;
+ n -= 1023;
+diff --git a/s_scalbnf.cpp b/s_scalbnf.cpp
+--- a/s_scalbnf.cpp
++++ b/s_scalbnf.cpp
+@@ -8,16 +4,17 @@
+ * SPDX-License-Identifier: MIT
+ *
+ * Please see https://git.musl-libc.org/cgit/musl/tree/COPYRIGHT
+ * for all contributors to musl.
+ */
+-#include <math.h>
+ #include <stdint.h>
+
++#include "math_private.h"
++
+ float scalbnf(float x, int n)
+ {
+ union {float f; uint32_t i;} u;
+ float_t y = x;
+
+ if (n > 127) {
+ y *= 0x1p127f;
+ n -= 127;
+diff --git a/s_sin.cpp b/s_sin.cpp
+--- a/s_sin.cpp
++++ b/s_sin.cpp
+@@ -41,17 +41,16 @@ __FBSDID("$FreeBSD$");
+ * trig(NaN) is that NaN;
+ *
+ * Accuracy:
+ * TRIG(x) returns trig(x) nearly rounded
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #define INLINE_REM_PIO2
+ #include "math_private.h"
+ #include "e_rem_pio2.c"
+
+ double
+ sin(double x)
+ {
+ double y[2],z=0.0;
+diff --git a/s_tan.cpp b/s_tan.cpp
+--- a/s_tan.cpp
++++ b/s_tan.cpp
+@@ -40,17 +40,16 @@ __FBSDID("$FreeBSD$");
+ * trig(NaN) is that NaN;
+ *
+ * Accuracy:
+ * TRIG(x) returns trig(x) nearly rounded
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #define INLINE_REM_PIO2
+ #include "math_private.h"
+ #include "e_rem_pio2.c"
+
+ double
+ tan(double x)
+ {
+ double y[2],z=0.0;
+diff --git a/s_tanh.cpp b/s_tanh.cpp
+--- a/s_tanh.cpp
++++ b/s_tanh.cpp
+@@ -34,17 +34,16 @@ __FBSDID("$FreeBSD$");
+ *
+ * Special cases:
+ * tanh(NaN) is NaN;
+ * only tanh(0)=0 is exact for finite argument.
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const volatile double tiny = 1.0e-300;
+ static const double one = 1.0, two = 2.0, huge = 1.0e300;
+
+ double
+ tanh(double x)
+ {
+diff --git a/s_trunc.cpp b/s_trunc.cpp
+--- a/s_trunc.cpp
++++ b/s_trunc.cpp
+@@ -19,17 +19,16 @@ __FBSDID("$FreeBSD$");
+ * Method:
+ * Bit twiddling.
+ * Exception:
+ * Inexact flag raised if x not equal to trunc(x).
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const double huge = 1.0e300;
+
+ double
+ trunc(double x)
+ {
+ int32_t i0,i1,j0;
+diff --git a/s_truncf.cpp b/s_truncf.cpp
+--- a/s_truncf.cpp
++++ b/s_truncf.cpp
+@@ -17,17 +17,16 @@ __FBSDID("$FreeBSD$");
+ * truncf(x)
+ * Return x rounded toward 0 to integral value
+ * Method:
+ * Bit twiddling.
+ * Exception:
+ * Inexact flag raised if x not equal to truncf(x).
+ */
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const float huge = 1.0e30F;
+
+ float
+ truncf(float x)
+ {
+ int32_t i0,j0;
+--- a/e_acosf.cpp 2022-12-11 21:15:50.058431053 -0500
++++ b/e_acosf.cpp_new 2022-12-11 21:13:10.532515637 -0500
+@@ -12,15 +12,14 @@
+ * is preserved.
+ * ====================================================
+ */
+
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const float
+ one = 1.0000000000e+00, /* 0x3F800000 */
+ pi = 3.1415925026e+00, /* 0x40490fda */
+ pio2_hi = 1.5707962513e+00; /* 0x3fc90fda */
+ static volatile float
+--- a/e_asinf.cpp 2022-12-11 21:15:50.486437124 -0500
++++ b/e_asinf.cpp_new 2022-12-11 21:13:07.876498396 -0500
+@@ -12,15 +12,14 @@
+ * is preserved.
+ * ====================================================
+ */
+
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const float
+ one = 1.0000000000e+00, /* 0x3F800000 */
+ huge = 1.000e+30,
+ /* coefficient for R(x^2) */
+ pS0 = 1.6666586697e-01,
+--- a/e_expf.cpp 2022-12-11 21:15:51.250447960 -0500
++++ b/e_expf.cpp_new 2022-12-11 21:13:05.212481099 -0500
+@@ -14,15 +14,14 @@
+ */
+
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const float
+ one = 1.0,
+ halF[2] = {0.5,-0.5,},
+ o_threshold= 8.8721679688e+01, /* 0x42b17180 */
+ u_threshold= -1.0397208405e+02, /* 0xc2cff1b5 */
+--- a/e_logf.cpp 2022-12-11 21:15:51.978458278 -0500
++++ b/e_logf.cpp_new 2022-12-11 21:13:02.644464418 -0500
+@@ -12,15 +12,14 @@
+ * is preserved.
+ * ====================================================
+ */
+
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const float
+ ln2_hi = 6.9313812256e-01, /* 0x3f317180 */
+ ln2_lo = 9.0580006145e-06, /* 0x3717f7d1 */
+ two25 = 3.355443200e+07, /* 0x4c000000 */
+ /* |(log(1+s)-log(1-s))/s - Lg(s)| < 2**-34.24 (~[-4.95e-11, 4.97e-11]). */
+--- a/e_powf.cpp 2022-12-11 21:15:52.242462016 -0500
++++ b/e_powf.cpp_new 2022-12-11 21:13:00.164448306 -0500
+@@ -12,15 +12,14 @@
+ * is preserved.
+ * ====================================================
+ */
+
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const float
+ bp[] = {1.0, 1.5,},
+ dp_h[] = { 0.0, 5.84960938e-01,}, /* 0x3f15c000 */
+ dp_l[] = { 0.0, 1.56322085e-06,}, /* 0x35d1cfdc */
+ zero = 0.0,
+--- a/e_rem_pio2f.cpp 2022-12-11 21:15:52.626467454 -0500
++++ b/e_rem_pio2f.cpp_new 2022-12-11 21:12:57.300429693 -0500
+@@ -22,15 +22,14 @@
+ * return the remainder of x rem pi/2 in *y
+ * use double precision for everything except passing x
+ * use __kernel_rem_pio2() for large x
+ */
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ /*
+ * invpio2: 53 bits of 2/pi
+ * pio2_1: first 25 bits of pi/2
+ * pio2_1t: pi/2 - pio2_1
+ */
+--- a/e_sqrtf.cpp 2022-12-11 21:15:52.922471645 -0500
++++ b/e_sqrtf.cpp_new 2022-12-11 21:12:54.596412113 -0500
+@@ -13,15 +13,14 @@
+ * ====================================================
+ */
+
+ #ifndef lint
+ static char rcsid[] = "$FreeBSD$";
+ #endif
+
+-#include "math.h"
+ #include "math_private.h"
+
+ #ifdef USE_BUILTIN_SQRTF
+ float
+ __ieee754_sqrtf(float x)
+ {
+ return (__builtin_sqrtf(x));
+--- a/k_cosf.cpp 2022-12-11 21:15:53.178475268 -0500
++++ b/k_cosf.cpp_new 2022-12-11 21:12:52.012395309 -0500
+@@ -15,15 +15,14 @@
+ */
+
+ #ifndef INLINE_KERNEL_COSDF
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+ #endif
+
+-#include "math.h"
+ #include "math_private.h"
+
+ /* |cos(x) - c(x)| < 2**-34.1 (~[-5.37e-11, 5.295e-11]). */
+ static const double
+ one = 1.0,
+ C0 = -0x1ffffffd0c5e81.0p-54, /* -0.499999997251031003120 */
+ C1 = 0x155553e1053a42.0p-57, /* 0.0416666233237390631894 */
+--- a/k_expf.cpp 2022-12-11 21:15:53.434478888 -0500
++++ b/k_expf.cpp_new 2022-12-11 21:12:49.012375792 -0500
+@@ -27,15 +27,14 @@
+ */
+
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+
+ #include <complex.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const uint32_t k = 235; /* constant for reduction */
+ static const float kln2 = 162.88958740F; /* k * ln2 */
+
+ /*
+ * See k_exp.c for details.
+--- a/k_sinf.cpp 2022-12-11 21:15:53.946486130 -0500
++++ b/k_sinf.cpp_new 2022-12-11 21:12:46.020356322 -0500
+@@ -15,15 +15,14 @@
+ */
+
+ #ifndef INLINE_KERNEL_SINDF
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+ #endif
+
+-#include "math.h"
+ #include "math_private.h"
+
+ /* |sin(x)/x - s(x)| < 2**-37.5 (~[-4.89e-12, 4.824e-12]). */
+ static const double
+ S1 = -0x15555554cbac77.0p-55, /* -0.166666666416265235595 */
+ S2 = 0x111110896efbb2.0p-59, /* 0.0083333293858894631756 */
+ S3 = -0x1a00f9e2cae774.0p-65, /* -0.000198393348360966317347 */
+--- a/k_tanf.cpp 2022-12-11 21:15:54.254490484 -0500
++++ b/k_tanf.cpp_new 2022-12-11 21:12:42.972336479 -0500
+@@ -14,15 +14,14 @@
+ */
+
+ #ifndef INLINE_KERNEL_TANDF
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+ #endif
+
+-#include "math.h"
+ #include "math_private.h"
+
+ /* |tan(x)/x - t(x)| < 2**-25.5 (~[-2e-08, 2e-08]). */
+ static const double
+ T[] = {
+ 0x15554d3418c99f.0p-54, /* 0.333331395030791399758 */
+ 0x1112fd38999f72.0p-55, /* 0.133392002712976742718 */
+--- a/s_atanf.cpp 2022-12-11 21:15:54.670496364 -0500
++++ b/s_atanf.cpp_new 2022-12-11 21:12:39.940316733 -0500
+@@ -12,15 +12,14 @@
+ * is preserved.
+ * ====================================================
+ */
+
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+
+-#include "math.h"
+ #include "math_private.h"
+
+ static const float atanhi[] = {
+ 4.6364760399e-01, /* atan(0.5)hi 0x3eed6338 */
+ 7.8539812565e-01, /* atan(1.0)hi 0x3f490fda */
+ 9.8279368877e-01, /* atan(1.5)hi 0x3f7b985e */
+ 1.5707962513e+00, /* atan(inf)hi 0x3fc90fda */
+--- a/s_cosf.cpp 2022-12-11 21:15:55.494508001 -0500
++++ b/s_cosf.cpp_new 2022-12-11 21:12:37.396300160 -0500
+@@ -15,15 +15,14 @@
+ */
+
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+
+ #include <float.h>
+
+-#include "math.h"
+ #define INLINE_KERNEL_COSDF
+ #define INLINE_KERNEL_SINDF
+ #define INLINE_REM_PIO2F
+ #include "math_private.h"
+ #include "e_rem_pio2f.c"
+ #include "k_cosf.c"
+ #include "k_sinf.c"
+--- a/s_exp2.cpp 2022-12-11 21:15:55.790512181 -0500
++++ b/s_exp2.cpp_new 2022-12-11 21:12:34.476281131 -0500
+@@ -27,15 +27,14 @@
+ */
+
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ #define TBLBITS 8
+ #define TBLSIZE (1 << TBLBITS)
+
+ static const double
+ redux = 0x1.8p52 / TBLSIZE,
+--- a/s_exp2f.cpp 2022-12-11 21:15:55.926514101 -0500
++++ b/s_exp2f.cpp_new 2022-12-11 21:12:30.644256149 -0500
+@@ -27,15 +27,14 @@
+ */
+
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+
+ #include <float.h>
+
+-#include "math.h"
+ #include "math_private.h"
+
+ #define TBLBITS 4
+ #define TBLSIZE (1 << TBLBITS)
+
+ static const float
+ redux = 0x1.8p23f / TBLSIZE,
+--- a/s_fabsf.cpp 2022-12-11 21:15:56.202517995 -0500
++++ b/s_fabsf.cpp_new 2022-12-11 21:12:27.244233973 -0500
+@@ -16,15 +16,14 @@
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+
+ /*
+ * fabsf(x) returns the absolute value of x.
+ */
+
+-#include "math.h"
+ #include "math_private.h"
+
+ float
+ fabsf(float x)
+ {
+ u_int32_t ix;
+ GET_FLOAT_WORD(ix,x);
+--- a/s_sinf.cpp 2022-12-11 21:15:57.618537960 -0500
++++ b/s_sinf.cpp_new 2022-12-11 21:12:24.684217270 -0500
+@@ -15,15 +15,14 @@
+ */
+
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+
+ #include <float.h>
+
+-#include "math.h"
+ #define INLINE_KERNEL_COSDF
+ #define INLINE_KERNEL_SINDF
+ #define INLINE_REM_PIO2F
+ #include "math_private.h"
+ #include "e_rem_pio2f.c"
+ #include "k_cosf.c"
+ #include "k_sinf.c"
+--- a/s_tanf.cpp 2022-12-11 21:15:57.930542357 -0500
++++ b/s_tanf.cpp_new 2022-12-11 21:12:20.420189437 -0500
+@@ -15,15 +15,14 @@
+ */
+
+ #include <sys/cdefs.h>
+ __FBSDID("$FreeBSD$");
+
+ #include <float.h>
+
+-#include "math.h"
+ #define INLINE_KERNEL_TANDF
+ #define INLINE_REM_PIO2F
+ #include "math_private.h"
+ #include "e_rem_pio2f.c"
+ #include "k_tanf.c"
+
+ /* Small multiples of pi/2 rounded to double precision. */
diff --git a/modules/fdlibm/patches/05_include_stdint_h_in_math_private_h.patch b/modules/fdlibm/patches/05_include_stdint_h_in_math_private_h.patch
new file mode 100644
index 0000000000..aec7d94738
--- /dev/null
+++ b/modules/fdlibm/patches/05_include_stdint_h_in_math_private_h.patch
@@ -0,0 +1,21 @@
+diff --git a/math_private.h b/math_private.h
+--- a/math_private.h
++++ b/math_private.h
+@@ -12,16 +12,17 @@
+ /*
+ * from: @(#)fdlibm.h 5.1 93/09/24
+ * $FreeBSD$
+ */
+
+ #ifndef _MATH_PRIVATE_H_
+ #define _MATH_PRIVATE_H_
+
++#include <stdint.h>
+ #include <sys/types.h>
+ #include <machine/endian.h>
+
+ #include "fdlibm.h"
+
+ /*
+ * The original fdlibm code used statements like:
+ * n0 = ((*(int*)&one)>>29)^1; * index of high word *
diff --git a/modules/fdlibm/patches/06_use_mfbt_endian_h_in_math_private_h.patch b/modules/fdlibm/patches/06_use_mfbt_endian_h_in_math_private_h.patch
new file mode 100644
index 0000000000..b4374fc948
--- /dev/null
+++ b/modules/fdlibm/patches/06_use_mfbt_endian_h_in_math_private_h.patch
@@ -0,0 +1,116 @@
+diff --git a/math_private.h b/math_private.h
+--- a/math_private.h
++++ b/math_private.h
+@@ -14,52 +14,38 @@
+ * $FreeBSD$
+ */
+
+ #ifndef _MATH_PRIVATE_H_
+ #define _MATH_PRIVATE_H_
+
+ #include <stdint.h>
+ #include <sys/types.h>
+-#include <machine/endian.h>
++
++#include "mozilla/EndianUtils.h"
+
+ #include "fdlibm.h"
+
+ /*
+ * The original fdlibm code used statements like:
+ * n0 = ((*(int*)&one)>>29)^1; * index of high word *
+ * ix0 = *(n0+(int*)&x); * high word of x *
+ * ix1 = *((1-n0)+(int*)&x); * low word of x *
+ * to dig two 32 bit words out of the 64 bit IEEE floating point
+ * value. That is non-ANSI, and, moreover, the gcc instruction
+ * scheduler gets it wrong. We instead use the following macros.
+ * Unlike the original code, we determine the endianness at compile
+ * time, not at run time; I don't see much benefit to selecting
+ * endianness at run time.
+ */
+
+-/*
+- * A union which permits us to convert between a double and two 32 bit
+- * ints.
+- */
+-
+-#ifdef __arm__
+-#if defined(__VFP_FP__) || defined(__ARM_EABI__)
+-#define IEEE_WORD_ORDER BYTE_ORDER
+-#else
+-#define IEEE_WORD_ORDER BIG_ENDIAN
+-#endif
+-#else /* __arm__ */
+-#define IEEE_WORD_ORDER BYTE_ORDER
+-#endif
+-
+ /* A union which permits us to convert between a long double and
+ four 32 bit ints. */
+
+-#if IEEE_WORD_ORDER == BIG_ENDIAN
++#if MOZ_BIG_ENDIAN()
+
+ typedef union
+ {
+ long double value;
+ struct {
+ u_int32_t mswhi;
+ u_int32_t mswlo;
+ u_int32_t lswhi;
+@@ -68,17 +54,17 @@ typedef union
+ struct {
+ u_int64_t msw;
+ u_int64_t lsw;
+ } parts64;
+ } ieee_quad_shape_type;
+
+ #endif
+
+-#if IEEE_WORD_ORDER == LITTLE_ENDIAN
++#if MOZ_LITTLE_ENDIAN()
+
+ typedef union
+ {
+ long double value;
+ struct {
+ u_int32_t lswlo;
+ u_int32_t lswhi;
+ u_int32_t mswlo;
+@@ -87,17 +73,17 @@ typedef union
+ struct {
+ u_int64_t lsw;
+ u_int64_t msw;
+ } parts64;
+ } ieee_quad_shape_type;
+
+ #endif
+
+-#if IEEE_WORD_ORDER == BIG_ENDIAN
++#if MOZ_BIG_ENDIAN()
+
+ typedef union
+ {
+ double value;
+ struct
+ {
+ u_int32_t msw;
+ u_int32_t lsw;
+@@ -105,17 +91,17 @@ typedef union
+ struct
+ {
+ u_int64_t w;
+ } xparts;
+ } ieee_double_shape_type;
+
+ #endif
+
+-#if IEEE_WORD_ORDER == LITTLE_ENDIAN
++#if MOZ_LITTLE_ENDIAN()
+
+ typedef union
+ {
+ double value;
+ struct
+ {
+ u_int32_t lsw;
+ u_int32_t msw;
diff --git a/modules/fdlibm/patches/07_add_fdlibm_namespace_to_functions_defined_and_used_in_fdlibm.patch b/modules/fdlibm/patches/07_add_fdlibm_namespace_to_functions_defined_and_used_in_fdlibm.patch
new file mode 100644
index 0000000000..37e0017cbc
--- /dev/null
+++ b/modules/fdlibm/patches/07_add_fdlibm_namespace_to_functions_defined_and_used_in_fdlibm.patch
@@ -0,0 +1,71 @@
+diff --git a/math_private.h b/math_private.h
+--- a/math_private.h
++++ b/math_private.h
+@@ -867,16 +867,67 @@ irintl(long double x)
+ #define __ieee754_j1f j1f
+ #define __ieee754_y0f y0f
+ #define __ieee754_y1f y1f
+ #define __ieee754_jnf jnf
+ #define __ieee754_ynf ynf
+ #define __ieee754_remainderf remainderf
+ #define __ieee754_scalbf scalbf
+
++#define acos fdlibm::acos
++#define acosf fdlibm::acosf
++#define asin fdlibm::asin
++#define asinf fdlibm::asinf
++#define atan fdlibm::atan
++#define atanf fdlibm::atanf
++#define atan2 fdlibm::atan2
++#define cos fdlibm::cos
++#define cosf fdlibm::cosf
++#define sin fdlibm::sin
++#define sinf fdlibm::sinf
++#define tan fdlibm::tan
++#define tanf fdlibm::tanf
++#define cosh fdlibm::cosh
++#define sinh fdlibm::sinh
++#define tanh fdlibm::tanh
++#define exp fdlibm::exp
++#define expf fdlibm::expf
++#define exp2 fdlibm::exp2
++#define exp2f fdlibm::exp2f
++#define log fdlibm::log
++#define logf fdlibm::logf
++#define log10 fdlibm::log10
++#define pow fdlibm::pow
++#define powf fdlibm::powf
++#define ceil fdlibm::ceil
++#define ceilf fdlibm::ceilf
++#define fabs fdlibm::fabs
++#define fabsf fdlibm::fabsf
++#define floor fdlibm::floor
++#define acosh fdlibm::acosh
++#define asinh fdlibm::asinh
++#define atanh fdlibm::atanh
++#define cbrt fdlibm::cbrt
++#define expm1 fdlibm::expm1
++#define hypot fdlibm::hypot
++#define log1p fdlibm::log1p
++#define log2 fdlibm::log2
++#define scalb fdlibm::scalb
++#define copysign fdlibm::copysign
++#define scalbn fdlibm::scalbn
++#define scalbnf fdlibm::scalbnf
++#define trunc fdlibm::trunc
++#define truncf fdlibm::truncf
++#define floorf fdlibm::floorf
++#define nearbyint fdlibm::nearbyint
++#define nearbyintf fdlibm::nearbyintf
++#define rint fdlibm::rint
++#define rintf fdlibm::rintf
++#define sqrtf fdlibm::sqrtf
++
+ /* fdlibm kernel function */
+ int __kernel_rem_pio2(double*,double*,int,int,int);
+
+ /* double precision kernel functions */
+ #ifndef INLINE_REM_PIO2
+ int __ieee754_rem_pio2(double,double*);
+ #endif
+ double __kernel_sin(double,double,int);
diff --git a/modules/fdlibm/patches/08_remove_weak_reference_macro.patch b/modules/fdlibm/patches/08_remove_weak_reference_macro.patch
new file mode 100644
index 0000000000..338ff92b60
--- /dev/null
+++ b/modules/fdlibm/patches/08_remove_weak_reference_macro.patch
@@ -0,0 +1,459 @@
+diff --git a/e_acos.cpp b/e_acos.cpp
+--- a/e_acos.cpp
++++ b/e_acos.cpp
+@@ -99,12 +99,8 @@ __ieee754_acos(double x)
+ c = (z-df*df)/(s+df);
+ p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));
+ q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4)));
+ r = p/q;
+ w = r*s+c;
+ return 2.0*(df+w);
+ }
+ }
+-
+-#if LDBL_MANT_DIG == 53
+-__weak_reference(acos, acosl);
+-#endif
+diff --git a/e_acosh.cpp b/e_acosh.cpp
+--- a/e_acosh.cpp
++++ b/e_acosh.cpp
+@@ -56,12 +56,8 @@ __ieee754_acosh(double x)
+ } else if (hx > 0x40000000) { /* 2**28 > x > 2 */
+ t=x*x;
+ return __ieee754_log(2.0*x-one/(x+sqrt(t-one)));
+ } else { /* 1<x<2 */
+ t = x-one;
+ return log1p(t+sqrt(2.0*t+t*t));
+ }
+ }
+-
+-#if LDBL_MANT_DIG == 53
+-__weak_reference(acosh, acoshl);
+-#endif
+diff --git a/e_asin.cpp b/e_asin.cpp
+--- a/e_asin.cpp
++++ b/e_asin.cpp
+@@ -105,12 +105,8 @@ __ieee754_asin(double x)
+ c = (t-w*w)/(s+w);
+ r = p/q;
+ p = 2.0*s*r-(pio2_lo-2.0*c);
+ q = pio4_hi-2.0*w;
+ t = pio4_hi-(p-q);
+ }
+ if(hx>0) return t; else return -t;
+ }
+-
+-#if LDBL_MANT_DIG == 53
+-__weak_reference(asin, asinl);
+-#endif
+diff --git a/e_atan2.cpp b/e_atan2.cpp
+--- a/e_atan2.cpp
++++ b/e_atan2.cpp
+@@ -117,12 +117,8 @@ __ieee754_atan2(double y, double x)
+ switch (m) {
+ case 0: return z ; /* atan(+,+) */
+ case 1: return -z ; /* atan(-,+) */
+ case 2: return pi-(z-pi_lo);/* atan(+,-) */
+ default: /* case 3 */
+ return (z-pi_lo)-pi;/* atan(-,-) */
+ }
+ }
+-
+-#if LDBL_MANT_DIG == 53
+-__weak_reference(atan2, atan2l);
+-#endif
+diff --git a/e_atanh.cpp b/e_atanh.cpp
+--- a/e_atanh.cpp
++++ b/e_atanh.cpp
+@@ -56,12 +56,8 @@ __ieee754_atanh(double x)
+ SET_HIGH_WORD(x,ix);
+ if(ix<0x3fe00000) { /* x < 0.5 */
+ t = x+x;
+ t = 0.5*log1p(t+t*x/(one-x));
+ } else
+ t = 0.5*log1p((x+x)/(one-x));
+ if(hx>=0) return t; else return -t;
+ }
+-
+-#if LDBL_MANT_DIG == 53
+-__weak_reference(atanh, atanhl);
+-#endif
+diff --git a/e_cosh.cpp b/e_cosh.cpp
+--- a/e_cosh.cpp
++++ b/e_cosh.cpp
+@@ -73,12 +73,8 @@ __ieee754_cosh(double x)
+
+ /* |x| in [log(maxdouble), overflowthresold] */
+ if (ix<=0x408633CE)
+ return __ldexp_exp(fabs(x), -1);
+
+ /* |x| > overflowthresold, cosh(x) overflow */
+ return huge*huge;
+ }
+-
+-#if (LDBL_MANT_DIG == 53)
+-__weak_reference(cosh, coshl);
+-#endif
+diff --git a/e_exp.cpp b/e_exp.cpp
+--- a/e_exp.cpp
++++ b/e_exp.cpp
+@@ -152,12 +152,8 @@ __ieee754_exp(double x) /* default IEEE double exp */
+ else y = one-((lo-(x*c)/(2.0-c))-hi);
+ if(k >= -1021) {
+ if (k==1024) return y*2.0*0x1p1023;
+ return y*twopk;
+ } else {
+ return y*twopk*twom1000;
+ }
+ }
+-
+-#if (LDBL_MANT_DIG == 53)
+-__weak_reference(exp, expl);
+-#endif
+diff --git a/e_hypot.cpp b/e_hypot.cpp
+--- a/e_hypot.cpp
++++ b/e_hypot.cpp
+@@ -117,12 +117,8 @@ __ieee754_hypot(double x, double y)
+ w = sqrt(t1*y1-(w*(-w)-(t1*y2+t2*b)));
+ }
+ if(k!=0) {
+ t1 = 0.0;
+ SET_HIGH_WORD(t1,(1023+k)<<20);
+ return t1*w;
+ } else return w;
+ }
+-
+-#if LDBL_MANT_DIG == 53
+-__weak_reference(hypot, hypotl);
+-#endif
+diff --git a/e_log.cpp b/e_log.cpp
+--- a/e_log.cpp
++++ b/e_log.cpp
+@@ -135,12 +135,8 @@ __ieee754_log(double x)
+ hfsq=0.5*f*f;
+ if(k==0) return f-(hfsq-s*(hfsq+R)); else
+ return dk*ln2_hi-((hfsq-(s*(hfsq+R)+dk*ln2_lo))-f);
+ } else {
+ if(k==0) return f-s*(f-R); else
+ return dk*ln2_hi-((s*(f-R)-dk*ln2_lo)-f);
+ }
+ }
+-
+-#if (LDBL_MANT_DIG == 53)
+-__weak_reference(log, logl);
+-#endif
+diff --git a/e_log10.cpp b/e_log10.cpp
+--- a/e_log10.cpp
++++ b/e_log10.cpp
+@@ -82,12 +82,8 @@ __ieee754_log10(double x)
+ * with some parallelism and it reduces the error for many args.
+ */
+ w = y2 + val_hi;
+ val_lo += (y2 - w) + val_hi;
+ val_hi = w;
+
+ return val_lo + val_hi;
+ }
+-
+-#if (LDBL_MANT_DIG == 53)
+-__weak_reference(log10, log10l);
+-#endif
+diff --git a/e_log2.cpp b/e_log2.cpp
+--- a/e_log2.cpp
++++ b/e_log2.cpp
+@@ -105,12 +105,8 @@ __ieee754_log2(double x)
+
+ /* spadd(val_hi, val_lo, y), except for not using double_t: */
+ w = y + val_hi;
+ val_lo += (y - w) + val_hi;
+ val_hi = w;
+
+ return val_lo + val_hi;
+ }
+-
+-#if (LDBL_MANT_DIG == 53)
+-__weak_reference(log2, log2l);
+-#endif
+diff --git a/e_pow.cpp b/e_pow.cpp
+--- a/e_pow.cpp
++++ b/e_pow.cpp
+@@ -302,12 +302,8 @@ __ieee754_pow(double x, double y)
+ r = (z*t1)/(t1-two)-(w+z*w);
+ z = one-(r-z);
+ GET_HIGH_WORD(j,z);
+ j += (n<<20);
+ if((j>>20)<=0) z = scalbn(z,n); /* subnormal output */
+ else SET_HIGH_WORD(z,j);
+ return s*z;
+ }
+-
+-#if (LDBL_MANT_DIG == 53)
+-__weak_reference(pow, powl);
+-#endif
+diff --git a/e_sinh.cpp b/e_sinh.cpp
+--- a/e_sinh.cpp
++++ b/e_sinh.cpp
+@@ -67,12 +67,8 @@ __ieee754_sinh(double x)
+
+ /* |x| in [log(maxdouble), overflowthresold] */
+ if (ix<=0x408633CE)
+ return h*2.0*__ldexp_exp(fabs(x), -1);
+
+ /* |x| > overflowthresold, sinh(x) overflow */
+ return x*shuge;
+ }
+-
+-#if (LDBL_MANT_DIG == 53)
+-__weak_reference(sinh, sinhl);
+-#endif
+diff --git a/s_asinh.cpp b/s_asinh.cpp
+--- a/s_asinh.cpp
++++ b/s_asinh.cpp
+@@ -50,12 +50,8 @@ asinh(double x)
+ t = fabs(x);
+ w = __ieee754_log(2.0*t+one/(__ieee754_sqrt(x*x+one)+t));
+ } else { /* 2.0 > |x| > 2**-28 */
+ t = x*x;
+ w =log1p(fabs(x)+t/(one+__ieee754_sqrt(one+t)));
+ }
+ if(hx>0) return w; else return -w;
+ }
+-
+-#if LDBL_MANT_DIG == 53
+-__weak_reference(asinh, asinhl);
+-#endif
+diff --git a/s_atan.cpp b/s_atan.cpp
+--- a/s_atan.cpp
++++ b/s_atan.cpp
+@@ -112,12 +112,8 @@ atan(double x)
+ s1 = z*(aT[0]+w*(aT[2]+w*(aT[4]+w*(aT[6]+w*(aT[8]+w*aT[10])))));
+ s2 = w*(aT[1]+w*(aT[3]+w*(aT[5]+w*(aT[7]+w*aT[9]))));
+ if (id<0) return x - x*(s1+s2);
+ else {
+ z = atanhi[id] - ((x*(s1+s2) - atanlo[id]) - x);
+ return (hx<0)? -z:z;
+ }
+ }
+-
+-#if LDBL_MANT_DIG == 53
+-__weak_reference(atan, atanl);
+-#endif
+diff --git a/s_cbrt.cpp b/s_cbrt.cpp
+--- a/s_cbrt.cpp
++++ b/s_cbrt.cpp
+@@ -106,12 +106,8 @@ cbrt(double x)
+ s=t*t; /* t*t is exact */
+ r=x/s; /* error <= 0.5 ulps; |r| < |t| */
+ w=t+t; /* t+t is exact */
+ r=(r-t)/(w+r); /* r-t is exact; w+r ~= 3*t */
+ t=t+t*r; /* error <= (0.5 + 0.5/3) * ulp */
+
+ return(t);
+ }
+-
+-#if (LDBL_MANT_DIG == 53)
+-__weak_reference(cbrt, cbrtl);
+-#endif
+diff --git a/s_ceil.cpp b/s_ceil.cpp
+--- a/s_ceil.cpp
++++ b/s_ceil.cpp
+@@ -65,12 +65,8 @@ ceil(double x)
+ }
+ }
+ i1 &= (~i);
+ }
+ }
+ INSERT_WORDS(x,i0,i1);
+ return x;
+ }
+-
+-#if LDBL_MANT_DIG == 53
+-__weak_reference(ceil, ceill);
+-#endif
+diff --git a/s_cos.cpp b/s_cos.cpp
+--- a/s_cos.cpp
++++ b/s_cos.cpp
+@@ -77,12 +77,8 @@ cos(double x)
+ case 0: return __kernel_cos(y[0],y[1]);
+ case 1: return -__kernel_sin(y[0],y[1],1);
+ case 2: return -__kernel_cos(y[0],y[1]);
+ default:
+ return __kernel_sin(y[0],y[1],1);
+ }
+ }
+ }
+-
+-#if (LDBL_MANT_DIG == 53)
+-__weak_reference(cos, cosl);
+-#endif
+diff --git a/s_expm1.cpp b/s_expm1.cpp
+--- a/s_expm1.cpp
++++ b/s_expm1.cpp
+@@ -210,12 +210,8 @@ expm1(double x)
+ SET_HIGH_WORD(t,((0x3ff-k)<<20)); /* 2^-k */
+ y = x-(e+t);
+ y += one;
+ y = y*twopk;
+ }
+ }
+ return y;
+ }
+-
+-#if (LDBL_MANT_DIG == 53)
+-__weak_reference(expm1, expm1l);
+-#endif
+diff --git a/s_floor.cpp b/s_floor.cpp
+--- a/s_floor.cpp
++++ b/s_floor.cpp
+@@ -66,12 +66,8 @@ floor(double x)
+ }
+ }
+ i1 &= (~i);
+ }
+ }
+ INSERT_WORDS(x,i0,i1);
+ return x;
+ }
+-
+-#if LDBL_MANT_DIG == 53
+-__weak_reference(floor, floorl);
+-#endif
+diff --git a/s_log1p.cpp b/s_log1p.cpp
+--- a/s_log1p.cpp
++++ b/s_log1p.cpp
+@@ -168,12 +168,8 @@ log1p(double x)
+ return k*ln2_hi-((R-(k*ln2_lo+c))-f);
+ }
+ s = f/(2.0+f);
+ z = s*s;
+ R = z*(Lp1+z*(Lp2+z*(Lp3+z*(Lp4+z*(Lp5+z*(Lp6+z*Lp7))))));
+ if(k==0) return f-(hfsq-s*(hfsq+R)); else
+ return k*ln2_hi-((hfsq-(s*(hfsq+R)+(k*ln2_lo+c)))-f);
+ }
+-
+-#if (LDBL_MANT_DIG == 53)
+-__weak_reference(log1p, log1pl);
+-#endif
+diff --git a/s_rint.cpp b/s_rint.cpp
+--- a/s_rint.cpp
++++ b/s_rint.cpp
+@@ -80,12 +80,8 @@ rint(double x)
+ if((i1&i)==0) return x; /* x is integral */
+ i>>=1;
+ if((i1&i)!=0) i1 = (i1&(~i))|((0x40000000)>>(j0-20));
+ }
+ INSERT_WORDS(x,i0,i1);
+ STRICT_ASSIGN(double,w,TWO52[sx]+x);
+ return w-TWO52[sx];
+ }
+-
+-#if (LDBL_MANT_DIG == 53)
+-__weak_reference(rint, rintl);
+-#endif
+diff --git a/s_scalbn.cpp b/s_scalbn.cpp
+--- a/s_scalbn.cpp
++++ b/s_scalbn.cpp
+@@ -36,13 +36,8 @@ double scalbn(double x, int n)
+ if (n < -1022)
+ n = -1022;
+ }
+ }
+ u.i = (uint64_t)(0x3ff+n)<<52;
+ x = y * u.f;
+ return x;
+ }
+-
+-#if (LDBL_MANT_DIG == 53) && !defined(scalbn)
+-__weak_reference(scalbn, ldexpl);
+-__weak_reference(scalbn, scalbnl);
+-#endif
+diff --git a/s_sin.cpp b/s_sin.cpp
+--- a/s_sin.cpp
++++ b/s_sin.cpp
+@@ -77,12 +77,8 @@ sin(double x)
+ case 0: return __kernel_sin(y[0],y[1],1);
+ case 1: return __kernel_cos(y[0],y[1]);
+ case 2: return -__kernel_sin(y[0],y[1],1);
+ default:
+ return -__kernel_cos(y[0],y[1]);
+ }
+ }
+ }
+-
+-#if (LDBL_MANT_DIG == 53)
+-__weak_reference(sin, sinl);
+-#endif
+diff --git a/s_tan.cpp b/s_tan.cpp
+--- a/s_tan.cpp
++++ b/s_tan.cpp
+@@ -71,12 +71,8 @@ tan(double x)
+
+ /* argument reduction needed */
+ else {
+ n = __ieee754_rem_pio2(x,y);
+ return __kernel_tan(y[0],y[1],1-((n&1)<<1)); /* 1 -- n even
+ -1 -- n odd */
+ }
+ }
+-
+-#if (LDBL_MANT_DIG == 53)
+-__weak_reference(tan, tanl);
+-#endif
+diff --git a/s_tanh.cpp b/s_tanh.cpp
+--- a/s_tanh.cpp
++++ b/s_tanh.cpp
+@@ -72,12 +72,8 @@ tanh(double x)
+ z= -t/(t+two);
+ }
+ /* |x| >= 22, return +-1 */
+ } else {
+ z = one - tiny; /* raise inexact flag */
+ }
+ return (jx>=0)? z: -z;
+ }
+-
+-#if (LDBL_MANT_DIG == 53)
+-__weak_reference(tanh, tanhl);
+-#endif
+diff --git a/s_trunc.cpp b/s_trunc.cpp
+--- a/s_trunc.cpp
++++ b/s_trunc.cpp
+@@ -55,12 +55,8 @@ trunc(double x)
+ i = ((u_int32_t)(0xffffffff))>>(j0-20);
+ if((i1&i)==0) return x; /* x is integral */
+ if(huge+x>0.0) /* raise inexact flag */
+ i1 &= (~i);
+ }
+ INSERT_WORDS(x,i0,i1);
+ return x;
+ }
+-
+-#if LDBL_MANT_DIG == 53
+-__weak_reference(trunc, truncl);
+-#endif
+--- a/s_exp2.cpp 2022-12-11 21:20:40.290068458 -0500
++++ b/s_exp2.cpp_new 2022-12-11 21:20:57.490259247 -0500
+@@ -389,10 +389,6 @@
+ return (r * 2.0 * 0x1p1023);
+ return (r * twopk);
+ } else {
+ return (r * twopkp1000 * twom1000);
+ }
+ }
+-
+-#if (LDBL_MANT_DIG == 53)
+-__weak_reference(exp2, exp2l);
+-#endif
+--- a/s_scalbnf.cpp
++++ b/s_scalbnf.cpp
+@@ -32,10 +34,8 @@
+ if (n < -126)
+ n = -126;
+ }
+ }
+ u.i = (uint32_t)(0x7f+n)<<23;
+ x = y * u.f;
+ return x;
+ }
+-
+-__strong_reference(scalbnf, ldexpf);
diff --git a/modules/fdlibm/patches/09_comment_out_rcsid_variable.patch b/modules/fdlibm/patches/09_comment_out_rcsid_variable.patch
new file mode 100644
index 0000000000..715c859956
--- /dev/null
+++ b/modules/fdlibm/patches/09_comment_out_rcsid_variable.patch
@@ -0,0 +1,1362 @@
+diff --git a/e_acos.cpp b/e_acos.cpp
+--- a/e_acos.cpp
++++ b/e_acos.cpp
+@@ -6,18 +6,18 @@
+ *
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* __ieee754_acos(x)
+ * Method :
+ * acos(x) = pi/2 - asin(x)
+ * acos(-x) = pi/2 + asin(x)
+ * For |x|<=0.5
+ * acos(x) = pi/2 - (x + x*x^2*R(x^2)) (see asin.c)
+ * For x>0.5
+diff --git a/e_acosh.cpp b/e_acosh.cpp
+--- a/e_acosh.cpp
++++ b/e_acosh.cpp
+@@ -7,18 +7,18 @@
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ *
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* __ieee754_acosh(x)
+ * Method :
+ * Based on
+ * acosh(x) = log [ x + sqrt(x*x-1) ]
+ * we have
+ * acosh(x) := log(x)+ln2, if x is large; else
+ * acosh(x) := log(2x-1/(sqrt(x*x-1)+x)) if x>2; else
+diff --git a/e_asin.cpp b/e_asin.cpp
+--- a/e_asin.cpp
++++ b/e_asin.cpp
+@@ -6,18 +6,18 @@
+ *
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* __ieee754_asin(x)
+ * Method :
+ * Since asin(x) = x + x^3/6 + x^5*3/40 + x^7*15/336 + ...
+ * we approximate asin(x) on [0,0.5] by
+ * asin(x) = x + x*x^2*R(x^2)
+ * where
+ * R(x^2) is a rational approximation of (asin(x)-x)/x^3
+diff --git a/e_atan2.cpp b/e_atan2.cpp
+--- a/e_atan2.cpp
++++ b/e_atan2.cpp
+@@ -7,18 +7,18 @@
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ *
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* __ieee754_atan2(y,x)
+ * Method :
+ * 1. Reduce y to positive by atan2(y,x)=-atan2(-y,x).
+ * 2. Reduce x to positive by (if x and y are unexceptional):
+ * ARG (x+iy) = arctan(y/x) ... if x > 0,
+ * ARG (x+iy) = pi - arctan[y/(-x)] ... if x < 0,
+ *
+diff --git a/e_atanh.cpp b/e_atanh.cpp
+--- a/e_atanh.cpp
++++ b/e_atanh.cpp
+@@ -7,18 +7,18 @@
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ *
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* __ieee754_atanh(x)
+ * Method :
+ * 1.Reduced x to positive by atanh(-x) = -atanh(x)
+ * 2.For x>=0.5
+ * 1 2x x
+ * atanh(x) = --- * log(1 + -------) = 0.5 * log1p(2 * --------)
+ * 2 1 - x 1 - x
+diff --git a/e_cosh.cpp b/e_cosh.cpp
+--- a/e_cosh.cpp
++++ b/e_cosh.cpp
+@@ -6,18 +6,18 @@
+ *
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* __ieee754_cosh(x)
+ * Method :
+ * mathematically cosh(x) if defined to be (exp(x)+exp(-x))/2
+ * 1. Replace x by |x| (cosh(x) = cosh(-x)).
+ * 2.
+ * [ exp(x) - 1 ]^2
+ * 0 <= x <= ln2/2 : cosh(x) := 1 + -------------------
+diff --git a/e_exp.cpp b/e_exp.cpp
+--- a/e_exp.cpp
++++ b/e_exp.cpp
+@@ -5,18 +5,18 @@
+ * Copyright (C) 2004 by Sun Microsystems, Inc. All rights reserved.
+ *
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* __ieee754_exp(x)
+ * Returns the exponential of x.
+ *
+ * Method
+ * 1. Argument reduction:
+ * Reduce x to an r so that |r| <= 0.5*ln2 ~ 0.34658.
+ * Given x, find r and integer k such that
+diff --git a/e_hypot.cpp b/e_hypot.cpp
+--- a/e_hypot.cpp
++++ b/e_hypot.cpp
+@@ -6,18 +6,18 @@
+ *
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* __ieee754_hypot(x,y)
+ *
+ * Method :
+ * If (assume round-to-nearest) z=x*x+y*y
+ * has error less than sqrt(2)/2 ulp, than
+ * sqrt(z) has error less than 1 ulp (exercise).
+ *
+diff --git a/e_log.cpp b/e_log.cpp
+--- a/e_log.cpp
++++ b/e_log.cpp
+@@ -6,18 +6,18 @@
+ *
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* __ieee754_log(x)
+ * Return the logrithm of x
+ *
+ * Method :
+ * 1. Argument Reduction: find k and f such that
+ * x = 2^k * (1+f),
+ * where sqrt(2)/2 < 1+f < sqrt(2) .
+diff --git a/e_log10.cpp b/e_log10.cpp
+--- a/e_log10.cpp
++++ b/e_log10.cpp
+@@ -6,32 +6,32 @@
+ *
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /*
+ * Return the base 10 logarithm of x. See e_log.c and k_log.h for most
+ * comments.
+ *
+ * log10(x) = (f - 0.5*f*f + k_log1p(f)) / ln10 + k * log10(2)
+ * in not-quite-routine extra precision.
+ */
+
+ #include <float.h>
+
+ #include "math_private.h"
+ #include "k_log.h"
+-
++
+ static const double
+ two54 = 1.80143985094819840000e+16, /* 0x43500000, 0x00000000 */
+ ivln10hi = 4.34294481878168880939e-01, /* 0x3fdbcb7b, 0x15200000 */
+ ivln10lo = 2.50829467116452752298e-11, /* 0x3dbb9438, 0xca9aadd5 */
+ log10_2hi = 3.01029995663611771306e-01, /* 0x3FD34413, 0x509F6000 */
+ log10_2lo = 3.69423907715893078616e-13; /* 0x3D59FEF3, 0x11F12B36 */
+
+ static const double zero = 0.0;
+diff --git a/e_log2.cpp b/e_log2.cpp
+--- a/e_log2.cpp
++++ b/e_log2.cpp
+@@ -6,18 +6,18 @@
+ *
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /*
+ * Return the base 2 logarithm of x. See e_log.c and k_log.h for most
+ * comments.
+ *
+ * This reduces x to {k, 1+f} exactly as in e_log.c, then calls the kernel,
+ * then does the combining and scaling steps
+ * log2(x) = (f - 0.5*f*f + k_log1p(f)) / ln2 + k
+diff --git a/e_pow.cpp b/e_pow.cpp
+--- a/e_pow.cpp
++++ b/e_pow.cpp
+@@ -4,18 +4,18 @@
+ * Copyright (C) 2004 by Sun Microsystems, Inc. All rights reserved.
+ *
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* __ieee754_pow(x,y) return x**y
+ *
+ * n
+ * Method: Let x = 2 * (1+f)
+ * 1. Compute and return log2(x) in two pieces:
+ * log2(x) = w1 + w2,
+ * where w1 has 53-24 = 29 bit trailing zeros.
+diff --git a/e_rem_pio2.cpp b/e_rem_pio2.cpp
+--- a/e_rem_pio2.cpp
++++ b/e_rem_pio2.cpp
+@@ -8,18 +8,18 @@
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ *
+ * Optimized by Bruce D. Evans.
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* __ieee754_rem_pio2(x,y)
+ *
+ * return the remainder of x rem pi/2 in y[0]+y[1]
+ * use __kernel_rem_pio2()
+ */
+
+ #include <float.h>
+diff --git a/e_sinh.cpp b/e_sinh.cpp
+--- a/e_sinh.cpp
++++ b/e_sinh.cpp
+@@ -6,18 +6,18 @@
+ *
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* __ieee754_sinh(x)
+ * Method :
+ * mathematically sinh(x) if defined to be (exp(x)-exp(-x))/2
+ * 1. Replace x by |x| (sinh(-x) = -sinh(x)).
+ * 2.
+ * E + E/(E+1)
+ * 0 <= x <= 22 : sinh(x) := --------------, E=expm1(x)
+diff --git a/k_cos.cpp b/k_cos.cpp
+--- a/k_cos.cpp
++++ b/k_cos.cpp
+@@ -6,18 +6,18 @@
+ *
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /*
+ * __kernel_cos( x, y )
+ * kernel cos function on [-pi/4, pi/4], pi/4 ~ 0.785398164
+ * Input x is assumed to be bounded by ~pi/4 in magnitude.
+ * Input y is the tail of x.
+ *
+ * Algorithm
+diff --git a/k_exp.cpp b/k_exp.cpp
+--- a/k_exp.cpp
++++ b/k_exp.cpp
+@@ -21,22 +21,22 @@
+ * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
+ * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
+ * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+ * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
+ * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
+ * SUCH DAMAGE.
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ #include <complex.h>
+
+-#include "math_private.h"
++ #include "math_private.h"
+
+ static const uint32_t k = 1799; /* constant for reduction */
+ static const double kln2 = 1246.97177782734161156; /* k * ln2 */
+
+ /*
+ * Compute exp(x), scaled to avoid spurious overflow. An exponent is
+ * returned separately in 'expt'.
+ *
+diff --git a/k_log.h b/k_log.h
+--- a/k_log.h
++++ b/k_log.h
+@@ -6,18 +6,18 @@
+ *
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /*
+ * k_log1p(f):
+ * Return log(1+f) - f for 1+f in ~[sqrt(2)/2, sqrt(2)].
+ *
+ * The following describes the overall strategy for computing
+ * logarithms in base e. The argument reduction and adding the final
+ * term of the polynomial are done by the caller for increased accuracy
+diff --git a/k_rem_pio2.cpp b/k_rem_pio2.cpp
+--- a/k_rem_pio2.cpp
++++ b/k_rem_pio2.cpp
+@@ -6,18 +6,18 @@
+ *
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /*
+ * __kernel_rem_pio2(x,y,e0,nx,prec)
+ * double x[],y[]; int e0,nx,prec;
+ *
+ * __kernel_rem_pio2 return the last three digits of N with
+ * y = x - N*pi/2
+ * so that |y| < pi/2.
+diff --git a/k_sin.cpp b/k_sin.cpp
+--- a/k_sin.cpp
++++ b/k_sin.cpp
+@@ -6,18 +6,18 @@
+ *
+ * Developed at SunSoft, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* __kernel_sin( x, y, iy)
+ * kernel sin function on ~[-pi/4, pi/4] (except on -0), pi/4 ~ 0.7854
+ * Input x is assumed to be bounded by ~pi/4 in magnitude.
+ * Input y is the tail of x.
+ * Input iy indicates whether y is 0. (if iy=0, y assume to be 0).
+ *
+ * Algorithm
+diff --git a/k_tan.cpp b/k_tan.cpp
+--- a/k_tan.cpp
++++ b/k_tan.cpp
+@@ -6,18 +6,18 @@
+ *
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+ /* INDENT OFF */
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* __kernel_tan( x, y, k )
+ * kernel tan function on ~[-pi/4, pi/4] (except on -0), pi/4 ~ 0.7854
+ * Input x is assumed to be bounded by ~pi/4 in magnitude.
+ * Input y is the tail of x.
+ * Input k indicates whether tan (if k = 1) or -1/tan (if k = -1) is returned.
+ *
+ * Algorithm
+diff --git a/s_asinh.cpp b/s_asinh.cpp
+--- a/s_asinh.cpp
++++ b/s_asinh.cpp
+@@ -5,18 +5,18 @@
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* asinh(x)
+ * Method :
+ * Based on
+ * asinh(x) = sign(x) * log [ |x| + sqrt(x*x+1) ]
+ * we have
+ * asinh(x) := x if 1+x*x=1,
+ * := sign(x)*(log(x)+ln2)) for large |x|, else
+diff --git a/s_atan.cpp b/s_atan.cpp
+--- a/s_atan.cpp
++++ b/s_atan.cpp
+@@ -5,18 +5,18 @@
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* atan(x)
+ * Method
+ * 1. Reduce x to positive by atan(x) = -atan(-x).
+ * 2. According to the integer k=4t+0.25 chopped, t=x, the argument
+ * is further reduced to one of the following intervals and the
+ * arctangent of t is evaluated by the corresponding formula:
+ *
+diff --git a/s_cbrt.cpp b/s_cbrt.cpp
+--- a/s_cbrt.cpp
++++ b/s_cbrt.cpp
+@@ -7,18 +7,18 @@
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ *
+ * Optimized by Bruce D. Evans.
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ #include <float.h>
+ #include "math_private.h"
+
+ /* cbrt(x)
+ * Return cube root of x
+ */
+ static const u_int32_t
+diff --git a/s_ceil.cpp b/s_ceil.cpp
+--- a/s_ceil.cpp
++++ b/s_ceil.cpp
+@@ -5,18 +5,18 @@
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /*
+ * ceil(x)
+ * Return x rounded toward -inf to integral value
+ * Method:
+ * Bit twiddling.
+ * Exception:
+ * Inexact flag raised if x not equal to ceil(x).
+diff --git a/s_ceilf.cpp b/s_ceilf.cpp
+--- a/s_ceilf.cpp
++++ b/s_ceilf.cpp
+@@ -8,18 +8,18 @@
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ #include "math_private.h"
+
+ static const float huge = 1.0e30;
+
+ float
+ ceilf(float x)
+ {
+diff --git a/s_copysign.cpp b/s_copysign.cpp
+--- a/s_copysign.cpp
++++ b/s_copysign.cpp
+@@ -5,18 +5,18 @@
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /*
+ * copysign(double x, double y)
+ * copysign(x,y) returns a value with the magnitude of x and
+ * with the sign bit of y.
+ */
+
+ #include "math_private.h"
+diff --git a/s_cos.cpp b/s_cos.cpp
+--- a/s_cos.cpp
++++ b/s_cos.cpp
+@@ -5,18 +5,18 @@
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* cos(x)
+ * Return cosine function of x.
+ *
+ * kernel function:
+ * __kernel_sin ... sine function on [-pi/4,pi/4]
+ * __kernel_cos ... cosine function on [-pi/4,pi/4]
+ * __ieee754_rem_pio2 ... argument reduction routine
+diff --git a/s_expm1.cpp b/s_expm1.cpp
+--- a/s_expm1.cpp
++++ b/s_expm1.cpp
+@@ -5,18 +5,18 @@
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* expm1(x)
+ * Returns exp(x)-1, the exponential of x minus 1.
+ *
+ * Method
+ * 1. Argument reduction:
+ * Given x, find r and integer k such that
+ *
+diff --git a/s_fabs.cpp b/s_fabs.cpp
+--- a/s_fabs.cpp
++++ b/s_fabs.cpp
+@@ -5,18 +5,18 @@
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /*
+ * fabs(x) returns the absolute value of x.
+ */
+
+ #include "math_private.h"
+
+ double
+diff --git a/s_floor.cpp b/s_floor.cpp
+--- a/s_floor.cpp
++++ b/s_floor.cpp
+@@ -5,18 +5,18 @@
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /*
+ * floor(x)
+ * Return x rounded toward -inf to integral value
+ * Method:
+ * Bit twiddling.
+ * Exception:
+ * Inexact flag raised if x not equal to floor(x).
+diff --git a/s_floorf.cpp b/s_floorf.cpp
+--- a/s_floorf.cpp
++++ b/s_floorf.cpp
+@@ -8,18 +8,18 @@
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /*
+ * floorf(x)
+ * Return x rounded toward -inf to integral value
+ * Method:
+ * Bit twiddling.
+ * Exception:
+ * Inexact flag raised if x not equal to floorf(x).
+diff --git a/s_log1p.cpp b/s_log1p.cpp
+--- a/s_log1p.cpp
++++ b/s_log1p.cpp
+@@ -5,18 +5,18 @@
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* double log1p(double x)
+ *
+ * Method :
+ * 1. Argument Reduction: find k and f such that
+ * 1+x = 2^k * (1+f),
+ * where sqrt(2)/2 < 1+f < sqrt(2) .
+ *
+diff --git a/s_nearbyint.cpp b/s_nearbyint.cpp
+--- a/s_nearbyint.cpp
++++ b/s_nearbyint.cpp
+@@ -21,18 +21,18 @@
+ * DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS
+ * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
+ * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+ * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
+ * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
+ * SUCH DAMAGE.
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ #include <fenv.h>
+
+ #include "math_private.h"
+
+ /*
+ * We save and restore the floating-point environment to avoid raising
+ * an inexact exception. We can get away with using fesetenv()
+diff --git a/s_rint.cpp b/s_rint.cpp
+--- a/s_rint.cpp
++++ b/s_rint.cpp
+@@ -5,18 +5,18 @@
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /*
+ * rint(x)
+ * Return x rounded to integral value according to the prevailing
+ * rounding mode.
+ * Method:
+ * Using floating addition.
+ * Exception:
+diff --git a/s_rintf.cpp b/s_rintf.cpp
+--- a/s_rintf.cpp
++++ b/s_rintf.cpp
+@@ -8,18 +8,18 @@
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ #include <float.h>
+ #include <stdint.h>
+
+ #include "math_private.h"
+
+ static const float
+ TWO23[2]={
+diff --git a/s_sin.cpp b/s_sin.cpp
+--- a/s_sin.cpp
++++ b/s_sin.cpp
+@@ -5,18 +5,18 @@
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* sin(x)
+ * Return sine function of x.
+ *
+ * kernel function:
+ * __kernel_sin ... sine function on [-pi/4,pi/4]
+ * __kernel_cos ... cose function on [-pi/4,pi/4]
+ * __ieee754_rem_pio2 ... argument reduction routine
+diff --git a/s_tan.cpp b/s_tan.cpp
+--- a/s_tan.cpp
++++ b/s_tan.cpp
+@@ -5,18 +5,18 @@
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* tan(x)
+ * Return tangent function of x.
+ *
+ * kernel function:
+ * __kernel_tan ... tangent function on [-pi/4,pi/4]
+ * __ieee754_rem_pio2 ... argument reduction routine
+ *
+diff --git a/s_tanh.cpp b/s_tanh.cpp
+--- a/s_tanh.cpp
++++ b/s_tanh.cpp
+@@ -5,18 +5,18 @@
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* Tanh(x)
+ * Return the Hyperbolic Tangent of x
+ *
+ * Method :
+ * x -x
+ * e - e
+ * 0. tanh(x) is defined to be -----------
+diff --git a/s_trunc.cpp b/s_trunc.cpp
+--- a/s_trunc.cpp
++++ b/s_trunc.cpp
+@@ -5,18 +5,18 @@
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /*
+ * trunc(x)
+ * Return x rounded toward 0 to integral value
+ * Method:
+ * Bit twiddling.
+ * Exception:
+ * Inexact flag raised if x not equal to trunc(x).
+diff --git a/s_truncf.cpp b/s_truncf.cpp
+--- a/s_truncf.cpp
++++ b/s_truncf.cpp
+@@ -5,18 +5,18 @@
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /*
+ * truncf(x)
+ * Return x rounded toward 0 to integral value
+ * Method:
+ * Bit twiddling.
+ * Exception:
+ * Inexact flag raised if x not equal to truncf(x).
+--- a/e_acosf.cpp 2022-12-11 21:20:40.290068458 -0500
++++ b/e_acosf.cpp_new 2022-12-11 21:23:31.863880382 -0500
+@@ -9,16 +9,16 @@
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ #include "math_private.h"
+
+ static const float
+ one = 1.0000000000e+00, /* 0x3F800000 */
+ pi = 3.1415925026e+00, /* 0x40490fda */
+ pio2_hi = 1.5707962513e+00; /* 0x3fc90fda */
+--- a/e_asinf.cpp 2022-12-11 21:20:40.290068458 -0500
++++ b/e_asinf.cpp_new 2022-12-11 21:23:31.875880502 -0500
+@@ -9,16 +9,16 @@
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ #include "math_private.h"
+
+ static const float
+ one = 1.0000000000e+00, /* 0x3F800000 */
+ huge = 1.000e+30,
+ /* coefficient for R(x^2) */
+--- a/e_expf.cpp 2022-12-11 21:20:40.290068458 -0500
++++ b/e_expf.cpp_new 2022-12-11 21:23:31.895880703 -0500
+@@ -9,16 +9,16 @@
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ #include <float.h>
+
+ #include "math_private.h"
+
+ static const float
+ one = 1.0,
+--- a/e_logf.cpp 2022-12-11 21:20:40.290068458 -0500
++++ b/e_logf.cpp_new 2022-12-11 21:23:31.935881104 -0500
+@@ -9,16 +9,16 @@
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ #include "math_private.h"
+
+ static const float
+ ln2_hi = 6.9313812256e-01, /* 0x3f317180 */
+ ln2_lo = 9.0580006145e-06, /* 0x3717f7d1 */
+ two25 = 3.355443200e+07, /* 0x4c000000 */
+--- a/e_powf.cpp 2022-12-11 21:20:40.290068458 -0500
++++ b/e_powf.cpp_new 2022-12-11 21:23:31.971881465 -0500
+@@ -9,16 +9,16 @@
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ #include "math_private.h"
+
+ static const float
+ bp[] = {1.0, 1.5,},
+ dp_h[] = { 0.0, 5.84960938e-01,}, /* 0x3f15c000 */
+ dp_l[] = { 0.0, 1.56322085e-06,}, /* 0x35d1cfdc */
+--- a/e_rem_pio2f.cpp 2022-12-11 21:20:40.290068458 -0500
++++ b/e_rem_pio2f.cpp_new 2022-12-11 21:23:32.015881905 -0500
+@@ -10,16 +10,16 @@
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /* __ieee754_rem_pio2f(x,y)
+ *
+ * return the remainder of x rem pi/2 in *y
+ * use double precision for everything except passing x
+ * use __kernel_rem_pio2() for large x
+ */
+--- b/k_cosf.cpp 2022-12-11 21:20:40.290068458 -0500
++++ b/k_cosf.cpp_new 2022-12-11 21:23:32.051882266 -0500
+@@ -11,16 +11,16 @@
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+ #ifndef INLINE_KERNEL_COSDF
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+ #endif
+
+ #include "math_private.h"
+
+ /* |cos(x) - c(x)| < 2**-34.1 (~[-5.37e-11, 5.295e-11]). */
+ static const double
+ one = 1.0,
+--- b/k_expf.cpp 2022-12-11 21:20:40.290068458 -0500
++++ b/k_expf.cpp_new 2022-12-11 21:23:32.087882627 -0500
+@@ -22,16 +22,16 @@
+ * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
+ * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+ * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
+ * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
+ * SUCH DAMAGE.
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ #include <complex.h>
+
+ #include "math_private.h"
+
+ static const uint32_t k = 235; /* constant for reduction */
+ static const float kln2 = 162.88958740F; /* k * ln2 */
+--- b/k_sinf.cpp 2022-12-11 21:20:40.290068458 -0500
++++ b/k_sinf.cpp_new 2022-12-11 21:23:32.119882948 -0500
+@@ -11,16 +11,16 @@
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+ #ifndef INLINE_KERNEL_SINDF
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+ #endif
+
+ #include "math_private.h"
+
+ /* |sin(x)/x - s(x)| < 2**-37.5 (~[-4.89e-12, 4.824e-12]). */
+ static const double
+ S1 = -0x15555554cbac77.0p-55, /* -0.166666666416265235595 */
+--- b/k_tanf.cpp 2022-12-11 21:20:40.290068458 -0500
++++ b/k_tanf.cpp_new 2022-12-11 21:23:32.155883308 -0500
+@@ -10,16 +10,16 @@
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+ #ifndef INLINE_KERNEL_TANDF
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+ #endif
+
+ #include "math_private.h"
+
+ /* |tan(x)/x - t(x)| < 2**-25.5 (~[-2e-08, 2e-08]). */
+ static const double
+ T[] = {
+--- a/s_atanf.cpp 2022-12-11 21:20:40.290068458 -0500
++++ b/s_atanf.cpp_new 2022-12-11 21:23:32.199883749 -0500
+@@ -9,16 +9,16 @@
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ #include "math_private.h"
+
+ static const float atanhi[] = {
+ 4.6364760399e-01, /* atan(0.5)hi 0x3eed6338 */
+ 7.8539812565e-01, /* atan(1.0)hi 0x3f490fda */
+ 9.8279368877e-01, /* atan(1.5)hi 0x3f7b985e */
+--- a/s_cosf.cpp 2022-12-11 21:20:40.290068458 -0500
++++ b/s_cosf.cpp_new 2022-12-11 21:23:32.235884110 -0500
+@@ -10,16 +10,16 @@
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ #include <float.h>
+
+ #define INLINE_KERNEL_COSDF
+ #define INLINE_KERNEL_SINDF
+ #define INLINE_REM_PIO2F
+ #include "math_private.h"
+--- a/s_exp2.cpp 2022-12-11 21:20:40.290068458 -0500
++++ b/s_exp2.cpp_new 2022-12-11 21:23:32.275884511 -0500
+@@ -22,16 +22,16 @@
+ * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
+ * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+ * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
+ * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
+ * SUCH DAMAGE.
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ #include <float.h>
+
+ #include "math_private.h"
+
+ #define TBLBITS 8
+ #define TBLSIZE (1 << TBLBITS)
+--- a/s_exp2f.cpp 2022-12-11 21:20:40.290068458 -0500
++++ b/s_exp2f.cpp_new 2022-12-11 21:23:32.311884872 -0500
+@@ -22,16 +22,16 @@
+ * OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
+ * HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
+ * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
+ * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
+ * SUCH DAMAGE.
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ #include <float.h>
+
+ #include "math_private.h"
+
+ #define TBLBITS 4
+ #define TBLSIZE (1 << TBLBITS)
+--- a/s_fabsf.cpp 2022-12-11 21:20:40.290068458 -0500
++++ b/s_fabsf.cpp_new 2022-12-11 21:23:32.347885233 -0500
+@@ -9,16 +9,16 @@
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ /*
+ * fabsf(x) returns the absolute value of x.
+ */
+
+ #include "math_private.h"
+
+--- a/s_sinf.cpp 2022-12-11 21:20:40.290068458 -0500
++++ b/s_sinf.cpp_new 2022-12-11 21:23:32.383885594 -0500
+@@ -10,16 +10,16 @@
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ #include <float.h>
+
+ #define INLINE_KERNEL_COSDF
+ #define INLINE_KERNEL_SINDF
+ #define INLINE_REM_PIO2F
+ #include "math_private.h"
+--- a/s_tanf.cpp 2022-12-11 21:20:40.290068458 -0500
++++ b/s_tanf.cpp_new 2022-12-11 21:23:32.419885954 -0500
+@@ -10,16 +10,16 @@
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#include <sys/cdefs.h>
+-__FBSDID("$FreeBSD$");
++//#include <sys/cdefs.h>
++//__FBSDID("$FreeBSD$");
+
+ #include <float.h>
+
+ #define INLINE_KERNEL_TANDF
+ #define INLINE_REM_PIO2F
+ #include "math_private.h"
+ #include "e_rem_pio2f.c"
+--- a/e_sqrtf.cpp 2022-12-13 14:30:13.113415455 -0500
++++ b/e_sqrtf.cpp 2022-12-13 14:13:15.302434550 -0500
+@@ -8,19 +8,19 @@
+ *
+ * Developed at SunPro, a Sun Microsystems, Inc. business.
+ * Permission to use, copy, modify, and distribute this
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+-#ifndef lint
+-static char rcsid[] = "$FreeBSD$";
+-#endif
++//#ifndef lint
++//static char rcsid[] = "$FreeBSD$";
++//#endif
+
+ #include "math_private.h"
+
+ #ifdef USE_BUILTIN_SQRTF
+ float
+ __ieee754_sqrtf(float x)
+ {
+ return (__builtin_sqrtf(x));
diff --git a/modules/fdlibm/patches/10_remove_unused_function_from_k_exp_cpp.patch b/modules/fdlibm/patches/10_remove_unused_function_from_k_exp_cpp.patch
new file mode 100644
index 0000000000..0cd417bcd2
--- /dev/null
+++ b/modules/fdlibm/patches/10_remove_unused_function_from_k_exp_cpp.patch
@@ -0,0 +1,107 @@
+diff --git a/k_exp.cpp b/k_exp.cpp
+--- a/k_exp.cpp
++++ b/k_exp.cpp
+@@ -24,18 +24,16 @@
+ * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
+ * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
+ * SUCH DAMAGE.
+ */
+
+ //#include <sys/cdefs.h>
+ //__FBSDID("$FreeBSD$");
+
+-#include <complex.h>
+-
+ #include "math_private.h"
+
+ static const uint32_t k = 1799; /* constant for reduction */
+ static const double kln2 = 1246.97177782734161156; /* k * ln2 */
+
+ /*
+ * Compute exp(x), scaled to avoid spurious overflow. An exponent is
+ * returned separately in 'expt'.
+@@ -78,33 +76,8 @@ __ldexp_exp(double x, int expt)
+ double exp_x, scale;
+ int ex_expt;
+
+ exp_x = __frexp_exp(x, &ex_expt);
+ expt += ex_expt;
+ INSERT_WORDS(scale, (0x3ff + expt) << 20, 0);
+ return (exp_x * scale);
+ }
+-
+-double complex
+-__ldexp_cexp(double complex z, int expt)
+-{
+- double c, exp_x, s, scale1, scale2, x, y;
+- int ex_expt, half_expt;
+-
+- x = creal(z);
+- y = cimag(z);
+- exp_x = __frexp_exp(x, &ex_expt);
+- expt += ex_expt;
+-
+- /*
+- * Arrange so that scale1 * scale2 == 2**expt. We use this to
+- * compensate for scalbn being horrendously slow.
+- */
+- half_expt = expt / 2;
+- INSERT_WORDS(scale1, (0x3ff + half_expt) << 20, 0);
+- half_expt = expt - half_expt;
+- INSERT_WORDS(scale2, (0x3ff + half_expt) << 20, 0);
+-
+- sincos(y, &s, &c);
+- return (CMPLX(c * exp_x * scale1 * scale2,
+- s * exp_x * scale1 * scale2));
+-}
+--- a/k_expf.cpp 2022-12-13 14:39:11.791594163 -0500
++++ b/k_expf.cpp 2022-12-13 14:39:18.631622916 -0500
+@@ -24,18 +24,16 @@
+ * LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
+ * OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF
+ * SUCH DAMAGE.
+ */
+
+ //#include <sys/cdefs.h>
+ //__FBSDID("$FreeBSD$");
+
+-#include <complex.h>
+-
+ #include "math_private.h"
+
+ static const uint32_t k = 235; /* constant for reduction */
+ static const float kln2 = 162.88958740F; /* k * ln2 */
+
+ /*
+ * See k_exp.c for details.
+ *
+@@ -61,29 +59,8 @@
+ float exp_x, scale;
+ int ex_expt;
+
+ exp_x = __frexp_expf(x, &ex_expt);
+ expt += ex_expt;
+ SET_FLOAT_WORD(scale, (0x7f + expt) << 23);
+ return (exp_x * scale);
+ }
+-
+-float complex
+-__ldexp_cexpf(float complex z, int expt)
+-{
+- float c, exp_x, s, scale1, scale2, x, y;
+- int ex_expt, half_expt;
+-
+- x = crealf(z);
+- y = cimagf(z);
+- exp_x = __frexp_expf(x, &ex_expt);
+- expt += ex_expt;
+-
+- half_expt = expt / 2;
+- SET_FLOAT_WORD(scale1, (0x7f + half_expt) << 23);
+- half_expt = expt - half_expt;
+- SET_FLOAT_WORD(scale2, (0x7f + half_expt) << 23);
+-
+- sincosf(y, &s, &c);
+- return (CMPLXF(c * exp_x * scale1 * scale2,
+- s * exp_x * scale1 * scale2));
+-}
diff --git a/modules/fdlibm/patches/11_include_cfloat_to_use_flt_eval_method.patch b/modules/fdlibm/patches/11_include_cfloat_to_use_flt_eval_method.patch
new file mode 100644
index 0000000000..a2b7d959a6
--- /dev/null
+++ b/modules/fdlibm/patches/11_include_cfloat_to_use_flt_eval_method.patch
@@ -0,0 +1,21 @@
+diff --git a/math_private.h b/math_private.h
+--- a/math_private.h
++++ b/math_private.h
+@@ -12,16 +12,17 @@
+ /*
+ * from: @(#)fdlibm.h 5.1 93/09/24
+ * $FreeBSD$
+ */
+
+ #ifndef _MATH_PRIVATE_H_
+ #define _MATH_PRIVATE_H_
+
++#include <cfloat>
+ #include <stdint.h>
+ #include <sys/types.h>
+
+ #include "mozilla/EndianUtils.h"
+
+ #include "fdlibm.h"
+
+ /*
diff --git a/modules/fdlibm/patches/12_define_u_int32_t_and_u_int64_t_on_windows.patch b/modules/fdlibm/patches/12_define_u_int32_t_and_u_int64_t_on_windows.patch
new file mode 100644
index 0000000000..895509c8d2
--- /dev/null
+++ b/modules/fdlibm/patches/12_define_u_int32_t_and_u_int64_t_on_windows.patch
@@ -0,0 +1,27 @@
+diff --git a/math_private.h b/math_private.h
+--- a/math_private.h
++++ b/math_private.h
+@@ -33,16 +33,23 @@
+ * to dig two 32 bit words out of the 64 bit IEEE floating point
+ * value. That is non-ANSI, and, moreover, the gcc instruction
+ * scheduler gets it wrong. We instead use the following macros.
+ * Unlike the original code, we determine the endianness at compile
+ * time, not at run time; I don't see much benefit to selecting
+ * endianness at run time.
+ */
+
++#ifndef u_int32_t
++#define u_int32_t uint32_t
++#endif
++#ifndef u_int64_t
++#define u_int64_t uint64_t
++#endif
++
+ /* A union which permits us to convert between a long double and
+ four 32 bit ints. */
+
+ #if MOZ_BIG_ENDIAN()
+
+ typedef union
+ {
+ long double value;
diff --git a/modules/fdlibm/patches/13_define_strict_assign_even_if_flt_eval_method_is_not_defined.patch b/modules/fdlibm/patches/13_define_strict_assign_even_if_flt_eval_method_is_not_defined.patch
new file mode 100644
index 0000000000..a0417602d5
--- /dev/null
+++ b/modules/fdlibm/patches/13_define_strict_assign_even_if_flt_eval_method_is_not_defined.patch
@@ -0,0 +1,31 @@
+diff --git a/math_private.h b/math_private.h
+--- a/math_private.h
++++ b/math_private.h
+@@ -323,16 +323,27 @@ do { \
+ if (sizeof(type) >= sizeof(long double)) \
+ (lval) = (rval); \
+ else { \
+ __lval = (rval); \
+ (lval) = __lval; \
+ } \
+ } while (0)
+ #endif
++#else
++#define STRICT_ASSIGN(type, lval, rval) do { \
++ volatile type __lval; \
++ \
++ if (sizeof(type) >= sizeof(long double)) \
++ (lval) = (rval); \
++ else { \
++ __lval = (rval); \
++ (lval) = __lval; \
++ } \
++} while (0)
+ #endif /* FLT_EVAL_METHOD */
+
+ /* Support switching the mode to FP_PE if necessary. */
+ #if defined(__i386__) && !defined(NO_FPSETPREC)
+ #define ENTERI() ENTERIT(long double)
+ #define ENTERIT(returntype) \
+ returntype __retval; \
+ fp_prec_t __oprec; \
diff --git a/modules/fdlibm/patches/14_do_not_use_hexadecimal_floating_point_number.patch b/modules/fdlibm/patches/14_do_not_use_hexadecimal_floating_point_number.patch
new file mode 100644
index 0000000000..151a7a3b11
--- /dev/null
+++ b/modules/fdlibm/patches/14_do_not_use_hexadecimal_floating_point_number.patch
@@ -0,0 +1,68 @@
+diff --git a/e_exp.cpp b/e_exp.cpp
+--- a/e_exp.cpp
++++ b/e_exp.cpp
+@@ -146,14 +146,17 @@ __ieee754_exp(double x) /* default IEEE double exp */
+ if(k >= -1021)
+ INSERT_WORDS(twopk,((u_int32_t)(0x3ff+k))<<20, 0);
+ else
+ INSERT_WORDS(twopk,((u_int32_t)(0x3ff+(k+1000)))<<20, 0);
+ c = x - t*(P1+t*(P2+t*(P3+t*(P4+t*P5))));
+ if(k==0) return one-((x*c)/(c-2.0)-x);
+ else y = one-((lo-(x*c)/(2.0-c))-hi);
+ if(k >= -1021) {
+- if (k==1024) return y*2.0*0x1p1023;
++ if (k==1024) {
++ double const_0x1p1023 = pow(2, 1023);
++ return y*2.0*const_0x1p1023;
++ }
+ return y*twopk;
+ } else {
+ return y*twopk*twom1000;
+ }
+ }
+diff --git a/s_expm1.cpp b/s_expm1.cpp
+--- a/s_expm1.cpp
++++ b/s_expm1.cpp
+@@ -192,17 +192,20 @@ expm1(double x)
+ e -= hxs;
+ if(k== -1) return 0.5*(x-e)-0.5;
+ if(k==1) {
+ if(x < -0.25) return -2.0*(e-(x+0.5));
+ else return one+2.0*(x-e);
+ }
+ if (k <= -2 || k>56) { /* suffice to return exp(x)-1 */
+ y = one-(e-x);
+- if (k == 1024) y = y*2.0*0x1p1023;
++ if (k == 1024) {
++ double const_0x1p1023 = pow(2, 1023);
++ y = y*2.0*const_0x1p1023;
++ }
+ else y = y*twopk;
+ return y-one;
+ }
+ t = one;
+ if(k<20) {
+ SET_HIGH_WORD(t,0x3ff00000 - (0x200000>>k)); /* t=1-2^-k */
+ y = t-(e-x);
+ y = y*twopk;
+--- a/s_exp2.cpp 2022-12-11 21:26:36.045643147 -0500
++++ b/s_exp2.cpp_new 2022-12-11 21:28:42.058769925 -0500
+@@ -381,14 +381,16 @@
+ INSERT_WORDS(twopk, 0x3ff00000 + k, 0);
+ else
+ INSERT_WORDS(twopkp1000, 0x3ff00000 + k + (1000 << 20), 0);
+ r = t + t * z * (P1 + z * (P2 + z * (P3 + z * (P4 + z * P5))));
+
+ /* Scale by 2**(k>>20). */
+ if(k >= -(1021 << 20)) {
+- if (k == 1024 << 20)
+- return (r * 2.0 * 0x1p1023);
++ if (k == 1024 << 20) {
++ double const_0x1p1023 = pow(2, 1023);
++ return (r * 2.0 * const_0x1p1023);
++ }
+ return (r * twopk);
+ } else {
+ return (r * twopkp1000 * twom1000);
+ }
+ }
diff --git a/modules/fdlibm/patches/15_remove_unused_rintl_function_from_s_nearbyint_cpp.patch b/modules/fdlibm/patches/15_remove_unused_rintl_function_from_s_nearbyint_cpp.patch
new file mode 100644
index 0000000000..8c4f45a75a
--- /dev/null
+++ b/modules/fdlibm/patches/15_remove_unused_rintl_function_from_s_nearbyint_cpp.patch
@@ -0,0 +1,13 @@
+diff --git a/s_nearbyint.cpp b/s_nearbyint.cpp
+--- a/s_nearbyint.cpp
++++ b/s_nearbyint.cpp
+@@ -54,9 +54,8 @@ fn(type x) \
+ fegetenv(&env); \
+ ret = rint(x); \
+ fesetenv(&env); \
+ return (ret); \
+ }
+
+ DECL(double, nearbyint, rint)
+ DECL(float, nearbyintf, rintf)
+-DECL(long double, nearbyintl, rintl)
diff --git a/modules/fdlibm/patches/16_use_safer_strict_assign_on_visual_studio.patch b/modules/fdlibm/patches/16_use_safer_strict_assign_on_visual_studio.patch
new file mode 100644
index 0000000000..9f834149c6
--- /dev/null
+++ b/modules/fdlibm/patches/16_use_safer_strict_assign_on_visual_studio.patch
@@ -0,0 +1,22 @@
+diff --git a/math_private.h b/math_private.h
+--- a/math_private.h
++++ b/math_private.h
+@@ -309,17 +309,17 @@ do { \
+ /* The above works on non-i386 too, but we use this to check v. */
+ #define LD80C(m, ex, v) { .e = (v), }
+ #endif
+
+ #ifdef FLT_EVAL_METHOD
+ /*
+ * Attempt to get strict C99 semantics for assignment with non-C99 compilers.
+ */
+-#if FLT_EVAL_METHOD == 0 || __GNUC__ == 0
++#if !defined(_MSC_VER) && (FLT_EVAL_METHOD == 0 || __GNUC__ == 0)
+ #define STRICT_ASSIGN(type, lval, rval) ((lval) = (rval))
+ #else
+ #define STRICT_ASSIGN(type, lval, rval) do { \
+ volatile type __lval; \
+ \
+ if (sizeof(type) >= sizeof(long double)) \
+ (lval) = (rval); \
+ else { \
diff --git a/modules/fdlibm/patches/17_exp_exact_result_for_positive_one.patch b/modules/fdlibm/patches/17_exp_exact_result_for_positive_one.patch
new file mode 100644
index 0000000000..ada59b23a2
--- /dev/null
+++ b/modules/fdlibm/patches/17_exp_exact_result_for_positive_one.patch
@@ -0,0 +1,40 @@
+diff --git a/e_exp.cpp b/e_exp.cpp
+--- a/e_exp.cpp
++++ b/e_exp.cpp
+@@ -91,16 +91,18 @@ ln2LO[2] ={ 1.90821492927058770002e-10, /* 0x3dea39ef, 0x35793c76 */
+ -1.90821492927058770002e-10,},/* 0xbdea39ef, 0x35793c76 */
+ invln2 = 1.44269504088896338700e+00, /* 0x3ff71547, 0x652b82fe */
+ P1 = 1.66666666666666019037e-01, /* 0x3FC55555, 0x5555553E */
+ P2 = -2.77777777770155933842e-03, /* 0xBF66C16C, 0x16BEBD93 */
+ P3 = 6.61375632143793436117e-05, /* 0x3F11566A, 0xAF25DE2C */
+ P4 = -1.65339022054652515390e-06, /* 0xBEBBBD41, 0xC5D26BF1 */
+ P5 = 4.13813679705723846039e-08; /* 0x3E663769, 0x72BEA4D0 */
+
++static const double E = 2.7182818284590452354; /* e */
++
+ static volatile double
+ huge = 1.0e+300,
+ twom1000= 9.33263618503218878990e-302; /* 2**-1000=0x01700000,0*/
+
+ double
+ __ieee754_exp(double x) /* default IEEE double exp */
+ {
+ double y,hi=0.0,lo=0.0,c,t,twopk;
+@@ -122,16 +124,17 @@ __ieee754_exp(double x) /* default IEEE double exp */
+ }
+ if(x > o_threshold) return huge*huge; /* overflow */
+ if(x < u_threshold) return twom1000*twom1000; /* underflow */
+ }
+
+ /* argument reduction */
+ if(hx > 0x3fd62e42) { /* if |x| > 0.5 ln2 */
+ if(hx < 0x3FF0A2B2) { /* and |x| < 1.5 ln2 */
++ if (x == 1.0) return E;
+ hi = x-ln2HI[xsb]; lo=ln2LO[xsb]; k = 1-xsb-xsb;
+ } else {
+ k = (int)(invln2*x+halF[xsb]);
+ t = k;
+ hi = x - t*ln2HI[0]; /* t*ln2HI is exact here */
+ lo = t*ln2LO[0];
+ }
+ STRICT_ASSIGN(double, x, hi - lo);
diff --git a/modules/fdlibm/patches/18_use_stdlib_sqrt.patch b/modules/fdlibm/patches/18_use_stdlib_sqrt.patch
new file mode 100644
index 0000000000..a4f9f3f003
--- /dev/null
+++ b/modules/fdlibm/patches/18_use_stdlib_sqrt.patch
@@ -0,0 +1,289 @@
+diff --git a/e_acos.cpp b/e_acos.cpp
+--- a/e_acos.cpp
++++ b/e_acos.cpp
+@@ -33,16 +33,17 @@
+ *
+ * Special cases:
+ * if x is NaN, return x itself;
+ * if |x|>1, return NaN with invalid signal.
+ *
+ * Function needed: sqrt
+ */
+
++#include <cmath>
+ #include <float.h>
+
+ #include "math_private.h"
+
+ static const double
+ one= 1.00000000000000000000e+00, /* 0x3FF00000, 0x00000000 */
+ pi = 3.14159265358979311600e+00, /* 0x400921FB, 0x54442D18 */
+ pio2_hi = 1.57079632679489655800e+00; /* 0x3FF921FB, 0x54442D18 */
+@@ -82,23 +83,23 @@ __ieee754_acos(double x)
+ p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));
+ q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4)));
+ r = p/q;
+ return pio2_hi - (x - (pio2_lo-x*r));
+ } else if (hx<0) { /* x < -0.5 */
+ z = (one+x)*0.5;
+ p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));
+ q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4)));
+- s = sqrt(z);
++ s = std::sqrt(z);
+ r = p/q;
+ w = r*s-pio2_lo;
+ return pi - 2.0*(s+w);
+ } else { /* x > 0.5 */
+ z = (one-x)*0.5;
+- s = sqrt(z);
++ s = std::sqrt(z);
+ df = s;
+ SET_LOW_WORD(df,0);
+ c = (z-df*df)/(s+df);
+ p = z*(pS0+z*(pS1+z*(pS2+z*(pS3+z*(pS4+z*pS5)))));
+ q = one+z*(qS1+z*(qS2+z*(qS3+z*qS4)));
+ r = p/q;
+ w = r*s+c;
+ return 2.0*(df+w);
+diff --git a/e_acosh.cpp b/e_acosh.cpp
+--- a/e_acosh.cpp
++++ b/e_acosh.cpp
+@@ -24,16 +24,17 @@
+ * acosh(x) := log(2x-1/(sqrt(x*x-1)+x)) if x>2; else
+ * acosh(x) := log1p(t+sqrt(2.0*t+t*t)); where t=x-1.
+ *
+ * Special cases:
+ * acosh(x) is NaN with signal if x<1.
+ * acosh(NaN) is NaN without signal.
+ */
+
++#include <cmath>
+ #include <float.h>
+
+ #include "math_private.h"
+
+ static const double
+ one = 1.0,
+ ln2 = 6.93147180559945286227e-01; /* 0x3FE62E42, 0xFEFA39EF */
+
+@@ -50,14 +51,14 @@ __ieee754_acosh(double x)
+ if(hx >=0x7ff00000) { /* x is inf of NaN */
+ return x+x;
+ } else
+ return __ieee754_log(x)+ln2; /* acosh(huge)=log(2x) */
+ } else if(((hx-0x3ff00000)|lx)==0) {
+ return 0.0; /* acosh(1) = 0 */
+ } else if (hx > 0x40000000) { /* 2**28 > x > 2 */
+ t=x*x;
+- return __ieee754_log(2.0*x-one/(x+sqrt(t-one)));
++ return __ieee754_log(2.0*x-one/(x+std::sqrt(t-one)));
+ } else { /* 1<x<2 */
+ t = x-one;
+- return log1p(t+sqrt(2.0*t+t*t));
++ return log1p(t+std::sqrt(2.0*t+t*t));
+ }
+ }
+diff --git a/e_asin.cpp b/e_asin.cpp
+--- a/e_asin.cpp
++++ b/e_asin.cpp
+@@ -39,16 +39,17 @@
+ * = pio4_hi+(pio4-2f)-(2s*z*R(z)-(pio2_lo+2c))
+ *
+ * Special cases:
+ * if x is NaN, return x itself;
+ * if |x|>1, return NaN with invalid signal.
+ *
+ */
+
++#include <cmath>
+ #include <float.h>
+
+ #include "math_private.h"
+
+ static const double
+ one = 1.00000000000000000000e+00, /* 0x3FF00000, 0x00000000 */
+ huge = 1.000e+300,
+ pio2_hi = 1.57079632679489655800e+00, /* 0x3FF921FB, 0x54442D18 */
+@@ -90,17 +91,17 @@ __ieee754_asin(double x)
+ w = p/q;
+ return x+x*w;
+ }
+ /* 1> |x|>= 0.5 */
+ w = one-fabs(x);
+ t = w*0.5;
+ p = t*(pS0+t*(pS1+t*(pS2+t*(pS3+t*(pS4+t*pS5)))));
+ q = one+t*(qS1+t*(qS2+t*(qS3+t*qS4)));
+- s = sqrt(t);
++ s = std::sqrt(t);
+ if(ix>=0x3FEF3333) { /* if |x| > 0.975 */
+ w = p/q;
+ t = pio2_hi-(2.0*(s+s*w)-pio2_lo);
+ } else {
+ w = s;
+ SET_LOW_WORD(w,0);
+ c = (t-w*w)/(s+w);
+ r = p/q;
+diff --git a/e_hypot.cpp b/e_hypot.cpp
+--- a/e_hypot.cpp
++++ b/e_hypot.cpp
+@@ -41,16 +41,17 @@
+ * hypot(x,y) is INF if x or y is +INF or -INF; else
+ * hypot(x,y) is NAN if x or y is NAN.
+ *
+ * Accuracy:
+ * hypot(x,y) returns sqrt(x^2+y^2) with error less
+ * than 1 ulps (units in the last place)
+ */
+
++#include <cmath>
+ #include <float.h>
+
+ #include "math_private.h"
+
+ double
+ __ieee754_hypot(double x, double y)
+ {
+ double a,b,t1,t2,y1,y2,w;
+@@ -100,25 +101,25 @@ __ieee754_hypot(double x, double y)
+ }
+ }
+ /* medium size a and b */
+ w = a-b;
+ if (w>b) {
+ t1 = 0;
+ SET_HIGH_WORD(t1,ha);
+ t2 = a-t1;
+- w = sqrt(t1*t1-(b*(-b)-t2*(a+t1)));
++ w = std::sqrt(t1*t1-(b*(-b)-t2*(a+t1)));
+ } else {
+ a = a+a;
+ y1 = 0;
+ SET_HIGH_WORD(y1,hb);
+ y2 = b - y1;
+ t1 = 0;
+ SET_HIGH_WORD(t1,ha+0x00100000);
+ t2 = a - t1;
+- w = sqrt(t1*y1-(w*(-w)-(t1*y2+t2*b)));
++ w = std::sqrt(t1*y1-(w*(-w)-(t1*y2+t2*b)));
+ }
+ if(k!=0) {
+ t1 = 0.0;
+ SET_HIGH_WORD(t1,(1023+k)<<20);
+ return t1*w;
+ } else return w;
+ }
+diff --git a/e_pow.cpp b/e_pow.cpp
+--- a/e_pow.cpp
++++ b/e_pow.cpp
+@@ -52,16 +52,17 @@
+ *
+ * Constants :
+ * The hexadecimal values are the intended ones for the following
+ * constants. The decimal values may be used, provided that the
+ * compiler will convert from decimal to binary accurately enough
+ * to produce the hexadecimal values shown.
+ */
+
++#include <cmath>
+ #include <float.h>
+ #include "math_private.h"
+
+ static const double
+ bp[] = {1.0, 1.5,},
+ dp_h[] = { 0.0, 5.84962487220764160156e-01,}, /* 0x3FE2B803, 0x40000000 */
+ dp_l[] = { 0.0, 1.35003920212974897128e-08,}, /* 0x3E4CFDEB, 0x43CFD006 */
+ zero = 0.0,
+@@ -151,17 +152,17 @@ __ieee754_pow(double x, double y)
+ return (hy<0)?-y: zero;
+ }
+ if(iy==0x3ff00000) { /* y is +-1 */
+ if(hy<0) return one/x; else return x;
+ }
+ if(hy==0x40000000) return x*x; /* y is 2 */
+ if(hy==0x3fe00000) { /* y is 0.5 */
+ if(hx>=0) /* x >= +0 */
+- return sqrt(x);
++ return std::sqrt(x);
+ }
+ }
+
+ ax = fabs(x);
+ /* special value of x */
+ if(lx==0) {
+ if(ix==0x7ff00000||ix==0||ix==0x3ff00000){
+ z = ax; /*x is +-0,+-inf,+-1*/
+diff --git a/s_asinh.cpp b/s_asinh.cpp
+--- a/s_asinh.cpp
++++ b/s_asinh.cpp
+@@ -19,16 +19,17 @@
+ * asinh(x) = sign(x) * log [ |x| + sqrt(x*x+1) ]
+ * we have
+ * asinh(x) := x if 1+x*x=1,
+ * := sign(x)*(log(x)+ln2)) for large |x|, else
+ * := sign(x)*log(2|x|+1/(|x|+sqrt(x*x+1))) if|x|>2, else
+ * := sign(x)*log1p(|x| + x^2/(1 + sqrt(1+x^2)))
+ */
+
++#include <cmath>
+ #include <float.h>
+
+ #include "math_private.h"
+
+ static const double
+ one = 1.00000000000000000000e+00, /* 0x3FF00000, 0x00000000 */
+ ln2 = 6.93147180559945286227e-01, /* 0x3FE62E42, 0xFEFA39EF */
+ huge= 1.00000000000000000000e+300;
+@@ -43,15 +44,15 @@ asinh(double x)
+ if(ix>=0x7ff00000) return x+x; /* x is inf or NaN */
+ if(ix< 0x3e300000) { /* |x|<2**-28 */
+ if(huge+x>one) return x; /* return x inexact except 0 */
+ }
+ if(ix>0x41b00000) { /* |x| > 2**28 */
+ w = __ieee754_log(fabs(x))+ln2;
+ } else if (ix>0x40000000) { /* 2**28 > |x| > 2.0 */
+ t = fabs(x);
+- w = __ieee754_log(2.0*t+one/(__ieee754_sqrt(x*x+one)+t));
++ w = __ieee754_log(2.0*t+one/(std::sqrt(x*x+one)+t));
+ } else { /* 2.0 > |x| > 2**-28 */
+ t = x*x;
+- w =log1p(fabs(x)+t/(one+__ieee754_sqrt(one+t)));
++ w =log1p(fabs(x)+t/(one+std::sqrt(one+t)));
+ }
+ if(hx>0) return w; else return -w;
+ }
+--- a/e_asinf.cpp 2022-12-13 14:45:17.953154257 -0500
++++ b/e_asinf.cpp 2022-12-13 14:45:03.425091710 -0500
+@@ -11,16 +11,18 @@
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+ //#include <sys/cdefs.h>
+ //__FBSDID("$FreeBSD$");
+
++#include <cmath>
++
+ #include "math_private.h"
+
+ static const float
+ one = 1.0000000000e+00, /* 0x3F800000 */
+ huge = 1.000e+30,
+ /* coefficient for R(x^2) */
+ pS0 = 1.6666586697e-01,
+ pS1 = -4.2743422091e-02,
+@@ -52,13 +54,13 @@
+ w = p/q;
+ return x+x*w;
+ }
+ /* 1> |x|>= 0.5 */
+ w = one-fabsf(x);
+ t = w*(float)0.5;
+ p = t*(pS0+t*(pS1+t*pS2));
+ q = one+t*qS1;
+- s = sqrt(t);
++ s = std::sqrt(t);
+ w = p/q;
+ t = pio2-2.0*(s+s*w);
+ if(hx>0) return t; else return -t;
+ }
diff --git a/modules/fdlibm/patches/19_remove_unneeded_round_to_integer_helpers.patch b/modules/fdlibm/patches/19_remove_unneeded_round_to_integer_helpers.patch
new file mode 100644
index 0000000000..8720d8c0ed
--- /dev/null
+++ b/modules/fdlibm/patches/19_remove_unneeded_round_to_integer_helpers.patch
@@ -0,0 +1,50 @@
+diff --git a/math_private.h b/math_private.h
+--- a/math_private.h
++++ b/math_private.h
+@@ -612,46 +612,16 @@ rnint(__double_t x)
+ * magic number would need to be variable. Assuming that the
+ * rounding precision is always the default is too fragile. This
+ * and many other complications will move when the default is
+ * changed to FP_PE.
+ */
+ return ((double)(x + 0x1.8p52) - 0x1.8p52);
+ }
+
+-static inline float
+-rnintf(__float_t x)
+-{
+- /*
+- * As for rnint(), except we could just call that to handle the
+- * extra precision case, usually without losing efficiency.
+- */
+- return ((float)(x + 0x1.8p23F) - 0x1.8p23F);
+-}
+-
+-#ifdef LDBL_MANT_DIG
+-/*
+- * The complications for extra precision are smaller for rnintl() since it
+- * can safely assume that the rounding precision has been increased from
+- * its default to FP_PE on x86. We don't exploit that here to get small
+- * optimizations from limiting the rangle to double. We just need it for
+- * the magic number to work with long doubles. ld128 callers should use
+- * rnint() instead of this if possible. ld80 callers should prefer
+- * rnintl() since for amd64 this avoids swapping the register set, while
+- * for i386 it makes no difference (assuming FP_PE), and for other arches
+- * it makes little difference.
+- */
+-static inline long double
+-rnintl(long double x)
+-{
+- return (x + __CONCAT(0x1.8p, LDBL_MANT_DIG) / 2 -
+- __CONCAT(0x1.8p, LDBL_MANT_DIG) / 2);
+-}
+-#endif /* LDBL_MANT_DIG */
+-
+ /*
+ * irint() and i64rint() give the same result as casting to their integer
+ * return type provided their arg is a floating point integer. They can
+ * sometimes be more efficient because no rounding is required.
+ */
+ #if defined(amd64) || defined(__i386__)
+ #define irint(x) \
+ (sizeof(x) == sizeof(float) && \
diff --git a/modules/fdlibm/patches/20_emulate_freebsd_internal_double_types.patch b/modules/fdlibm/patches/20_emulate_freebsd_internal_double_types.patch
new file mode 100644
index 0000000000..d22073391e
--- /dev/null
+++ b/modules/fdlibm/patches/20_emulate_freebsd_internal_double_types.patch
@@ -0,0 +1,28 @@
+diff --git a/math_private.h b/math_private.h
+--- a/math_private.h
++++ b/math_private.h
+@@ -20,16 +20,24 @@
+ #include <cfloat>
+ #include <stdint.h>
+ #include <sys/types.h>
+
+ #include "mozilla/EndianUtils.h"
+
+ #include "fdlibm.h"
+
++/*
++ * Emulate FreeBSD internal double types.
++ * Adapted from https://github.com/freebsd/freebsd-src/search?q=__double_t
++ */
++
++typedef double __double_t;
++typedef __double_t double_t;
++
+ /*
+ * The original fdlibm code used statements like:
+ * n0 = ((*(int*)&one)>>29)^1; * index of high word *
+ * ix0 = *(n0+(int*)&x); * high word of x *
+ * ix1 = *((1-n0)+(int*)&x); * low word of x *
+ * to dig two 32 bit words out of the 64 bit IEEE floating point
+ * value. That is non-ANSI, and, moreover, the gcc instruction
+ * scheduler gets it wrong. We instead use the following macros.
diff --git a/modules/fdlibm/patches/21_rem_pio_fixups.patch b/modules/fdlibm/patches/21_rem_pio_fixups.patch
new file mode 100644
index 0000000000..8961db62c3
--- /dev/null
+++ b/modules/fdlibm/patches/21_rem_pio_fixups.patch
@@ -0,0 +1,226 @@
+diff --git a/e_rem_pio2.cpp b/e_rem_pio2.cpp
+--- a/e_rem_pio2.cpp
++++ b/e_rem_pio2.cpp
+@@ -43,17 +43,17 @@ invpio2 = 6.36619772367581382433e-01, /* 0x3FE45F30, 0x6DC9C883 */
+ pio2_1 = 1.57079632673412561417e+00, /* 0x3FF921FB, 0x54400000 */
+ pio2_1t = 6.07710050650619224932e-11, /* 0x3DD0B461, 0x1A626331 */
+ pio2_2 = 6.07710050630396597660e-11, /* 0x3DD0B461, 0x1A600000 */
+ pio2_2t = 2.02226624879595063154e-21, /* 0x3BA3198A, 0x2E037073 */
+ pio2_3 = 2.02226624871116645580e-21, /* 0x3BA3198A, 0x2E000000 */
+ pio2_3t = 8.47842766036889956997e-32; /* 0x397B839A, 0x252049C1 */
+
+ #ifdef INLINE_REM_PIO2
+-static __inline __always_inline
++static inline
+ #endif
+ int
+ __ieee754_rem_pio2(double x, double *y)
+ {
+ double z,w,t,r,fn;
+ double tx[3],ty[2];
+ int32_t e0,i,j,nx,n,ix,hx;
+ u_int32_t low;
+diff --git a/k_rem_pio2.cpp b/k_rem_pio2.cpp
+--- a/k_rem_pio2.cpp
++++ b/k_rem_pio2.cpp
+@@ -305,17 +305,18 @@ __kernel_rem_pio2(double *x, double *y, int e0, int nx, int prec)
+ q0 = e0-24*(jv+1);
+
+ /* set up f[0] to f[jx+jk] where f[jx+jk] = ipio2[jv+jk] */
+ j = jv-jx; m = jx+jk;
+ for(i=0;i<=m;i++,j++) f[i] = (j<0)? zero : (double) ipio2[j];
+
+ /* compute q[0],q[1],...q[jk] */
+ for (i=0;i<=jk;i++) {
+- for(j=0,fw=0.0;j<=jx;j++) fw += x[j]*f[jx+i-j]; q[i] = fw;
++ for(j=0,fw=0.0;j<=jx;j++) fw += x[j]*f[jx+i-j];
++ q[i] = fw;
+ }
+
+ jz = jk;
+ recompute:
+ /* distill q[] into iq[] reversingly */
+ for(i=0,j=jz,z=q[jz];j>0;i++,j--) {
+ fw = (double)((int32_t)(twon24* z));
+ iq[i] = (int32_t)(z-two24*fw);
+diff --git a/s_cos.cpp b/s_cos.cpp
+--- a/s_cos.cpp
++++ b/s_cos.cpp
+@@ -43,17 +43,17 @@
+ * Accuracy:
+ * TRIG(x) returns trig(x) nearly rounded
+ */
+
+ #include <float.h>
+
+ #define INLINE_REM_PIO2
+ #include "math_private.h"
+-#include "e_rem_pio2.c"
++#include "e_rem_pio2.cpp"
+
+ double
+ cos(double x)
+ {
+ double y[2],z=0.0;
+ int32_t n, ix;
+
+ /* High word of x. */
+diff --git a/s_sin.cpp b/s_sin.cpp
+--- a/s_sin.cpp
++++ b/s_sin.cpp
+@@ -43,17 +43,17 @@
+ * Accuracy:
+ * TRIG(x) returns trig(x) nearly rounded
+ */
+
+ #include <float.h>
+
+ #define INLINE_REM_PIO2
+ #include "math_private.h"
+-#include "e_rem_pio2.c"
++#include "e_rem_pio2.cpp"
+
+ double
+ sin(double x)
+ {
+ double y[2],z=0.0;
+ int32_t n, ix;
+
+ /* High word of x. */
+diff --git a/s_tan.cpp b/s_tan.cpp
+--- a/s_tan.cpp
++++ b/s_tan.cpp
+@@ -42,17 +42,17 @@
+ * Accuracy:
+ * TRIG(x) returns trig(x) nearly rounded
+ */
+
+ #include <float.h>
+
+ #define INLINE_REM_PIO2
+ #include "math_private.h"
+-#include "e_rem_pio2.c"
++#include "e_rem_pio2.cpp"
+
+ double
+ tan(double x)
+ {
+ double y[2],z=0.0;
+ int32_t n, ix;
+
+ /* High word of x. */
+--- a/s_cosf.cpp 2022-12-11 21:26:36.041643110 -0500
++++ b/s_cosf.cpp 2022-12-11 21:34:08.197486828 -0500
+@@ -18,19 +18,19 @@
+ //__FBSDID("$FreeBSD$");
+
+ #include <float.h>
+
+ #define INLINE_KERNEL_COSDF
+ #define INLINE_KERNEL_SINDF
+ #define INLINE_REM_PIO2F
+ #include "math_private.h"
+-#include "e_rem_pio2f.c"
+-#include "k_cosf.c"
+-#include "k_sinf.c"
++#include "e_rem_pio2f.cpp"
++#include "k_cosf.cpp"
++#include "k_sinf.cpp"
+
+ /* Small multiples of pi/2 rounded to double precision. */
+ static const double
+ c1pio2 = 1*M_PI_2, /* 0x3FF921FB, 0x54442D18 */
+ c2pio2 = 2*M_PI_2, /* 0x400921FB, 0x54442D18 */
+ c3pio2 = 3*M_PI_2, /* 0x4012D97C, 0x7F3321D2 */
+ c4pio2 = 4*M_PI_2; /* 0x401921FB, 0x54442D18 */
+
+--- a/s_sinf.cpp 2022-12-11 21:26:36.045643147 -0500
++++ b/s_sinf.cpp 2022-12-11 21:33:52.741362926 -0500
+@@ -18,19 +18,19 @@
+ //__FBSDID("$FreeBSD$");
+
+ #include <float.h>
+
+ #define INLINE_KERNEL_COSDF
+ #define INLINE_KERNEL_SINDF
+ #define INLINE_REM_PIO2F
+ #include "math_private.h"
+-#include "e_rem_pio2f.c"
+-#include "k_cosf.c"
+-#include "k_sinf.c"
++#include "e_rem_pio2f.cpp"
++#include "k_cosf.cpp"
++#include "k_sinf.cpp"
+
+ /* Small multiples of pi/2 rounded to double precision. */
+ static const double
+ s1pio2 = 1*M_PI_2, /* 0x3FF921FB, 0x54442D18 */
+ s2pio2 = 2*M_PI_2, /* 0x400921FB, 0x54442D18 */
+ s3pio2 = 3*M_PI_2, /* 0x4012D97C, 0x7F3321D2 */
+ s4pio2 = 4*M_PI_2; /* 0x401921FB, 0x54442D18 */
+
+--- a/s_tanf.cpp 2022-12-11 21:26:36.045643147 -0500
++++ b/s_tanf.cpp 2022-12-11 21:33:45.269302886 -0500
+@@ -17,18 +17,18 @@
+ //#include <sys/cdefs.h>
+ //__FBSDID("$FreeBSD$");
+
+ #include <float.h>
+
+ #define INLINE_KERNEL_TANDF
+ #define INLINE_REM_PIO2F
+ #include "math_private.h"
+-#include "e_rem_pio2f.c"
+-#include "k_tanf.c"
++#include "e_rem_pio2f.cpp"
++#include "k_tanf.cpp"
+
+ /* Small multiples of pi/2 rounded to double precision. */
+ static const double
+ t1pio2 = 1*M_PI_2, /* 0x3FF921FB, 0x54442D18 */
+ t2pio2 = 2*M_PI_2, /* 0x400921FB, 0x54442D18 */
+ t3pio2 = 3*M_PI_2, /* 0x4012D97C, 0x7F3321D2 */
+ t4pio2 = 4*M_PI_2; /* 0x401921FB, 0x54442D18 */
+
+--- a/e_rem_pio2f.cpp 2022-12-13 14:50:12.582776365 -0500
++++ b/e_rem_pio2f.cpp 2022-12-13 14:51:05.003200166 -0500
+@@ -35,17 +35,17 @@
+ */
+
+ static const double
+ invpio2 = 6.36619772367581382433e-01, /* 0x3FE45F30, 0x6DC9C883 */
+ pio2_1 = 1.57079631090164184570e+00, /* 0x3FF921FB, 0x50000000 */
+ pio2_1t = 1.58932547735281966916e-08; /* 0x3E5110b4, 0x611A6263 */
+
+ #ifdef INLINE_REM_PIO2F
+-static __inline __always_inline
++static inline
+ #endif
+ int
+ __ieee754_rem_pio2f(float x, double *y)
+ {
+ double w,r,fn;
+ double tx[1],ty[1];
+ float z;
+ int32_t e0,n,ix,hx;
+--- a/e_powf.cpp
++++ b/e_powf.cpp
+@@ -11,16 +11,18 @@
+ * software is freely granted, provided that this notice
+ * is preserved.
+ * ====================================================
+ */
+
+ //#include <sys/cdefs.h>
+ //__FBSDID("$FreeBSD$");
+
++#include "s_scalbnf.cpp"
++
+ #include "math_private.h"
+
+ static const float
+ bp[] = {1.0, 1.5,},
+ dp_h[] = { 0.0, 5.84960938e-01,}, /* 0x3f15c000 */
+ dp_l[] = { 0.0, 1.56322085e-06,}, /* 0x35d1cfdc */
+ zero = 0.0,
+ half = 0.5,
diff --git a/modules/fdlibm/patches/22_add_float_types.patch b/modules/fdlibm/patches/22_add_float_types.patch
new file mode 100644
index 0000000000..39b991a486
--- /dev/null
+++ b/modules/fdlibm/patches/22_add_float_types.patch
@@ -0,0 +1,66 @@
+diff --git a/math_private.h b/math_private.h
+--- a/math_private.h
++++ b/math_private.h
+@@ -27,16 +27,17 @@
+
+ /*
+ * Emulate FreeBSD internal double types.
+ * Adapted from https://github.com/freebsd/freebsd-src/search?q=__double_t
+ */
+
+ typedef double __double_t;
+ typedef __double_t double_t;
++typedef float __float_t;
+
+ /*
+ * The original fdlibm code used statements like:
+ * n0 = ((*(int*)&one)>>29)^1; * index of high word *
+ * ix0 = *(n0+(int*)&x); * high word of x *
+ * ix1 = *((1-n0)+(int*)&x); * low word of x *
+ * to dig two 32 bit words out of the 64 bit IEEE floating point
+ * value. That is non-ANSI, and, moreover, the gcc instruction
+diff --git a/e_rem_pio2f.cpp b/e_rem_pio2f.cpp
+--- a/e_rem_pio2f.cpp
++++ b/e_rem_pio2f.cpp
+@@ -49,17 +49,17 @@ int
+ double tx[1],ty[1];
+ float z;
+ int32_t e0,n,ix,hx;
+
+ GET_FLOAT_WORD(hx,x);
+ ix = hx&0x7fffffff;
+ /* 33+53 bit pi is good enough for medium size */
+ if(ix<0x4dc90fdb) { /* |x| ~< 2^28*(pi/2), medium size */
+- fn = rnint((float_t)x*invpio2);
++ fn = rnint((float)x*invpio2);
+ n = irint(fn);
+ r = x-fn*pio2_1;
+ w = fn*pio2_1t;
+ *y = r-w;
+ return n;
+ }
+ /*
+ * all other (large) arguments
+diff --git a/s_scalbnf.cpp b/s_scalbnf.cpp
+--- a/s_scalbnf.cpp
++++ b/s_scalbnf.cpp
+@@ -7,18 +7,18 @@
+ * for all contributors to musl.
+ */
+ #include <stdint.h>
+
+ #include "math_private.h"
+
+ float scalbnf(float x, int n)
+ {
+ union {float f; uint32_t i;} u;
+- float_t y = x;
++ float y = x;
+
+ if (n > 127) {
+ y *= 0x1p127f;
+ n -= 127;
+ if (n > 127) {
+ y *= 0x1p127f;
+ n -= 127;
+ if (n > 127)