| author | |
| committer | |
| log | eab93d3574887512ac86c5ae47805b0184277fd8 |
| tree | d521a46d54f8f5e6a1c11a5657960b1348ec989c |
| parent | 814b1e9a5835902f143bafbb98c7cdf56f1d1f19 |
These symbols are already provided by compiler_rt11 files changed, 0 insertions(+), 591 deletions(-)
lib/libc/mingw/math/log10f.c deleted-11| ... | @@ -1,11 +0,0 @@ | ||
| 1 | /** | ||
| 2 | * This file has no copyright assigned and is placed in the Public Domain. | ||
| 3 | * This file is part of the mingw-w64 runtime package. | ||
| 4 | * No warranty is given; refer to the file DISCLAIMER.PD within this package. | ||
| 5 | */ | ||
| 6 | #include <math.h> | ||
| 7 | |||
| 8 | float log10f(float _X) | ||
| 9 | { | ||
| 10 | return ((float)log10((double)_X)); | ||
| 11 | } | ||
lib/libc/mingw/math/logf.c deleted-11| ... | @@ -1,11 +0,0 @@ | ||
| 1 | /** | ||
| 2 | * This file has no copyright assigned and is placed in the Public Domain. | ||
| 3 | * This file is part of the mingw-w64 runtime package. | ||
| 4 | * No warranty is given; refer to the file DISCLAIMER.PD within this package. | ||
| 5 | */ | ||
| 6 | #include <math.h> | ||
| 7 | |||
| 8 | float logf(float _X) | ||
| 9 | { | ||
| 10 | return ((float)log((double)_X)); | ||
| 11 | } | ||
lib/libc/musl/src/math/log.c deleted-112| ... | @@ -1,112 +0,0 @@ | ||
| 1 | /* | ||
| 2 | * Double-precision log(x) function. | ||
| 3 | * | ||
| 4 | * Copyright (c) 2018, Arm Limited. | ||
| 5 | * SPDX-License-Identifier: MIT | ||
| 6 | */ | ||
| 7 | |||
| 8 | #include <math.h> | ||
| 9 | #include <stdint.h> | ||
| 10 | #include "libm.h" | ||
| 11 | #include "log_data.h" | ||
| 12 | |||
| 13 | #define T __log_data.tab | ||
| 14 | #define T2 __log_data.tab2 | ||
| 15 | #define B __log_data.poly1 | ||
| 16 | #define A __log_data.poly | ||
| 17 | #define Ln2hi __log_data.ln2hi | ||
| 18 | #define Ln2lo __log_data.ln2lo | ||
| 19 | #define N (1 << LOG_TABLE_BITS) | ||
| 20 | #define OFF 0x3fe6000000000000 | ||
| 21 | |||
| 22 | /* Top 16 bits of a double. */ | ||
| 23 | static inline uint32_t top16(double x) | ||
| 24 | { | ||
| 25 | 	return asuint64(x) >> 48; | ||
| 26 | } | ||
| 27 | |||
| 28 | double log(double x) | ||
| 29 | { | ||
| 30 | 	double_t w, z, r, r2, r3, y, invc, logc, kd, hi, lo; | ||
| 31 | 	uint64_t ix, iz, tmp; | ||
| 32 | 	uint32_t top; | ||
| 33 | 	int k, i; | ||
| 34 | |||
| 35 | 	ix = asuint64(x); | ||
| 36 | 	top = top16(x); | ||
| 37 | #define LO asuint64(1.0 - 0x1p-4) | ||
| 38 | #define HI asuint64(1.0 + 0x1.09p-4) | ||
| 39 | 	if (predict_false(ix - LO < HI - LO)) { | ||
| 40 | 		/* Handle close to 1.0 inputs separately. */ | ||
| 41 | 		/* Fix sign of zero with downward rounding when x==1. */ | ||
| 42 | 		if (WANT_ROUNDING && predict_false(ix == asuint64(1.0))) | ||
| 43 | 			return 0; | ||
| 44 | 		r = x - 1.0; | ||
| 45 | 		r2 = r * r; | ||
| 46 | 		r3 = r * r2; | ||
| 47 | 		y = r3 * | ||
| 48 | 		 (B[1] + r * B[2] + r2 * B[3] + | ||
| 49 | 		 r3 * (B[4] + r * B[5] + r2 * B[6] + | ||
| 50 | 			 r3 * (B[7] + r * B[8] + r2 * B[9] + r3 * B[10]))); | ||
| 51 | 		/* Worst-case error is around 0.507 ULP. */ | ||
| 52 | 		w = r * 0x1p27; | ||
| 53 | 		double_t rhi = r + w - w; | ||
| 54 | 		double_t rlo = r - rhi; | ||
| 55 | 		w = rhi * rhi * B[0]; /* B[0] == -0.5. */ | ||
| 56 | 		hi = r + w; | ||
| 57 | 		lo = r - hi + w; | ||
| 58 | 		lo += B[0] * rlo * (rhi + r); | ||
| 59 | 		y += lo; | ||
| 60 | 		y += hi; | ||
| 61 | 		return eval_as_double(y); | ||
| 62 | 	} | ||
| 63 | 	if (predict_false(top - 0x0010 >= 0x7ff0 - 0x0010)) { | ||
| 64 | 		/* x < 0x1p-1022 or inf or nan. */ | ||
| 65 | 		if (ix * 2 == 0) | ||
| 66 | 			return __math_divzero(1); | ||
| 67 | 		if (ix == asuint64(INFINITY)) /* log(inf) == inf. */ | ||
| 68 | 			return x; | ||
| 69 | 		if ((top & 0x8000) || (top & 0x7ff0) == 0x7ff0) | ||
| 70 | 			return __math_invalid(x); | ||
| 71 | 		/* x is subnormal, normalize it. */ | ||
| 72 | 		ix = asuint64(x * 0x1p52); | ||
| 73 | 		ix -= 52ULL << 52; | ||
| 74 | 	} | ||
| 75 | |||
| 76 | 	/* x = 2^k z; where z is in range [OFF,2*OFF) and exact. | ||
| 77 | 	 The range is split into N subintervals. | ||
| 78 | 	 The ith subinterval contains z and c is near its center. */ | ||
| 79 | 	tmp = ix - OFF; | ||
| 80 | 	i = (tmp >> (52 - LOG_TABLE_BITS)) % N; | ||
| 81 | 	k = (int64_t)tmp >> 52; /* arithmetic shift */ | ||
| 82 | 	iz = ix - (tmp & 0xfffULL << 52); | ||
| 83 | 	invc = T[i].invc; | ||
| 84 | 	logc = T[i].logc; | ||
| 85 | 	z = asdouble(iz); | ||
| 86 | |||
| 87 | 	/* log(x) = log1p(z/c-1) + log(c) + k*Ln2. */ | ||
| 88 | 	/* r ~= z/c - 1, |r| < 1/(2*N). */ | ||
| 89 | #if __FP_FAST_FMA | ||
| 90 | 	/* rounding error: 0x1p-55/N. */ | ||
| 91 | 	r = __builtin_fma(z, invc, -1.0); | ||
| 92 | #else | ||
| 93 | 	/* rounding error: 0x1p-55/N + 0x1p-66. */ | ||
| 94 | 	r = (z - T2[i].chi - T2[i].clo) * invc; | ||
| 95 | #endif | ||
| 96 | 	kd = (double_t)k; | ||
| 97 | |||
| 98 | 	/* hi + lo = r + log(c) + k*Ln2. */ | ||
| 99 | 	w = kd * Ln2hi + logc; | ||
| 100 | 	hi = w + r; | ||
| 101 | 	lo = w - hi + r + kd * Ln2lo; | ||
| 102 | |||
| 103 | 	/* log(x) = lo + (log1p(r) - r) + hi. */ | ||
| 104 | 	r2 = r * r; /* rounding error: 0x1p-54/N^2. */ | ||
| 105 | 	/* Worst case error if |y| > 0x1p-5: | ||
| 106 | 	 0.5 + 4.13/N + abs-poly-error*2^57 ULP (+ 0.002 ULP without fma) | ||
| 107 | 	 Worst case error if |y| > 0x1p-4: | ||
| 108 | 	 0.5 + 2.06/N + abs-poly-error*2^56 ULP (+ 0.001 ULP without fma). */ | ||
| 109 | 	y = lo + r2 * A[0] + | ||
| 110 | 	 r * r2 * (A[1] + r * A[2] + r2 * (A[3] + r * A[4])) + hi; | ||
| 111 | 	return eval_as_double(y); | ||
| 112 | } | ||
lib/libc/musl/src/math/log10.c deleted-101| ... | @@ -1,101 +0,0 @@ | ||
| 1 | /* origin: FreeBSD /usr/src/lib/msun/src/e_log10.c */ | ||
| 2 | /* | ||
| 3 | * ==================================================== | ||
| 4 | * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. | ||
| 5 | * | ||
| 6 | * Developed at SunSoft, a Sun Microsystems, Inc. business. | ||
| 7 | * Permission to use, copy, modify, and distribute this | ||
| 8 | * software is freely granted, provided that this notice | ||
| 9 | * is preserved. | ||
| 10 | * ==================================================== | ||
| 11 | */ | ||
| 12 | /* | ||
| 13 | * Return the base 10 logarithm of x. See log.c for most comments. | ||
| 14 | * | ||
| 15 | * Reduce x to 2^k (1+f) and calculate r = log(1+f) - f + f*f/2 | ||
| 16 | * as in log.c, then combine and scale in extra precision: | ||
| 17 | * log10(x) = (f - f*f/2 + r)/log(10) + k*log10(2) | ||
| 18 | */ | ||
| 19 | |||
| 20 | #include <math.h> | ||
| 21 | #include <stdint.h> | ||
| 22 | |||
| 23 | static const double | ||
| 24 | ivln10hi = 4.34294481878168880939e-01, /* 0x3fdbcb7b, 0x15200000 */ | ||
| 25 | ivln10lo = 2.50829467116452752298e-11, /* 0x3dbb9438, 0xca9aadd5 */ | ||
| 26 | log10_2hi = 3.01029995663611771306e-01, /* 0x3FD34413, 0x509F6000 */ | ||
| 27 | log10_2lo = 3.69423907715893078616e-13, /* 0x3D59FEF3, 0x11F12B36 */ | ||
| 28 | Lg1 = 6.666666666666735130e-01, /* 3FE55555 55555593 */ | ||
| 29 | Lg2 = 3.999999999940941908e-01, /* 3FD99999 9997FA04 */ | ||
| 30 | Lg3 = 2.857142874366239149e-01, /* 3FD24924 94229359 */ | ||
| 31 | Lg4 = 2.222219843214978396e-01, /* 3FCC71C5 1D8E78AF */ | ||
| 32 | Lg5 = 1.818357216161805012e-01, /* 3FC74664 96CB03DE */ | ||
| 33 | Lg6 = 1.531383769920937332e-01, /* 3FC39A09 D078C69F */ | ||
| 34 | Lg7 = 1.479819860511658591e-01; /* 3FC2F112 DF3E5244 */ | ||
| 35 | |||
| 36 | double log10(double x) | ||
| 37 | { | ||
| 38 | 	union {double f; uint64_t i;} u = {x}; | ||
| 39 | 	double_t hfsq,f,s,z,R,w,t1,t2,dk,y,hi,lo,val_hi,val_lo; | ||
| 40 | 	uint32_t hx; | ||
| 41 | 	int k; | ||
| 42 | |||
| 43 | 	hx = u.i>>32; | ||
| 44 | 	k = 0; | ||
| 45 | 	if (hx < 0x00100000 || hx>>31) { | ||
| 46 | 		if (u.i<<1 == 0) | ||
| 47 | 			return -1/(x*x); /* log(+-0)=-inf */ | ||
| 48 | 		if (hx>>31) | ||
| 49 | 			return (x-x)/0.0; /* log(-#) = NaN */ | ||
| 50 | 		/* subnormal number, scale x up */ | ||
| 51 | 		k -= 54; | ||
| 52 | 		x *= 0x1p54; | ||
| 53 | 		u.f = x; | ||
| 54 | 		hx = u.i>>32; | ||
| 55 | 	} else if (hx >= 0x7ff00000) { | ||
| 56 | 		return x; | ||
| 57 | 	} else if (hx == 0x3ff00000 && u.i<<32 == 0) | ||
| 58 | 		return 0; | ||
| 59 | |||
| 60 | 	/* reduce x into [sqrt(2)/2, sqrt(2)] */ | ||
| 61 | 	hx += 0x3ff00000 - 0x3fe6a09e; | ||
| 62 | 	k += (int)(hx>>20) - 0x3ff; | ||
| 63 | 	hx = (hx&0x000fffff) + 0x3fe6a09e; | ||
| 64 | 	u.i = (uint64_t)hx<<32 | (u.i&0xffffffff); | ||
| 65 | 	x = u.f; | ||
| 66 | |||
| 67 | 	f = x - 1.0; | ||
| 68 | 	hfsq = 0.5*f*f; | ||
| 69 | 	s = f/(2.0+f); | ||
| 70 | 	z = s*s; | ||
| 71 | 	w = z*z; | ||
| 72 | 	t1 = w*(Lg2+w*(Lg4+w*Lg6)); | ||
| 73 | 	t2 = z*(Lg1+w*(Lg3+w*(Lg5+w*Lg7))); | ||
| 74 | 	R = t2 + t1; | ||
| 75 | |||
| 76 | 	/* See log2.c for details. */ | ||
| 77 | 	/* hi+lo = f - hfsq + s*(hfsq+R) ~ log(1+f) */ | ||
| 78 | 	hi = f - hfsq; | ||
| 79 | 	u.f = hi; | ||
| 80 | 	u.i &= (uint64_t)-1<<32; | ||
| 81 | 	hi = u.f; | ||
| 82 | 	lo = f - hi - hfsq + s*(hfsq+R); | ||
| 83 | |||
| 84 | 	/* val_hi+val_lo ~ log10(1+f) + k*log10(2) */ | ||
| 85 | 	val_hi = hi*ivln10hi; | ||
| 86 | 	dk = k; | ||
| 87 | 	y = dk*log10_2hi; | ||
| 88 | 	val_lo = dk*log10_2lo + (lo+hi)*ivln10lo + lo*ivln10hi; | ||
| 89 | |||
| 90 | 	/* | ||
| 91 | 	 * Extra precision in for adding y is not strictly needed | ||
| 92 | 	 * since there is no very large cancellation near x = sqrt(2) or | ||
| 93 | 	 * x = 1/sqrt(2), but we do it anyway since it costs little on CPUs | ||
| 94 | 	 * with some parallelism and it reduces the error for many args. | ||
| 95 | 	 */ | ||
| 96 | 	w = y + val_hi; | ||
| 97 | 	val_lo += (y - w) + val_hi; | ||
| 98 | 	val_hi = w; | ||
| 99 | |||
| 100 | 	return val_lo + val_hi; | ||
| 101 | } | ||
lib/libc/musl/src/math/log10f.c deleted-77| ... | @@ -1,77 +0,0 @@ | ||
| 1 | /* origin: FreeBSD /usr/src/lib/msun/src/e_log10f.c */ | ||
| 2 | /* | ||
| 3 | * ==================================================== | ||
| 4 | * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. | ||
| 5 | * | ||
| 6 | * Developed at SunPro, a Sun Microsystems, Inc. business. | ||
| 7 | * Permission to use, copy, modify, and distribute this | ||
| 8 | * software is freely granted, provided that this notice | ||
| 9 | * is preserved. | ||
| 10 | * ==================================================== | ||
| 11 | */ | ||
| 12 | /* | ||
| 13 | * See comments in log10.c. | ||
| 14 | */ | ||
| 15 | |||
| 16 | #include <math.h> | ||
| 17 | #include <stdint.h> | ||
| 18 | |||
| 19 | static const float | ||
| 20 | ivln10hi = 4.3432617188e-01, /* 0x3ede6000 */ | ||
| 21 | ivln10lo = -3.1689971365e-05, /* 0xb804ead9 */ | ||
| 22 | log10_2hi = 3.0102920532e-01, /* 0x3e9a2080 */ | ||
| 23 | log10_2lo = 7.9034151668e-07, /* 0x355427db */ | ||
| 24 | /* |(log(1+s)-log(1-s))/s - Lg(s)| < 2**-34.24 (~[-4.95e-11, 4.97e-11]). */ | ||
| 25 | Lg1 = 0xaaaaaa.0p-24, /* 0.66666662693 */ | ||
| 26 | Lg2 = 0xccce13.0p-25, /* 0.40000972152 */ | ||
| 27 | Lg3 = 0x91e9ee.0p-25, /* 0.28498786688 */ | ||
| 28 | Lg4 = 0xf89e26.0p-26; /* 0.24279078841 */ | ||
| 29 | |||
| 30 | float log10f(float x) | ||
| 31 | { | ||
| 32 | 	union {float f; uint32_t i;} u = {x}; | ||
| 33 | 	float_t hfsq,f,s,z,R,w,t1,t2,dk,hi,lo; | ||
| 34 | 	uint32_t ix; | ||
| 35 | 	int k; | ||
| 36 | |||
| 37 | 	ix = u.i; | ||
| 38 | 	k = 0; | ||
| 39 | 	if (ix < 0x00800000 || ix>>31) { /* x < 2**-126 */ | ||
| 40 | 		if (ix<<1 == 0) | ||
| 41 | 			return -1/(x*x); /* log(+-0)=-inf */ | ||
| 42 | 		if (ix>>31) | ||
| 43 | 			return (x-x)/0.0f; /* log(-#) = NaN */ | ||
| 44 | 		/* subnormal number, scale up x */ | ||
| 45 | 		k -= 25; | ||
| 46 | 		x *= 0x1p25f; | ||
| 47 | 		u.f = x; | ||
| 48 | 		ix = u.i; | ||
| 49 | 	} else if (ix >= 0x7f800000) { | ||
| 50 | 		return x; | ||
| 51 | 	} else if (ix == 0x3f800000) | ||
| 52 | 		return 0; | ||
| 53 | |||
| 54 | 	/* reduce x into [sqrt(2)/2, sqrt(2)] */ | ||
| 55 | 	ix += 0x3f800000 - 0x3f3504f3; | ||
| 56 | 	k += (int)(ix>>23) - 0x7f; | ||
| 57 | 	ix = (ix&0x007fffff) + 0x3f3504f3; | ||
| 58 | 	u.i = ix; | ||
| 59 | 	x = u.f; | ||
| 60 | |||
| 61 | 	f = x - 1.0f; | ||
| 62 | 	s = f/(2.0f + f); | ||
| 63 | 	z = s*s; | ||
| 64 | 	w = z*z; | ||
| 65 | 	t1= w*(Lg2+w*Lg4); | ||
| 66 | 	t2= z*(Lg1+w*Lg3); | ||
| 67 | 	R = t2 + t1; | ||
| 68 | 	hfsq = 0.5f*f*f; | ||
| 69 | |||
| 70 | 	hi = f - hfsq; | ||
| 71 | 	u.f = hi; | ||
| 72 | 	u.i &= 0xfffff000; | ||
| 73 | 	hi = u.f; | ||
| 74 | 	lo = f - hi - hfsq + s*(hfsq+R); | ||
| 75 | 	dk = k; | ||
| 76 | 	return dk*log10_2lo + (lo+hi)*ivln10lo + lo*ivln10hi + hi*ivln10hi + dk*log10_2hi; | ||
| 77 | } | ||
lib/libc/musl/src/math/log2.c deleted-122| ... | @@ -1,122 +0,0 @@ | ||
| 1 | /* | ||
| 2 | * Double-precision log2(x) function. | ||
| 3 | * | ||
| 4 | * Copyright (c) 2018, Arm Limited. | ||
| 5 | * SPDX-License-Identifier: MIT | ||
| 6 | */ | ||
| 7 | |||
| 8 | #include <math.h> | ||
| 9 | #include <stdint.h> | ||
| 10 | #include "libm.h" | ||
| 11 | #include "log2_data.h" | ||
| 12 | |||
| 13 | #define T __log2_data.tab | ||
| 14 | #define T2 __log2_data.tab2 | ||
| 15 | #define B __log2_data.poly1 | ||
| 16 | #define A __log2_data.poly | ||
| 17 | #define InvLn2hi __log2_data.invln2hi | ||
| 18 | #define InvLn2lo __log2_data.invln2lo | ||
| 19 | #define N (1 << LOG2_TABLE_BITS) | ||
| 20 | #define OFF 0x3fe6000000000000 | ||
| 21 | |||
| 22 | /* Top 16 bits of a double. */ | ||
| 23 | static inline uint32_t top16(double x) | ||
| 24 | { | ||
| 25 | 	return asuint64(x) >> 48; | ||
| 26 | } | ||
| 27 | |||
| 28 | double log2(double x) | ||
| 29 | { | ||
| 30 | 	double_t z, r, r2, r4, y, invc, logc, kd, hi, lo, t1, t2, t3, p; | ||
| 31 | 	uint64_t ix, iz, tmp; | ||
| 32 | 	uint32_t top; | ||
| 33 | 	int k, i; | ||
| 34 | |||
| 35 | 	ix = asuint64(x); | ||
| 36 | 	top = top16(x); | ||
| 37 | #define LO asuint64(1.0 - 0x1.5b51p-5) | ||
| 38 | #define HI asuint64(1.0 + 0x1.6ab2p-5) | ||
| 39 | 	if (predict_false(ix - LO < HI - LO)) { | ||
| 40 | 		/* Handle close to 1.0 inputs separately. */ | ||
| 41 | 		/* Fix sign of zero with downward rounding when x==1. */ | ||
| 42 | 		if (WANT_ROUNDING && predict_false(ix == asuint64(1.0))) | ||
| 43 | 			return 0; | ||
| 44 | 		r = x - 1.0; | ||
| 45 | #if __FP_FAST_FMA | ||
| 46 | 		hi = r * InvLn2hi; | ||
| 47 | 		lo = r * InvLn2lo + __builtin_fma(r, InvLn2hi, -hi); | ||
| 48 | #else | ||
| 49 | 		double_t rhi, rlo; | ||
| 50 | 		rhi = asdouble(asuint64(r) & -1ULL << 32); | ||
| 51 | 		rlo = r - rhi; | ||
| 52 | 		hi = rhi * InvLn2hi; | ||
| 53 | 		lo = rlo * InvLn2hi + r * InvLn2lo; | ||
| 54 | #endif | ||
| 55 | 		r2 = r * r; /* rounding error: 0x1p-62. */ | ||
| 56 | 		r4 = r2 * r2; | ||
| 57 | 		/* Worst-case error is less than 0.54 ULP (0.55 ULP without fma). */ | ||
| 58 | 		p = r2 * (B[0] + r * B[1]); | ||
| 59 | 		y = hi + p; | ||
| 60 | 		lo += hi - y + p; | ||
| 61 | 		lo += r4 * (B[2] + r * B[3] + r2 * (B[4] + r * B[5]) + | ||
| 62 | 			 r4 * (B[6] + r * B[7] + r2 * (B[8] + r * B[9]))); | ||
| 63 | 		y += lo; | ||
| 64 | 		return eval_as_double(y); | ||
| 65 | 	} | ||
| 66 | 	if (predict_false(top - 0x0010 >= 0x7ff0 - 0x0010)) { | ||
| 67 | 		/* x < 0x1p-1022 or inf or nan. */ | ||
| 68 | 		if (ix * 2 == 0) | ||
| 69 | 			return __math_divzero(1); | ||
| 70 | 		if (ix == asuint64(INFINITY)) /* log(inf) == inf. */ | ||
| 71 | 			return x; | ||
| 72 | 		if ((top & 0x8000) || (top & 0x7ff0) == 0x7ff0) | ||
| 73 | 			return __math_invalid(x); | ||
| 74 | 		/* x is subnormal, normalize it. */ | ||
| 75 | 		ix = asuint64(x * 0x1p52); | ||
| 76 | 		ix -= 52ULL << 52; | ||
| 77 | 	} | ||
| 78 | |||
| 79 | 	/* x = 2^k z; where z is in range [OFF,2*OFF) and exact. | ||
| 80 | 	 The range is split into N subintervals. | ||
| 81 | 	 The ith subinterval contains z and c is near its center. */ | ||
| 82 | 	tmp = ix - OFF; | ||
| 83 | 	i = (tmp >> (52 - LOG2_TABLE_BITS)) % N; | ||
| 84 | 	k = (int64_t)tmp >> 52; /* arithmetic shift */ | ||
| 85 | 	iz = ix - (tmp & 0xfffULL << 52); | ||
| 86 | 	invc = T[i].invc; | ||
| 87 | 	logc = T[i].logc; | ||
| 88 | 	z = asdouble(iz); | ||
| 89 | 	kd = (double_t)k; | ||
| 90 | |||
| 91 | 	/* log2(x) = log2(z/c) + log2(c) + k. */ | ||
| 92 | 	/* r ~= z/c - 1, |r| < 1/(2*N). */ | ||
| 93 | #if __FP_FAST_FMA | ||
| 94 | 	/* rounding error: 0x1p-55/N. */ | ||
| 95 | 	r = __builtin_fma(z, invc, -1.0); | ||
| 96 | 	t1 = r * InvLn2hi; | ||
| 97 | 	t2 = r * InvLn2lo + __builtin_fma(r, InvLn2hi, -t1); | ||
| 98 | #else | ||
| 99 | 	double_t rhi, rlo; | ||
| 100 | 	/* rounding error: 0x1p-55/N + 0x1p-65. */ | ||
| 101 | 	r = (z - T2[i].chi - T2[i].clo) * invc; | ||
| 102 | 	rhi = asdouble(asuint64(r) & -1ULL << 32); | ||
| 103 | 	rlo = r - rhi; | ||
| 104 | 	t1 = rhi * InvLn2hi; | ||
| 105 | 	t2 = rlo * InvLn2hi + r * InvLn2lo; | ||
| 106 | #endif | ||
| 107 | |||
| 108 | 	/* hi + lo = r/ln2 + log2(c) + k. */ | ||
| 109 | 	t3 = kd + logc; | ||
| 110 | 	hi = t3 + t1; | ||
| 111 | 	lo = t3 - hi + t1 + t2; | ||
| 112 | |||
| 113 | 	/* log2(r+1) = r/ln2 + r^2*poly(r). */ | ||
| 114 | 	/* Evaluation is optimized assuming superscalar pipelined execution. */ | ||
| 115 | 	r2 = r * r; /* rounding error: 0x1p-54/N^2. */ | ||
| 116 | 	r4 = r2 * r2; | ||
| 117 | 	/* Worst-case error if |y| > 0x1p-4: 0.547 ULP (0.550 ULP without fma). | ||
| 118 | 	 ~ 0.5 + 2/N/ln2 + abs-poly-error*0x1p56 ULP (+ 0.003 ULP without fma). */ | ||
| 119 | 	p = A[0] + r * A[1] + r2 * (A[2] + r * A[3]) + r4 * (A[4] + r * A[5]); | ||
| 120 | 	y = lo + r2 * p + hi; | ||
| 121 | 	return eval_as_double(y); | ||
| 122 | } | ||
lib/libc/musl/src/math/log2f.c deleted-72| ... | @@ -1,72 +0,0 @@ | ||
| 1 | /* | ||
| 2 | * Single-precision log2 function. | ||
| 3 | * | ||
| 4 | * Copyright (c) 2017-2018, Arm Limited. | ||
| 5 | * SPDX-License-Identifier: MIT | ||
| 6 | */ | ||
| 7 | |||
| 8 | #include <math.h> | ||
| 9 | #include <stdint.h> | ||
| 10 | #include "libm.h" | ||
| 11 | #include "log2f_data.h" | ||
| 12 | |||
| 13 | /* | ||
| 14 | LOG2F_TABLE_BITS = 4 | ||
| 15 | LOG2F_POLY_ORDER = 4 | ||
| 16 | |||
| 17 | ULP error: 0.752 (nearest rounding.) | ||
| 18 | Relative error: 1.9 * 2^-26 (before rounding.) | ||
| 19 | */ | ||
| 20 | |||
| 21 | #define N (1 << LOG2F_TABLE_BITS) | ||
| 22 | #define T __log2f_data.tab | ||
| 23 | #define A __log2f_data.poly | ||
| 24 | #define OFF 0x3f330000 | ||
| 25 | |||
| 26 | float log2f(float x) | ||
| 27 | { | ||
| 28 | 	double_t z, r, r2, p, y, y0, invc, logc; | ||
| 29 | 	uint32_t ix, iz, top, tmp; | ||
| 30 | 	int k, i; | ||
| 31 | |||
| 32 | 	ix = asuint(x); | ||
| 33 | 	/* Fix sign of zero with downward rounding when x==1. */ | ||
| 34 | 	if (WANT_ROUNDING && predict_false(ix == 0x3f800000)) | ||
| 35 | 		return 0; | ||
| 36 | 	if (predict_false(ix - 0x00800000 >= 0x7f800000 - 0x00800000)) { | ||
| 37 | 		/* x < 0x1p-126 or inf or nan. */ | ||
| 38 | 		if (ix * 2 == 0) | ||
| 39 | 			return __math_divzerof(1); | ||
| 40 | 		if (ix == 0x7f800000) /* log2(inf) == inf. */ | ||
| 41 | 			return x; | ||
| 42 | 		if ((ix & 0x80000000) || ix * 2 >= 0xff000000) | ||
| 43 | 			return __math_invalidf(x); | ||
| 44 | 		/* x is subnormal, normalize it. */ | ||
| 45 | 		ix = asuint(x * 0x1p23f); | ||
| 46 | 		ix -= 23 << 23; | ||
| 47 | 	} | ||
| 48 | |||
| 49 | 	/* x = 2^k z; where z is in range [OFF,2*OFF] and exact. | ||
| 50 | 	 The range is split into N subintervals. | ||
| 51 | 	 The ith subinterval contains z and c is near its center. */ | ||
| 52 | 	tmp = ix - OFF; | ||
| 53 | 	i = (tmp >> (23 - LOG2F_TABLE_BITS)) % N; | ||
| 54 | 	top = tmp & 0xff800000; | ||
| 55 | 	iz = ix - top; | ||
| 56 | 	k = (int32_t)tmp >> 23; /* arithmetic shift */ | ||
| 57 | 	invc = T[i].invc; | ||
| 58 | 	logc = T[i].logc; | ||
| 59 | 	z = (double_t)asfloat(iz); | ||
| 60 | |||
| 61 | 	/* log2(x) = log1p(z/c-1)/ln2 + log2(c) + k */ | ||
| 62 | 	r = z * invc - 1; | ||
| 63 | 	y0 = logc + (double_t)k; | ||
| 64 | |||
| 65 | 	/* Pipelined polynomial evaluation to approximate log1p(r)/ln2. */ | ||
| 66 | 	r2 = r * r; | ||
| 67 | 	y = A[1] * r + A[2]; | ||
| 68 | 	y = A[0] * r2 + y; | ||
| 69 | 	p = A[3] * r + y0; | ||
| 70 | 	y = y * r2 + p; | ||
| 71 | 	return eval_as_float(y); | ||
| 72 | } | ||
lib/libc/musl/src/math/logf.c deleted-71| ... | @@ -1,71 +0,0 @@ | ||
| 1 | /* | ||
| 2 | * Single-precision log function. | ||
| 3 | * | ||
| 4 | * Copyright (c) 2017-2018, Arm Limited. | ||
| 5 | * SPDX-License-Identifier: MIT | ||
| 6 | */ | ||
| 7 | |||
| 8 | #include <math.h> | ||
| 9 | #include <stdint.h> | ||
| 10 | #include "libm.h" | ||
| 11 | #include "logf_data.h" | ||
| 12 | |||
| 13 | /* | ||
| 14 | LOGF_TABLE_BITS = 4 | ||
| 15 | LOGF_POLY_ORDER = 4 | ||
| 16 | |||
| 17 | ULP error: 0.818 (nearest rounding.) | ||
| 18 | Relative error: 1.957 * 2^-26 (before rounding.) | ||
| 19 | */ | ||
| 20 | |||
| 21 | #define T __logf_data.tab | ||
| 22 | #define A __logf_data.poly | ||
| 23 | #define Ln2 __logf_data.ln2 | ||
| 24 | #define N (1 << LOGF_TABLE_BITS) | ||
| 25 | #define OFF 0x3f330000 | ||
| 26 | |||
| 27 | float logf(float x) | ||
| 28 | { | ||
| 29 | 	double_t z, r, r2, y, y0, invc, logc; | ||
| 30 | 	uint32_t ix, iz, tmp; | ||
| 31 | 	int k, i; | ||
| 32 | |||
| 33 | 	ix = asuint(x); | ||
| 34 | 	/* Fix sign of zero with downward rounding when x==1. */ | ||
| 35 | 	if (WANT_ROUNDING && predict_false(ix == 0x3f800000)) | ||
| 36 | 		return 0; | ||
| 37 | 	if (predict_false(ix - 0x00800000 >= 0x7f800000 - 0x00800000)) { | ||
| 38 | 		/* x < 0x1p-126 or inf or nan. */ | ||
| 39 | 		if (ix * 2 == 0) | ||
| 40 | 			return __math_divzerof(1); | ||
| 41 | 		if (ix == 0x7f800000) /* log(inf) == inf. */ | ||
| 42 | 			return x; | ||
| 43 | 		if ((ix & 0x80000000) || ix * 2 >= 0xff000000) | ||
| 44 | 			return __math_invalidf(x); | ||
| 45 | 		/* x is subnormal, normalize it. */ | ||
| 46 | 		ix = asuint(x * 0x1p23f); | ||
| 47 | 		ix -= 23 << 23; | ||
| 48 | 	} | ||
| 49 | |||
| 50 | 	/* x = 2^k z; where z is in range [OFF,2*OFF] and exact. | ||
| 51 | 	 The range is split into N subintervals. | ||
| 52 | 	 The ith subinterval contains z and c is near its center. */ | ||
| 53 | 	tmp = ix - OFF; | ||
| 54 | 	i = (tmp >> (23 - LOGF_TABLE_BITS)) % N; | ||
| 55 | 	k = (int32_t)tmp >> 23; /* arithmetic shift */ | ||
| 56 | 	iz = ix - (tmp & 0xff800000); | ||
| 57 | 	invc = T[i].invc; | ||
| 58 | 	logc = T[i].logc; | ||
| 59 | 	z = (double_t)asfloat(iz); | ||
| 60 | |||
| 61 | 	/* log(x) = log1p(z/c-1) + log(c) + k*Ln2 */ | ||
| 62 | 	r = z * invc - 1; | ||
| 63 | 	y0 = logc + (double_t)k * Ln2; | ||
| 64 | |||
| 65 | 	/* Pipelined polynomial evaluation to approximate log1p(r). */ | ||
| 66 | 	r2 = r * r; | ||
| 67 | 	y = A[1] * r + A[2]; | ||
| 68 | 	y = A[0] * r2 + y; | ||
| 69 | 	y = y * r2 + (y0 + r); | ||
| 70 | 	return eval_as_float(y); | ||
| 71 | } | ||
src/libs/mingw.zig-2| ... | @@ -997,8 +997,6 @@ const mingw32_x86_src = [_][]const u8{ | ... | @@ -997,8 +997,6 @@ const mingw32_x86_src = [_][]const u8{ |
| 997 | const mingw32_x86_32_src = [_][]const u8{ | 997 | const mingw32_x86_32_src = [_][]const u8{ |
| 998 | // ucrtbase | 998 | // ucrtbase |
| 999 | "math" ++ path.sep_str ++ "coshf.c", | 999 | "math" ++ path.sep_str ++ "coshf.c", |
| 1000 | "math" ++ path.sep_str ++ "log10f.c", | ||
| 1001 | "math" ++ path.sep_str ++ "logf.c", | ||
| 1002 | "math" ++ path.sep_str ++ "modff.c", | 1000 | "math" ++ path.sep_str ++ "modff.c", |
| 1003 | "math" ++ path.sep_str ++ "powf.c", | 1001 | "math" ++ path.sep_str ++ "powf.c", |
| 1004 | "math" ++ path.sep_str ++ "sinhf.c", | 1002 | "math" ++ path.sep_str ++ "sinhf.c", |
src/libs/musl.zig-6| ... | @@ -1016,23 +1016,17 @@ const src_files = [_][]const u8{ | ... | @@ -1016,23 +1016,17 @@ const src_files = [_][]const u8{ |
| 1016 | "musl/src/math/llround.c", | 1016 | "musl/src/math/llround.c", |
| 1017 | "musl/src/math/llroundf.c", | 1017 | "musl/src/math/llroundf.c", |
| 1018 | "musl/src/math/llroundl.c", | 1018 | "musl/src/math/llroundl.c", |
| 1019 | "musl/src/math/log10.c", | ||
| 1020 | "musl/src/math/log10f.c", | ||
| 1021 | "musl/src/math/log10l.c", | 1019 | "musl/src/math/log10l.c", |
| 1022 | "musl/src/math/log1p.c", | 1020 | "musl/src/math/log1p.c", |
| 1023 | "musl/src/math/log1pf.c", | 1021 | "musl/src/math/log1pf.c", |
| 1024 | "musl/src/math/log1pl.c", | 1022 | "musl/src/math/log1pl.c", |
| 1025 | "musl/src/math/log2.c", | ||
| 1026 | "musl/src/math/log2_data.c", | 1023 | "musl/src/math/log2_data.c", |
| 1027 | "musl/src/math/log2f.c", | ||
| 1028 | "musl/src/math/log2f_data.c", | 1024 | "musl/src/math/log2f_data.c", |
| 1029 | "musl/src/math/log2l.c", | 1025 | "musl/src/math/log2l.c", |
| 1030 | "musl/src/math/logb.c", | 1026 | "musl/src/math/logb.c", |
| 1031 | "musl/src/math/logbf.c", | 1027 | "musl/src/math/logbf.c", |
| 1032 | "musl/src/math/logbl.c", | 1028 | "musl/src/math/logbl.c", |
| 1033 | "musl/src/math/log.c", | ||
| 1034 | "musl/src/math/log_data.c", | 1029 | "musl/src/math/log_data.c", |
| 1035 | "musl/src/math/logf.c", | ||
| 1036 | "musl/src/math/logf_data.c", | 1030 | "musl/src/math/logf_data.c", |
| 1037 | "musl/src/math/logl.c", | 1031 | "musl/src/math/logl.c", |
| 1038 | "musl/src/math/lrint.c", | 1032 | "musl/src/math/lrint.c", |
src/libs/wasi_libc.zig-6| ... | @@ -785,23 +785,17 @@ const libc_top_half_src_files = [_][]const u8{ | ... | @@ -785,23 +785,17 @@ const libc_top_half_src_files = [_][]const u8{ |
| 785 | "musl/src/math/llround.c", | 785 | "musl/src/math/llround.c", |
| 786 | "musl/src/math/llroundf.c", | 786 | "musl/src/math/llroundf.c", |
| 787 | "musl/src/math/llroundl.c", | 787 | "musl/src/math/llroundl.c", |
| 788 | "musl/src/math/log10.c", | ||
| 789 | "musl/src/math/log10f.c", | ||
| 790 | "musl/src/math/log10l.c", | 788 | "musl/src/math/log10l.c", |
| 791 | "musl/src/math/log1p.c", | 789 | "musl/src/math/log1p.c", |
| 792 | "musl/src/math/log1pf.c", | 790 | "musl/src/math/log1pf.c", |
| 793 | "musl/src/math/log1pl.c", | 791 | "musl/src/math/log1pl.c", |
| 794 | "musl/src/math/log2.c", | ||
| 795 | "musl/src/math/log2_data.c", | 792 | "musl/src/math/log2_data.c", |
| 796 | "musl/src/math/log2f.c", | ||
| 797 | "musl/src/math/log2f_data.c", | 793 | "musl/src/math/log2f_data.c", |
| 798 | "musl/src/math/log2l.c", | 794 | "musl/src/math/log2l.c", |
| 799 | "musl/src/math/logb.c", | 795 | "musl/src/math/logb.c", |
| 800 | "musl/src/math/logbf.c", | 796 | "musl/src/math/logbf.c", |
| 801 | "musl/src/math/logbl.c", | 797 | "musl/src/math/logbl.c", |
| 802 | "musl/src/math/log.c", | ||
| 803 | "musl/src/math/log_data.c", | 798 | "musl/src/math/log_data.c", |
| 804 | "musl/src/math/logf.c", | ||
| 805 | "musl/src/math/logf_data.c", | 799 | "musl/src/math/logf_data.c", |
| 806 | "musl/src/math/logl.c", | 800 | "musl/src/math/logl.c", |
| 807 | "musl/src/math/lrint.c", | 801 | "musl/src/math/lrint.c", |