Adding upstream version 115.7.0esr.upstream/115.7.0esr

Signed-off-by: Daniel Baumann <daniel.baumann@progress-linux.org>

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)