Remove the msun source files not mentioned in a .bp file.
You'll note suspicious items here, like "why s_sincosf.c but not the other two s_sincos*.c files?", but they seem to be legit cases where arm-optimized-routines has optimized one variant (in this case sincosf()) but not the others (sincos() and sincosl()). The long double variants in particular receiving little attention.
But it's annoying we don't really know/can't prove we don't have problems here. My kingdom for a way to opt in to strict ODR violation checking!
Test: treehugger
Change-Id: I62148316e54501de1db10d10859fb5369efa6041
diff --git a/libm/NOTICE b/libm/NOTICE
index 424725e..3c0e783 100644
--- a/libm/NOTICE
+++ b/libm/NOTICE
@@ -140,15 +140,6 @@
-------------------------------------------------------------------
====================================================
-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.
-
--------------------------------------------------------------------
-
-====================================================
Copyright 2004 Sun Microsystems, Inc. All Rights Reserved.
Permission to use, copy, modify, and distribute this
@@ -423,32 +414,6 @@
-------------------------------------------------------------------
-Copyright (c) 2007 David Schultz
-All rights reserved.
-
-Redistribution and use in source and binary forms, with or without
-modification, are permitted provided that the following conditions
-are met:
-1. Redistributions of source code must retain the above copyright
- notice, this list of conditions and the following disclaimer.
-2. Redistributions in binary form must reproduce the above copyright
- notice, this list of conditions and the following disclaimer in the
- documentation and/or other materials provided with the distribution.
-
-THIS SOFTWARE IS PROVIDED BY AUTHOR AND CONTRIBUTORS ``AS IS'' AND
-ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
-IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
-ARE DISCLAIMED. IN NO EVENT SHALL AUTHOR OR CONTRIBUTORS BE LIABLE
-FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
-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.
-
--------------------------------------------------------------------
-
Copyright (c) 2007 The NetBSD Foundation, Inc.
All rights reserved.
@@ -1273,33 +1238,3 @@
-------------------------------------------------------------------
-SPDX-License-Identifier: BSD-3-Clause
-
-Copyright (c) 2003 Dag-Erling Smørgrav
-All rights reserved.
-
-Redistribution and use in source and binary forms, with or without
-modification, are permitted provided that the following conditions
-are met:
-1. Redistributions of source code must retain the above copyright
- notice, this list of conditions and the following disclaimer
- in this position and unchanged.
-2. Redistributions in binary form must reproduce the above copyright
- notice, this list of conditions and the following disclaimer in the
- documentation and/or other materials provided with the distribution.
-3. The name of the author may not be used to endorse or promote products
- derived from this software without specific prior written permission.
-
-THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR
-IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
-OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
-IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT,
-INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL 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.
-
--------------------------------------------------------------------
-
diff --git a/libm/upstream-freebsd/lib/msun/src/e_exp.c b/libm/upstream-freebsd/lib/msun/src/e_exp.c
deleted file mode 100644
index 59da392..0000000
--- a/libm/upstream-freebsd/lib/msun/src/e_exp.c
+++ /dev/null
@@ -1,164 +0,0 @@
-
-/* @(#)e_exp.c 1.6 04/04/22 */
-/*
- * ====================================================
- * 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$");
-
-/* 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
- *
- * x = k*ln2 + r, |r| <= 0.5*ln2.
- *
- * Here r will be represented as r = hi-lo for better
- * accuracy.
- *
- * 2. Approximation of exp(r) by a special rational function on
- * the interval [0,0.34658]:
- * Write
- * R(r**2) = r*(exp(r)+1)/(exp(r)-1) = 2 + r*r/6 - r**4/360 + ...
- * We use a special Remes algorithm on [0,0.34658] to generate
- * a polynomial of degree 5 to approximate R. The maximum error
- * of this polynomial approximation is bounded by 2**-59. In
- * other words,
- * R(z) ~ 2.0 + P1*z + P2*z**2 + P3*z**3 + P4*z**4 + P5*z**5
- * (where z=r*r, and the values of P1 to P5 are listed below)
- * and
- * | 5 | -59
- * | 2.0+P1*z+...+P5*z - R(z) | <= 2
- * | |
- * The computation of exp(r) thus becomes
- * 2*r
- * exp(r) = 1 + -------
- * R - r
- * r*R1(r)
- * = 1 + r + ----------- (for better accuracy)
- * 2 - R1(r)
- * where
- * 2 4 10
- * R1(r) = r - (P1*r + P2*r + ... + P5*r ).
- *
- * 3. Scale back to obtain exp(x):
- * From step 1, we have
- * exp(x) = 2^k * exp(r)
- *
- * Special cases:
- * exp(INF) is INF, exp(NaN) is NaN;
- * exp(-INF) is 0, and
- * for finite argument, only exp(0)=1 is exact.
- *
- * Accuracy:
- * according to an error analysis, the error is always less than
- * 1 ulp (unit in the last place).
- *
- * Misc. info.
- * For IEEE double
- * if x > 7.09782712893383973096e+02 then exp(x) overflow
- * if x < -7.45133219101941108420e+02 then exp(x) underflow
- *
- * 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
-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 */
- -6.93147180369123816490e-01,},/* 0xbfe62e42, 0xfee00000 */
-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 volatile double
-huge = 1.0e+300,
-twom1000= 9.33263618503218878990e-302; /* 2**-1000=0x01700000,0*/
-
-double
-exp(double x) /* default IEEE double exp */
-{
- double y,hi=0.0,lo=0.0,c,t,twopk;
- int32_t k=0,xsb;
- u_int32_t hx;
-
- GET_HIGH_WORD(hx,x);
- xsb = (hx>>31)&1; /* sign bit of x */
- hx &= 0x7fffffff; /* high word of |x| */
-
- /* filter out non-finite argument */
- if(hx >= 0x40862E42) { /* if |x|>=709.78... */
- if(hx>=0x7ff00000) {
- u_int32_t lx;
- GET_LOW_WORD(lx,x);
- if(((hx&0xfffff)|lx)!=0)
- return x+x; /* NaN */
- else return (xsb==0)? x:0.0; /* exp(+-inf)={inf,0} */
- }
- 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 */
- 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);
- }
- else if(hx < 0x3e300000) { /* when |x|<2**-28 */
- if(huge+x>one) return one+x;/* trigger inexact */
- }
- else k = 0;
-
- /* x is now in primary range */
- t = x*x;
- 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;
- return y*twopk;
- } else {
- return y*twopk*twom1000;
- }
-}
-
-#if (LDBL_MANT_DIG == 53)
-__weak_reference(exp, expl);
-#endif
diff --git a/libm/upstream-freebsd/lib/msun/src/e_expf.c b/libm/upstream-freebsd/lib/msun/src/e_expf.c
deleted file mode 100644
index 620d341..0000000
--- a/libm/upstream-freebsd/lib/msun/src/e_expf.c
+++ /dev/null
@@ -1,98 +0,0 @@
-/* e_expf.c -- float version of e_exp.c.
- * Conversion to float by Ian Lance Taylor, Cygnus Support, ian@cygnus.com.
- */
-
-/*
- * ====================================================
- * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
- *
- * 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 <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 */
-ln2HI[2] ={ 6.9314575195e-01, /* 0x3f317200 */
- -6.9314575195e-01,}, /* 0xbf317200 */
-ln2LO[2] ={ 1.4286067653e-06, /* 0x35bfbe8e */
- -1.4286067653e-06,}, /* 0xb5bfbe8e */
-invln2 = 1.4426950216e+00, /* 0x3fb8aa3b */
-/*
- * Domain [-0.34568, 0.34568], range ~[-4.278e-9, 4.447e-9]:
- * |x*(exp(x)+1)/(exp(x)-1) - p(x)| < 2**-27.74
- */
-P1 = 1.6666625440e-1, /* 0xaaaa8f.0p-26 */
-P2 = -2.7667332906e-3; /* -0xb55215.0p-32 */
-
-static volatile float
-huge = 1.0e+30,
-twom100 = 7.8886090522e-31; /* 2**-100=0x0d800000 */
-
-float
-expf(float x)
-{
- float y,hi=0.0,lo=0.0,c,t,twopk;
- int32_t k=0,xsb;
- u_int32_t hx;
-
- GET_FLOAT_WORD(hx,x);
- xsb = (hx>>31)&1; /* sign bit of x */
- hx &= 0x7fffffff; /* high word of |x| */
-
- /* filter out non-finite argument */
- if(hx >= 0x42b17218) { /* if |x|>=88.721... */
- if(hx>0x7f800000)
- return x+x; /* NaN */
- if(hx==0x7f800000)
- return (xsb==0)? x:0.0; /* exp(+-inf)={inf,0} */
- if(x > o_threshold) return huge*huge; /* overflow */
- if(x < u_threshold) return twom100*twom100; /* underflow */
- }
-
- /* argument reduction */
- if(hx > 0x3eb17218) { /* if |x| > 0.5 ln2 */
- if(hx < 0x3F851592) { /* and |x| < 1.5 ln2 */
- hi = x-ln2HI[xsb]; lo=ln2LO[xsb]; k = 1-xsb-xsb;
- } else {
- k = invln2*x+halF[xsb];
- t = k;
- hi = x - t*ln2HI[0]; /* t*ln2HI is exact here */
- lo = t*ln2LO[0];
- }
- STRICT_ASSIGN(float, x, hi - lo);
- }
- else if(hx < 0x39000000) { /* when |x|<2**-14 */
- if(huge+x>one) return one+x;/* trigger inexact */
- }
- else k = 0;
-
- /* x is now in primary range */
- t = x*x;
- if(k >= -125)
- SET_FLOAT_WORD(twopk,((u_int32_t)(0x7f+k))<<23);
- else
- SET_FLOAT_WORD(twopk,((u_int32_t)(0x7f+(k+100)))<<23);
- c = x - t*(P1+t*P2);
- if(k==0) return one-((x*c)/(c-(float)2.0)-x);
- else y = one-((lo-(x*c)/((float)2.0-c))-hi);
- if(k >= -125) {
- if(k==128) return y*2.0F*0x1p127F;
- return y*twopk;
- } else {
- return y*twopk*twom100;
- }
-}
diff --git a/libm/upstream-freebsd/lib/msun/src/e_log.c b/libm/upstream-freebsd/lib/msun/src/e_log.c
deleted file mode 100644
index 03ce820..0000000
--- a/libm/upstream-freebsd/lib/msun/src/e_log.c
+++ /dev/null
@@ -1,147 +0,0 @@
-
-/* @(#)e_log.c 1.3 95/01/18 */
-/*
- * ====================================================
- * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
- *
- * 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$");
-
-/* 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) .
- *
- * 2. Approximation of log(1+f).
- * Let s = f/(2+f) ; based on log(1+f) = log(1+s) - log(1-s)
- * = 2s + 2/3 s**3 + 2/5 s**5 + .....,
- * = 2s + s*R
- * We use a special Reme algorithm on [0,0.1716] to generate
- * a polynomial of degree 14 to approximate R The maximum error
- * of this polynomial approximation is bounded by 2**-58.45. In
- * other words,
- * 2 4 6 8 10 12 14
- * R(z) ~ Lg1*s +Lg2*s +Lg3*s +Lg4*s +Lg5*s +Lg6*s +Lg7*s
- * (the values of Lg1 to Lg7 are listed in the program)
- * and
- * | 2 14 | -58.45
- * | Lg1*s +...+Lg7*s - R(z) | <= 2
- * | |
- * Note that 2s = f - s*f = f - hfsq + s*hfsq, where hfsq = f*f/2.
- * In order to guarantee error in log below 1ulp, we compute log
- * by
- * log(1+f) = f - s*(f - R) (if f is not too large)
- * log(1+f) = f - (hfsq - s*(hfsq+R)). (better accuracy)
- *
- * 3. Finally, log(x) = k*ln2 + log(1+f).
- * = k*ln2_hi+(f-(hfsq-(s*(hfsq+R)+k*ln2_lo)))
- * Here ln2 is split into two floating point number:
- * ln2_hi + ln2_lo,
- * where n*ln2_hi is always exact for |n| < 2000.
- *
- * Special cases:
- * log(x) is NaN with signal if x < 0 (including -INF) ;
- * log(+INF) is +INF; log(0) is -INF with signal;
- * log(NaN) is that NaN with no signal.
- *
- * Accuracy:
- * according to an error analysis, the error is always less than
- * 1 ulp (unit in the last place).
- *
- * 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
-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 */
-Lg3 = 2.857142874366239149e-01, /* 3FD24924 94229359 */
-Lg4 = 2.222219843214978396e-01, /* 3FCC71C5 1D8E78AF */
-Lg5 = 1.818357216161805012e-01, /* 3FC74664 96CB03DE */
-Lg6 = 1.531383769920937332e-01, /* 3FC39A09 D078C69F */
-Lg7 = 1.479819860511658591e-01; /* 3FC2F112 DF3E5244 */
-
-static const double zero = 0.0;
-static volatile double vzero = 0.0;
-
-double
-log(double x)
-{
- double hfsq,f,s,z,R,w,t1,t2,dk;
- int32_t k,hx,i,j;
- u_int32_t lx;
-
- EXTRACT_WORDS(hx,lx,x);
-
- k=0;
- if (hx < 0x00100000) { /* x < 2**-1022 */
- if (((hx&0x7fffffff)|lx)==0)
- return -two54/vzero; /* log(+-0)=-inf */
- if (hx<0) return (x-x)/zero; /* log(-#) = NaN */
- k -= 54; x *= two54; /* subnormal number, scale up x */
- GET_HIGH_WORD(hx,x);
- }
- if (hx >= 0x7ff00000) return x+x;
- k += (hx>>20)-1023;
- hx &= 0x000fffff;
- i = (hx+0x95f64)&0x100000;
- SET_HIGH_WORD(x,hx|(i^0x3ff00000)); /* normalize x or x/2 */
- k += (i>>20);
- f = x-1.0;
- if((0x000fffff&(2+hx))<3) { /* -2**-20 <= f < 2**-20 */
- if(f==zero) {
- if(k==0) {
- return zero;
- } else {
- dk=(double)k;
- return dk*ln2_hi+dk*ln2_lo;
- }
- }
- R = f*f*(0.5-0.33333333333333333*f);
- if(k==0) return f-R; else {dk=(double)k;
- return dk*ln2_hi-((R-dk*ln2_lo)-f);}
- }
- s = f/(2.0+f);
- dk = (double)k;
- z = s*s;
- i = hx-0x6147a;
- w = z*z;
- j = 0x6b851-hx;
- t1= w*(Lg2+w*(Lg4+w*Lg6));
- t2= z*(Lg1+w*(Lg3+w*(Lg5+w*Lg7)));
- i |= j;
- R = t2+t1;
- if(i>0) {
- 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/libm/upstream-freebsd/lib/msun/src/e_log2.c b/libm/upstream-freebsd/lib/msun/src/e_log2.c
deleted file mode 100644
index 10b1c00..0000000
--- a/libm/upstream-freebsd/lib/msun/src/e_log2.c
+++ /dev/null
@@ -1,117 +0,0 @@
-
-/* @(#)e_log10.c 1.3 95/01/18 */
-/*
- * ====================================================
- * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
- *
- * 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$");
-
-/*
- * 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
- * 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 */
-
-static const double zero = 0.0;
-static volatile double vzero = 0.0;
-
-double
-log2(double x)
-{
- double f,hfsq,hi,lo,r,val_hi,val_lo,w,y;
- int32_t i,k,hx;
- u_int32_t lx;
-
- EXTRACT_WORDS(hx,lx,x);
-
- k=0;
- if (hx < 0x00100000) { /* x < 2**-1022 */
- if (((hx&0x7fffffff)|lx)==0)
- return -two54/vzero; /* log(+-0)=-inf */
- if (hx<0) return (x-x)/zero; /* log(-#) = NaN */
- k -= 54; x *= two54; /* subnormal number, scale up x */
- GET_HIGH_WORD(hx,x);
- }
- if (hx >= 0x7ff00000) return x+x;
- if (hx == 0x3ff00000 && lx == 0)
- return zero; /* log(1) = +0 */
- k += (hx>>20)-1023;
- hx &= 0x000fffff;
- i = (hx+0x95f64)&0x100000;
- SET_HIGH_WORD(x,hx|(i^0x3ff00000)); /* normalize x or x/2 */
- k += (i>>20);
- y = (double)k;
- f = x - 1.0;
- hfsq = 0.5*f*f;
- r = k_log1p(f);
-
- /*
- * f-hfsq must (for args near 1) be evaluated in extra precision
- * to avoid a large cancellation when x is near sqrt(2) or 1/sqrt(2).
- * This is fairly efficient since f-hfsq only depends on f, so can
- * be evaluated in parallel with R. Not combining hfsq with R also
- * keeps R small (though not as small as a true `lo' term would be),
- * so that extra precision is not needed for terms involving R.
- *
- * Compiler bugs involving extra precision used to break Dekker's
- * theorem for spitting f-hfsq as hi+lo, unless double_t was used
- * or the multi-precision calculations were avoided when double_t
- * has extra precision. These problems are now automatically
- * avoided as a side effect of the optimization of combining the
- * Dekker splitting step with the clear-low-bits step.
- *
- * y must (for args near sqrt(2) and 1/sqrt(2)) be added in extra
- * precision to avoid a very large cancellation when x is very near
- * these values. Unlike the above cancellations, this problem is
- * specific to base 2. It is strange that adding +-1 is so much
- * harder than adding +-ln2 or +-log10_2.
- *
- * This uses Dekker's theorem to normalize y+val_hi, so the
- * compiler bugs are back in some configurations, sigh. And I
- * don't want to used double_t to avoid them, since that gives a
- * pessimization and the support for avoiding the pessimization
- * is not yet available.
- *
- * The multi-precision calculations for the multiplications are
- * routine.
- */
- hi = f - hfsq;
- SET_LOW_WORD(hi,0);
- lo = (f - hi) - hfsq + r;
- val_hi = hi*ivln2hi;
- val_lo = (lo+hi)*ivln2lo + lo*ivln2hi;
-
- /* 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/libm/upstream-freebsd/lib/msun/src/e_log2f.c b/libm/upstream-freebsd/lib/msun/src/e_log2f.c
deleted file mode 100644
index 956f33a..0000000
--- a/libm/upstream-freebsd/lib/msun/src/e_log2f.c
+++ /dev/null
@@ -1,82 +0,0 @@
-/*
- * ====================================================
- * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
- *
- * 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$");
-
-/*
- * Float version of e_log2.c. See the latter for most comments.
- */
-
-#include "math.h"
-#include "math_private.h"
-#include "k_logf.h"
-
-static const float
-two25 = 3.3554432000e+07, /* 0x4c000000 */
-ivln2hi = 1.4428710938e+00, /* 0x3fb8b000 */
-ivln2lo = -1.7605285393e-04; /* 0xb9389ad4 */
-
-static const float zero = 0.0;
-static volatile float vzero = 0.0;
-
-float
-log2f(float x)
-{
- float f,hfsq,hi,lo,r,y;
- int32_t i,k,hx;
-
- GET_FLOAT_WORD(hx,x);
-
- k=0;
- if (hx < 0x00800000) { /* x < 2**-126 */
- if ((hx&0x7fffffff)==0)
- return -two25/vzero; /* log(+-0)=-inf */
- if (hx<0) return (x-x)/zero; /* log(-#) = NaN */
- k -= 25; x *= two25; /* subnormal number, scale up x */
- GET_FLOAT_WORD(hx,x);
- }
- if (hx >= 0x7f800000) return x+x;
- if (hx == 0x3f800000)
- return zero; /* log(1) = +0 */
- k += (hx>>23)-127;
- hx &= 0x007fffff;
- i = (hx+(0x4afb0d))&0x800000;
- SET_FLOAT_WORD(x,hx|(i^0x3f800000)); /* normalize x or x/2 */
- k += (i>>23);
- y = (float)k;
- f = x - (float)1.0;
- hfsq = (float)0.5*f*f;
- r = k_log1pf(f);
-
- /*
- * We no longer need to avoid falling into the multi-precision
- * calculations due to compiler bugs breaking Dekker's theorem.
- * Keep avoiding this as an optimization. See e_log2.c for more
- * details (some details are here only because the optimization
- * is not yet available in double precision).
- *
- * Another compiler bug turned up. With gcc on i386,
- * (ivln2lo + ivln2hi) would be evaluated in float precision
- * despite runtime evaluations using double precision. So we
- * must cast one of its terms to float_t. This makes the whole
- * expression have type float_t, so return is forced to waste
- * time clobbering its extra precision.
- */
- if (sizeof(float_t) > sizeof(float))
- return (r - hfsq + f) * ((float_t)ivln2lo + ivln2hi) + y;
-
- hi = f - hfsq;
- GET_FLOAT_WORD(hx,hi);
- SET_FLOAT_WORD(hi,hx&0xfffff000);
- lo = (f - hi) - hfsq + r;
- return (lo+hi)*ivln2lo + lo*ivln2hi + hi*ivln2hi + y;
-}
diff --git a/libm/upstream-freebsd/lib/msun/src/e_logf.c b/libm/upstream-freebsd/lib/msun/src/e_logf.c
deleted file mode 100644
index 68a4d5d..0000000
--- a/libm/upstream-freebsd/lib/msun/src/e_logf.c
+++ /dev/null
@@ -1,89 +0,0 @@
-/* e_logf.c -- float version of e_log.c.
- * Conversion to float by Ian Lance Taylor, Cygnus Support, ian@cygnus.com.
- */
-
-/*
- * ====================================================
- * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
- *
- * 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 "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]). */
-Lg1 = 0xaaaaaa.0p-24, /* 0.66666662693 */
-Lg2 = 0xccce13.0p-25, /* 0.40000972152 */
-Lg3 = 0x91e9ee.0p-25, /* 0.28498786688 */
-Lg4 = 0xf89e26.0p-26; /* 0.24279078841 */
-
-static const float zero = 0.0;
-static volatile float vzero = 0.0;
-
-float
-logf(float x)
-{
- float hfsq,f,s,z,R,w,t1,t2,dk;
- int32_t k,ix,i,j;
-
- GET_FLOAT_WORD(ix,x);
-
- k=0;
- if (ix < 0x00800000) { /* x < 2**-126 */
- if ((ix&0x7fffffff)==0)
- return -two25/vzero; /* log(+-0)=-inf */
- if (ix<0) return (x-x)/zero; /* log(-#) = NaN */
- k -= 25; x *= two25; /* subnormal number, scale up x */
- GET_FLOAT_WORD(ix,x);
- }
- if (ix >= 0x7f800000) return x+x;
- k += (ix>>23)-127;
- ix &= 0x007fffff;
- i = (ix+(0x95f64<<3))&0x800000;
- SET_FLOAT_WORD(x,ix|(i^0x3f800000)); /* normalize x or x/2 */
- k += (i>>23);
- f = x-(float)1.0;
- if((0x007fffff&(0x8000+ix))<0xc000) { /* -2**-9 <= f < 2**-9 */
- if(f==zero) {
- if(k==0) {
- return zero;
- } else {
- dk=(float)k;
- return dk*ln2_hi+dk*ln2_lo;
- }
- }
- R = f*f*((float)0.5-(float)0.33333333333333333*f);
- if(k==0) return f-R; else {dk=(float)k;
- return dk*ln2_hi-((R-dk*ln2_lo)-f);}
- }
- s = f/((float)2.0+f);
- dk = (float)k;
- z = s*s;
- i = ix-(0x6147a<<3);
- w = z*z;
- j = (0x6b851<<3)-ix;
- t1= w*(Lg2+w*Lg4);
- t2= z*(Lg1+w*Lg3);
- i |= j;
- R = t2+t1;
- if(i>0) {
- hfsq=(float)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);
- }
-}
diff --git a/libm/upstream-freebsd/lib/msun/src/e_pow.c b/libm/upstream-freebsd/lib/msun/src/e_pow.c
deleted file mode 100644
index adc64c9..0000000
--- a/libm/upstream-freebsd/lib/msun/src/e_pow.c
+++ /dev/null
@@ -1,314 +0,0 @@
-/* @(#)e_pow.c 1.5 04/04/22 SMI */
-/*
- * ====================================================
- * 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$");
-
-/* 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.
- * 2. Perform y*log2(x) = n+y' by simulating multi-precision
- * arithmetic, where |y'|<=0.5.
- * 3. Return x**y = 2**n*exp(y'*log2)
- *
- * Special cases:
- * 1. (anything) ** 0 is 1
- * 2. (anything) ** 1 is itself
- * 3. (anything) ** NAN is NAN except 1 ** NAN = 1
- * 4. NAN ** (anything except 0) is NAN
- * 5. +-(|x| > 1) ** +INF is +INF
- * 6. +-(|x| > 1) ** -INF is +0
- * 7. +-(|x| < 1) ** +INF is +0
- * 8. +-(|x| < 1) ** -INF is +INF
- * 9. +-1 ** +-INF is 1
- * 10. +0 ** (+anything except 0, NAN) is +0
- * 11. -0 ** (+anything except 0, NAN, odd integer) is +0
- * 12. +0 ** (-anything except 0, NAN) is +INF
- * 13. -0 ** (-anything except 0, NAN, odd integer) is +INF
- * 14. -0 ** (odd integer) = -( +0 ** (odd integer) )
- * 15. +INF ** (+anything except 0,NAN) is +INF
- * 16. +INF ** (-anything except 0,NAN) is +0
- * 17. -INF ** (anything) = -0 ** (-anything)
- * 18. (-anything) ** (integer) is (-1)**(integer)*(+anything**integer)
- * 19. (-anything except 0 and inf) ** (non-integer) is NAN
- *
- * Accuracy:
- * pow(x,y) returns x**y nearly rounded. In particular
- * pow(integer,integer)
- * always returns the correct integer provided it is
- * representable.
- *
- * 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,
-qrtr = 0.25,
-thrd = 3.3333333333333331e-01, /* 0x3fd55555, 0x55555555 */
-one = 1.0,
-two = 2.0,
-two53 = 9007199254740992.0, /* 0x43400000, 0x00000000 */
-huge = 1.0e300,
-tiny = 1.0e-300,
- /* poly coefs for (3/2)*(log(x)-2s-2/3*s**3 */
-L1 = 5.99999999999994648725e-01, /* 0x3FE33333, 0x33333303 */
-L2 = 4.28571428578550184252e-01, /* 0x3FDB6DB6, 0xDB6FABFF */
-L3 = 3.33333329818377432918e-01, /* 0x3FD55555, 0x518F264D */
-L4 = 2.72728123808534006489e-01, /* 0x3FD17460, 0xA91D4101 */
-L5 = 2.30660745775561754067e-01, /* 0x3FCD864A, 0x93C9DB65 */
-L6 = 2.06975017800338417784e-01, /* 0x3FCA7E28, 0x4A454EEF */
-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 */
-lg2 = 6.93147180559945286227e-01, /* 0x3FE62E42, 0xFEFA39EF */
-lg2_h = 6.93147182464599609375e-01, /* 0x3FE62E43, 0x00000000 */
-lg2_l = -1.90465429995776804525e-09, /* 0xBE205C61, 0x0CA86C39 */
-ovt = 8.0085662595372944372e-0017, /* -(1024-log2(ovfl+.5ulp)) */
-cp = 9.61796693925975554329e-01, /* 0x3FEEC709, 0xDC3A03FD =2/(3ln2) */
-cp_h = 9.61796700954437255859e-01, /* 0x3FEEC709, 0xE0000000 =(float)cp */
-cp_l = -7.02846165095275826516e-09, /* 0xBE3E2FE0, 0x145B01F5 =tail of cp_h*/
-ivln2 = 1.44269504088896338700e+00, /* 0x3FF71547, 0x652B82FE =1/ln2 */
-ivln2_h = 1.44269502162933349609e+00, /* 0x3FF71547, 0x60000000 =24b 1/ln2*/
-ivln2_l = 1.92596299112661746887e-08; /* 0x3E54AE0B, 0xF85DDF44 =1/ln2 tail*/
-
-double
-pow(double x, double y)
-{
- double z,ax,z_h,z_l,p_h,p_l;
- double y1,t1,t2,r,s,t,u,v,w;
- int32_t i,j,k,yisint,n;
- int32_t hx,hy,ix,iy;
- u_int32_t lx,ly;
-
- EXTRACT_WORDS(hx,lx,x);
- EXTRACT_WORDS(hy,ly,y);
- ix = hx&0x7fffffff; iy = hy&0x7fffffff;
-
- /* y==zero: x**0 = 1 */
- if((iy|ly)==0) return one;
-
- /* x==1: 1**y = 1, even if y is NaN */
- if (hx==0x3ff00000 && lx == 0) return one;
-
- /* y!=zero: result is NaN if either arg is NaN */
- if(ix > 0x7ff00000 || ((ix==0x7ff00000)&&(lx!=0)) ||
- iy > 0x7ff00000 || ((iy==0x7ff00000)&&(ly!=0)))
- return nan_mix(x, y);
-
- /* determine if y is an odd int when x < 0
- * yisint = 0 ... y is not an integer
- * yisint = 1 ... y is an odd int
- * yisint = 2 ... y is an even int
- */
- yisint = 0;
- if(hx<0) {
- if(iy>=0x43400000) yisint = 2; /* even integer y */
- else if(iy>=0x3ff00000) {
- k = (iy>>20)-0x3ff; /* exponent */
- if(k>20) {
- j = ly>>(52-k);
- if(((u_int32_t)j<<(52-k))==ly) yisint = 2-(j&1);
- } else if(ly==0) {
- j = iy>>(20-k);
- if((j<<(20-k))==iy) yisint = 2-(j&1);
- }
- }
- }
-
- /* special value of y */
- if(ly==0) {
- if (iy==0x7ff00000) { /* y is +-inf */
- if(((ix-0x3ff00000)|lx)==0)
- return one; /* (-1)**+-inf is 1 */
- else if (ix >= 0x3ff00000)/* (|x|>1)**+-inf = inf,0 */
- return (hy>=0)? y: zero;
- else /* (|x|<1)**-,+inf = inf,0 */
- 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);
- }
- }
-
- ax = fabs(x);
- /* special value of x */
- if(lx==0) {
- if(ix==0x7ff00000||ix==0||ix==0x3ff00000){
- z = ax; /*x is +-0,+-inf,+-1*/
- if(hy<0) z = one/z; /* z = (1/|x|) */
- if(hx<0) {
- if(((ix-0x3ff00000)|yisint)==0) {
- z = (z-z)/(z-z); /* (-1)**non-int is NaN */
- } else if(yisint==1)
- z = -z; /* (x<0)**odd = -(|x|**odd) */
- }
- return z;
- }
- }
-
- /* CYGNUS LOCAL + fdlibm-5.3 fix: This used to be
- n = (hx>>31)+1;
- but ANSI C says a right shift of a signed negative quantity is
- implementation defined. */
- n = ((u_int32_t)hx>>31)-1;
-
- /* (x<0)**(non-int) is NaN */
- if((n|yisint)==0) return (x-x)/(x-x);
-
- s = one; /* s (sign of result -ve**odd) = -1 else = 1 */
- if((n|(yisint-1))==0) s = -one;/* (-ve)**(odd int) */
-
- /* |y| is huge */
- if(iy>0x41e00000) { /* if |y| > 2**31 */
- if(iy>0x43f00000){ /* if |y| > 2**64, must o/uflow */
- if(ix<=0x3fefffff) return (hy<0)? huge*huge:tiny*tiny;
- if(ix>=0x3ff00000) return (hy>0)? huge*huge:tiny*tiny;
- }
- /* over/underflow if x is not close to one */
- if(ix<0x3fefffff) return (hy<0)? s*huge*huge:s*tiny*tiny;
- if(ix>0x3ff00000) return (hy>0)? s*huge*huge:s*tiny*tiny;
- /* now |1-x| is tiny <= 2**-20, suffice to compute
- log(x) by x-x^2/2+x^3/3-x^4/4 */
- t = ax-one; /* t has 20 trailing zeros */
- w = (t*t)*(half-t*(thrd-t*qrtr));
- u = ivln2_h*t; /* ivln2_h has 21 sig. bits */
- v = t*ivln2_l-w*ivln2;
- t1 = u+v;
- SET_LOW_WORD(t1,0);
- t2 = v-(t1-u);
- } else {
- double ss,s2,s_h,s_l,t_h,t_l;
- n = 0;
- /* take care subnormal number */
- if(ix<0x00100000)
- {ax *= two53; n -= 53; GET_HIGH_WORD(ix,ax); }
- n += ((ix)>>20)-0x3ff;
- j = ix&0x000fffff;
- /* determine interval */
- ix = j|0x3ff00000; /* normalize ix */
- if(j<=0x3988E) k=0; /* |x|<sqrt(3/2) */
- else if(j<0xBB67A) k=1; /* |x|<sqrt(3) */
- else {k=0;n+=1;ix -= 0x00100000;}
- SET_HIGH_WORD(ax,ix);
-
- /* compute ss = s_h+s_l = (x-1)/(x+1) or (x-1.5)/(x+1.5) */
- u = ax-bp[k]; /* bp[0]=1.0, bp[1]=1.5 */
- v = one/(ax+bp[k]);
- ss = u*v;
- s_h = ss;
- SET_LOW_WORD(s_h,0);
- /* t_h=ax+bp[k] High */
- t_h = zero;
- SET_HIGH_WORD(t_h,((ix>>1)|0x20000000)+0x00080000+(k<<18));
- t_l = ax - (t_h-bp[k]);
- s_l = v*((u-s_h*t_h)-s_h*t_l);
- /* compute log(ax) */
- s2 = ss*ss;
- r = s2*s2*(L1+s2*(L2+s2*(L3+s2*(L4+s2*(L5+s2*L6)))));
- r += s_l*(s_h+ss);
- s2 = s_h*s_h;
- t_h = 3+s2+r;
- SET_LOW_WORD(t_h,0);
- t_l = r-((t_h-3)-s2);
- /* u+v = ss*(1+...) */
- u = s_h*t_h;
- v = s_l*t_h+t_l*ss;
- /* 2/(3log2)*(ss+...) */
- p_h = u+v;
- SET_LOW_WORD(p_h,0);
- p_l = v-(p_h-u);
- z_h = cp_h*p_h; /* cp_h+cp_l = 2/(3*log2) */
- z_l = cp_l*p_h+p_l*cp+dp_l[k];
- /* log2(ax) = (ss+..)*2/(3*log2) = n + dp_h + z_h + z_l */
- t = n;
- t1 = (((z_h+z_l)+dp_h[k])+t);
- SET_LOW_WORD(t1,0);
- t2 = z_l-(((t1-t)-dp_h[k])-z_h);
- }
-
- /* split up y into y1+y2 and compute (y1+y2)*(t1+t2) */
- y1 = y;
- SET_LOW_WORD(y1,0);
- p_l = (y-y1)*t1+y*t2;
- p_h = y1*t1;
- z = p_l+p_h;
- EXTRACT_WORDS(j,i,z);
- if (j>=0x40900000) { /* z >= 1024 */
- if(((j-0x40900000)|i)!=0) /* if z > 1024 */
- return s*huge*huge; /* overflow */
- else {
- if(p_l+ovt>z-p_h) return s*huge*huge; /* overflow */
- }
- } else if((j&0x7fffffff)>=0x4090cc00 ) { /* z <= -1075 */
- if(((j-0xc090cc00)|i)!=0) /* z < -1075 */
- return s*tiny*tiny; /* underflow */
- else {
- if(p_l<=z-p_h) return s*tiny*tiny; /* underflow */
- }
- }
- /*
- * compute 2**(p_h+p_l)
- */
- i = j&0x7fffffff;
- k = (i>>20)-0x3ff;
- n = 0;
- if(i>0x3fe00000) { /* if |z| > 0.5, set n = [z+0.5] */
- n = j+(0x00100000>>(k+1));
- k = ((n&0x7fffffff)>>20)-0x3ff; /* new k for n */
- t = zero;
- SET_HIGH_WORD(t,n&~(0x000fffff>>k));
- n = ((n&0x000fffff)|0x00100000)>>(20-k);
- if(j<0) n = -n;
- p_h -= t;
- }
- t = p_l+p_h;
- SET_LOW_WORD(t,0);
- u = t*lg2_h;
- v = (p_l-(t-p_h))*lg2+t*lg2_l;
- z = u+v;
- w = v-(z-u);
- t = z*z;
- t1 = z - t*(P1+t*(P2+t*(P3+t*(P4+t*P5))));
- 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/libm/upstream-freebsd/lib/msun/src/e_powf.c b/libm/upstream-freebsd/lib/msun/src/e_powf.c
deleted file mode 100644
index f5a2c70..0000000
--- a/libm/upstream-freebsd/lib/msun/src/e_powf.c
+++ /dev/null
@@ -1,252 +0,0 @@
-/* e_powf.c -- float version of e_pow.c.
- * Conversion to float by Ian Lance Taylor, Cygnus Support, ian@cygnus.com.
- */
-
-/*
- * ====================================================
- * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
- *
- * 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 "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,
-half = 0.5,
-qrtr = 0.25,
-thrd = 3.33333343e-01, /* 0x3eaaaaab */
-one = 1.0,
-two = 2.0,
-two24 = 16777216.0, /* 0x4b800000 */
-huge = 1.0e30,
-tiny = 1.0e-30,
- /* poly coefs for (3/2)*(log(x)-2s-2/3*s**3 */
-L1 = 6.0000002384e-01, /* 0x3f19999a */
-L2 = 4.2857143283e-01, /* 0x3edb6db7 */
-L3 = 3.3333334327e-01, /* 0x3eaaaaab */
-L4 = 2.7272811532e-01, /* 0x3e8ba305 */
-L5 = 2.3066075146e-01, /* 0x3e6c3255 */
-L6 = 2.0697501302e-01, /* 0x3e53f142 */
-P1 = 1.6666667163e-01, /* 0x3e2aaaab */
-P2 = -2.7777778450e-03, /* 0xbb360b61 */
-P3 = 6.6137559770e-05, /* 0x388ab355 */
-P4 = -1.6533901999e-06, /* 0xb5ddea0e */
-P5 = 4.1381369442e-08, /* 0x3331bb4c */
-lg2 = 6.9314718246e-01, /* 0x3f317218 */
-lg2_h = 6.93145752e-01, /* 0x3f317200 */
-lg2_l = 1.42860654e-06, /* 0x35bfbe8c */
-ovt = 4.2995665694e-08, /* -(128-log2(ovfl+.5ulp)) */
-cp = 9.6179670095e-01, /* 0x3f76384f =2/(3ln2) */
-cp_h = 9.6191406250e-01, /* 0x3f764000 =12b cp */
-cp_l = -1.1736857402e-04, /* 0xb8f623c6 =tail of cp_h */
-ivln2 = 1.4426950216e+00, /* 0x3fb8aa3b =1/ln2 */
-ivln2_h = 1.4426879883e+00, /* 0x3fb8aa00 =16b 1/ln2*/
-ivln2_l = 7.0526075433e-06; /* 0x36eca570 =1/ln2 tail*/
-
-float
-powf(float x, float y)
-{
- float z,ax,z_h,z_l,p_h,p_l;
- float y1,t1,t2,r,s,sn,t,u,v,w;
- int32_t i,j,k,yisint,n;
- int32_t hx,hy,ix,iy,is;
-
- GET_FLOAT_WORD(hx,x);
- GET_FLOAT_WORD(hy,y);
- ix = hx&0x7fffffff; iy = hy&0x7fffffff;
-
- /* y==zero: x**0 = 1 */
- if(iy==0) return one;
-
- /* x==1: 1**y = 1, even if y is NaN */
- if (hx==0x3f800000) return one;
-
- /* y!=zero: result is NaN if either arg is NaN */
- if(ix > 0x7f800000 ||
- iy > 0x7f800000)
- return nan_mix(x, y);
-
- /* determine if y is an odd int when x < 0
- * yisint = 0 ... y is not an integer
- * yisint = 1 ... y is an odd int
- * yisint = 2 ... y is an even int
- */
- yisint = 0;
- if(hx<0) {
- if(iy>=0x4b800000) yisint = 2; /* even integer y */
- else if(iy>=0x3f800000) {
- k = (iy>>23)-0x7f; /* exponent */
- j = iy>>(23-k);
- if((j<<(23-k))==iy) yisint = 2-(j&1);
- }
- }
-
- /* special value of y */
- if (iy==0x7f800000) { /* y is +-inf */
- if (ix==0x3f800000)
- return one; /* (-1)**+-inf is NaN */
- else if (ix > 0x3f800000)/* (|x|>1)**+-inf = inf,0 */
- return (hy>=0)? y: zero;
- else /* (|x|<1)**-,+inf = inf,0 */
- return (hy<0)?-y: zero;
- }
- if(iy==0x3f800000) { /* y is +-1 */
- if(hy<0) return one/x; else return x;
- }
- if(hy==0x40000000) return x*x; /* y is 2 */
- if(hy==0x3f000000) { /* y is 0.5 */
- if(hx>=0) /* x >= +0 */
- return sqrtf(x);
- }
-
- ax = fabsf(x);
- /* special value of x */
- if(ix==0x7f800000||ix==0||ix==0x3f800000){
- z = ax; /*x is +-0,+-inf,+-1*/
- if(hy<0) z = one/z; /* z = (1/|x|) */
- if(hx<0) {
- if(((ix-0x3f800000)|yisint)==0) {
- z = (z-z)/(z-z); /* (-1)**non-int is NaN */
- } else if(yisint==1)
- z = -z; /* (x<0)**odd = -(|x|**odd) */
- }
- return z;
- }
-
- n = ((u_int32_t)hx>>31)-1;
-
- /* (x<0)**(non-int) is NaN */
- if((n|yisint)==0) return (x-x)/(x-x);
-
- sn = one; /* s (sign of result -ve**odd) = -1 else = 1 */
- if((n|(yisint-1))==0) sn = -one;/* (-ve)**(odd int) */
-
- /* |y| is huge */
- if(iy>0x4d000000) { /* if |y| > 2**27 */
- /* over/underflow if x is not close to one */
- if(ix<0x3f7ffff6) return (hy<0)? sn*huge*huge:sn*tiny*tiny;
- if(ix>0x3f800007) return (hy>0)? sn*huge*huge:sn*tiny*tiny;
- /* now |1-x| is tiny <= 2**-20, suffice to compute
- log(x) by x-x^2/2+x^3/3-x^4/4 */
- t = ax-1; /* t has 20 trailing zeros */
- w = (t*t)*(half-t*(thrd-t*qrtr));
- u = ivln2_h*t; /* ivln2_h has 16 sig. bits */
- v = t*ivln2_l-w*ivln2;
- t1 = u+v;
- GET_FLOAT_WORD(is,t1);
- SET_FLOAT_WORD(t1,is&0xfffff000);
- t2 = v-(t1-u);
- } else {
- float s2,s_h,s_l,t_h,t_l;
- n = 0;
- /* take care subnormal number */
- if(ix<0x00800000)
- {ax *= two24; n -= 24; GET_FLOAT_WORD(ix,ax); }
- n += ((ix)>>23)-0x7f;
- j = ix&0x007fffff;
- /* determine interval */
- ix = j|0x3f800000; /* normalize ix */
- if(j<=0x1cc471) k=0; /* |x|<sqrt(3/2) */
- else if(j<0x5db3d7) k=1; /* |x|<sqrt(3) */
- else {k=0;n+=1;ix -= 0x00800000;}
- SET_FLOAT_WORD(ax,ix);
-
- /* compute s = s_h+s_l = (x-1)/(x+1) or (x-1.5)/(x+1.5) */
- u = ax-bp[k]; /* bp[0]=1.0, bp[1]=1.5 */
- v = one/(ax+bp[k]);
- s = u*v;
- s_h = s;
- GET_FLOAT_WORD(is,s_h);
- SET_FLOAT_WORD(s_h,is&0xfffff000);
- /* t_h=ax+bp[k] High */
- is = ((ix>>1)&0xfffff000)|0x20000000;
- SET_FLOAT_WORD(t_h,is+0x00400000+(k<<21));
- t_l = ax - (t_h-bp[k]);
- s_l = v*((u-s_h*t_h)-s_h*t_l);
- /* compute log(ax) */
- s2 = s*s;
- r = s2*s2*(L1+s2*(L2+s2*(L3+s2*(L4+s2*(L5+s2*L6)))));
- r += s_l*(s_h+s);
- s2 = s_h*s_h;
- t_h = 3+s2+r;
- GET_FLOAT_WORD(is,t_h);
- SET_FLOAT_WORD(t_h,is&0xfffff000);
- t_l = r-((t_h-3)-s2);
- /* u+v = s*(1+...) */
- u = s_h*t_h;
- v = s_l*t_h+t_l*s;
- /* 2/(3log2)*(s+...) */
- p_h = u+v;
- GET_FLOAT_WORD(is,p_h);
- SET_FLOAT_WORD(p_h,is&0xfffff000);
- p_l = v-(p_h-u);
- z_h = cp_h*p_h; /* cp_h+cp_l = 2/(3*log2) */
- z_l = cp_l*p_h+p_l*cp+dp_l[k];
- /* log2(ax) = (s+..)*2/(3*log2) = n + dp_h + z_h + z_l */
- t = n;
- t1 = (((z_h+z_l)+dp_h[k])+t);
- GET_FLOAT_WORD(is,t1);
- SET_FLOAT_WORD(t1,is&0xfffff000);
- t2 = z_l-(((t1-t)-dp_h[k])-z_h);
- }
-
- /* split up y into y1+y2 and compute (y1+y2)*(t1+t2) */
- GET_FLOAT_WORD(is,y);
- SET_FLOAT_WORD(y1,is&0xfffff000);
- p_l = (y-y1)*t1+y*t2;
- p_h = y1*t1;
- z = p_l+p_h;
- GET_FLOAT_WORD(j,z);
- if (j>0x43000000) /* if z > 128 */
- return sn*huge*huge; /* overflow */
- else if (j==0x43000000) { /* if z == 128 */
- if(p_l+ovt>z-p_h) return sn*huge*huge; /* overflow */
- }
- else if ((j&0x7fffffff)>0x43160000) /* z <= -150 */
- return sn*tiny*tiny; /* underflow */
- else if (j==0xc3160000){ /* z == -150 */
- if(p_l<=z-p_h) return sn*tiny*tiny; /* underflow */
- }
- /*
- * compute 2**(p_h+p_l)
- */
- i = j&0x7fffffff;
- k = (i>>23)-0x7f;
- n = 0;
- if(i>0x3f000000) { /* if |z| > 0.5, set n = [z+0.5] */
- n = j+(0x00800000>>(k+1));
- k = ((n&0x7fffffff)>>23)-0x7f; /* new k for n */
- SET_FLOAT_WORD(t,n&~(0x007fffff>>k));
- n = ((n&0x007fffff)|0x00800000)>>(23-k);
- if(j<0) n = -n;
- p_h -= t;
- }
- t = p_l+p_h;
- GET_FLOAT_WORD(is,t);
- SET_FLOAT_WORD(t,is&0xffff8000);
- u = t*lg2_h;
- v = (p_l-(t-p_h))*lg2+t*lg2_l;
- z = u+v;
- w = v-(z-u);
- t = z*z;
- t1 = z - t*(P1+t*(P2+t*(P3+t*(P4+t*P5))));
- r = (z*t1)/(t1-two)-(w+z*w);
- z = one-(r-z);
- GET_FLOAT_WORD(j,z);
- j += (n<<23);
- if((j>>23)<=0) z = scalbnf(z,n); /* subnormal output */
- else SET_FLOAT_WORD(z,j);
- return sn*z;
-}
diff --git a/libm/upstream-freebsd/lib/msun/src/s_copysign.c b/libm/upstream-freebsd/lib/msun/src/s_copysign.c
deleted file mode 100644
index a5f3870..0000000
--- a/libm/upstream-freebsd/lib/msun/src/s_copysign.c
+++ /dev/null
@@ -1,33 +0,0 @@
-/* @(#)s_copysign.c 5.1 93/09/24 */
-/*
- * ====================================================
- * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
- *
- * 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$");
-
-/*
- * 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);
- SET_HIGH_WORD(x,(hx&0x7fffffff)|(hy&0x80000000));
- return x;
-}
diff --git a/libm/upstream-freebsd/lib/msun/src/s_copysignf.c b/libm/upstream-freebsd/lib/msun/src/s_copysignf.c
deleted file mode 100644
index 05ca1e3..0000000
--- a/libm/upstream-freebsd/lib/msun/src/s_copysignf.c
+++ /dev/null
@@ -1,36 +0,0 @@
-/* s_copysignf.c -- float version of s_copysign.c.
- * Conversion to float by Ian Lance Taylor, Cygnus Support, ian@cygnus.com.
- */
-
-/*
- * ====================================================
- * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
- *
- * 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$");
-
-/*
- * copysignf(float x, float y)
- * copysignf(x,y) returns a value with the magnitude of x and
- * with the sign bit of y.
- */
-
-#include "math.h"
-#include "math_private.h"
-
-float
-copysignf(float x, float y)
-{
- u_int32_t ix,iy;
- GET_FLOAT_WORD(ix,x);
- GET_FLOAT_WORD(iy,y);
- SET_FLOAT_WORD(x,(ix&0x7fffffff)|(iy&0x80000000));
- return x;
-}
diff --git a/libm/upstream-freebsd/lib/msun/src/s_copysignl.c b/libm/upstream-freebsd/lib/msun/src/s_copysignl.c
deleted file mode 100644
index 666a2ba..0000000
--- a/libm/upstream-freebsd/lib/msun/src/s_copysignl.c
+++ /dev/null
@@ -1,44 +0,0 @@
-/*-
- * SPDX-License-Identifier: BSD-2-Clause
- *
- * Copyright (c) 2004 Stefan Farfeleder
- * All rights reserved.
- *
- * Redistribution and use in source and binary forms, with or without
- * modification, are permitted provided that the following conditions
- * are met:
- * 1. Redistributions of source code must retain the above copyright
- * notice, this list of conditions and the following disclaimer.
- * 2. Redistributions in binary form must reproduce the above copyright
- * notice, this list of conditions and the following disclaimer in the
- * documentation and/or other materials provided with the distribution.
- *
- * THIS SOFTWARE IS PROVIDED BY AUTHOR AND CONTRIBUTORS ``AS IS'' AND
- * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
- * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
- * ARE DISCLAIMED. IN NO EVENT SHALL AUTHOR OR CONTRIBUTORS BE LIABLE
- * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
- * 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.
- *
- * $FreeBSD$
- */
-
-#include <math.h>
-
-#include "fpmath.h"
-
-long double
-copysignl(long double x, long double y)
-{
- union IEEEl2bits ux, uy;
-
- ux.e = x;
- uy.e = y;
- ux.bits.sign = uy.bits.sign;
- return (ux.e);
-}
diff --git a/libm/upstream-freebsd/lib/msun/src/s_cosf.c b/libm/upstream-freebsd/lib/msun/src/s_cosf.c
deleted file mode 100644
index b701fd2..0000000
--- a/libm/upstream-freebsd/lib/msun/src/s_cosf.c
+++ /dev/null
@@ -1,87 +0,0 @@
-/* s_cosf.c -- float version of s_cos.c.
- * Conversion to float by Ian Lance Taylor, Cygnus Support, ian@cygnus.com.
- * Optimized by Bruce D. Evans.
- */
-
-/*
- * ====================================================
- * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
- *
- * 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 <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"
-
-/* 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 */
-
-float
-cosf(float x)
-{
- double y;
- int32_t n, hx, ix;
-
- GET_FLOAT_WORD(hx,x);
- ix = hx & 0x7fffffff;
-
- if(ix <= 0x3f490fda) { /* |x| ~<= pi/4 */
- if(ix<0x39800000) /* |x| < 2**-12 */
- if(((int)x)==0) return 1.0; /* 1 with inexact if x != 0 */
- return __kernel_cosdf(x);
- }
- if(ix<=0x407b53d1) { /* |x| ~<= 5*pi/4 */
- if(ix>0x4016cbe3) /* |x| ~> 3*pi/4 */
- return -__kernel_cosdf(x + (hx > 0 ? -c2pio2 : c2pio2));
- else {
- if(hx>0)
- return __kernel_sindf(c1pio2 - x);
- else
- return __kernel_sindf(x + c1pio2);
- }
- }
- if(ix<=0x40e231d5) { /* |x| ~<= 9*pi/4 */
- if(ix>0x40afeddf) /* |x| ~> 7*pi/4 */
- return __kernel_cosdf(x + (hx > 0 ? -c4pio2 : c4pio2));
- else {
- if(hx>0)
- return __kernel_sindf(x - c3pio2);
- else
- return __kernel_sindf(-c3pio2 - x);
- }
- }
-
- /* cos(Inf or NaN) is NaN */
- else if (ix>=0x7f800000) return x-x;
-
- /* general argument reduction needed */
- else {
- n = __ieee754_rem_pio2f(x,&y);
- switch(n&3) {
- case 0: return __kernel_cosdf(y);
- case 1: return __kernel_sindf(-y);
- case 2: return -__kernel_cosdf(y);
- default:
- return __kernel_sindf(y);
- }
- }
-}
diff --git a/libm/upstream-freebsd/lib/msun/src/s_exp2.c b/libm/upstream-freebsd/lib/msun/src/s_exp2.c
deleted file mode 100644
index 1dd9673..0000000
--- a/libm/upstream-freebsd/lib/msun/src/s_exp2.c
+++ /dev/null
@@ -1,399 +0,0 @@
-/*-
- * SPDX-License-Identifier: BSD-2-Clause
- *
- * Copyright (c) 2005 David Schultz <das@FreeBSD.ORG>
- * All rights reserved.
- *
- * Redistribution and use in source and binary forms, with or without
- * modification, are permitted provided that the following conditions
- * are met:
- * 1. Redistributions of source code must retain the above copyright
- * notice, this list of conditions and the following disclaimer.
- * 2. Redistributions in binary form must reproduce the above copyright
- * notice, this list of conditions and the following disclaimer in the
- * documentation and/or other materials provided with the distribution.
- *
- * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND
- * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
- * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
- * ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE
- * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
- * 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 <float.h>
-
-#include "math.h"
-#include "math_private.h"
-
-#define TBLBITS 8
-#define TBLSIZE (1 << TBLBITS)
-
-static const double
- redux = 0x1.8p52 / TBLSIZE,
- P1 = 0x1.62e42fefa39efp-1,
- P2 = 0x1.ebfbdff82c575p-3,
- P3 = 0x1.c6b08d704a0a6p-5,
- P4 = 0x1.3b2ab88f70400p-7,
- P5 = 0x1.5d88003875c74p-10;
-
-static volatile double
- huge = 0x1p1000,
- twom1000 = 0x1p-1000;
-
-static const double tbl[TBLSIZE * 2] = {
-/* exp2(z + eps) eps */
- 0x1.6a09e667f3d5dp-1, 0x1.9880p-44,
- 0x1.6b052fa751744p-1, 0x1.8000p-50,
- 0x1.6c012750bd9fep-1, -0x1.8780p-45,
- 0x1.6cfdcddd476bfp-1, 0x1.ec00p-46,
- 0x1.6dfb23c651a29p-1, -0x1.8000p-50,
- 0x1.6ef9298593ae3p-1, -0x1.c000p-52,
- 0x1.6ff7df9519386p-1, -0x1.fd80p-45,
- 0x1.70f7466f42da3p-1, -0x1.c880p-45,
- 0x1.71f75e8ec5fc3p-1, 0x1.3c00p-46,
- 0x1.72f8286eacf05p-1, -0x1.8300p-44,
- 0x1.73f9a48a58152p-1, -0x1.0c00p-47,
- 0x1.74fbd35d7ccfcp-1, 0x1.f880p-45,
- 0x1.75feb564267f1p-1, 0x1.3e00p-47,
- 0x1.77024b1ab6d48p-1, -0x1.7d00p-45,
- 0x1.780694fde5d38p-1, -0x1.d000p-50,
- 0x1.790b938ac1d00p-1, 0x1.3000p-49,
- 0x1.7a11473eb0178p-1, -0x1.d000p-49,
- 0x1.7b17b0976d060p-1, 0x1.0400p-45,
- 0x1.7c1ed0130c133p-1, 0x1.0000p-53,
- 0x1.7d26a62ff8636p-1, -0x1.6900p-45,
- 0x1.7e2f336cf4e3bp-1, -0x1.2e00p-47,
- 0x1.7f3878491c3e8p-1, -0x1.4580p-45,
- 0x1.80427543e1b4ep-1, 0x1.3000p-44,
- 0x1.814d2add1071ap-1, 0x1.f000p-47,
- 0x1.82589994ccd7ep-1, -0x1.1c00p-45,
- 0x1.8364c1eb942d0p-1, 0x1.9d00p-45,
- 0x1.8471a4623cab5p-1, 0x1.7100p-43,
- 0x1.857f4179f5bbcp-1, 0x1.2600p-45,
- 0x1.868d99b4491afp-1, -0x1.2c40p-44,
- 0x1.879cad931a395p-1, -0x1.3000p-45,
- 0x1.88ac7d98a65b8p-1, -0x1.a800p-45,
- 0x1.89bd0a4785800p-1, -0x1.d000p-49,
- 0x1.8ace5422aa223p-1, 0x1.3280p-44,
- 0x1.8be05bad619fap-1, 0x1.2b40p-43,
- 0x1.8cf3216b54383p-1, -0x1.ed00p-45,
- 0x1.8e06a5e08664cp-1, -0x1.0500p-45,
- 0x1.8f1ae99157807p-1, 0x1.8280p-45,
- 0x1.902fed0282c0ep-1, -0x1.cb00p-46,
- 0x1.9145b0b91ff96p-1, -0x1.5e00p-47,
- 0x1.925c353aa2ff9p-1, 0x1.5400p-48,
- 0x1.93737b0cdc64ap-1, 0x1.7200p-46,
- 0x1.948b82b5f98aep-1, -0x1.9000p-47,
- 0x1.95a44cbc852cbp-1, 0x1.5680p-45,
- 0x1.96bdd9a766f21p-1, -0x1.6d00p-44,
- 0x1.97d829fde4e2ap-1, -0x1.1000p-47,
- 0x1.98f33e47a23a3p-1, 0x1.d000p-45,
- 0x1.9a0f170ca0604p-1, -0x1.8a40p-44,
- 0x1.9b2bb4d53ff89p-1, 0x1.55c0p-44,
- 0x1.9c49182a3f15bp-1, 0x1.6b80p-45,
- 0x1.9d674194bb8c5p-1, -0x1.c000p-49,
- 0x1.9e86319e3238ep-1, 0x1.7d00p-46,
- 0x1.9fa5e8d07f302p-1, 0x1.6400p-46,
- 0x1.a0c667b5de54dp-1, -0x1.5000p-48,
- 0x1.a1e7aed8eb8f6p-1, 0x1.9e00p-47,
- 0x1.a309bec4a2e27p-1, 0x1.ad80p-45,
- 0x1.a42c980460a5dp-1, -0x1.af00p-46,
- 0x1.a5503b23e259bp-1, 0x1.b600p-47,
- 0x1.a674a8af46213p-1, 0x1.8880p-44,
- 0x1.a799e1330b3a7p-1, 0x1.1200p-46,
- 0x1.a8bfe53c12e8dp-1, 0x1.6c00p-47,
- 0x1.a9e6b5579fcd2p-1, -0x1.9b80p-45,
- 0x1.ab0e521356fb8p-1, 0x1.b700p-45,
- 0x1.ac36bbfd3f381p-1, 0x1.9000p-50,
- 0x1.ad5ff3a3c2780p-1, 0x1.4000p-49,
- 0x1.ae89f995ad2a3p-1, -0x1.c900p-45,
- 0x1.afb4ce622f367p-1, 0x1.6500p-46,
- 0x1.b0e07298db790p-1, 0x1.fd40p-45,
- 0x1.b20ce6c9a89a9p-1, 0x1.2700p-46,
- 0x1.b33a2b84f1a4bp-1, 0x1.d470p-43,
- 0x1.b468415b747e7p-1, -0x1.8380p-44,
- 0x1.b59728de5593ap-1, 0x1.8000p-54,
- 0x1.b6c6e29f1c56ap-1, 0x1.ad00p-47,
- 0x1.b7f76f2fb5e50p-1, 0x1.e800p-50,
- 0x1.b928cf22749b2p-1, -0x1.4c00p-47,
- 0x1.ba5b030a10603p-1, -0x1.d700p-47,
- 0x1.bb8e0b79a6f66p-1, 0x1.d900p-47,
- 0x1.bcc1e904bc1ffp-1, 0x1.2a00p-47,
- 0x1.bdf69c3f3a16fp-1, -0x1.f780p-46,
- 0x1.bf2c25bd71db8p-1, -0x1.0a00p-46,
- 0x1.c06286141b2e9p-1, -0x1.1400p-46,
- 0x1.c199bdd8552e0p-1, 0x1.be00p-47,
- 0x1.c2d1cd9fa64eep-1, -0x1.9400p-47,
- 0x1.c40ab5fffd02fp-1, -0x1.ed00p-47,
- 0x1.c544778fafd15p-1, 0x1.9660p-44,
- 0x1.c67f12e57d0cbp-1, -0x1.a100p-46,
- 0x1.c7ba88988c1b6p-1, -0x1.8458p-42,
- 0x1.c8f6d9406e733p-1, -0x1.a480p-46,
- 0x1.ca3405751c4dfp-1, 0x1.b000p-51,
- 0x1.cb720dcef9094p-1, 0x1.1400p-47,
- 0x1.ccb0f2e6d1689p-1, 0x1.0200p-48,
- 0x1.cdf0b555dc412p-1, 0x1.3600p-48,
- 0x1.cf3155b5bab3bp-1, -0x1.6900p-47,
- 0x1.d072d4a0789bcp-1, 0x1.9a00p-47,
- 0x1.d1b532b08c8fap-1, -0x1.5e00p-46,
- 0x1.d2f87080d8a85p-1, 0x1.d280p-46,
- 0x1.d43c8eacaa203p-1, 0x1.1a00p-47,
- 0x1.d5818dcfba491p-1, 0x1.f000p-50,
- 0x1.d6c76e862e6a1p-1, -0x1.3a00p-47,
- 0x1.d80e316c9834ep-1, -0x1.cd80p-47,
- 0x1.d955d71ff6090p-1, 0x1.4c00p-48,
- 0x1.da9e603db32aep-1, 0x1.f900p-48,
- 0x1.dbe7cd63a8325p-1, 0x1.9800p-49,
- 0x1.dd321f301b445p-1, -0x1.5200p-48,
- 0x1.de7d5641c05bfp-1, -0x1.d700p-46,
- 0x1.dfc97337b9aecp-1, -0x1.6140p-46,
- 0x1.e11676b197d5ep-1, 0x1.b480p-47,
- 0x1.e264614f5a3e7p-1, 0x1.0ce0p-43,
- 0x1.e3b333b16ee5cp-1, 0x1.c680p-47,
- 0x1.e502ee78b3fb4p-1, -0x1.9300p-47,
- 0x1.e653924676d68p-1, -0x1.5000p-49,
- 0x1.e7a51fbc74c44p-1, -0x1.7f80p-47,
- 0x1.e8f7977cdb726p-1, -0x1.3700p-48,
- 0x1.ea4afa2a490e8p-1, 0x1.5d00p-49,
- 0x1.eb9f4867ccae4p-1, 0x1.61a0p-46,
- 0x1.ecf482d8e680dp-1, 0x1.5500p-48,
- 0x1.ee4aaa2188514p-1, 0x1.6400p-51,
- 0x1.efa1bee615a13p-1, -0x1.e800p-49,
- 0x1.f0f9c1cb64106p-1, -0x1.a880p-48,
- 0x1.f252b376bb963p-1, -0x1.c900p-45,
- 0x1.f3ac948dd7275p-1, 0x1.a000p-53,
- 0x1.f50765b6e4524p-1, -0x1.4f00p-48,
- 0x1.f6632798844fdp-1, 0x1.a800p-51,
- 0x1.f7bfdad9cbe38p-1, 0x1.abc0p-48,
- 0x1.f91d802243c82p-1, -0x1.4600p-50,
- 0x1.fa7c1819e908ep-1, -0x1.b0c0p-47,
- 0x1.fbdba3692d511p-1, -0x1.0e00p-51,
- 0x1.fd3c22b8f7194p-1, -0x1.0de8p-46,
- 0x1.fe9d96b2a23eep-1, 0x1.e430p-49,
- 0x1.0000000000000p+0, 0x0.0000p+0,
- 0x1.00b1afa5abcbep+0, -0x1.3400p-52,
- 0x1.0163da9fb3303p+0, -0x1.2170p-46,
- 0x1.02168143b0282p+0, 0x1.a400p-52,
- 0x1.02c9a3e77806cp+0, 0x1.f980p-49,
- 0x1.037d42e11bbcap+0, -0x1.7400p-51,
- 0x1.04315e86e7f89p+0, 0x1.8300p-50,
- 0x1.04e5f72f65467p+0, -0x1.a3f0p-46,
- 0x1.059b0d315855ap+0, -0x1.2840p-47,
- 0x1.0650a0e3c1f95p+0, 0x1.1600p-48,
- 0x1.0706b29ddf71ap+0, 0x1.5240p-46,
- 0x1.07bd42b72a82dp+0, -0x1.9a00p-49,
- 0x1.0874518759bd0p+0, 0x1.6400p-49,
- 0x1.092bdf66607c8p+0, -0x1.0780p-47,
- 0x1.09e3ecac6f383p+0, -0x1.8000p-54,
- 0x1.0a9c79b1f3930p+0, 0x1.fa00p-48,
- 0x1.0b5586cf988fcp+0, -0x1.ac80p-48,
- 0x1.0c0f145e46c8ap+0, 0x1.9c00p-50,
- 0x1.0cc922b724816p+0, 0x1.5200p-47,
- 0x1.0d83b23395dd8p+0, -0x1.ad00p-48,
- 0x1.0e3ec32d3d1f3p+0, 0x1.bac0p-46,
- 0x1.0efa55fdfa9a6p+0, -0x1.4e80p-47,
- 0x1.0fb66affed2f0p+0, -0x1.d300p-47,
- 0x1.1073028d7234bp+0, 0x1.1500p-48,
- 0x1.11301d0125b5bp+0, 0x1.c000p-49,
- 0x1.11edbab5e2af9p+0, 0x1.6bc0p-46,
- 0x1.12abdc06c31d5p+0, 0x1.8400p-49,
- 0x1.136a814f2047dp+0, -0x1.ed00p-47,
- 0x1.1429aaea92de9p+0, 0x1.8e00p-49,
- 0x1.14e95934f3138p+0, 0x1.b400p-49,
- 0x1.15a98c8a58e71p+0, 0x1.5300p-47,
- 0x1.166a45471c3dfp+0, 0x1.3380p-47,
- 0x1.172b83c7d5211p+0, 0x1.8d40p-45,
- 0x1.17ed48695bb9fp+0, -0x1.5d00p-47,
- 0x1.18af9388c8d93p+0, -0x1.c880p-46,
- 0x1.1972658375d66p+0, 0x1.1f00p-46,
- 0x1.1a35beb6fcba7p+0, 0x1.0480p-46,
- 0x1.1af99f81387e3p+0, -0x1.7390p-43,
- 0x1.1bbe084045d54p+0, 0x1.4e40p-45,
- 0x1.1c82f95281c43p+0, -0x1.a200p-47,
- 0x1.1d4873168b9b2p+0, 0x1.3800p-49,
- 0x1.1e0e75eb44031p+0, 0x1.ac00p-49,
- 0x1.1ed5022fcd938p+0, 0x1.1900p-47,
- 0x1.1f9c18438cdf7p+0, -0x1.b780p-46,
- 0x1.2063b88628d8fp+0, 0x1.d940p-45,
- 0x1.212be3578a81ep+0, 0x1.8000p-50,
- 0x1.21f49917ddd41p+0, 0x1.b340p-45,
- 0x1.22bdda2791323p+0, 0x1.9f80p-46,
- 0x1.2387a6e7561e7p+0, -0x1.9c80p-46,
- 0x1.2451ffb821427p+0, 0x1.2300p-47,
- 0x1.251ce4fb2a602p+0, -0x1.3480p-46,
- 0x1.25e85711eceb0p+0, 0x1.2700p-46,
- 0x1.26b4565e27d16p+0, 0x1.1d00p-46,
- 0x1.2780e341de00fp+0, 0x1.1ee0p-44,
- 0x1.284dfe1f5633ep+0, -0x1.4c00p-46,
- 0x1.291ba7591bb30p+0, -0x1.3d80p-46,
- 0x1.29e9df51fdf09p+0, 0x1.8b00p-47,
- 0x1.2ab8a66d10e9bp+0, -0x1.27c0p-45,
- 0x1.2b87fd0dada3ap+0, 0x1.a340p-45,
- 0x1.2c57e39771af9p+0, -0x1.0800p-46,
- 0x1.2d285a6e402d9p+0, -0x1.ed00p-47,
- 0x1.2df961f641579p+0, -0x1.4200p-48,
- 0x1.2ecafa93e2ecfp+0, -0x1.4980p-45,
- 0x1.2f9d24abd8822p+0, -0x1.6300p-46,
- 0x1.306fe0a31b625p+0, -0x1.2360p-44,
- 0x1.31432edeea50bp+0, -0x1.0df8p-40,
- 0x1.32170fc4cd7b8p+0, -0x1.2480p-45,
- 0x1.32eb83ba8e9a2p+0, -0x1.5980p-45,
- 0x1.33c08b2641766p+0, 0x1.ed00p-46,
- 0x1.3496266e3fa27p+0, -0x1.c000p-50,
- 0x1.356c55f929f0fp+0, -0x1.0d80p-44,
- 0x1.36431a2de88b9p+0, 0x1.2c80p-45,
- 0x1.371a7373aaa39p+0, 0x1.0600p-45,
- 0x1.37f26231e74fep+0, -0x1.6600p-46,
- 0x1.38cae6d05d838p+0, -0x1.ae00p-47,
- 0x1.39a401b713ec3p+0, -0x1.4720p-43,
- 0x1.3a7db34e5a020p+0, 0x1.8200p-47,
- 0x1.3b57fbfec6e95p+0, 0x1.e800p-44,
- 0x1.3c32dc313a8f2p+0, 0x1.f800p-49,
- 0x1.3d0e544ede122p+0, -0x1.7a00p-46,
- 0x1.3dea64c1234bbp+0, 0x1.6300p-45,
- 0x1.3ec70df1c4eccp+0, -0x1.8a60p-43,
- 0x1.3fa4504ac7e8cp+0, -0x1.cdc0p-44,
- 0x1.40822c367a0bbp+0, 0x1.5b80p-45,
- 0x1.4160a21f72e95p+0, 0x1.ec00p-46,
- 0x1.423fb27094646p+0, -0x1.3600p-46,
- 0x1.431f5d950a920p+0, 0x1.3980p-45,
- 0x1.43ffa3f84b9ebp+0, 0x1.a000p-48,
- 0x1.44e0860618919p+0, -0x1.6c00p-48,
- 0x1.45c2042a7d201p+0, -0x1.bc00p-47,
- 0x1.46a41ed1d0016p+0, -0x1.2800p-46,
- 0x1.4786d668b3326p+0, 0x1.0e00p-44,
- 0x1.486a2b5c13c00p+0, -0x1.d400p-45,
- 0x1.494e1e192af04p+0, 0x1.c200p-47,
- 0x1.4a32af0d7d372p+0, -0x1.e500p-46,
- 0x1.4b17dea6db801p+0, 0x1.7800p-47,
- 0x1.4bfdad53629e1p+0, -0x1.3800p-46,
- 0x1.4ce41b817c132p+0, 0x1.0800p-47,
- 0x1.4dcb299fddddbp+0, 0x1.c700p-45,
- 0x1.4eb2d81d8ab96p+0, -0x1.ce00p-46,
- 0x1.4f9b2769d2d02p+0, 0x1.9200p-46,
- 0x1.508417f4531c1p+0, -0x1.8c00p-47,
- 0x1.516daa2cf662ap+0, -0x1.a000p-48,
- 0x1.5257de83f51eap+0, 0x1.a080p-43,
- 0x1.5342b569d4edap+0, -0x1.6d80p-45,
- 0x1.542e2f4f6ac1ap+0, -0x1.2440p-44,
- 0x1.551a4ca5d94dbp+0, 0x1.83c0p-43,
- 0x1.56070dde9116bp+0, 0x1.4b00p-45,
- 0x1.56f4736b529dep+0, 0x1.15a0p-43,
- 0x1.57e27dbe2c40ep+0, -0x1.9e00p-45,
- 0x1.58d12d497c76fp+0, -0x1.3080p-45,
- 0x1.59c0827ff0b4cp+0, 0x1.dec0p-43,
- 0x1.5ab07dd485427p+0, -0x1.4000p-51,
- 0x1.5ba11fba87af4p+0, 0x1.0080p-44,
- 0x1.5c9268a59460bp+0, -0x1.6c80p-45,
- 0x1.5d84590998e3fp+0, 0x1.69a0p-43,
- 0x1.5e76f15ad20e1p+0, -0x1.b400p-46,
- 0x1.5f6a320dcebcap+0, 0x1.7700p-46,
- 0x1.605e1b976dcb8p+0, 0x1.6f80p-45,
- 0x1.6152ae6cdf715p+0, 0x1.1000p-47,
- 0x1.6247eb03a5531p+0, -0x1.5d00p-46,
- 0x1.633dd1d1929b5p+0, -0x1.2d00p-46,
- 0x1.6434634ccc313p+0, -0x1.a800p-49,
- 0x1.652b9febc8efap+0, -0x1.8600p-45,
- 0x1.6623882553397p+0, 0x1.1fe0p-40,
- 0x1.671c1c708328ep+0, -0x1.7200p-44,
- 0x1.68155d44ca97ep+0, 0x1.6800p-49,
- 0x1.690f4b19e9471p+0, -0x1.9780p-45,
-};
-
-/*
- * exp2(x): compute the base 2 exponential of x
- *
- * Accuracy: Peak error < 0.503 ulp for normalized results.
- *
- * Method: (accurate tables)
- *
- * Reduce x:
- * x = 2**k + y, for integer k and |y| <= 1/2.
- * Thus we have exp2(x) = 2**k * exp2(y).
- *
- * Reduce y:
- * y = i/TBLSIZE + z - eps[i] for integer i near y * TBLSIZE.
- * Thus we have exp2(y) = exp2(i/TBLSIZE) * exp2(z - eps[i]),
- * with |z - eps[i]| <= 2**-9 + 2**-39 for the table used.
- *
- * We compute exp2(i/TBLSIZE) via table lookup and exp2(z - eps[i]) via
- * a degree-5 minimax polynomial with maximum error under 1.3 * 2**-61.
- * The values in exp2t[] and eps[] are chosen such that
- * exp2t[i] = exp2(i/TBLSIZE + eps[i]), and eps[i] is a small offset such
- * that exp2t[i] is accurate to 2**-64.
- *
- * Note that the range of i is +-TBLSIZE/2, so we actually index the tables
- * by i0 = i + TBLSIZE/2. For cache efficiency, exp2t[] and eps[] are
- * virtual tables, interleaved in the real table tbl[].
- *
- * This method is due to Gal, with many details due to Gal and Bachelis:
- *
- * Gal, S. and Bachelis, B. An Accurate Elementary Mathematical Library
- * for the IEEE Floating Point Standard. TOMS 17(1), 26-46 (1991).
- */
-double
-exp2(double x)
-{
- double r, t, twopk, twopkp1000, z;
- uint32_t hx, ix, lx, i0;
- int k;
-
- /* Filter out exceptional cases. */
- GET_HIGH_WORD(hx,x);
- ix = hx & 0x7fffffff; /* high word of |x| */
- if(ix >= 0x40900000) { /* |x| >= 1024 */
- if(ix >= 0x7ff00000) {
- GET_LOW_WORD(lx,x);
- if(((ix & 0xfffff) | lx) != 0 || (hx & 0x80000000) == 0)
- return (x + x); /* x is NaN or +Inf */
- else
- return (0.0); /* x is -Inf */
- }
- if(x >= 0x1.0p10)
- return (huge * huge); /* overflow */
- if(x <= -0x1.0ccp10)
- return (twom1000 * twom1000); /* underflow */
- } else if (ix < 0x3c900000) { /* |x| < 0x1p-54 */
- return (1.0 + x);
- }
-
- /* Reduce x, computing z, i0, and k. */
- STRICT_ASSIGN(double, t, x + redux);
- GET_LOW_WORD(i0, t);
- i0 += TBLSIZE / 2;
- k = (i0 >> TBLBITS) << 20;
- i0 = (i0 & (TBLSIZE - 1)) << 1;
- t -= redux;
- z = x - t;
-
- /* Compute r = exp2(y) = exp2t[i0] * p(z - eps[i]). */
- t = tbl[i0]; /* exp2t[i0] */
- z -= tbl[i0 + 1]; /* eps[i0] */
- if (k >= -(1021 << 20))
- 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);
- return (r * twopk);
- } else {
- return (r * twopkp1000 * twom1000);
- }
-}
-
-#if (LDBL_MANT_DIG == 53)
-__weak_reference(exp2, exp2l);
-#endif
diff --git a/libm/upstream-freebsd/lib/msun/src/s_exp2f.c b/libm/upstream-freebsd/lib/msun/src/s_exp2f.c
deleted file mode 100644
index c5b4c8e..0000000
--- a/libm/upstream-freebsd/lib/msun/src/s_exp2f.c
+++ /dev/null
@@ -1,139 +0,0 @@
-/*-
- * SPDX-License-Identifier: BSD-2-Clause
- *
- * Copyright (c) 2005 David Schultz <das@FreeBSD.ORG>
- * All rights reserved.
- *
- * Redistribution and use in source and binary forms, with or without
- * modification, are permitted provided that the following conditions
- * are met:
- * 1. Redistributions of source code must retain the above copyright
- * notice, this list of conditions and the following disclaimer.
- * 2. Redistributions in binary form must reproduce the above copyright
- * notice, this list of conditions and the following disclaimer in the
- * documentation and/or other materials provided with the distribution.
- *
- * THIS SOFTWARE IS PROVIDED BY THE AUTHOR AND CONTRIBUTORS ``AS IS'' AND
- * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
- * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
- * ARE DISCLAIMED. IN NO EVENT SHALL THE AUTHOR OR CONTRIBUTORS BE LIABLE
- * FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
- * 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 <float.h>
-
-#include "math.h"
-#include "math_private.h"
-
-#define TBLBITS 4
-#define TBLSIZE (1 << TBLBITS)
-
-static const float
- redux = 0x1.8p23f / TBLSIZE,
- P1 = 0x1.62e430p-1f,
- P2 = 0x1.ebfbe0p-3f,
- P3 = 0x1.c6b348p-5f,
- P4 = 0x1.3b2c9cp-7f;
-
-static volatile float
- huge = 0x1p100f,
- twom100 = 0x1p-100f;
-
-static const double exp2ft[TBLSIZE] = {
- 0x1.6a09e667f3bcdp-1,
- 0x1.7a11473eb0187p-1,
- 0x1.8ace5422aa0dbp-1,
- 0x1.9c49182a3f090p-1,
- 0x1.ae89f995ad3adp-1,
- 0x1.c199bdd85529cp-1,
- 0x1.d5818dcfba487p-1,
- 0x1.ea4afa2a490dap-1,
- 0x1.0000000000000p+0,
- 0x1.0b5586cf9890fp+0,
- 0x1.172b83c7d517bp+0,
- 0x1.2387a6e756238p+0,
- 0x1.306fe0a31b715p+0,
- 0x1.3dea64c123422p+0,
- 0x1.4bfdad5362a27p+0,
- 0x1.5ab07dd485429p+0,
-};
-
-/*
- * exp2f(x): compute the base 2 exponential of x
- *
- * Accuracy: Peak error < 0.501 ulp; location of peak: -0.030110927.
- *
- * Method: (equally-spaced tables)
- *
- * Reduce x:
- * x = 2**k + y, for integer k and |y| <= 1/2.
- * Thus we have exp2f(x) = 2**k * exp2(y).
- *
- * Reduce y:
- * y = i/TBLSIZE + z for integer i near y * TBLSIZE.
- * Thus we have exp2(y) = exp2(i/TBLSIZE) * exp2(z),
- * with |z| <= 2**-(TBLSIZE+1).
- *
- * We compute exp2(i/TBLSIZE) via table lookup and exp2(z) via a
- * degree-4 minimax polynomial with maximum error under 1.4 * 2**-33.
- * Using double precision for everything except the reduction makes
- * roundoff error insignificant and simplifies the scaling step.
- *
- * This method is due to Tang, but I do not use his suggested parameters:
- *
- * Tang, P. Table-driven Implementation of the Exponential Function
- * in IEEE Floating-Point Arithmetic. TOMS 15(2), 144-157 (1989).
- */
-float
-exp2f(float x)
-{
- double tv, twopk, u, z;
- float t;
- uint32_t hx, ix, i0;
- int32_t k;
-
- /* Filter out exceptional cases. */
- GET_FLOAT_WORD(hx, x);
- ix = hx & 0x7fffffff; /* high word of |x| */
- if(ix >= 0x43000000) { /* |x| >= 128 */
- if(ix >= 0x7f800000) {
- if ((ix & 0x7fffff) != 0 || (hx & 0x80000000) == 0)
- return (x + x); /* x is NaN or +Inf */
- else
- return (0.0); /* x is -Inf */
- }
- if(x >= 0x1.0p7f)
- return (huge * huge); /* overflow */
- if(x <= -0x1.2cp7f)
- return (twom100 * twom100); /* underflow */
- } else if (ix <= 0x33000000) { /* |x| <= 0x1p-25 */
- return (1.0f + x);
- }
-
- /* Reduce x, computing z, i0, and k. */
- STRICT_ASSIGN(float, t, x + redux);
- GET_FLOAT_WORD(i0, t);
- i0 += TBLSIZE / 2;
- k = (i0 >> TBLBITS) << 20;
- i0 &= TBLSIZE - 1;
- t -= redux;
- z = x - t;
- INSERT_WORDS(twopk, 0x3ff00000 + k, 0);
-
- /* Compute r = exp2(y) = exp2ft[i0] * p(z). */
- tv = exp2ft[i0];
- u = tv * z;
- tv = tv + u * (P1 + z * P2) + u * (z * z) * (P3 + z * P4);
-
- /* Scale by 2**(k>>20). */
- return (tv * twopk);
-}
diff --git a/libm/upstream-freebsd/lib/msun/src/s_fabsl.c b/libm/upstream-freebsd/lib/msun/src/s_fabsl.c
deleted file mode 100644
index 5076d8a..0000000
--- a/libm/upstream-freebsd/lib/msun/src/s_fabsl.c
+++ /dev/null
@@ -1,45 +0,0 @@
-/*-
- * SPDX-License-Identifier: BSD-3-Clause
- *
- * Copyright (c) 2003 Dag-Erling Smørgrav
- * All rights reserved.
- *
- * Redistribution and use in source and binary forms, with or without
- * modification, are permitted provided that the following conditions
- * are met:
- * 1. Redistributions of source code must retain the above copyright
- * notice, this list of conditions and the following disclaimer
- * in this position and unchanged.
- * 2. Redistributions in binary form must reproduce the above copyright
- * notice, this list of conditions and the following disclaimer in the
- * documentation and/or other materials provided with the distribution.
- * 3. The name of the author may not be used to endorse or promote products
- * derived from this software without specific prior written permission.
- *
- * THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR
- * IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
- * OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED.
- * IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT, INDIRECT,
- * INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL 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.
- *
- * $FreeBSD$
- */
-
-#include <math.h>
-
-#include "fpmath.h"
-
-long double
-fabsl(long double x)
-{
- union IEEEl2bits u;
-
- u.e = x;
- u.bits.sign = 0;
- return (u.e);
-}
diff --git a/libm/upstream-freebsd/lib/msun/src/s_sincosf.c b/libm/upstream-freebsd/lib/msun/src/s_sincosf.c
deleted file mode 100644
index 755ff05..0000000
--- a/libm/upstream-freebsd/lib/msun/src/s_sincosf.c
+++ /dev/null
@@ -1,126 +0,0 @@
-/*-
- * ====================================================
- * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
- *
- * 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.
- * ====================================================
- */
-
-/* s_sincosf.c -- float version of s_sincos.c.
- * Conversion to float by Ian Lance Taylor, Cygnus Support, ian@cygnus.com.
- * Optimized by Bruce D. Evans.
- * Merged s_sinf.c and s_cosf.c by Steven G. Kargl.
- */
-
-#include <sys/cdefs.h>
-__FBSDID("$FreeBSD$");
-
-#include <float.h>
-
-#include "math.h"
-#define INLINE_REM_PIO2F
-#include "math_private.h"
-#include "e_rem_pio2f.c"
-#include "k_sincosf.h"
-
-/* Small multiples of pi/2 rounded to double precision. */
-static const double
-p1pio2 = 1*M_PI_2, /* 0x3FF921FB, 0x54442D18 */
-p2pio2 = 2*M_PI_2, /* 0x400921FB, 0x54442D18 */
-p3pio2 = 3*M_PI_2, /* 0x4012D97C, 0x7F3321D2 */
-p4pio2 = 4*M_PI_2; /* 0x401921FB, 0x54442D18 */
-
-void
-sincosf(float x, float *sn, float *cs)
-{
- float c, s;
- double y;
- int32_t n, hx, ix;
-
- GET_FLOAT_WORD(hx, x);
- ix = hx & 0x7fffffff;
-
- if (ix <= 0x3f490fda) { /* |x| ~<= pi/4 */
- if (ix < 0x39800000) { /* |x| < 2**-12 */
- if ((int)x == 0) {
- *sn = x; /* x with inexact if x != 0 */
- *cs = 1;
- return;
- }
- }
- __kernel_sincosdf(x, sn, cs);
- return;
- }
-
- if (ix <= 0x407b53d1) { /* |x| ~<= 5*pi/4 */
- if (ix <= 0x4016cbe3) { /* |x| ~<= 3pi/4 */
- if (hx > 0) {
- __kernel_sincosdf(x - p1pio2, cs, sn);
- *cs = -*cs;
- } else {
- __kernel_sincosdf(x + p1pio2, cs, sn);
- *sn = -*sn;
- }
- } else {
- if (hx > 0)
- __kernel_sincosdf(x - p2pio2, sn, cs);
- else
- __kernel_sincosdf(x + p2pio2, sn, cs);
- *sn = -*sn;
- *cs = -*cs;
- }
- return;
- }
-
- if (ix <= 0x40e231d5) { /* |x| ~<= 9*pi/4 */
- if (ix <= 0x40afeddf) { /* |x| ~<= 7*pi/4 */
- if (hx > 0) {
- __kernel_sincosdf(x - p3pio2, cs, sn);
- *sn = -*sn;
- } else {
- __kernel_sincosdf(x + p3pio2, cs, sn);
- *cs = -*cs;
- }
- } else {
- if (hx > 0)
- __kernel_sincosdf(x - p4pio2, sn, cs);
- else
- __kernel_sincosdf(x + p4pio2, sn, cs);
- }
- return;
- }
-
- /* If x = Inf or NaN, then sin(x) = NaN and cos(x) = NaN. */
- if (ix >= 0x7f800000) {
- *sn = x - x;
- *cs = x - x;
- return;
- }
-
- /* Argument reduction. */
- n = __ieee754_rem_pio2f(x, &y);
- __kernel_sincosdf(y, &s, &c);
-
- switch(n & 3) {
- case 0:
- *sn = s;
- *cs = c;
- break;
- case 1:
- *sn = c;
- *cs = -s;
- break;
- case 2:
- *sn = -s;
- *cs = -c;
- break;
- default:
- *sn = -c;
- *cs = s;
- }
-}
-
-
diff --git a/libm/upstream-freebsd/lib/msun/src/s_sinf.c b/libm/upstream-freebsd/lib/msun/src/s_sinf.c
deleted file mode 100644
index 41b5dc1..0000000
--- a/libm/upstream-freebsd/lib/msun/src/s_sinf.c
+++ /dev/null
@@ -1,85 +0,0 @@
-/* s_sinf.c -- float version of s_sin.c.
- * Conversion to float by Ian Lance Taylor, Cygnus Support, ian@cygnus.com.
- * Optimized by Bruce D. Evans.
- */
-
-/*
- * ====================================================
- * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
- *
- * 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 <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"
-
-/* 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 */
-
-float
-sinf(float x)
-{
- double y;
- int32_t n, hx, ix;
-
- GET_FLOAT_WORD(hx,x);
- ix = hx & 0x7fffffff;
-
- if(ix <= 0x3f490fda) { /* |x| ~<= pi/4 */
- if(ix<0x39800000) /* |x| < 2**-12 */
- if(((int)x)==0) return x; /* x with inexact if x != 0 */
- return __kernel_sindf(x);
- }
- if(ix<=0x407b53d1) { /* |x| ~<= 5*pi/4 */
- if(ix<=0x4016cbe3) { /* |x| ~<= 3pi/4 */
- if(hx>0)
- return __kernel_cosdf(x - s1pio2);
- else
- return -__kernel_cosdf(x + s1pio2);
- } else
- return __kernel_sindf((hx > 0 ? s2pio2 : -s2pio2) - x);
- }
- if(ix<=0x40e231d5) { /* |x| ~<= 9*pi/4 */
- if(ix<=0x40afeddf) { /* |x| ~<= 7*pi/4 */
- if(hx>0)
- return -__kernel_cosdf(x - s3pio2);
- else
- return __kernel_cosdf(x + s3pio2);
- } else
- return __kernel_sindf(x + (hx > 0 ? -s4pio2 : s4pio2));
- }
-
- /* sin(Inf or NaN) is NaN */
- else if (ix>=0x7f800000) return x-x;
-
- /* general argument reduction needed */
- else {
- n = __ieee754_rem_pio2f(x,&y);
- switch(n&3) {
- case 0: return __kernel_sindf(y);
- case 1: return __kernel_cosdf(y);
- case 2: return __kernel_sindf(-y);
- default:
- return -__kernel_cosdf(y);
- }
- }
-}