FreeBSD Bugzilla – Attachment 229010 Details for
Bug 218514
[LIBM] implementations of sinpi[fl], cospi[fl], and tanpi[fl]
Home
|
New
|
Browse
|
Search
|
[?]
|
Reports
|
Help
|
New Account
|
Log In
Remember
[x]
|
Forgot Password
Login:
[x]
[patch]
new new updated patch
cospi_20211024.diff (text/plain), 71.06 KB, created by
Steve Kargl
on 2021-10-25 03:11:17 UTC
(
hide
)
Description:
new new updated patch
Filename:
MIME Type:
Creator:
Steve Kargl
Created:
2021-10-25 03:11:17 UTC
Size:
71.06 KB
patch
obsolete
>diff --git a/lib/msun/Makefile b/lib/msun/Makefile >index 7107aad56aa..5164aef6020 100644 >--- a/lib/msun/Makefile >+++ b/lib/msun/Makefile >@@ -138,6 +138,11 @@ COMMON_SRCS+= catrig.c catrigf.c \ > s_cproj.c s_cprojf.c s_creal.c s_crealf.c s_creall.c \ > s_csinh.c s_csinhf.c s_ctanh.c s_ctanhf.c > >+# IEEE-754 2008 and ISO/IEC TS 18661-4 half-cycle trignometric functions >+COMMON_SRCS+= s_cospi.c s_cospif.c s_cospil.c \ >+ s_sinpi.c s_sinpif.c s_sinpil.c \ >+ s_tanpi.c s_tanpif.c s_tanpil.c >+ > # FreeBSD's C library supplies these functions: > #COMMON_SRCS+= s_fabs.c s_frexp.c s_isnan.c s_ldexp.c s_modf.c > >@@ -154,7 +159,8 @@ INCS+= fenv.h math.h > > MAN= acos.3 acosh.3 asin.3 asinh.3 atan.3 atan2.3 atanh.3 \ > ceil.3 cacos.3 ccos.3 ccosh.3 cexp.3 \ >- cimag.3 clog.3 copysign.3 cos.3 cosh.3 cpow.3 csqrt.3 erf.3 \ >+ cimag.3 clog.3 copysign.3 cos.3 cosh.3 cospi.3 \ >+ cpow.3 csqrt.3 erf.3 \ > exp.3 fabs.3 fdim.3 \ > feclearexcept.3 feenableexcept.3 fegetenv.3 \ > fegetround.3 fenv.3 floor.3 \ >@@ -162,7 +168,7 @@ MAN= acos.3 acosh.3 asin.3 asinh.3 atan.3 atan2.3 atanh.3 \ > lgamma.3 log.3 lrint.3 lround.3 math.3 nan.3 \ > nextafter.3 remainder.3 rint.3 \ > round.3 scalbn.3 signbit.3 sin.3 sincos.3 \ >- sinh.3 sqrt.3 tan.3 tanh.3 trunc.3 \ >+ sinh.3 sinpi.3 sqrt.3 tan.3 tanh.3 tanpi.3 trunc.3 \ > complex.3 > > MLINKS+=acos.3 acosf.3 acos.3 acosl.3 >@@ -192,6 +198,7 @@ MLINKS+=clog.3 clogf.3 clog.3 clogl.3 > MLINKS+=copysign.3 copysignf.3 copysign.3 copysignl.3 > MLINKS+=cos.3 cosf.3 cos.3 cosl.3 > MLINKS+=cosh.3 coshf.3 cosh.3 coshl.3 >+MLINKS+=cospi.3 cospif.3 cospi.3 cospil.3 > MLINKS+=cpow.3 cpowf.3 cpow.3 cpowl.3 > MLINKS+=csqrt.3 csqrtf.3 csqrt.3 csqrtl.3 > MLINKS+=erf.3 erfc.3 erf.3 erff.3 erf.3 erfcf.3 erf.3 erfl.3 erf.3 erfcl.3 >@@ -244,10 +251,12 @@ MLINKS+=scalbn.3 scalbnf.3 scalbn.3 scalbnl.3 > MLINKS+=sin.3 sinf.3 sin.3 sinl.3 > MLINKS+=sincos.3 sincosf.3 sin.3 sincosl.3 > MLINKS+=sinh.3 sinhf.3 sinh.3 sinhl.3 >+MLINKS+=sinpi.3 sinpif.3 sinpi.3 sinpil.3 > MLINKS+=sqrt.3 cbrt.3 sqrt.3 cbrtf.3 sqrt.3 cbrtl.3 sqrt.3 sqrtf.3 \ > sqrt.3 sqrtl.3 > MLINKS+=tan.3 tanf.3 tan.3 tanl.3 > MLINKS+=tanh.3 tanhf.3 tanh.3 tanhl.3 >+MLINKS+=tanpi.3 tanpif.3 tanpi.3 tanpil.3 > MLINKS+=trunc.3 truncf.3 trunc.3 truncl.3 > > .include <src.opts.mk> >diff --git a/lib/msun/Symbol.map b/lib/msun/Symbol.map >index 51deded732f..317c0de77f9 100644 >--- a/lib/msun/Symbol.map >+++ b/lib/msun/Symbol.map >@@ -304,3 +304,16 @@ FBSD_1.5 { > sincosf; > sincosl; > }; >+ >+/* First added in 14.0-CURRENT */ >+FBSD_1.7{ >+ cospi; >+ cospif; >+ cospil; >+ sinpi; >+ sinpif; >+ sinpil; >+ tanpi; >+ tanpif; >+ tanpil; >+}; >diff --git a/lib/msun/ld128/k_cospil.h b/lib/msun/ld128/k_cospil.h >new file mode 100644 >index 00000000000..592f1922953 >--- /dev/null >+++ b/lib/msun/ld128/k_cospil.h >@@ -0,0 +1,42 @@ >+/*- >+ * Copyright (c) 2017 Steven G. Kargl >+ * 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 unmodified, 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 ``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. >+ */ >+ >+/* >+ * See ../src/k_cospi.c for implementation details. >+ */ >+ >+static inline long double >+__kernel_cospil(long double x) >+{ >+ long double hi, lo; >+ >+ hi = (double)x; >+ lo = x - hi; >+ lo = lo * (pi_lo + pi_hi) + hi * pi_lo; >+ hi *= pi_hi; >+ _2sumF(hi, lo); >+ return (__kernel_cosl(hi, lo)); >+} >diff --git a/lib/msun/ld128/k_sinpil.h b/lib/msun/ld128/k_sinpil.h >new file mode 100644 >index 00000000000..fa4e7d6336d >--- /dev/null >+++ b/lib/msun/ld128/k_sinpil.h >@@ -0,0 +1,42 @@ >+/*- >+ * Copyright (c) 2017 Steven G. Kargl >+ * 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 unmodified, 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 ``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. >+ */ >+ >+/* >+ * See ../src/k_sinpi.c for implementation details. >+ */ >+ >+static inline long double >+__kernel_sinpil(long double x) >+{ >+ long double hi, lo; >+ >+ hi = (double)x; >+ lo = x - hi; >+ lo = lo * (pi_lo + pi_hi) + hi * pi_lo; >+ hi *= pi_hi; >+ _2sumF(hi, lo); >+ return (__kernel_sinl(hi, lo, 1)); >+} >diff --git a/lib/msun/ld128/s_cospil.c b/lib/msun/ld128/s_cospil.c >new file mode 100644 >index 00000000000..b4bc50bb4d8 >--- /dev/null >+++ b/lib/msun/ld128/s_cospil.c >@@ -0,0 +1,109 @@ >+/*- >+ * Copyright (c) 2017 Steven G. Kargl >+ * 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 unmodified, 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 ``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. >+ */ >+ >+/* >+ * See ../src/s_cospi.c for implementation details. >+ * >+ * FIXME: This has not been compiled nor has it been tested for accuracy. >+ * FIXME: This should use bit twiddling. >+ */ >+ >+#include "math.h" >+#include "math_private.h" >+ >+/* >+ * pi_hi contains the leading 56 bits of a 169 bit approximation for pi. >+ */ >+static const long double >+pi_hi = 3.14159265358979322702026593105983920e+00L, >+pi_lo = 1.14423774522196636802434264184180742e-17L; >+ >+#include "k_cospil.h" >+#include "k_sinpil.h" >+ >+volatile static const double vzero = 0; >+ >+long double >+cospil(long double x) >+{ >+ long double ax, c, xf; >+ uint32_t ix; >+ >+ ax = fabsl(x); >+ >+ if (ax < 1) { >+ if (ax < 0.25) { >+ if (ax < 0x1p-60) { >+ if ((int)x == 0) >+ return (1); >+ } >+ return (__kernel_cospil(ax)); >+ } >+ >+ if (ax < 0.5) >+ c = __kernel_sinpil(0.5 - ax); >+ else if (ax < 0.75) { >+ if (ax == 0.5) >+ return (0); >+ c = -__kernel_sinpil(ax - 0.5); >+ } else >+ c = -__kernel_cospil(1 - ax); >+ return (c); >+ } >+ >+ if (ax < 0x1p112) { >+ xf = floorl(ax); >+ ax -= xf; >+ if (x < 0.5) { >+ if (x < 0.25) >+ c = ax == 0 ? 1 : __kernel_cospil(ax); >+ else >+ c = __kernel_sinpil(0.5 - ax); >+ } else { >+ if (x < 0.75) { >+ if (ax == 0.5) >+ return (0); >+ c = -__kernel_sinpil(ax - 0.5); >+ } else >+ c = -__kernel_cospil(1 - ax); >+ } >+ >+ if (xf > 0x1p50) >+ xf -= 0x1p50; >+ if (xf > 0x1p30) >+ xf -= 0x1p30; >+ ix = (uint32_t)xf; >+ return (ix & 1 ? -c : c); >+ } >+ >+ if (isinf(x) || isnan(x)) >+ return (vzero / vzero); >+ >+ /* >+ * |x| >= 0x1p112 is always an even integer, so return 1. >+ */ >+ return (1); >+} >diff --git a/lib/msun/ld128/s_sinpil.c b/lib/msun/ld128/s_sinpil.c >new file mode 100644 >index 00000000000..56ede66dc32 >--- /dev/null >+++ b/lib/msun/ld128/s_sinpil.c >@@ -0,0 +1,118 @@ >+/*- >+ * Copyright (c) 2017 Steven G. Kargl >+ * 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 unmodified, 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 ``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. >+ */ >+ >+/* >+ * See ../src/s_sinpi.c for implementation details. >+ * >+ * FIXME: This has not been compiled nor has it been tested for accuracy. >+ * FIXME: This should use bit twiddling. >+ */ >+ >+#include "math.h" >+#include "math_private.h" >+ >+/* >+ * pi_hi contains the leading 56 bits of a 169 bit approximation for pi. >+ */ >+static const long double >+pi_hi = 3.14159265358979322702026593105983920e+00L, >+pi_lo = 1.14423774522196636802434264184180742e-17L; >+ >+#include "k_cospil.h" >+#include "k_sinpil.h" >+ >+volatile static const double vzero = 0; >+ >+long double >+sinpil(long double x) >+{ >+ long double ax, hi, lo, s, xf, xhi, xlo; >+ uint32_t ix; >+ >+ ax = fabsl(x); >+ >+ if (ax < 1) { >+ if (ax < 0.25) { >+ if (ax < 0x1p-60) { >+ if (x == 0) >+ return (x); >+ hi = (double)x; >+ hi *= 0x1p113L >+ lo = x * 0x1p113L - hi; >+ s = (pi_lo + pi_hi) * lo + pi_lo * lo + >+ pi_hi * hi; >+ return (s * 0x1p-113L); >+ } >+ >+ s = __kernel_sinpil(ax); >+ return (copysignl(s, x)); >+ } >+ >+ if (ax < 0.5) >+ s = __kernel_cospil(0.5 - ax); >+ else if (ax < 0.75) >+ s = __kernel_cospil(ax - 0.5); >+ else >+ s = __kernel_sinpil(1 - ax); >+ return (copysignl(s, x)); >+ } >+ >+ if (ax < 0x1p112) { >+ xf = floorl(ax); >+ ax -= xf; >+ if (ax == 0) { >+ s = 0; >+ } else { >+ if (ax < 0.5) { >+ if (ax <= 0.25) >+ s = __kernel_sinpil(ax); >+ else >+ s = __kernel_cospil(0.5 - ax); >+ } else { >+ if (ax < 0.75) >+ s = __kernel_cospil(ax - 0.5); >+ else >+ s = __kernel_sinpil(1 - ax); >+ } >+ >+ if (xf > 0x1p50) >+ xf -= 0x1p50; >+ if (xf > 0x1p30) >+ xf -= 0x1p30; >+ ix = (uint32_t)xf; >+ if (ix & 1) s = -s; >+ } >+ return (copysignl(s, x)); >+ } >+ >+ if (isinf(x) || isnan(x)) >+ return (vzero / vzero); >+ >+ /* >+ * |x| >= 0x1p112 is always an integer, so return +-0. >+ */ >+ return (copysignl(0, x)); >+} >diff --git a/lib/msun/ld128/s_tanpil.c b/lib/msun/ld128/s_tanpil.c >new file mode 100644 >index 00000000000..33a61cf3115 >--- /dev/null >+++ b/lib/msun/ld128/s_tanpil.c >@@ -0,0 +1,120 @@ >+/*- >+ * Copyright (c) 2017 Steven G. Kargl >+ * 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 unmodified, 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 ``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. >+ */ >+ >+/* >+ * See ../src/s_tanpi.c for implementation details. >+ * >+ * FIXME: This has not been compiled nor has it been tested for accuracy. >+ * FIXME: This should use bit twiddling. >+ */ >+ >+#include "math.h" >+#include "math_private.h" >+ >+/* >+ * pi_hi contains the leading 56 bits of a 169 bit approximation for pi. >+ */ >+static const long double >+pi_hi = 3.14159265358979322702026593105983920e+00L, >+pi_lo = 1.14423774522196636802434264184180742e-17L; >+ >+static inline long double >+__kernel_tanpi(long double x) >+{ >+ long double hi, lo, t; >+ >+ if (x < 0.25) { >+ hi = (double)x; >+ lo = x - hi; >+ lo = lo * (pi_lo + pi_hi) + hi * pi_lo; >+ hi *= pi_hi; >+ _2sumF(hi, lo); >+ t = __kernel_tanl(hi, lo, -1); >+ } else if (x > 0.25) { >+ x = 0.5 - x; >+ hi = (double)x; >+ lo = x - hi; >+ lo = lo * (pi_lo + pi_hi) + hi * pi_lo; >+ hi *= pi_hi; >+ _2sumF(hi, lo); >+ t = - __kernel_tanl(hi, lo, 1); >+ } else >+ t = 1; >+ >+ return (t); >+} >+ >+volatile static const double vzero = 0; >+ >+long double >+tanpil(long double x) >+{ >+ long double ax, hi, lo, xf; >+ uint32_t ix; >+ >+ ax = fabsl(ax); >+ >+ if (ax < 1) { >+ if (ax < 0.5) { >+ if (ax < 0x1p-60) { >+ if (x == 0) >+ return (x); >+ hi = (double)x; >+ hi *= 0x1p113L >+ lo = x * 0x1p113L - hi; >+ t = (pi_lo + pi_hi) * lo + pi_lo * lo + >+ pi_hi * hi; >+ return (t * 0x1p-113L); >+ } >+ t = __kernel_tanpil(ax); >+ } else if (ax == 0.5) >+ return ((ax - ax) / (ax - ax)); >+ else >+ t = -__kernel_tanpil(1 - ax); >+ return (copysignl(t, x)); >+ } >+ >+ if (ix < 0x1p112) { >+ xf = floorl(ax); >+ ax -= xf; >+ if (ax < 0.5) >+ t = ax == 0 ? 0 : __kernel_tanpil(ax); >+ else if (ax == 0.5) >+ return ((ax - ax) / (ax - ax)); >+ else >+ t = -__kernel_tanpil(1 - ax); >+ return (copysignl(t, x)); >+ } >+ >+ /* x = +-inf or nan. */ >+ if (isinf(x) || isnan(x)) >+ return (vzero / vzero); >+ >+ /* >+ * |x| >= 0x1p53 is always an integer, so return +-0. >+ */ >+ return (copysignl(0, x)); >+} >diff --git a/lib/msun/ld80/k_cospil.h b/lib/msun/ld80/k_cospil.h >new file mode 100644 >index 00000000000..6e13ef02aea >--- /dev/null >+++ b/lib/msun/ld80/k_cospil.h >@@ -0,0 +1,42 @@ >+/*- >+ * Copyright (c) 2017 Steven G. Kargl >+ * 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 unmodified, 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 ``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. >+ */ >+ >+/* >+ * See ../src/k_cospi.c for implementation details. >+ */ >+ >+static inline long double >+__kernel_cospil(long double x) >+{ >+ long double hi, lo; >+ >+ hi = (float)x; >+ lo = x - hi; >+ lo = lo * (pi_lo + pi_hi) + hi * pi_lo; >+ hi *= pi_hi; >+ _2sumF(hi, lo); >+ return (__kernel_cosl(hi, lo)); >+} >diff --git a/lib/msun/ld80/k_sinpil.h b/lib/msun/ld80/k_sinpil.h >new file mode 100644 >index 00000000000..00241b932e9 >--- /dev/null >+++ b/lib/msun/ld80/k_sinpil.h >@@ -0,0 +1,42 @@ >+/*- >+ * Copyright (c) 2017 Steven G. Kargl >+ * 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 unmodified, 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 ``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. >+ */ >+ >+/* >+ * See ../src/k_sinpi.c for implementation details. >+ */ >+ >+static inline long double >+__kernel_sinpil(long double x) >+{ >+ long double hi, lo; >+ >+ hi = (float)x; >+ lo = x - hi; >+ lo = lo * (pi_lo + pi_hi) + hi * pi_lo; >+ hi *= pi_hi; >+ _2sumF(hi, lo); >+ return (__kernel_sinl(hi, lo, 1)); >+} >diff --git a/lib/msun/ld80/s_cospil.c b/lib/msun/ld80/s_cospil.c >new file mode 100644 >index 00000000000..199479e9eaf >--- /dev/null >+++ b/lib/msun/ld80/s_cospil.c >@@ -0,0 +1,129 @@ >+/*- >+ * Copyright (c) 2017 Steven G. Kargl >+ * 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 unmodified, 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 ``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. >+ */ >+ >+/* >+ * See ../src/s_cospi.c for implementation details. >+ */ >+ >+#ifdef __i386__ >+#include <ieeefp.h> >+#endif >+#include <stdint.h> >+ >+#include "fpmath.h" >+#include "math.h" >+#include "math_private.h" >+ >+static const double >+pi_hi = 3.1415926814079285e+00, /* 0x400921fb 0x58000000 */ >+pi_lo =-2.7818135228334233e-08; /* 0xbe5dde97 0x3dcb3b3a */ >+ >+#include "k_cospil.h" >+#include "k_sinpil.h" >+ >+volatile static const double vzero = 0; >+ >+long double >+cospil(long double x) >+{ >+ long double ax, c; >+ uint64_t lx, m; >+ uint32_t j0; >+ uint16_t hx, ix; >+ >+ EXTRACT_LDBL80_WORDS(hx, lx, x); >+ ix = hx & 0x7fff; >+ INSERT_LDBL80_WORDS(ax, ix, lx); >+ >+ ENTERI(); >+ >+ if (ix < 0x3fff) { /* |x| < 1 */ >+ if (ix < 0x3ffd) { /* |x| < 0.25 */ >+ if (ix < 0x3fdd) { /* |x| < 0x1p-34 */ >+ if ((int)x == 0) >+ RETURNI(1); >+ } >+ RETURNI(__kernel_cospil(ax)); >+ } >+ >+ if (ix < 0x3ffe) /* |x| < 0.5 */ >+ c = __kernel_sinpil(0.5 - ax); >+ else if (lx < 0xc000000000000000ull) { /* |x| < 0.75 */ >+ if (ax == 0.5) >+ RETURNI(0); >+ c = -__kernel_sinpil(ax - 0.5); >+ } else >+ c = -__kernel_cospil(1 - ax); >+ RETURNI(c); >+ } >+ >+ if (ix < 0x403e) { /* 1 <= |x| < 0x1p63 */ >+ /* Determine integer part of ax. */ >+ j0 = ix - 0x3fff + 1; >+ if (j0 < 32) { >+ lx = (lx >> 32) << 32; >+ lx &= ~(((lx << 32)-1) >> j0); >+ } else { >+ m = (uint64_t)-1 >> (j0 + 1); >+ if (lx & m) lx &= ~m; >+ } >+ INSERT_LDBL80_WORDS(x, ix, lx); >+ >+ ax -= x; >+ EXTRACT_LDBL80_WORDS(ix, lx, ax); >+ >+ if (ix < 0x3ffe) { /* |x| < 0.5 */ >+ if (ix < 0x3ffd) /* |x| < 0.25 */ >+ c = ix == 0 ? 1 : __kernel_cospil(ax); >+ else >+ c = __kernel_sinpil(0.5 - ax); >+ >+ } else { >+ if (lx < 0xc000000000000000ull) { /* |x| < 0.75 */ >+ if (ax == 0.5) >+ RETURNI(0); >+ c = -__kernel_sinpil(ax - 0.5); >+ } else >+ c = -__kernel_cospil(1 - ax); >+ } >+ >+ if (j0 > 40) >+ x -= 0x1p40; >+ if (j0 > 30) >+ x -= 0x1p30; >+ j0 = (uint32_t)x; >+ >+ RETURNI(j0 & 1 ? -c : c); >+ } >+ >+ if (ix >= 0x7fff) >+ RETURNI(vzero / vzero); >+ >+ /* >+ * |x| >= 0x1p63 is always an even integer, so return 1. >+ */ >+ RETURNI(1); >+} >diff --git a/lib/msun/ld80/s_sinpil.c b/lib/msun/ld80/s_sinpil.c >new file mode 100644 >index 00000000000..4cefa92352e >--- /dev/null >+++ b/lib/msun/ld80/s_sinpil.c >@@ -0,0 +1,140 @@ >+/*- >+ * Copyright (c) 2017 Steven G. Kargl >+ * 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 unmodified, 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 ``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. >+ */ >+ >+/* >+ * See ../src/s_sinpi.c for implementation details. >+ */ >+ >+#ifdef __i386__ >+#include <ieeefp.h> >+#endif >+#include <stdint.h> >+ >+#include "fpmath.h" >+#include "math.h" >+#include "math_private.h" >+ >+static const union IEEEl2bits >+pi_hi_u = LD80C(0xc90fdaa200000000, 1, 3.14159265346825122833e+00L), >+pi_lo_u = LD80C(0x85a308d313198a2e, -33, 1.21542010130123852029e-10L); >+#define pi_hi (pi_hi_u.e) >+#define pi_lo (pi_lo_u.e) >+ >+#include "k_cospil.h" >+#include "k_sinpil.h" >+ >+volatile static const double vzero = 0; >+ >+long double >+sinpil(long double x) >+{ >+ long double ax, hi, lo, s; >+ uint64_t lx, m; >+ uint32_t j0; >+ uint16_t hx, ix; >+ >+ EXTRACT_LDBL80_WORDS(hx, lx, x); >+ ix = hx & 0x7fff; >+ INSERT_LDBL80_WORDS(ax, ix, lx); >+ >+ ENTERI(); >+ >+ if (ix < 0x3fff) { /* |x| < 1 */ >+ if (ix < 0x3ffd) { /* |x| < 0.25 */ >+ if (ix < 0x3fdd) { /* |x| < 0x1p-34 */ >+ if (x == 0) >+ RETURNI(x); >+ INSERT_LDBL80_WORDS(hi, hx, >+ lx & 0xffffffff00000000ull); >+ hi *= 0x1p63L; >+ lo = x * 0x1p63L - hi; >+ s = (pi_lo + pi_hi) * lo + pi_lo * hi + >+ pi_hi * hi; >+ RETURNI(s * 0x1p-63L); >+ } >+ s = __kernel_sinpil(ax); >+ RETURNI((hx & 0x8000) ? -s : s); >+ } >+ >+ if (ix < 0x3ffe) /* |x| < 0.5 */ >+ s = __kernel_cospil(0.5 - ax); >+ else if (lx < 0xc000000000000000ull) /* |x| < 0.75 */ >+ s = __kernel_cospil(ax - 0.5); >+ else >+ s = __kernel_sinpil(1 - ax); >+ RETURNI((hx & 0x8000) ? -s : s); >+ } >+ >+ if (ix < 0x403e) { /* 1 <= |x| < 0x1p63 */ >+ /* Determine integer part of ax. */ >+ j0 = ix - 0x3fff + 1; >+ if (j0 < 32) { >+ lx = (lx >> 32) << 32; >+ lx &= ~(((lx << 32)-1) >> j0); >+ } else { >+ m = (uint64_t)-1 >> (j0 + 1); >+ if (lx & m) lx &= ~m; >+ } >+ INSERT_LDBL80_WORDS(x, ix, lx); >+ >+ ax -= x; >+ EXTRACT_LDBL80_WORDS(ix, lx, ax); >+ >+ if (ix == 0) { >+ s = 0; >+ } else { >+ if (ix < 0x3ffe) { /* |x| < 0.5 */ >+ if (ix < 0x3ffd) /* |x| < 0.25 */ >+ s = __kernel_sinpil(ax); >+ else >+ s = __kernel_cospil(0.5 - ax); >+ } else { >+ /* |x| < 0.75 */ >+ if (lx < 0xc000000000000000ull) >+ s = __kernel_cospil(ax - 0.5); >+ else >+ s = __kernel_sinpil(1 - ax); >+ } >+ >+ if (j0 > 40) >+ x -= 0x1p40; >+ if (j0 > 30) >+ x -= 0x1p30; >+ j0 = (uint32_t)x; >+ if (j0 & 1) s = -s; >+ } >+ RETURNI((hx & 0x8000) ? -s : s); >+ } >+ >+ /* x = +-inf or nan. */ >+ if (ix >= 0x7fff) >+ RETURNI(vzero / vzero); >+ >+ /* >+ * |x| >= 0x1p63 is always an integer, so return +-0. >+ */ >+ RETURNI(copysignl(0, x)); >+} >diff --git a/lib/msun/ld80/s_tanpil.c b/lib/msun/ld80/s_tanpil.c >new file mode 100644 >index 00000000000..02451e56202 >--- /dev/null >+++ b/lib/msun/ld80/s_tanpil.c >@@ -0,0 +1,139 @@ >+/*- >+ * Copyright (c) 2017 Steven G. Kargl >+ * 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 unmodified, 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 ``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. >+ */ >+ >+/* >+ * See ../src/s_tanpi.c for implementation details. >+ */ >+ >+#ifdef __i386__ >+#include <ieeefp.h> >+#endif >+#include <stdint.h> >+ >+#include "fpmath.h" >+#include "math.h" >+#include "math_private.h" >+ >+static const double >+pi_hi = 3.1415926814079285e+00, /* 0x400921fb 0x58000000 */ >+pi_lo = -2.7818135228334233e-08; /* 0xbe5dde97 0x3dcb3b3a */ >+ >+static inline long double >+__kernel_tanpil(long double x) >+{ >+ long double hi, lo, t; >+ >+ if (x < 0.25) { >+ hi = (float)x; >+ lo = x - hi; >+ lo = lo * (pi_lo + pi_hi) + hi * pi_lo; >+ hi *= pi_hi; >+ _2sumF(hi, lo); >+ t = __kernel_tanl(hi, lo, -1); >+ } else if (x > 0.25) { >+ x = 0.5 - x; >+ hi = (float)x; >+ lo = x - hi; >+ lo = lo * (pi_lo + pi_hi) + hi * pi_lo; >+ hi *= pi_hi; >+ _2sumF(hi, lo); >+ t = - __kernel_tanl(hi, lo, 1); >+ } else >+ t = 1; >+ >+ return (t); >+} >+ >+volatile static const double vzero = 0; >+ >+long double >+tanpil(long double x) >+{ >+ long double ax, hi, lo, t; >+ uint64_t lx, m; >+ uint32_t j0; >+ uint16_t hx, ix; >+ >+ EXTRACT_LDBL80_WORDS(hx, lx, x); >+ ix = hx & 0x7fff; >+ INSERT_LDBL80_WORDS(ax, ix, lx); >+ >+ ENTERI(); >+ >+ if (ix < 0x3fff) { /* |x| < 1 */ >+ if (ix < 0x3ffe) { /* |x| < 0.5 */ >+ if (ix < 0x3fdd) { /* |x| < 0x1p-34 */ >+ if (x == 0) >+ RETURNI(x); >+ INSERT_LDBL80_WORDS(hi, hx, >+ lx & 0xffffffff00000000ull); >+ hi *= 0x1p63L; >+ lo = x * 0x1p63L - hi; >+ t = (pi_lo + pi_hi) * lo + pi_lo * hi + >+ pi_hi * hi; >+ RETURNI(t * 0x1p-63L); >+ } >+ t = __kernel_tanpil(ax); >+ } else if (ax == 0.5) >+ RETURNI((ax - ax) / (ax - ax)); >+ else >+ t = -__kernel_tanpil(1 - ax); >+ RETURNI((hx & 0x8000) ? -t : t); >+ } >+ >+ if (ix < 0x403e) { /* 1 <= |x| < 0x1p63 */ >+ /* Determine integer part of ax. */ >+ j0 = ix - 0x3fff + 1; >+ if (j0 < 32) { >+ lx = (lx >> 32) << 32; >+ lx &= ~(((lx << 32)-1) >> j0); >+ } else { >+ m = (uint64_t)-1 >> (j0 + 1); >+ if (lx & m) lx &= ~m; >+ } >+ INSERT_LDBL80_WORDS(x, ix, lx); >+ >+ ax -= x; >+ EXTRACT_LDBL80_WORDS(ix, lx, ax); >+ >+ if (ix < 0x3ffe) /* |x| < 0.5 */ >+ t = ax == 0 ? 0 : __kernel_tanpil(ax); >+ else if (ax == 0.5) >+ RETURNI((ax - ax) / (ax - ax)); >+ else >+ t = -__kernel_tanpil(1 - ax); >+ RETURNI((hx & 0x8000) ? -t : t); >+ } >+ >+ /* x = +-inf or nan. */ >+ if (ix >= 0x7fff) >+ RETURNI(vzero / vzero); >+ >+ /* >+ * |x| >= 0x1p63 is always an integer, so return +-0. >+ */ >+ RETURNI(copysignl(0, x)); >+} >diff --git a/lib/msun/man/cospi.3 b/lib/msun/man/cospi.3 >new file mode 100644 >index 00000000000..b76e1858ed3 >--- /dev/null >+++ b/lib/msun/man/cospi.3 >@@ -0,0 +1,111 @@ >+.\" Copyright (c) 2017 Steven G. Kargl <kargl@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 REGENTS 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 REGENTS 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$ >+.\" >+.Dd April 1, 2017 >+.Dt COSPI 3 >+.Os >+.Sh NAME >+.Nm cospi , >+.Nm cospif , >+.Nm cospil >+.Nd half\(encycle cosine functions >+.Sh LIBRARY >+.Lb libm >+.Sh SYNOPSIS >+.In math.h >+.Ft double >+.Fn cospi "double x" >+.Ft float >+.Fn cospif "float x" >+.Ft long double >+.Fn cospil "long double x" >+.Sh DESCRIPTION >+The >+.Fn cospi , >+.Fn cospif , >+and >+.Fn cospil >+functions compute the cosine of >+.Fa "\(*p \(mu x" . >+and measure angles in half-cycles. >+.Sh RETURN VALUES >+The >+.Fn cospi , >+.Fn cospif , >+and >+.Fn cospil >+functions returns >+.Fn cos "\(*p \(mu x" . >+If \*(Bax\*(Ba \*(Ge 2^(p - 1) >+where p is the floating\(enpoint precision of >+.Ar x , >+then the returned value is 1 and it has no significance. >+.Sh SPECIAL VALUES >+.Bl -tag >+.It >+.Fn cospi \*(Pm0 >+returns 1. >+.It >+.Fn cospi \*(Pmn/2 >+returns 0 for positive integers >+.Ar n . >+.It >+.Fn cospi n >+returns 1 for even integers >+.Ar n . >+.It >+.Fn cospi n >+returns \-1 for odd integers >+.Ar n . >+.It >+.Fn cospi \*(Pm\(if >+return an \*(Na and raises an FE_INVALID exception. >+.It >+.Fn cospi \*(Na >+return an \*(Na and raises an FE_INVALID exception. >+.El >+.Sh SEE ALSO >+.Xr cos 3 , >+.Xr fenv 3 , >+.Xr math 3 , >+.Xr sin 3 , >+.Xr sinpi 3 , >+.Xr tan 3 , >+.Xr tanpi 3 >+.Sh AUTHORS >+The half\(encycle trignometric functions were written by >+.An Steven G. Kargl Aq Mt kargl@FreeBSD.org . >+.Sh STANDARDS >+These functions conform to >+IEEE Std 754\(tm\(en2008 , >+\(dqIEEE Standard for Floating-Point Arithmetic\(dq >+and to >+ISO/IEC TS 18661-4 , >+\(dqInformation technology \(em Programming languages, their environments, >+and system software interfaces \(em Floating\(enpoint extensions for >+C\(dq \(em Part 4: Supplementary functions. >+ >+ >diff --git a/lib/msun/man/sinpi.3 b/lib/msun/man/sinpi.3 >new file mode 100644 >index 00000000000..4ae162bd1c7 >--- /dev/null >+++ b/lib/msun/man/sinpi.3 >@@ -0,0 +1,102 @@ >+.\" Copyright (c) 2017 Steven G. Kargl <kargl@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 REGENTS 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 REGENTS 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$ >+.\" >+.Dd April 1, 2017 >+.Dt SINPI 3 >+.Os >+.Sh NAME >+.Nm sinpi , >+.Nm sinpif , >+.Nm sinpil >+.Nd half\(encycle sine functions >+.Sh LIBRARY >+.Lb libm >+.Sh SYNOPSIS >+.In math.h >+.Ft double >+.Fn sinpi "double x" >+.Ft float >+.Fn sinpif "float x" >+.Ft long double >+.Fn sinpil "long double x" >+.Sh DESCRIPTION >+The >+.Fn sinpi , >+.Fn sinpif , >+and >+.Fn sinpil >+functions compute the sine of >+.Fa "\(*p \(mu x" . >+and measure angles in half-cycles. >+.Sh RETURN VALUES >+The >+.Fn sinpi , >+.Fn sinpif , >+and >+.Fn sinpil >+functions returns >+.Fn sin "\(*p \(mu x" . >+If \*(Bax\*(Ba \*(Ge 2^(p - 1) >+where p is the floating\(enpoint precision of >+.Ar x , >+then the returned value is \*(Pm0 and it has no significance. >+.Sh SPECIAL VALUES >+.Bl -tag >+.It >+.Fn sinpi \*(Pm0 >+returns \*(Pm0. >+.It >+.Fn sinpi \*(Pmn >+returns \*(Pm0 for positive integers >+.Ar n . >+.It >+.Fn sinpi \*(Pm\(if >+return an \*(Na and raises an FE_INVALID exception. >+.It >+.Fn sinpi \*(Na >+return an \*(Na and raises an FE_INVALID exception. >+.El >+.Sh SEE ALSO >+.Xr cos 3 , >+.Xr cospi 3 , >+.Xr fenv 3 , >+.Xr math 3 , >+.Xr sin 3 , >+.Xr tan 3 , >+.Xr tanpi 3 >+.Sh AUTHORS >+The half\(encycle trignometric functions were written by >+.An Steven G. Kargl Aq Mt kargl@FreeBSD.org . >+.Sh STANDARDS >+These functions conform to >+IEEE Std 754\(tm\(en2008 , >+\(dqIEEE Standard for Floating-Point Arithmetic\(dq >+and to >+ISO/IEC TS 18661-4 , >+\(dqInformation technology \(em Programming languages, their environments, >+and system software interfaces \(em Floating\(enpoint extensions for >+C\(dq \(em Part 4: Supplementary functions. >+ >diff --git a/lib/msun/man/tanpi.3 b/lib/msun/man/tanpi.3 >new file mode 100644 >index 00000000000..9cbbf4efaf2 >--- /dev/null >+++ b/lib/msun/man/tanpi.3 >@@ -0,0 +1,106 @@ >+.\" Copyright (c) 2017 Steven G. Kargl <kargl@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 REGENTS 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 REGENTS 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$ >+.\" >+.Dd April 1, 2017 >+.Dt TANPI 3 >+.Os >+.Sh NAME >+.Nm tanpi , >+.Nm tanpif , >+.Nm tanpil >+.Nd half\(encycle tangent functions >+.Sh LIBRARY >+.Lb libm >+.Sh SYNOPSIS >+.In math.h >+.Ft double >+.Fn tanpi "double x" >+.Ft float >+.Fn tanpif "float x" >+.Ft long double >+.Fn tanpil "long double x" >+.Sh DESCRIPTION >+The >+.Fn tanpi , >+.Fn tanpif , >+and >+.Fn tanpil >+functions compute the tangent of >+.Fa "\(*p \(mu x" >+and measure angles in half-cycles. >+.Sh RETURN VALUES >+The >+.Fn tanpi , >+.Fn tanpif , >+and >+.Fn tanpil >+functions returns >+.Fn tan "\(*p \(mu x" . >+If \*(Bax\*(Ba \*(Ge 2^(p - 1) >+where p is the floating\(enpoint precision of >+.Ar x , >+then the returned value is \*(Pm0 and it has no significance. >+.Sh SPECIAL VALUES >+.Bl -tag >+.It >+.Fn tanpi \*(Pm0 >+returns \*(Pm0. >+.It >+.Fn tanpi \*(Pmn >+returns \*(Pm0 for positive integers >+.Ar n . >+.It >+.Fn tanpi \*(Pmn/2 >+returns \*(Na for n > 0 and raises an FE_INVALID exception. >+.It >+.Fn tanpi \*(Pm\(if >+return an \*(Na and raises an FE_INVALID exception. >+.It >+.Fn tanpi \*(Na >+return an \*(Na and raises an FE_INVALID exception. >+.El >+.Sh SEE ALSO >+.Xr cos 3 , >+.Xr cospi 3 , >+.Xr fenv 3 , >+.Xr math 3 , >+.Xr sin 3 , >+.Xr sinpi 3 , >+.Xr tan 3 , >+.Sh AUTHORS >+The half\(encycle trignometric functions were written by >+.An Steven G. Kargl Aq Mt kargl@FreeBSD.org . >+.Sh STANDARDS >+These functions conform to >+IEEE Std 754\(tm\(en2008 , >+\(dqIEEE Standard for Floating-Point Arithmetic\(dq >+and to >+ISO/IEC TS 18661-4 , >+\(dqInformation technology \(em Programming languages, their environments, >+and system software interfaces \(em Floating\(enpoint extensions for >+C\(dq \(em Part 4: Supplementary functions. >+ >+ >diff --git a/lib/msun/src/k_cospi.h b/lib/msun/src/k_cospi.h >new file mode 100644 >index 00000000000..32985be809a >--- /dev/null >+++ b/lib/msun/src/k_cospi.h >@@ -0,0 +1,44 @@ >+/*- >+ * Copyright (c) 2017 Steven G. Kargl >+ * 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 unmodified, 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 ``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. >+ */ >+ >+/* >+ * The basic kernel for x in [0,0.25]. To use the kernel for cos(x), the >+ * argument to __kernel_cospi() must be multiplied by pi. >+ */ >+ >+static inline double >+__kernel_cospi(double x) >+{ >+ double_t hi, lo; >+ >+ hi = (float)x; >+ lo = x - hi; >+ lo = lo * (pi_lo + pi_hi) + hi * pi_lo; >+ hi *= pi_hi; >+ _2sumF(hi, lo); >+ return (__kernel_cos(hi, lo)); >+} >+ >diff --git a/lib/msun/src/k_sinpi.h b/lib/msun/src/k_sinpi.h >new file mode 100644 >index 00000000000..66152ba9931 >--- /dev/null >+++ b/lib/msun/src/k_sinpi.h >@@ -0,0 +1,43 @@ >+/*- >+ * Copyright (c) 2017 Steven G. Kargl >+ * 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 unmodified, 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 ``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. >+ */ >+ >+/* >+ * The basic kernel for x in [0,0.25]. To use the kernel for sin(x), the >+ * argument to __kernel_sinpi() must be multiplied by pi. >+ */ >+ >+static inline double >+__kernel_sinpi(double x) >+{ >+ double_t hi, lo; >+ >+ hi = (float)x; >+ lo = x - hi; >+ lo = lo * (pi_lo + pi_hi) + hi * pi_lo; >+ hi *= pi_hi; >+ _2sumF(hi, lo); >+ return (__kernel_sin(hi, lo, 1)); >+} >diff --git a/lib/msun/src/math.h b/lib/msun/src/math.h >index 8156094313e..98601f825bb 100644 >--- a/lib/msun/src/math.h >+++ b/lib/msun/src/math.h >@@ -508,6 +508,15 @@ long double lgammal_r(long double, int *); > void sincos(double, double *, double *); > void sincosf(float, float *, float *); > void sincosl(long double, long double *, long double *); >+double cospi(double); >+float cospif(float); >+long double cospil(long double); >+double sinpi(double); >+float sinpif(float); >+long double sinpil(long double); >+double tanpi(double); >+float tanpif(float); >+long double tanpil(long double); > #endif > > __END_DECLS >diff --git a/lib/msun/src/s_cospi.c b/lib/msun/src/s_cospi.c >new file mode 100644 >index 00000000000..92a5f467c97 >--- /dev/null >+++ b/lib/msun/src/s_cospi.c >@@ -0,0 +1,151 @@ >+/*- >+ * Copyright (c) 2017 Steven G. Kargl >+ * 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 unmodified, 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 ``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. >+ */ >+ >+/** >+ * cospi(x) computes cos(pi*x) without multiplication by pi (almost). First, >+ * note that cospi(-x) = cospi(x), so the algorithm considers only |x|. The >+ * method used depends on the magnitude of x. >+ * >+ * 1. For small |x|, cospi(x) = 1 with FE_INEXACT raised where a sloppy >+ * threshold is used. The threshold is |x| < 0x1pN with N = -(P/2+M). >+ * P is the precision of the floating-point type and M = 2 to 4. >+ * >+ * 2. For |x| < 1, argument reduction is not required and sinpi(x) is >+ * computed by calling a kernel that leverages the kernels for sin(x) >+ * ans cos(x). See k_sinpi.c and k_cospi.c for details. >+ * >+ * 3. For 1 <= |x| < 0x1p(P-1), argument reduction is required where >+ * |x| = j0 + r with j0 an integer and the remainder r satisfies >+ * 0 <= r < 1. With the given domain, a simplified inline floor(x) >+ * is used. Also, note the following identity >+ * >+ * cospi(x) = cos(pi*(j0+r)) >+ * = cos(pi*j0) * cos(pi*r) - sin(pi*j0) * sin(pi*r) >+ * = cos(pi*j0) * cos(pi*r) >+ * = +-cospi(r) >+ * >+ * If j0 is even, then cos(pi*j0) = 1. If j0 is odd, then cos(pi*j0) = -1. >+ * cospi(r) is then computed via an appropriate kernel. >+ * >+ * 4. For |x| >= 0x1p(P-1), |x| is integral and cospi(x) = 1. >+ * >+ * 5. Special cases: >+ * >+ * cospi(+-0) = 1. >+ * cospi(n.5) = 0 for n an integer. >+ * cospi(+-inf) = nan. Raises the "invalid" floating-point exception. >+ * cospi(nan) = nan. Raises the "invalid" floating-point exception. >+ */ >+ >+#include "math.h" >+#include "math_private.h" >+ >+static const double >+pi_hi = 3.1415926814079285e+00, /* 0x400921fb 0x58000000 */ >+pi_lo =-2.7818135228334233e-08; /* 0xbe5dde97 0x3dcb3b3a */ >+ >+#include "k_cospi.h" >+#include "k_sinpi.h" >+ >+volatile static const double vzero = 0; >+ >+double >+cospi(double x) >+{ >+ double ax, c; >+ uint32_t hx, ix, j0, lx; >+ >+ EXTRACT_WORDS(hx, lx, x); >+ ix = hx & 0x7fffffff; >+ INSERT_WORDS(ax, ix, lx); >+ >+ if (ix < 0x3ff00000) { /* |x| < 1 */ >+ if (ix < 0x3fd00000) { /* |x| < 0.25 */ >+ if (ix < 0x3e200000) { /* |x| < 0x1p-29 */ >+ if ((int)ax == 0) >+ return (1); >+ } >+ return (__kernel_cospi(ax)); >+ } >+ >+ if (ix < 0x3fe00000) /* |x| < 0.5 */ >+ c = __kernel_sinpi(0.5 - ax); >+ else if (ix < 0x3fe80000){ /* |x| < 0.75 */ >+ if (ax == 0.5) >+ return (0); >+ c = -__kernel_sinpi(ax - 0.5); >+ } else >+ c = -__kernel_cospi(1 - ax); >+ return (c); >+ } >+ >+ if (ix < 0x43300000) { /* 1 <= |x| < 0x1p52 */ >+ /* Determine integer part of ax. */ >+ j0 = ((ix >> 20) & 0x7ff) - 0x3ff; >+ if (j0 < 20) { >+ ix &= ~(0x000fffff >> j0); >+ lx = 0; >+ } else { >+ lx &= ~((uint32_t)0xffffffff >> (j0 - 20)); >+ } >+ INSERT_WORDS(x, ix, lx); >+ >+ ax -= x; >+ EXTRACT_WORDS(ix, lx, ax); >+ >+ >+ if (ix < 0x3fe00000) { /* |x| < 0.5 */ >+ if (ix < 0x3fd00000) /* |x| < 0.25 */ >+ c = ix == 0 ? 1 : __kernel_cospi(ax); >+ else >+ c = __kernel_sinpi(0.5 - ax); >+ } else { >+ if (ix < 0x3fe80000) { /* |x| < 0.75 */ >+ if (ax == 0.5) >+ return (0); >+ c = -__kernel_sinpi(ax - 0.5); >+ } else >+ c = -__kernel_cospi(1 - ax); >+ } >+ >+ if (j0 > 30) >+ x -= 0x1p30; >+ j0 = (uint32_t)x; >+ return (j0 & 1 ? -c : c); >+ } >+ >+ if (ix >= 0x7f800000) >+ return (vzero / vzero); >+ >+ /* >+ * |x| >= 0x1p52 is always an even integer, so return 1. >+ */ >+ return (1); >+} >+ >+#if LDBL_MANT_DIG == 53 >+__weak_reference(cospi, cospil); >+#endif >diff --git a/lib/msun/src/s_cospif.c b/lib/msun/src/s_cospif.c >new file mode 100644 >index 00000000000..4dd881395ba >--- /dev/null >+++ b/lib/msun/src/s_cospif.c >@@ -0,0 +1,109 @@ >+/*- >+ * Copyright (c) 2017 Steven G. Kargl >+ * 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 unmodified, 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 ``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. >+ */ >+ >+/* >+ * See ../src/s_cospi.c for implementation details. >+ */ >+#define INLINE_KERNEL_SINDF >+#define INLINE_KERNEL_COSDF >+ >+#include "math.h" >+#include "math_private.h" >+#include "k_cosf.c" >+#include "k_sinf.c" >+ >+#define __kernel_cospif(x) (__kernel_cosdf(M_PI * (x))) >+#define __kernel_sinpif(x) (__kernel_sindf(M_PI * (x))) >+ >+volatile static const float vzero = 0; >+ >+float >+cospif(float x) >+{ >+ float ax, c; >+ uint32_t ix, j0; >+ >+ GET_FLOAT_WORD(ix, x); >+ ix = ix & 0x7fffffff; >+ SET_FLOAT_WORD(ax, ix); >+ >+ if (ix < 0x3f800000) { /* |x| < 1 */ >+ if (ix < 0x3e800000) { /* |x| < 0.25 */ >+ if (ix < 0x38800000) { /* |x| < 0x1p-14 */ >+ /* Raise inexact iff != 0. */ >+ if ((int)ax == 0) >+ return (1); >+ } >+ return (__kernel_cospif(ax)); >+ } >+ >+ if (ix < 0x3f000000) /* |x| < 0.5 */ >+ c = __kernel_sinpif(0.5F - ax); >+ else if (ix < 0x3f400000) { /* |x| < 0.75 */ >+ if (ix == 0x3f000000) >+ return (0); >+ c = -__kernel_sinpif(ax - 0.5F); >+ } else >+ c = -__kernel_cospif(1 - ax); >+ return (c); >+ } >+ >+ if (ix < 0x4b000000) { /* 1 <= |x| < 0x1p23 */ >+ /* Determine integer part of ax. */ >+ j0 = ((ix >> 23) & 0xff) - 0x7f; >+ ix &= ~(0x007fffff >> j0); >+ SET_FLOAT_WORD(x, ix); >+ >+ ax -= x; >+ GET_FLOAT_WORD(ix, ax); >+ >+ if (ix < 0x3f000000) { /* |x| < 0.5 */ >+ if (ix < 0x3e800000) /* |x| < 0.25 */ >+ c = ix == 0 ? 1 : __kernel_cospif(ax); >+ else >+ c = __kernel_sinpif(0.5F - ax); >+ } else { >+ if (ix < 0x3f400000) { /* |x| < 0.75 */ >+ if (ix == 0x3f000000) >+ return (0); >+ c = -__kernel_sinpif(ax - 0.5F); >+ } else >+ c = -__kernel_cospif(1 - ax); >+ } >+ >+ j0 = (uint32_t)x; >+ return (j0 & 1 ? -c : c); >+ } >+ >+ /* x = +-inf or nan. */ >+ if (ix >= 0x7f800000) >+ return (vzero / vzero); >+ >+ /* >+ * |x| >= 0x1p23 is always an even integer, so return 1. >+ */ >+ return (1); >+} >diff --git a/lib/msun/src/s_sinpi.c b/lib/msun/src/s_sinpi.c >new file mode 100644 >index 00000000000..b9731112a7e >--- /dev/null >+++ b/lib/msun/src/s_sinpi.c >@@ -0,0 +1,168 @@ >+/*- >+ * Copyright (c) 2017 Steven G. Kargl >+ * 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 unmodified, 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 ``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. >+ */ >+ >+/** >+ * sinpi(x) computes sin(pi*x) without multiplication by pi (almost). First, >+ * note that sinpi(-x) = -sinpi(x), so the algorithm considers only |x| and >+ * includes reflection symmetry by considering the sign of x on output. The >+ * method used depends on the magnitude of x. >+ * >+ * 1. For small |x|, sinpi(x) = pi * x where a sloppy threshold is used. The >+ * threshold is |x| < 0x1pN with N = -(P/2+M). P is the precision of the >+ * floating-point type and M = 2 to 4. To achieve high accuracy, pi is >+ * decomposed into high and low parts with the high part containing a >+ * number of trailing zero bits. x is also split into high and low parts. >+ * >+ * 2. For |x| < 1, argument reduction is not required and sinpi(x) is >+ * computed by calling a kernel that leverages the kernels for sin(x) >+ * ans cos(x). See k_sinpi.c and k_cospi.c for details. >+ * >+ * 3. For 1 <= |x| < 0x1p(P-1), argument reduction is required where >+ * |x| = j0 + r with j0 an integer and the remainder r satisfies >+ * 0 <= r < 1. With the given domain, a simplified inline floor(x) >+ * is used. Also, note the following identity >+ * >+ * sinpi(x) = sin(pi*(j0+r)) >+ * = sin(pi*j0) * cos(pi*r) + cos(pi*j0) * sin(pi*r) >+ * = cos(pi*j0) * sin(pi*r) >+ * = +-sinpi(r) >+ * >+ * If j0 is even, then cos(pi*j0) = 1. If j0 is odd, then cos(pi*j0) = -1. >+ * sinpi(r) is then computed via an appropriate kernel. >+ * >+ * 4. For |x| >= 0x1p(P-1), |x| is integral and sinpi(x) = copysign(0,x). >+ * >+ * 5. Special cases: >+ * >+ * sinpi(+-0) = +-0 >+ * sinpi(+-n) = +-0, for positive integers n. >+ * sinpi(+-inf) = nan. Raises the "invalid" floating-point exception. >+ * sinpi(nan) = nan. Raises the "invalid" floating-point exception. >+ */ >+ >+#include "math.h" >+#include "math_private.h" >+ >+static const double >+pi_hi = 3.1415926814079285e+00, /* 0x400921fb 0x58000000 */ >+pi_lo =-2.7818135228334233e-08; /* 0xbe5dde97 0x3dcb3b3a */ >+ >+#include "k_cospi.h" >+#include "k_sinpi.h" >+ >+volatile static const double vzero = 0; >+ >+double >+sinpi(double x) >+{ >+ double ax, hi, lo, s; >+ uint32_t hx, ix, j0, lx; >+ >+ EXTRACT_WORDS(hx, lx, x); >+ ix = hx & 0x7fffffff; >+ INSERT_WORDS(ax, ix, lx); >+ >+ if (ix < 0x3ff00000) { /* |x| < 1 */ >+ if (ix < 0x3fd00000) { /* |x| < 0.25 */ >+ if (ix < 0x3e200000) { /* |x| < 0x1p-29 */ >+ if (x == 0) >+ return (x); >+ /* >+ * To avoid issues with subnormal values, >+ * scale the computation and rescale on >+ * return. >+ */ >+ INSERT_WORDS(hi, hx, 0); >+ hi *= 0x1p53; >+ lo = x * 0x1p53 - hi; >+ s = (pi_lo + pi_hi) * lo + pi_lo * hi + >+ pi_hi * hi; >+ return (s * 0x1p-53); >+ } >+ >+ s = __kernel_sinpi(ax); >+ return ((hx & 0x80000000) ? -s : s); >+ } >+ >+ if (ix < 0x3fe00000) /* |x| < 0.5 */ >+ s = __kernel_cospi(0.5 - ax); >+ else if (ix < 0x3fe80000) /* |x| < 0.75 */ >+ s = __kernel_cospi(ax - 0.5); >+ else >+ s = __kernel_sinpi(1 - ax); >+ return ((hx & 0x80000000) ? -s : s); >+ } >+ >+ if (ix < 0x43300000) { /* 1 <= |x| < 0x1p52 */ >+ /* Determine integer part of ax. */ >+ j0 = ((ix >> 20) & 0x7ff) - 0x3ff; >+ if (j0 < 20) { >+ ix &= ~(0x000fffff >> j0); >+ lx = 0; >+ } else { >+ lx &= ~((uint32_t)0xffffffff >> (j0 - 20)); >+ } >+ INSERT_WORDS(x, ix, lx); >+ >+ ax -= x; >+ EXTRACT_WORDS(ix, lx, ax); >+ >+ if (ix == 0) >+ s = 0; >+ else { >+ if (ix < 0x3fe00000) { /* |x| < 0.5 */ >+ if (ix < 0x3fd00000) /* |x| < 0.25 */ >+ s = __kernel_sinpi(ax); >+ else >+ s = __kernel_cospi(0.5 - ax); >+ } else { >+ if (ix < 0x3fe80000) /* |x| < 0.75 */ >+ s = __kernel_cospi(ax - 0.5); >+ else >+ s = __kernel_sinpi(1 - ax); >+ } >+ >+ if (j0 > 30) >+ x -= 0x1p30; >+ j0 = (uint32_t)x; >+ if (j0 & 1) s = -s; >+ } >+ >+ return ((hx & 0x80000000) ? -s : s); >+ } >+ >+ if (ix >= 0x7f800000) >+ return (vzero / vzero); >+ >+ /* >+ * |x| >= 0x1p52 is always an integer, so return +-0. >+ */ >+ return (copysign(0, x)); >+} >+ >+#if LDBL_MANT_DIG == 53 >+__weak_reference(sinpi, sinpil); >+#endif >diff --git a/lib/msun/src/s_sinpif.c b/lib/msun/src/s_sinpif.c >new file mode 100644 >index 00000000000..c9f76f8a235 >--- /dev/null >+++ b/lib/msun/src/s_sinpif.c >@@ -0,0 +1,122 @@ >+/*- >+ * Copyright (c) 2017 Steven G. Kargl >+ * 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 unmodified, 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 ``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. >+ */ >+ >+/* >+ * See ../src/s_sinpi.c for implementation details. >+ */ >+ >+#define INLINE_KERNEL_SINDF >+#define INLINE_KERNEL_COSDF >+ >+#include "math.h" >+#include "math_private.h" >+#include "k_cosf.c" >+#include "k_sinf.c" >+ >+#define __kernel_cospif(x) (__kernel_cosdf(M_PI * (x))) >+#define __kernel_sinpif(x) (__kernel_sindf(M_PI * (x))) >+ >+static const float >+pi_hi = 3.14160156e+00F, /* 0x40491000 */ >+pi_lo = -8.90890988e-06F; /* 0xb715777a */ >+ >+volatile static const float vzero = 0; >+ >+float >+sinpif(float x) >+{ >+ float ax, hi, lo, s; >+ uint32_t hx, ix, j0; >+ >+ GET_FLOAT_WORD(hx, x); >+ ix = hx & 0x7fffffff; >+ SET_FLOAT_WORD(ax, ix); >+ >+ if (ix < 0x3f800000) { /* |x| < 1 */ >+ if (ix < 0x3e800000) { /* |x| < 0.25 */ >+ if (ix < 0x38800000) { /* |x| < 0x1p-14 */ >+ if (x == 0) >+ return (x); >+ SET_FLOAT_WORD(hi, hx & 0xffff0000); >+ hi *= 0x1p23F; >+ lo = x * 0x1p23F - hi; >+ s = (pi_lo + pi_hi) * lo + pi_lo * hi + >+ pi_hi * hi; >+ return (s * 0x1p-23F); >+ } >+ >+ s = __kernel_sinpif(ax); >+ return ((hx & 0x80000000) ? -s : s); >+ } >+ >+ if (ix < 0x3f000000) /* |x| < 0.5 */ >+ s = __kernel_cospif(0.5F - ax); >+ else if (ix < 0x3f400000) /* |x| < 0.75 */ >+ s = __kernel_cospif(ax - 0.5F); >+ else >+ s = __kernel_sinpif(1 - ax); >+ return ((hx & 0x80000000) ? -s : s); >+ } >+ >+ if (ix < 0x4b000000) { /* 1 <= |x| < 0x1p23 */ >+ /* Determine integer part of ax. */ >+ j0 = ((ix >> 23) & 0xff) - 0x7f; >+ ix &= ~(0x007fffff >> j0); >+ SET_FLOAT_WORD(x, ix); >+ >+ ax -= x; >+ GET_FLOAT_WORD(ix, ax); >+ >+ if (ix == 0) >+ s = 0; >+ else { >+ if (ix < 0x3f000000) { /* |x| < 0.5 */ >+ if (ix < 0x3e800000) /* |x| < 0.25 */ >+ s = __kernel_sinpif(ax); >+ else >+ s = __kernel_cospif(0.5F - ax); >+ } else { >+ if (ix < 0x3f400000) /* |x| < 0.75 */ >+ s = __kernel_cospif(ax - 0.5F); >+ else >+ s = __kernel_sinpif(1 - ax); >+ } >+ >+ j0 = (uint32_t)x; >+ s = (j0 & 1) ? -s : s; >+ } >+ return ((hx & 0x80000000) ? -s : s); >+ } >+ >+ /* x = +-inf or nan. */ >+ if (ix >= 0x7f800000) >+ return (vzero / vzero); >+ >+ /* >+ * |x| >= 0x1p23 is always an integer, so return +-0. >+ */ >+ return (copysignf(0, x)); >+} >diff --git a/lib/msun/src/s_tanpi.c b/lib/msun/src/s_tanpi.c >new file mode 100644 >index 00000000000..e01917c94c1 >--- /dev/null >+++ b/lib/msun/src/s_tanpi.c >@@ -0,0 +1,176 @@ >+/*- >+ * Copyright (c) 2017 Steven G. Kargl >+ * 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 unmodified, 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 ``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. >+ */ >+ >+/** >+ * tanpi(x) computes tan(pi*x) without multiplication by pi (almost). First, >+ * note that tanpi(-x) = -tanpi(x), so the algorithm considers only |x| and >+ * includes reflection symmetry by considering the sign of x on output. The >+ * method used depends on the magnitude of x. >+ * >+ * 1. For small |x|, tanpi(x) = pi * x where a sloppy threshold is used. The >+ * threshold is |x| < 0x1pN with N = -(P/2+M). P is the precision of the >+ * floating-point type and M = 2 to 4. To achieve high accuracy, pi is >+ * decomposed into high and low parts with the high part containing a >+ * number of trailing zero bits. x is also split into high and low parts. >+ * >+ * 2. For |x| < 1, argument reduction is not required and tanpi(x) is >+ * computed by a direct call to a kernel, which uses the kernel for >+ * tan(x). See below. >+ * >+ * 3. For 1 <= |x| < 0x1p(P-1), argument reduction is required where >+ * |x| = j0 + r with j0 an integer and the remainder r satisfies >+ * 0 <= r < 1. With the given domain, a simplified inline floor(x) >+ * is used. Also, note the following identity >+ * >+ * tan(pi*j0) + tan(pi*r) >+ * tanpi(x) = tan(pi*(j0+r)) = ---------------------------- = tanpi(r) >+ * 1 - tan(pi*j0) * tan(pi*r) >+ * >+ * So, after argument reduction, the kernel is again invoked. >+ * >+ * 4. For |x| >= 0x1p(P-1), |x| is integral and tanpi(x) = copysign(0,x). >+ * >+ * 5. Special cases: >+ * >+ * tanpi(+-0) = +-0 >+ * tanpi(+-n) = +-0, for positive integers n. >+ * tanpi(+-n+1/4) = +-1, for positive integers n. >+ * tanpi(+-n+1/2) = NaN, for positive integers n. >+ * tanpi(+-inf) = NaN. Raises the "invalid" floating-point exception. >+ * tanpi(nan) = NaN. Raises the "invalid" floating-point exception. >+ */ >+ >+#include "math.h" >+#include "math_private.h" >+ >+static const double >+pi_hi = 3.1415926814079285e+00, /* 0x400921fb 0x58000000 */ >+pi_lo = -2.7818135228334233e-08; /* 0xbe5dde97 0x3dcb3b3a */ >+ >+/* >+ * The kernel for tanpi(x) multiplies x by an 80-bit approximation of >+ * pi, where the hi and lo parts are used with with kernel for tan(x). >+ */ >+static inline double >+__kernel_tanpi(double x) >+{ >+ double_t hi, lo, t; >+ >+ if (x < 0.25) { >+ hi = (float)x; >+ lo = x - hi; >+ lo = lo * (pi_lo + pi_hi) + hi * pi_lo; >+ hi *= pi_hi; >+ _2sumF(hi, lo); >+ t = __kernel_tan(hi, lo, 1); >+ } else if (x > 0.25) { >+ x = 0.5 - x; >+ hi = (float)x; >+ lo = x - hi; >+ lo = lo * (pi_lo + pi_hi) + hi * pi_lo; >+ hi *= pi_hi; >+ _2sumF(hi, lo); >+ t = - __kernel_tan(hi, lo, -1); >+ } else >+ t = 1; >+ >+ return (t); >+} >+ >+volatile static const double vzero = 0; >+ >+double >+tanpi(double x) >+{ >+ double ax, hi, lo, t; >+ uint32_t hx, ix, j0, lx; >+ >+ EXTRACT_WORDS(hx, lx, x); >+ ix = hx & 0x7fffffff; >+ INSERT_WORDS(ax, ix, lx); >+ >+ if (ix < 0x3ff00000) { /* |x| < 1 */ >+ if (ix < 0x3fe00000) { /* |x| < 0.5 */ >+ if (ix < 0x3e200000) { /* |x| < 0x1p-29 */ >+ if (x == 0) >+ return (x); >+ /* >+ * To avoid issues with subnormal values, >+ * scale the computation and rescale on >+ * return. >+ */ >+ INSERT_WORDS(hi, hx, 0); >+ hi *= 0x1p53; >+ lo = x * 0x1p53 - hi; >+ t = (pi_lo + pi_hi) * lo + pi_lo * hi + >+ pi_hi * hi; >+ return (t * 0x1p-53); >+ } >+ t = __kernel_tanpi(ax); >+ } else if (ax == 0.5) >+ return ((ax - ax) / (ax - ax)); >+ else >+ t = - __kernel_tanpi(1 - ax); >+ return ((hx & 0x80000000) ? -t : t); >+ } >+ >+ if (ix < 0x43300000) { /* 1 <= |x| < 0x1p52 */ >+ /* Determine integer part of ax. */ >+ j0 = ((ix >> 20) & 0x7ff) - 0x3ff; >+ if (j0 < 20) { >+ ix &= ~(0x000fffff >> j0); >+ lx = 0; >+ } else { >+ lx &= ~(((uint32_t)(0xffffffff)) >> (j0 - 20)); >+ } >+ INSERT_WORDS(x,ix,lx); >+ >+ ax -= x; >+ EXTRACT_WORDS(ix, lx, ax); >+ >+ if (ix < 0x3fe00000) /* |x| < 0.5 */ >+ t = ax == 0 ? 0 : __kernel_tanpi(ax); >+ else if (ax == 0.5) >+ return ((ax - ax) / (ax - ax)); >+ else >+ t = - __kernel_tanpi(1 - ax); >+ >+ return ((hx & 0x80000000) ? -t : t); >+ } >+ >+ /* x = +-inf or nan. */ >+ if (ix >= 0x7f800000) >+ return (vzero / vzero); >+ >+ /* >+ * |x| >= 0x1p52 is always an integer, so return +-0. >+ */ >+ return (copysign(0, x)); >+} >+ >+#if LDBL_MANT_DIG == 53 >+__weak_reference(tanpi, tanpil); >+#endif >diff --git a/lib/msun/src/s_tanpif.c b/lib/msun/src/s_tanpif.c >new file mode 100644 >index 00000000000..6d4b627d1cf >--- /dev/null >+++ b/lib/msun/src/s_tanpif.c >@@ -0,0 +1,114 @@ >+/*- >+ * Copyright (c) 2017 Steven G. Kargl >+ * 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 unmodified, 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 ``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. >+ */ >+ >+/* >+ * See ../src/s_tanpi.c for implementation details. >+ */ >+ >+#define INLINE_KERNEL_TANDF >+ >+#include "math.h" >+#include "math_private.h" >+#include "k_tanf.c" >+ >+static const float >+pi_hi = 3.14160156e+00F, /* 0x40491000 */ >+pi_lo = -8.90890988e-06F; /* 0xb715777a */ >+ >+static inline float >+__kernel_tanpif(float x) >+{ >+ float t; >+ >+ if (x < 0.25F) >+ t = __kernel_tandf(M_PI * x, 1); >+ else if (x > 0.25F) >+ t = -__kernel_tandf(M_PI * (0.5 - x), -1); >+ else >+ t = 1; >+ >+ return (t); >+} >+ >+volatile static const float vzero = 0; >+ >+float >+tanpif(float x) >+{ >+ float ax, hi, lo, t; >+ uint32_t hx, ix, j0; >+ >+ GET_FLOAT_WORD(hx, x); >+ ix = hx & 0x7fffffff; >+ SET_FLOAT_WORD(ax, ix); >+ >+ if (ix < 0x3f800000) { /* |x| < 1 */ >+ if (ix < 0x3f000000) { /* |x| < 0.5 */ >+ if (ix < 0x38800000) { /* |x| < 0x1p-14 */ >+ if (ix == 0) >+ return (x); >+ SET_FLOAT_WORD(hi, hx & 0xffff0000); >+ hi *= 0x1p23F; >+ lo = x * 0x1p23F - hi; >+ t = (pi_lo + pi_hi) * lo + pi_lo * hi + >+ pi_hi * hi; >+ return (t * 0x1p-23F); >+ } >+ t = __kernel_tanpif(ax); >+ } else if (ix == 0x3f000000) >+ return ((ax - ax) / (ax - ax)); >+ else >+ t = - __kernel_tanpif(1 - ax); >+ return ((hx & 0x80000000) ? -t : t); >+ } >+ >+ if (ix < 0x4b000000) { /* 1 <= |x| < 0x1p23 */ >+ /* Determine integer part of ax. */ >+ j0 = ((ix >> 23) & 0xff) - 0x7f; >+ ix &= ~(0x007fffff >> j0); >+ SET_FLOAT_WORD(x, ix); >+ >+ ax -= x; >+ GET_FLOAT_WORD(ix, ax); >+ >+ if (ix < 0x3f000000) /* |x| < 0.5 */ >+ t = ix == 0 ? 0 : __kernel_tanpif(ax); >+ else if (ix == 0x3f000000) >+ return ((ax - ax) / (ax - ax)); >+ else >+ t = - __kernel_tanpif(1 - ax); >+ return ((hx & 0x80000000) ? -t : t); >+ } >+ >+ /* x = +-inf or nan. */ >+ if (ix >= 0x7f800000) >+ return (vzero / vzero); >+ >+ /* >+ * |x| >= 0x1p23 is always an integer, so return +-0. >+ */ >+ return (copysignf(0, x)); >+}
You cannot view the attachment while viewing its details because your browser does not support IFRAMEs.
View the attachment on a separate page
.
View Attachment As Diff
View Attachment As Raw
Actions:
View
|
Diff
Attachments on
bug 218514
:
181622
|
182139
|
182648
|
229008
| 229010 |
229037
|
229047
|
229051
|
229108
|
229180
|
229215