This is the mail archive of the
libc-alpha@sourceware.org
mailing list for the glibc project.
sinl/cosl for ldbl-96
- From: "Joseph S. Myers" <joseph at codesourcery dot com>
- To: libc-alpha at sourceware dot org
- Date: Fri, 16 Mar 2012 00:52:28 +0000 (UTC)
- Subject: sinl/cosl for ldbl-96
This patch fixes the problems with range reduction for sinl / cosl /
sincosl (ldbl-96, bug 13851) by providing a C implementation and using
it instead of the x87 instructions.
In general the C implementation is based on Jakub's ldbl-128
implementation, and I did not try to optimize the polynomials used for
the lower precision of ldbl-96 (it's probably possible to use fewer
terms without any significant effect on accuracy); I used plain
floating point in ldbl-96 in place of various bit-manipulations of the
floating-point representation done in ldbl-128. I did regenerate the
table of precomputed values (t_sincosl.c) with the attached program
(having verified that a variant of it for 113-bit precision appeared
to replicate the values in ldbl-128/t_sincosl.c).
It turned out that the general range reduction code, not previously
used for x86, needed a __FLT_EVAL_METHOD__ != 0 fix in the case prec
== 2 (the value to be used for 64-bit mantissa) similar to that
already present for prec == 3 (113-bit mantissa).
In the course of preparing this I found what appears to be a bug in
the ldbl-128 and ldbl-128ibm implementations of k_sinl.c, for which
I'll send a patch separately.
Tested x86 and x86_64 and ULPs updated based on that testing.
2012-03-16 Joseph Myers <joseph@codesourcery.com>
[BZ #13851]
* sysdeps/ieee754/ldbl-96/e_rem_pio2l.c: New file.
* sysdeps/ieee754/ldbl-96/k_cosl.c: Likewise.
* sysdeps/ieee754/ldbl-96/k_sinl.c: Likewise.
* sysdeps/ieee754/ldbl-96/t_sincosl.c: Likewise.
* sysdeps/ieee754/ldbl-96/s_cosl.c (__cosl): Correct test for
infinite argument.
* sysdeps/ieee754/ldbl-96/s_sinl.c (__sinl): Likewise.
* sysdeps/ieee754/dbl-64/k_rem_pio2.c: Handle __FLT_EVAL_METHOD__
!= 0 for prec == 2.
* sysdeps/i386/fpu/e_rem_pio2l.c: Remove.
* sysdeps/i386/fpu/k_rem_pio2.c: Likewise.
* sysdeps/i386/fpu/s_cosl.S: Likewise.
* sysdeps/i386/fpu/s_sincosl.S: Likewise.
* sysdeps/i386/fpu/s_sinl.S: Likewise.
* sysdeps/x86_64/fpu/e_rem_pio2l.c: Likewise.
* sysdeps/x86_64/fpu/k_cosl.c: Likewise.
* sysdeps/x86_64/fpu/k_sinl.c: Likewise.
* sysdeps/x86_64/fpu/s_cosl.S: Likewise.
* sysdeps/x86_64/fpu/s_sincosl.S: Likewise.
* sysdeps/x86_64/fpu/s_sinl.S: Likewise.
* math/libm-test.inc (cos_test): Add more tests and enable more
tests for long double.
(sin_test): Likewise.
(sincos_test): Likewise.
* sysdeps/i386/fpu/libm-test-ulps: Update.
* sysdeps/x86_64/fpu/libm-test-ulps: Likewise.
diff --git a/math/libm-test.inc b/math/libm-test.inc
index 89d0eb1..fb82926 100644
--- a/math/libm-test.inc
+++ b/math/libm-test.inc
@@ -2112,16 +2112,20 @@ cos_test (void)
TEST_f_f (cos, 0.75L, 0.731688868873820886311838753000084544L);
-#ifndef TEST_LDOUBLE
- /* Enable for long double once x86 and x86-64 implementations are fixed. */
TEST_f_f (cos, 0x1p65, 0.99888622066058013610642172179340364209972L);
TEST_f_f (cos, -0x1p65, 0.99888622066058013610642172179340364209972L);
-#endif
#ifdef TEST_DOUBLE
TEST_f_f (cos, 0.80190127184058835, 0.69534156199418473);
- TEST_f_f (cos, 1e22, 0.5232147853951389454975944733847);
- TEST_f_f (cos, 0x1p1023, -0.8263698346141479945007856808117);
+#endif
+
+#ifndef TEST_FLOAT
+ TEST_f_f (cos, 1e22, 0.5232147853951389454975944733847094921409L);
+ TEST_f_f (cos, 0x1p1023, -0.826369834614147994500785680811743734805L);
+#endif
+
+#if defined TEST_LDOUBLE && LDBL_MAX_EXP >= 16384
+ TEST_f_f (cos, 0x1p16383L, 0.9210843909921906206874509522505756251609L);
#endif
END (cos);
@@ -6395,18 +6399,21 @@ sin_test (void)
TEST_f_f (sin, -M_PI_2l, -1);
TEST_f_f (sin, 0.75L, 0.681638760023334166733241952779893935L);
-#ifndef TEST_LDOUBLE
-
- /* Enable for long double once x86 and x86-64 implementations are fixed. */
TEST_f_f (sin, 0x1p65, -0.047183876212354673805106149805700013943218L);
TEST_f_f (sin, -0x1p65, 0.047183876212354673805106149805700013943218L);
-#endif
#ifdef TEST_DOUBLE
TEST_f_f (sin, 0.80190127184058835, 0.71867942238767868);
TEST_f_f (sin, 2.522464e-1, 2.4957989804940911e-1);
- TEST_f_f (sin, 1e22, -0.8522008497671888017727058937530);
- TEST_f_f (sin, 0x1p1023, 0.5631277798508840248814522055909);
+#endif
+
+#ifndef TEST_FLOAT
+ TEST_f_f (sin, 1e22, -0.8522008497671888017727058937530293682618L);
+ TEST_f_f (sin, 0x1p1023, 0.5631277798508840134529434079444683477104L);
+#endif
+
+#if defined TEST_LDOUBLE && LDBL_MAX_EXP >= 16384
+ TEST_f_f (sin, 0x1p16383L, 0.3893629985894208126948115852610595405563L);
#endif
END (sin);
@@ -6576,16 +6583,20 @@ sincos_test (void)
TEST_extra (sincos, M_PI_6l*2.0, 0.86602540378443864676372317075293616L, 0.5);
TEST_extra (sincos, 0.75L, 0.681638760023334166733241952779893935L, 0.731688868873820886311838753000084544L);
-#ifndef TEST_LDOUBLE
- /* Enable for long double once x86 and x86-64 implementations are fixed. */
TEST_extra (sincos, 0x1p65, -0.047183876212354673805106149805700013943218L, 0.99888622066058013610642172179340364209972L);
TEST_extra (sincos, -0x1p65, 0.047183876212354673805106149805700013943218L, 0.99888622066058013610642172179340364209972L);
-#endif
#ifdef TEST_DOUBLE
TEST_extra (sincos, 0.80190127184058835, 0.71867942238767868, 0.69534156199418473);
- TEST_extra (sincos, 1e22, -0.8522008497671888017727058937530, 0.5232147853951389454975944733847);
- TEST_extra (sincos, 0x1p1023, 0.5631277798508840248814522055909, -0.8263698346141479945007856808117);
+#endif
+
+#ifndef TEST_FLOAT
+ TEST_extra (sincos, 1e22, -0.8522008497671888017727058937530293682618L, 0.5232147853951389454975944733847094921409L);
+ TEST_extra (sincos, 0x1p1023, 0.5631277798508840134529434079444683477104L, -0.826369834614147994500785680811743734805L);
+#endif
+
+#if defined TEST_LDOUBLE && LDBL_MAX_EXP >= 16384
+ TEST_extra (sincos, 0x1p16383L, 0.3893629985894208126948115852610595405563L, 0.9210843909921906206874509522505756251609L);
#endif
END (sincos);
diff --git a/sysdeps/i386/fpu/e_rem_pio2l.c b/sysdeps/i386/fpu/e_rem_pio2l.c
deleted file mode 100644
index 1347b04..0000000
--- a/sysdeps/i386/fpu/e_rem_pio2l.c
+++ /dev/null
@@ -1,3 +0,0 @@
-/* Empty. This file is only meant to avoid compiling the file with the
- same name in the libm-ieee754 directory. The code is not used since
- there is an assembler version for all users of this file. */
diff --git a/sysdeps/i386/fpu/k_rem_pio2.c b/sysdeps/i386/fpu/k_rem_pio2.c
deleted file mode 100644
index 1347b04..0000000
--- a/sysdeps/i386/fpu/k_rem_pio2.c
+++ /dev/null
@@ -1,3 +0,0 @@
-/* Empty. This file is only meant to avoid compiling the file with the
- same name in the libm-ieee754 directory. The code is not used since
- there is an assembler version for all users of this file. */
diff --git a/sysdeps/i386/fpu/libm-test-ulps b/sysdeps/i386/fpu/libm-test-ulps
index 76d25d7..ac9b777 100644
--- a/sysdeps/i386/fpu/libm-test-ulps
+++ b/sysdeps/i386/fpu/libm-test-ulps
@@ -569,12 +569,18 @@ double: 1
float: 1
idouble: 1
ifloat: 1
+ildouble: 1
+ldouble: 1
Test "cos_downward (2) == -0.4161468365471423869975682295007621897660":
float: 1
ifloat: 1
+ildouble: 1
+ldouble: 1
Test "cos_downward (3) == -0.9899924966004454572715727947312613023937":
double: 1
idouble: 1
+ildouble: 1
+ldouble: 1
Test "cos_downward (4) == -0.6536436208636119146391681830977503814241":
float: 1
ifloat: 1
@@ -649,6 +655,8 @@ double: 1
float: 1
idouble: 1
ifloat: 1
+ildouble: 1
+ldouble: 1
# cos_upward
Test "cos_upward (1) == 0.5403023058681397174009366074429766037323":
@@ -672,6 +680,8 @@ ldouble: 1
Test "cos_upward (4) == -0.6536436208636119146391681830977503814241":
double: 1
idouble: 1
+ildouble: 1
+ldouble: 1
Test "cos_upward (5) == 0.2836621854632262644666391715135573083344":
double: 1
idouble: 1
@@ -692,6 +702,11 @@ double: 1
float: 1
idouble: 1
ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "cos_upward (9) == -0.9111302618846769883682947111811653112463":
+ildouble: 1
+ldouble: 1
# cosh
Test "cosh (0.75) == 1.29468328467684468784170818539018176":
@@ -1014,8 +1029,8 @@ double: 1
float: 2
idouble: 1
ifloat: 2
-ildouble: 1
-ldouble: 1
+ildouble: 2
+ldouble: 2
Test "j0 (10.0) == -0.245935764451348335197760862485328754":
double: 3
float: 1
@@ -1031,8 +1046,8 @@ double: 1
float: 2
idouble: 1
ifloat: 2
-ildouble: 1
-ldouble: 1
+ildouble: 2
+ldouble: 2
Test "j0 (8.0) == 0.171650807137553906090869407851972001":
float: 1
ifloat: 1
@@ -1065,8 +1080,8 @@ double: 1
float: 2
idouble: 1
ifloat: 2
-ildouble: 1
-ldouble: 1
+ildouble: 2
+ldouble: 2
Test "jn (0, 10.0) == -0.245935764451348335197760862485328754":
double: 3
float: 1
@@ -1082,8 +1097,8 @@ double: 1
float: 2
idouble: 1
ifloat: 2
-ildouble: 1
-ldouble: 1
+ildouble: 2
+ldouble: 2
Test "jn (0, 8.0) == 0.171650807137553906090869407851972001":
float: 1
ifloat: 1
@@ -1343,6 +1358,8 @@ ldouble: 1
Test "sin_tonearest (1) == 0.8414709848078965066525023216302989996226":
float: 1
ifloat: 1
+ildouble: 1
+ldouble: 1
Test "sin_tonearest (10) == -0.5440211108893698134047476618513772816836":
ildouble: 1
ldouble: 1
@@ -1369,6 +1386,8 @@ ldouble: 1
Test "sin_towardzero (4) == -0.7568024953079282513726390945118290941359":
float: 1
ifloat: 1
+ildouble: 1
+ldouble: 1
Test "sin_towardzero (5) == -0.9589242746631384688931544061559939733525":
float: 1
ifloat: 1
@@ -1397,9 +1416,13 @@ double: 1
float: 1
idouble: 1
ifloat: 1
+ildouble: 1
+ldouble: 1
Test "sin_upward (10) == -0.5440211108893698134047476618513772816836":
float: 1
ifloat: 1
+ildouble: 1
+ldouble: 1
Test "sin_upward (2) == 0.9092974268256816953960198659117448427023":
float: 1
ifloat: 1
@@ -1410,12 +1433,18 @@ double: 1
float: 1
idouble: 1
ifloat: 1
+ildouble: 1
+ldouble: 1
Test "sin_upward (4) == -0.7568024953079282513726390945118290941359":
float: 1
ifloat: 1
+ildouble: 2
+ldouble: 2
Test "sin_upward (5) == -0.9589242746631384688931544061559939733525":
float: 1
ifloat: 1
+ildouble: 1
+ldouble: 1
Test "sin_upward (6) == -0.2794154981989258728115554466118947596280":
ildouble: 1
ldouble: 1
@@ -2209,8 +2238,8 @@ double: 3
float: 2
idouble: 3
ifloat: 2
-ildouble: 1
-ldouble: 1
+ildouble: 2
+ldouble: 2
Function: "j1":
double: 2
@@ -2301,8 +2330,8 @@ double: 1
float: 1
idouble: 1
ifloat: 1
-ildouble: 1
-ldouble: 1
+ildouble: 2
+ldouble: 2
Function: "sincos":
double: 1
diff --git a/sysdeps/i386/fpu/s_cosl.S b/sysdeps/i386/fpu/s_cosl.S
deleted file mode 100644
index 27dd74f..0000000
--- a/sysdeps/i386/fpu/s_cosl.S
+++ /dev/null
@@ -1,54 +0,0 @@
-/*
- * Written by J.T. Conklin <jtc@netbsd.org>.
- * Public domain.
- *
- * Adapted for `long double' by Ulrich Drepper <drepper@cygnus.com>.
- * Fixed errno handling by Ulrich Drepper <drepper@redhat.com>.
- */
-
-#define __need_Emath
-#include <bits/errno.h>
-#include <machine/asm.h>
-
-ENTRY(__cosl)
- fldt 4(%esp)
- fxam
- fstsw %ax
- movb $0x45, %dh
- andb %ah, %dh
- cmpb $0x05, %dh
- je 3f
-4: fcos
- fnstsw %ax
- testl $0x400,%eax
- jnz 1f
- ret
- .align ALIGNARG(4)
-1: fldpi
- fadd %st(0)
- fxch %st(1)
-2: fprem1
- fnstsw %ax
- testl $0x400,%eax
- jnz 2b
- fstp %st(1)
- fcos
- ret
-3:
-#ifdef PIC
- pushl %ebx
- cfi_adjust_cfa_offset (4)
- cfi_rel_offset (ebx, 0)
- LOAD_PIC_REG (bx)
- call __errno_location@PLT
- movl $EDOM, (%eax)
- popl %ebx
- cfi_adjust_cfa_offset (-4)
- cfi_restore (ebx)
-#else
- call __errno_location@PLT
- movl $EDOM, (%eax)
-#endif
- jmp 4b
-END (__cosl)
-weak_alias (__cosl, cosl)
diff --git a/sysdeps/i386/fpu/s_sincosl.S b/sysdeps/i386/fpu/s_sincosl.S
deleted file mode 100644
index 42a9625..0000000
--- a/sysdeps/i386/fpu/s_sincosl.S
+++ /dev/null
@@ -1,64 +0,0 @@
-/* Compute sine and cosine of argument.
- Copyright (C) 1997, 2000 Free Software Foundation, Inc.
- This file is part of the GNU C Library.
- Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997.
-
- The GNU C Library is free software; you can redistribute it and/or
- modify it under the terms of the GNU Lesser General Public
- License as published by the Free Software Foundation; either
- version 2.1 of the License, or (at your option) any later version.
-
- The GNU C Library is distributed in the hope that it will be useful,
- but WITHOUT ANY WARRANTY; without even the implied warranty of
- MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
- Lesser General Public License for more details.
-
- You should have received a copy of the GNU Lesser General Public
- License along with the GNU C Library; if not, see
- <http://www.gnu.org/licenses/>. */
-
-#include <machine/asm.h>
-#include "bp-sym.h"
-#include "bp-asm.h"
-
-#define PARMS LINKAGE /* no space for saved regs */
-#define ANGLE PARMS
-#define SINP ANGLE+12
-#define COSP SINP+PTR_SIZE
-
- .text
-ENTRY (BP_SYM (__sincosl))
- ENTER
-
- fldt ANGLE(%esp)
- fsincos
- movl SINP(%esp), %ecx
- CHECK_BOUNDS_BOTH_WIDE (%ecx, SINP(%esp), $12)
- movl COSP(%esp), %edx
- CHECK_BOUNDS_BOTH_WIDE (%edx, COSP(%esp), $12)
- fnstsw %ax
- testl $0x400,%eax
- jnz 1f
- fstpt (%edx)
- fstpt (%ecx)
-
- LEAVE
- ret
-
- .align ALIGNARG(4)
-1: fldpi
- fadd %st(0)
- fxch %st(1)
-2: fprem1
- fnstsw %ax
- testl $0x400,%eax
- jnz 2b
- fstp %st(1)
- fsincos
- fstpt (%edx)
- fstpt (%ecx)
-
- LEAVE
- ret
-END (BP_SYM (__sincosl))
-weak_alias (BP_SYM (__sincosl), BP_SYM (sincosl))
diff --git a/sysdeps/i386/fpu/s_sinl.S b/sysdeps/i386/fpu/s_sinl.S
deleted file mode 100644
index 68c4f99..0000000
--- a/sysdeps/i386/fpu/s_sinl.S
+++ /dev/null
@@ -1,54 +0,0 @@
-/*
- * Written by J.T. Conklin <jtc@netbsd.org>.
- * Public domain.
- *
- * Adapted for `long double' by Ulrich Drepper <drepper@cygnus.com>.
- * Fixed errno handling by Ulrich Drepper <drepper@redhat.com>.
- */
-
-#define __need_Emath
-#include <bits/errno.h>
-#include <machine/asm.h>
-
-ENTRY(__sinl)
- fldt 4(%esp)
- fxam
- fstsw %ax
- movb $0x45, %dh
- andb %ah, %dh
- cmpb $0x05, %dh
- je 3f
-4: fsin
- fnstsw %ax
- testl $0x400,%eax
- jnz 1f
- ret
- .align ALIGNARG(4)
-1: fldpi
- fadd %st(0)
- fxch %st(1)
-2: fprem1
- fnstsw %ax
- testl $0x400,%eax
- jnz 2b
- fstp %st(1)
- fsin
- ret
-3:
-#ifdef PIC
- pushl %ebx
- cfi_adjust_cfa_offset (4)
- cfi_rel_offset (ebx, 0)
- LOAD_PIC_REG (bx)
- call __errno_location@PLT
- movl $EDOM, (%eax)
- popl %ebx
- cfi_adjust_cfa_offset (-4)
- cfi_restore (ebx)
-#else
- call __errno_location@PLT
- movl $EDOM, (%eax)
-#endif
- jmp 4b
-END (__sinl)
-weak_alias (__sinl, sinl)
diff --git a/sysdeps/ieee754/dbl-64/k_rem_pio2.c b/sysdeps/ieee754/dbl-64/k_rem_pio2.c
index a1e0c6d..fcf956a 100644
--- a/sysdeps/ieee754/dbl-64/k_rem_pio2.c
+++ b/sysdeps/ieee754/dbl-64/k_rem_pio2.c
@@ -273,13 +273,16 @@ recompute:
y[0] = (ih==0)? fw: -fw;
break;
case 1:
- case 2:
- fw = 0.0;
- for (i=jz;i>=0;i--) fw += fq[i];
- y[0] = (ih==0)? fw: -fw;
- fw = fq[0]-fw;
- for (i=1;i<=jz;i++) fw += fq[i];
- y[1] = (ih==0)? fw: -fw;
+ case 2:;
+#if __FLT_EVAL_METHOD__ != 0
+ volatile
+#endif
+ double fv = 0.0;
+ for (i=jz;i>=0;i--) fv += fq[i];
+ y[0] = (ih==0)? fv: -fv;
+ fv = fq[0]-fv;
+ for (i=1;i<=jz;i++) fv += fq[i];
+ y[1] = (ih==0)? fv: -fv;
break;
case 3: /* painful */
for (i=jz;i>0;i--) {
diff --git a/sysdeps/ieee754/ldbl-96/e_rem_pio2l.c b/sysdeps/ieee754/ldbl-96/e_rem_pio2l.c
new file mode 100644
index 0000000..c35f2fe
--- /dev/null
+++ b/sysdeps/ieee754/ldbl-96/e_rem_pio2l.c
@@ -0,0 +1,236 @@
+/* Extended-precision floating point argument reduction.
+ Copyright (C) 1999-2012 Free Software Foundation, Inc.
+ This file is part of the GNU C Library.
+ Based on quad-precision code by Jakub Jelinek <jj@ultra.linux.cz>
+
+ The GNU C Library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ The GNU C Library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with the GNU C Library; if not, see
+ <http://www.gnu.org/licenses/>. */
+
+#include <math.h>
+#include <math_private.h>
+
+/* Table of constants for 2/pi, 5628 hexadecimal digits of 2/pi. */
+static const int32_t two_over_pi[] = {
+0xa2f983, 0x6e4e44, 0x1529fc, 0x2757d1, 0xf534dd, 0xc0db62,
+0x95993c, 0x439041, 0xfe5163, 0xabdebb, 0xc561b7, 0x246e3a,
+0x424dd2, 0xe00649, 0x2eea09, 0xd1921c, 0xfe1deb, 0x1cb129,
+0xa73ee8, 0x8235f5, 0x2ebb44, 0x84e99c, 0x7026b4, 0x5f7e41,
+0x3991d6, 0x398353, 0x39f49c, 0x845f8b, 0xbdf928, 0x3b1ff8,
+0x97ffde, 0x05980f, 0xef2f11, 0x8b5a0a, 0x6d1f6d, 0x367ecf,
+0x27cb09, 0xb74f46, 0x3f669e, 0x5fea2d, 0x7527ba, 0xc7ebe5,
+0xf17b3d, 0x0739f7, 0x8a5292, 0xea6bfb, 0x5fb11f, 0x8d5d08,
+0x560330, 0x46fc7b, 0x6babf0, 0xcfbc20, 0x9af436, 0x1da9e3,
+0x91615e, 0xe61b08, 0x659985, 0x5f14a0, 0x68408d, 0xffd880,
+0x4d7327, 0x310606, 0x1556ca, 0x73a8c9, 0x60e27b, 0xc08c6b,
+0x47c419, 0xc367cd, 0xdce809, 0x2a8359, 0xc4768b, 0x961ca6,
+0xddaf44, 0xd15719, 0x053ea5, 0xff0705, 0x3f7e33, 0xe832c2,
+0xde4f98, 0x327dbb, 0xc33d26, 0xef6b1e, 0x5ef89f, 0x3a1f35,
+0xcaf27f, 0x1d87f1, 0x21907c, 0x7c246a, 0xfa6ed5, 0x772d30,
+0x433b15, 0xc614b5, 0x9d19c3, 0xc2c4ad, 0x414d2c, 0x5d000c,
+0x467d86, 0x2d71e3, 0x9ac69b, 0x006233, 0x7cd2b4, 0x97a7b4,
+0xd55537, 0xf63ed7, 0x1810a3, 0xfc764d, 0x2a9d64, 0xabd770,
+0xf87c63, 0x57b07a, 0xe71517, 0x5649c0, 0xd9d63b, 0x3884a7,
+0xcb2324, 0x778ad6, 0x23545a, 0xb91f00, 0x1b0af1, 0xdfce19,
+0xff319f, 0x6a1e66, 0x615799, 0x47fbac, 0xd87f7e, 0xb76522,
+0x89e832, 0x60bfe6, 0xcdc4ef, 0x09366c, 0xd43f5d, 0xd7de16,
+0xde3b58, 0x929bde, 0x2822d2, 0xe88628, 0x4d58e2, 0x32cac6,
+0x16e308, 0xcb7de0, 0x50c017, 0xa71df3, 0x5be018, 0x34132e,
+0x621283, 0x014883, 0x5b8ef5, 0x7fb0ad, 0xf2e91e, 0x434a48,
+0xd36710, 0xd8ddaa, 0x425fae, 0xce616a, 0xa4280a, 0xb499d3,
+0xf2a606, 0x7f775c, 0x83c2a3, 0x883c61, 0x78738a, 0x5a8caf,
+0xbdd76f, 0x63a62d, 0xcbbff4, 0xef818d, 0x67c126, 0x45ca55,
+0x36d9ca, 0xd2a828, 0x8d61c2, 0x77c912, 0x142604, 0x9b4612,
+0xc459c4, 0x44c5c8, 0x91b24d, 0xf31700, 0xad43d4, 0xe54929,
+0x10d5fd, 0xfcbe00, 0xcc941e, 0xeece70, 0xf53e13, 0x80f1ec,
+0xc3e7b3, 0x28f8c7, 0x940593, 0x3e71c1, 0xb3092e, 0xf3450b,
+0x9c1288, 0x7b20ab, 0x9fb52e, 0xc29247, 0x2f327b, 0x6d550c,
+0x90a772, 0x1fe76b, 0x96cb31, 0x4a1679, 0xe27941, 0x89dff4,
+0x9794e8, 0x84e6e2, 0x973199, 0x6bed88, 0x365f5f, 0x0efdbb,
+0xb49a48, 0x6ca467, 0x427271, 0x325d8d, 0xb8159f, 0x09e5bc,
+0x25318d, 0x3974f7, 0x1c0530, 0x010c0d, 0x68084b, 0x58ee2c,
+0x90aa47, 0x02e774, 0x24d6bd, 0xa67df7, 0x72486e, 0xef169f,
+0xa6948e, 0xf691b4, 0x5153d1, 0xf20acf, 0x339820, 0x7e4bf5,
+0x6863b2, 0x5f3edd, 0x035d40, 0x7f8985, 0x295255, 0xc06437,
+0x10d86d, 0x324832, 0x754c5b, 0xd4714e, 0x6e5445, 0xc1090b,
+0x69f52a, 0xd56614, 0x9d0727, 0x50045d, 0xdb3bb4, 0xc576ea,
+0x17f987, 0x7d6b49, 0xba271d, 0x296996, 0xacccc6, 0x5414ad,
+0x6ae290, 0x89d988, 0x50722c, 0xbea404, 0x940777, 0x7030f3,
+0x27fc00, 0xa871ea, 0x49c266, 0x3de064, 0x83dd97, 0x973fa3,
+0xfd9443, 0x8c860d, 0xde4131, 0x9d3992, 0x8c70dd, 0xe7b717,
+0x3bdf08, 0x2b3715, 0xa0805c, 0x93805a, 0x921110, 0xd8e80f,
+0xaf806c, 0x4bffdb, 0x0f9038, 0x761859, 0x15a562, 0xbbcb61,
+0xb989c7, 0xbd4010, 0x04f2d2, 0x277549, 0xf6b6eb, 0xbb22db,
+0xaa140a, 0x2f2689, 0x768364, 0x333b09, 0x1a940e, 0xaa3a51,
+0xc2a31d, 0xaeedaf, 0x12265c, 0x4dc26d, 0x9c7a2d, 0x9756c0,
+0x833f03, 0xf6f009, 0x8c402b, 0x99316d, 0x07b439, 0x15200c,
+0x5bc3d8, 0xc492f5, 0x4badc6, 0xa5ca4e, 0xcd37a7, 0x36a9e6,
+0x9492ab, 0x6842dd, 0xde6319, 0xef8c76, 0x528b68, 0x37dbfc,
+0xaba1ae, 0x3115df, 0xa1ae00, 0xdafb0c, 0x664d64, 0xb705ed,
+0x306529, 0xbf5657, 0x3aff47, 0xb9f96a, 0xf3be75, 0xdf9328,
+0x3080ab, 0xf68c66, 0x15cb04, 0x0622fa, 0x1de4d9, 0xa4b33d,
+0x8f1b57, 0x09cd36, 0xe9424e, 0xa4be13, 0xb52333, 0x1aaaf0,
+0xa8654f, 0xa5c1d2, 0x0f3f0b, 0xcd785b, 0x76f923, 0x048b7b,
+0x721789, 0x53a6c6, 0xe26e6f, 0x00ebef, 0x584a9b, 0xb7dac4,
+0xba66aa, 0xcfcf76, 0x1d02d1, 0x2df1b1, 0xc1998c, 0x77adc3,
+0xda4886, 0xa05df7, 0xf480c6, 0x2ff0ac, 0x9aecdd, 0xbc5c3f,
+0x6dded0, 0x1fc790, 0xb6db2a, 0x3a25a3, 0x9aaf00, 0x9353ad,
+0x0457b6, 0xb42d29, 0x7e804b, 0xa707da, 0x0eaa76, 0xa1597b,
+0x2a1216, 0x2db7dc, 0xfde5fa, 0xfedb89, 0xfdbe89, 0x6c76e4,
+0xfca906, 0x70803e, 0x156e85, 0xff87fd, 0x073e28, 0x336761,
+0x86182a, 0xeabd4d, 0xafe7b3, 0x6e6d8f, 0x396795, 0x5bbf31,
+0x48d784, 0x16df30, 0x432dc7, 0x356125, 0xce70c9, 0xb8cb30,
+0xfd6cbf, 0xa200a4, 0xe46c05, 0xa0dd5a, 0x476f21, 0xd21262,
+0x845cb9, 0x496170, 0xe0566b, 0x015299, 0x375550, 0xb7d51e,
+0xc4f133, 0x5f6e13, 0xe4305d, 0xa92e85, 0xc3b21d, 0x3632a1,
+0xa4b708, 0xd4b1ea, 0x21f716, 0xe4698f, 0x77ff27, 0x80030c,
+0x2d408d, 0xa0cd4f, 0x99a520, 0xd3a2b3, 0x0a5d2f, 0x42f9b4,
+0xcbda11, 0xd0be7d, 0xc1db9b, 0xbd17ab, 0x81a2ca, 0x5c6a08,
+0x17552e, 0x550027, 0xf0147f, 0x8607e1, 0x640b14, 0x8d4196,
+0xdebe87, 0x2afdda, 0xb6256b, 0x34897b, 0xfef305, 0x9ebfb9,
+0x4f6a68, 0xa82a4a, 0x5ac44f, 0xbcf82d, 0x985ad7, 0x95c7f4,
+0x8d4d0d, 0xa63a20, 0x5f57a4, 0xb13f14, 0x953880, 0x0120cc,
+0x86dd71, 0xb6dec9, 0xf560bf, 0x11654d, 0x6b0701, 0xacb08c,
+0xd0c0b2, 0x485551, 0x0efb1e, 0xc37295, 0x3b06a3, 0x3540c0,
+0x7bdc06, 0xcc45e0, 0xfa294e, 0xc8cad6, 0x41f3e8, 0xde647c,
+0xd8649b, 0x31bed9, 0xc397a4, 0xd45877, 0xc5e369, 0x13daf0,
+0x3c3aba, 0x461846, 0x5f7555, 0xf5bdd2, 0xc6926e, 0x5d2eac,
+0xed440e, 0x423e1c, 0x87c461, 0xe9fd29, 0xf3d6e7, 0xca7c22,
+0x35916f, 0xc5e008, 0x8dd7ff, 0xe26a6e, 0xc6fdb0, 0xc10893,
+0x745d7c, 0xb2ad6b, 0x9d6ecd, 0x7b723e, 0x6a11c6, 0xa9cff7,
+0xdf7329, 0xbac9b5, 0x5100b7, 0x0db2e2, 0x24ba74, 0x607de5,
+0x8ad874, 0x2c150d, 0x0c1881, 0x94667e, 0x162901, 0x767a9f,
+0xbefdfd, 0xef4556, 0x367ed9, 0x13d9ec, 0xb9ba8b, 0xfc97c4,
+0x27a831, 0xc36ef1, 0x36c594, 0x56a8d8, 0xb5a8b4, 0x0ecccf,
+0x2d8912, 0x34576f, 0x89562c, 0xe3ce99, 0xb920d6, 0xaa5e6b,
+0x9c2a3e, 0xcc5f11, 0x4a0bfd, 0xfbf4e1, 0x6d3b8e, 0x2c86e2,
+0x84d4e9, 0xa9b4fc, 0xd1eeef, 0xc9352e, 0x61392f, 0x442138,
+0xc8d91b, 0x0afc81, 0x6a4afb, 0xd81c2f, 0x84b453, 0x8c994e,
+0xcc2254, 0xdc552a, 0xd6c6c0, 0x96190b, 0xb8701a, 0x649569,
+0x605a26, 0xee523f, 0x0f117f, 0x11b5f4, 0xf5cbfc, 0x2dbc34,
+0xeebc34, 0xcc5de8, 0x605edd, 0x9b8e67, 0xef3392, 0xb817c9,
+0x9b5861, 0xbc57e1, 0xc68351, 0x103ed8, 0x4871dd, 0xdd1c2d,
+0xa118af, 0x462c21, 0xd7f359, 0x987ad9, 0xc0549e, 0xfa864f,
+0xfc0656, 0xae79e5, 0x362289, 0x22ad38, 0xdc9367, 0xaae855,
+0x382682, 0x9be7ca, 0xa40d51, 0xb13399, 0x0ed7a9, 0x480569,
+0xf0b265, 0xa7887f, 0x974c88, 0x36d1f9, 0xb39221, 0x4a827b,
+0x21cf98, 0xdc9f40, 0x5547dc, 0x3a74e1, 0x42eb67, 0xdf9dfe,
+0x5fd45e, 0xa4677b, 0x7aacba, 0xa2f655, 0x23882b, 0x55ba41,
+0x086e59, 0x862a21, 0x834739, 0xe6e389, 0xd49ee5, 0x40fb49,
+0xe956ff, 0xca0f1c, 0x8a59c5, 0x2bfa94, 0xc5c1d3, 0xcfc50f,
+0xae5adb, 0x86c547, 0x624385, 0x3b8621, 0x94792c, 0x876110,
+0x7b4c2a, 0x1a2c80, 0x12bf43, 0x902688, 0x893c78, 0xe4c4a8,
+0x7bdbe5, 0xc23ac4, 0xeaf426, 0x8a67f7, 0xbf920d, 0x2ba365,
+0xb1933d, 0x0b7cbd, 0xdc51a4, 0x63dd27, 0xdde169, 0x19949a,
+0x9529a8, 0x28ce68, 0xb4ed09, 0x209f44, 0xca984e, 0x638270,
+0x237c7e, 0x32b90f, 0x8ef5a7, 0xe75614, 0x08f121, 0x2a9db5,
+0x4d7e6f, 0x5119a5, 0xabf9b5, 0xd6df82, 0x61dd96, 0x023616,
+0x9f3ac4, 0xa1a283, 0x6ded72, 0x7a8d39, 0xa9b882, 0x5c326b,
+0x5b2746, 0xed3400, 0x7700d2, 0x55f4fc, 0x4d5901, 0x8071e0,
+0xe13f89, 0xb295f3, 0x64a8f1, 0xaea74b, 0x38fc4c, 0xeab2bb,
+0x47270b, 0xabc3a7, 0x34ba60, 0x52dd34, 0xf8563a, 0xeb7e8a,
+0x31bb36, 0x5895b7, 0x47f7a9, 0x94c3aa, 0xd39225, 0x1e7f3e,
+0xd8974e, 0xbba94f, 0xd8ae01, 0xe661b4, 0x393d8e, 0xa523aa,
+0x33068e, 0x1633b5, 0x3bb188, 0x1d3a9d, 0x4013d0, 0xcc1be5,
+0xf862e7, 0x3bf28f, 0x39b5bf, 0x0bc235, 0x22747e, 0xa247c0,
+0xd52d1f, 0x19add3, 0x9094df, 0x9311d0, 0xb42b25, 0x496db2,
+0xe264b2, 0x5ef135, 0x3bc6a4, 0x1a4ad0, 0xaac92e, 0x64e886,
+0x573091, 0x982cfb, 0x311b1a, 0x08728b, 0xbdcee1, 0x60e142,
+0xeb641d, 0xd0bba3, 0xe559d4, 0x597b8c, 0x2a4483, 0xf332ba,
+0xf84867, 0x2c8d1b, 0x2fa9b0, 0x50f3dd, 0xf9f573, 0xdb61b4,
+0xfe233e, 0x6c41a6, 0xeea318, 0x775a26, 0xbc5e5c, 0xcea708,
+0x94dc57, 0xe20196, 0xf1e839, 0xbe4851, 0x5d2d2f, 0x4e9555,
+0xd96ec2, 0xe7d755, 0x6304e0, 0xc02e0e, 0xfc40a0, 0xbbf9b3,
+0x7125a7, 0x222dfb, 0xf619d8, 0x838c1c, 0x6619e6, 0xb20d55,
+0xbb5137, 0x79e809, 0xaf9149, 0x0d73de, 0x0b0da5, 0xce7f58,
+0xac1934, 0x724667, 0x7a1a13, 0x9e26bc, 0x4555e7, 0x585cb5,
+0x711d14, 0x486991, 0x480d60, 0x56adab, 0xd62f64, 0x96ee0c,
+0x212ff3, 0x5d6d88, 0xa67684, 0x95651e, 0xab9e0a, 0x4ddefe,
+0x571010, 0x836a39, 0xf8ea31, 0x9e381d, 0xeac8b1, 0xcac96b,
+0x37f21e, 0xd505e9, 0x984743, 0x9fc56c, 0x0331b7, 0x3b8bf8,
+0x86e56a, 0x8dc343, 0x6230e7, 0x93cfd5, 0x6a8f2d, 0x733005,
+0x1af021, 0xa09fcb, 0x7415a1, 0xd56b23, 0x6ff725, 0x2f4bc7,
+0xb8a591, 0x7fac59, 0x5c55de, 0x212c38, 0xb13296, 0x5cff50,
+0x366262, 0xfa7b16, 0xf4d9a6, 0x2acfe7, 0xf07403, 0xd4d604,
+0x6fd916, 0x31b1bf, 0xcbb450, 0x5bd7c8, 0x0ce194, 0x6bd643,
+0x4fd91c, 0xdf4543, 0x5f3453, 0xe2b5aa, 0xc9aec8, 0x131485,
+0xf9d2bf, 0xbadb9e, 0x76f5b9, 0xaf15cf, 0xca3182, 0x14b56d,
+0xe9fe4d, 0x50fc35, 0xf5aed5, 0xa2d0c1, 0xc96057, 0x192eb6,
+0xe91d92, 0x07d144, 0xaea3c6, 0x343566, 0x26d5b4, 0x3161e2,
+0x37f1a2, 0x209eff, 0x958e23, 0x493798, 0x35f4a6, 0x4bdc02,
+0xc2be13, 0xbe80a0, 0x0b72a3, 0x115c5f, 0x1e1bd1, 0x0db4d3,
+0x869e85, 0x96976b, 0x2ac91f, 0x8a26c2, 0x3070f0, 0x041412,
+0xfc9fa5, 0xf72a38, 0x9c6878, 0xe2aa76, 0x50cfe1, 0x559274,
+0x934e38, 0x0a92f7, 0x5533f0, 0xa63db4, 0x399971, 0xe2b755,
+0xa98a7c, 0x008f19, 0xac54d2, 0x2ea0b4, 0xf5f3e0, 0x60c849,
+0xffd269, 0xae52ce, 0x7a5fdd, 0xe9ce06, 0xfb0ae8, 0xa50cce,
+0xea9d3e, 0x3766dd, 0xb834f5, 0x0da090, 0x846f88, 0x4ae3d5,
+0x099a03, 0x2eae2d, 0xfcb40a, 0xfb9b33, 0xe281dd, 0x1b16ba,
+0xd8c0af, 0xd96b97, 0xb52dc9, 0x9c277f, 0x5951d5, 0x21ccd6,
+0xb6496b, 0x584562, 0xb3baf2, 0xa1a5c4, 0x7ca2cf, 0xa9b93d,
+0x7b7b89, 0x483d38,
+};
+
+int32_t
+__ieee754_rem_pio2l (long double x, long double *y)
+{
+ double tx[3], ty[3];
+ int32_t se, j0;
+ u_int32_t i0, i1;
+ int sx;
+ int n, exp;
+
+ GET_LDOUBLE_WORDS (se, i0, i1, x);
+ sx = (se >> 15) & 1;
+ j0 = (se & 0x7fff) - 0x3fff;
+
+ if (j0 < -1)
+ {
+ /* |x| < pi/4. */
+ y[0] = x;
+ y[1] = 0;
+ return 0;
+ }
+
+ if (j0 >= 0x8000)
+ {
+ /* x is infinite or NaN. */
+ y[0] = x - x;
+ y[1] = y[0];
+ return 0;
+ }
+
+ /* Split the 64 bits of the mantissa into three 24-bit integers
+ stored in a double array. */
+ exp = j0 - 23;
+ tx[0] = (double) (i0 >> 8);
+ tx[1] = (double) (((i0 << 16) | (i1 >> 16)) & 0xffffff);
+ tx[2] = (double) ((i1 << 8) & 0xffffff);
+
+ n = __kernel_rem_pio2 (tx, ty, exp, 3, 2, two_over_pi);
+
+ /* The result is now stored in two double values, we need to convert
+ it into two long double values. */
+ if (sx == 0)
+ {
+ y[0] = (long double) ty[0] + (long double) ty[1];
+ y[1] = ty[1] - (y[0] - ty[0]);
+ return n;
+ }
+ else
+ {
+ y[0] = -((long double) ty[0] + (long double) ty[1]);
+ y[1] = -ty[1] - (y[0] + ty[0]);
+ return -n;
+ }
+}
diff --git a/sysdeps/ieee754/ldbl-96/k_cosl.c b/sysdeps/ieee754/ldbl-96/k_cosl.c
new file mode 100644
index 0000000..d684f34
--- /dev/null
+++ b/sysdeps/ieee754/ldbl-96/k_cosl.c
@@ -0,0 +1,119 @@
+/* Extended-precision floating point cosine on <-pi/4,pi/4>.
+ Copyright (C) 1999-2012 Free Software Foundation, Inc.
+ This file is part of the GNU C Library.
+ Based on quad-precision cosine by Jakub Jelinek <jj@ultra.linux.cz>
+
+ The GNU C Library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ The GNU C Library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with the GNU C Library; if not, see
+ <http://www.gnu.org/licenses/>. */
+
+#include <math.h>
+#include <math_private.h>
+
+/* The polynomials have not been optimized for extended-precision and
+ may contain more terms than needed. */
+
+static const long double c[] = {
+#define ONE c[0]
+ 1.00000000000000000000000000000000000E+00L,
+
+/* cos x ~ ONE + x^2 ( SCOS1 + SCOS2 * x^2 + ... + SCOS4 * x^6 + SCOS5 * x^8 )
+ x in <0,1/256> */
+#define SCOS1 c[1]
+#define SCOS2 c[2]
+#define SCOS3 c[3]
+#define SCOS4 c[4]
+#define SCOS5 c[5]
+-5.00000000000000000000000000000000000E-01L,
+ 4.16666666666666666666666666556146073E-02L,
+-1.38888888888888888888309442601939728E-03L,
+ 2.48015873015862382987049502531095061E-05L,
+-2.75573112601362126593516899592158083E-07L,
+
+/* cos x ~ ONE + x^2 ( COS1 + COS2 * x^2 + ... + COS7 * x^12 + COS8 * x^14 )
+ x in <0,0.1484375> */
+#define COS1 c[6]
+#define COS2 c[7]
+#define COS3 c[8]
+#define COS4 c[9]
+#define COS5 c[10]
+#define COS6 c[11]
+#define COS7 c[12]
+#define COS8 c[13]
+-4.99999999999999999999999999999999759E-01L,
+ 4.16666666666666666666666666651287795E-02L,
+-1.38888888888888888888888742314300284E-03L,
+ 2.48015873015873015867694002851118210E-05L,
+-2.75573192239858811636614709689300351E-07L,
+ 2.08767569877762248667431926878073669E-09L,
+-1.14707451049343817400420280514614892E-11L,
+ 4.77810092804389587579843296923533297E-14L,
+
+/* sin x ~ ONE * x + x^3 ( SSIN1 + SSIN2 * x^2 + ... + SSIN4 * x^6 + SSIN5 * x^8 )
+ x in <0,1/256> */
+#define SSIN1 c[14]
+#define SSIN2 c[15]
+#define SSIN3 c[16]
+#define SSIN4 c[17]
+#define SSIN5 c[18]
+-1.66666666666666666666666666666666659E-01L,
+ 8.33333333333333333333333333146298442E-03L,
+-1.98412698412698412697726277416810661E-04L,
+ 2.75573192239848624174178393552189149E-06L,
+-2.50521016467996193495359189395805639E-08L,
+};
+
+#define SINCOSL_COS_HI 0
+#define SINCOSL_COS_LO 1
+#define SINCOSL_SIN_HI 2
+#define SINCOSL_SIN_LO 3
+extern const long double __sincosl_table[];
+
+long double
+__kernel_cosl(long double x, long double y)
+{
+ long double h, l, z, sin_l, cos_l_m1;
+ int index;
+
+ x = fabsl (x);
+ if (x < 0.1484375L)
+ {
+ /* Argument is small enough to approximate it by a Chebyshev
+ polynomial of degree 16. */
+ if (x < 0x1p-33L)
+ if (!((int)x)) return ONE; /* generate inexact */
+ z = x * x;
+ return ONE + (z*(COS1+z*(COS2+z*(COS3+z*(COS4+
+ z*(COS5+z*(COS6+z*(COS7+z*COS8))))))));
+ }
+ else
+ {
+ /* So that we don't have to use too large polynomial, we find
+ l and h such that x = l + h, where fabsl(l) <= 1.0/256 with 83
+ possible values for h. We look up cosl(h) and sinl(h) in
+ pre-computed tables, compute cosl(l) and sinl(l) using a
+ Chebyshev polynomial of degree 10(11) and compute
+ cosl(h+l) = cosl(h)cosl(l) - sinl(h)sinl(l). */
+ index = (int) (128 * (x - (0.1484375L - 1.0L / 256.0L)));
+ h = 0.1484375L + index / 128.0;
+ index *= 4;
+ l = y - (h - x);
+ z = l * l;
+ sin_l = l*(ONE+z*(SSIN1+z*(SSIN2+z*(SSIN3+z*(SSIN4+z*SSIN5)))));
+ cos_l_m1 = z*(SCOS1+z*(SCOS2+z*(SCOS3+z*(SCOS4+z*SCOS5))));
+ return __sincosl_table [index + SINCOSL_COS_HI]
+ + (__sincosl_table [index + SINCOSL_COS_LO]
+ - (__sincosl_table [index + SINCOSL_SIN_HI] * sin_l
+ - __sincosl_table [index + SINCOSL_COS_HI] * cos_l_m1));
+ }
+}
diff --git a/sysdeps/ieee754/ldbl-96/k_sinl.c b/sysdeps/ieee754/ldbl-96/k_sinl.c
new file mode 100644
index 0000000..feb24d9
--- /dev/null
+++ b/sysdeps/ieee754/ldbl-96/k_sinl.c
@@ -0,0 +1,126 @@
+/* Quad-precision floating point sine on <-pi/4,pi/4>.
+ Copyright (C) 1999-2012 Free Software Foundation, Inc.
+ This file is part of the GNU C Library.
+ Based on quad-precision sine by Jakub Jelinek <jj@ultra.linux.cz>
+
+ The GNU C Library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ The GNU C Library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with the GNU C Library; if not, see
+ <http://www.gnu.org/licenses/>. */
+
+/* The polynomials have not been optimized for extended-precision and
+ may contain more terms than needed. */
+
+#include <math.h>
+#include <math_private.h>
+
+/* The polynomials have not been optimized for extended-precision and
+ may contain more terms than needed. */
+
+static const long double c[] = {
+#define ONE c[0]
+ 1.00000000000000000000000000000000000E+00L,
+
+/* cos x ~ ONE + x^2 ( SCOS1 + SCOS2 * x^2 + ... + SCOS4 * x^6 + SCOS5 * x^8 )
+ x in <0,1/256> */
+#define SCOS1 c[1]
+#define SCOS2 c[2]
+#define SCOS3 c[3]
+#define SCOS4 c[4]
+#define SCOS5 c[5]
+-5.00000000000000000000000000000000000E-01L,
+ 4.16666666666666666666666666556146073E-02L,
+-1.38888888888888888888309442601939728E-03L,
+ 2.48015873015862382987049502531095061E-05L,
+-2.75573112601362126593516899592158083E-07L,
+
+/* sin x ~ ONE * x + x^3 ( SIN1 + SIN2 * x^2 + ... + SIN7 * x^12 + SIN8 * x^14 )
+ x in <0,0.1484375> */
+#define SIN1 c[6]
+#define SIN2 c[7]
+#define SIN3 c[8]
+#define SIN4 c[9]
+#define SIN5 c[10]
+#define SIN6 c[11]
+#define SIN7 c[12]
+#define SIN8 c[13]
+-1.66666666666666666666666666666666538e-01L,
+ 8.33333333333333333333333333307532934e-03L,
+-1.98412698412698412698412534478712057e-04L,
+ 2.75573192239858906520896496653095890e-06L,
+-2.50521083854417116999224301266655662e-08L,
+ 1.60590438367608957516841576404938118e-10L,
+-7.64716343504264506714019494041582610e-13L,
+ 2.81068754939739570236322404393398135e-15L,
+
+/* sin x ~ ONE * x + x^3 ( SSIN1 + SSIN2 * x^2 + ... + SSIN4 * x^6 + SSIN5 * x^8 )
+ x in <0,1/256> */
+#define SSIN1 c[14]
+#define SSIN2 c[15]
+#define SSIN3 c[16]
+#define SSIN4 c[17]
+#define SSIN5 c[18]
+-1.66666666666666666666666666666666659E-01L,
+ 8.33333333333333333333333333146298442E-03L,
+-1.98412698412698412697726277416810661E-04L,
+ 2.75573192239848624174178393552189149E-06L,
+-2.50521016467996193495359189395805639E-08L,
+};
+
+#define SINCOSL_COS_HI 0
+#define SINCOSL_COS_LO 1
+#define SINCOSL_SIN_HI 2
+#define SINCOSL_SIN_LO 3
+extern const long double __sincosl_table[];
+
+long double
+__kernel_sinl(long double x, long double y, int iy)
+{
+ long double absx, h, l, z, sin_l, cos_l_m1;
+ int index;
+
+ absx = fabsl (x);
+ if (absx < 0.1484375L)
+ {
+ /* Argument is small enough to approximate it by a Chebyshev
+ polynomial of degree 17. */
+ if (absx < 0x1p-33L)
+ if (!((int)x)) return x; /* generate inexact */
+ z = x * x;
+ return x + (x * (z*(SIN1+z*(SIN2+z*(SIN3+z*(SIN4+
+ z*(SIN5+z*(SIN6+z*(SIN7+z*SIN8)))))))));
+ }
+ else
+ {
+ /* So that we don't have to use too large polynomial, we find
+ l and h such that x = l + h, where fabsl(l) <= 1.0/256 with 83
+ possible values for h. We look up cosl(h) and sinl(h) in
+ pre-computed tables, compute cosl(l) and sinl(l) using a
+ Chebyshev polynomial of degree 10(11) and compute
+ sinl(h+l) = sinl(h)cosl(l) + cosl(h)sinl(l). */
+ index = (int) (128 * (absx - (0.1484375L - 1.0L / 256.0L)));
+ h = 0.1484375L + index / 128.0;
+ index *= 4;
+ if (iy)
+ l = (x < 0 ? -y : y) - (h - absx);
+ else
+ l = absx - h;
+ z = l * l;
+ sin_l = l*(ONE+z*(SSIN1+z*(SSIN2+z*(SSIN3+z*(SSIN4+z*SSIN5)))));
+ cos_l_m1 = z*(SCOS1+z*(SCOS2+z*(SCOS3+z*(SCOS4+z*SCOS5))));
+ z = __sincosl_table [index + SINCOSL_SIN_HI]
+ + (__sincosl_table [index + SINCOSL_SIN_LO]
+ + (__sincosl_table [index + SINCOSL_SIN_HI] * cos_l_m1)
+ + (__sincosl_table [index + SINCOSL_COS_HI] * sin_l));
+ return (x < 0) ? -z : z;
+ }
+}
diff --git a/sysdeps/ieee754/ldbl-96/s_cosl.c b/sysdeps/ieee754/ldbl-96/s_cosl.c
index 9d5606c..8b0b7d3 100644
--- a/sysdeps/ieee754/ldbl-96/s_cosl.c
+++ b/sysdeps/ieee754/ldbl-96/s_cosl.c
@@ -68,7 +68,7 @@ long double __cosl(long double x)
/* cos(Inf or NaN) is NaN */
else if (se==0x7fff) {
- if ((i0 | i1) == 0)
+ if (i1 == 0 && i0 == 0x80000000)
__set_errno (EDOM);
return x-x;
}
diff --git a/sysdeps/ieee754/ldbl-96/s_sinl.c b/sysdeps/ieee754/ldbl-96/s_sinl.c
index f15eb74..11e1899 100644
--- a/sysdeps/ieee754/ldbl-96/s_sinl.c
+++ b/sysdeps/ieee754/ldbl-96/s_sinl.c
@@ -68,7 +68,7 @@ long double __sinl(long double x)
/* sin(Inf or NaN) is NaN */
else if (se==0x7fff) {
- if ((i0 | i1) == 0)
+ if (i1 == 0 && i0 == 0x80000000)
__set_errno (EDOM);
return x-x;
}
diff --git a/sysdeps/ieee754/ldbl-96/t_sincosl.c b/sysdeps/ieee754/ldbl-96/t_sincosl.c
new file mode 100644
index 0000000..1ab9c2c
--- /dev/null
+++ b/sysdeps/ieee754/ldbl-96/t_sincosl.c
@@ -0,0 +1,444 @@
+/* Extended-precision floating point sine and cosine tables.
+ Copyright (C) 1999-2012 Free Software Foundation, Inc.
+ This file is part of the GNU C Library.
+ Based on quad-precision tables by Jakub Jelinek <jj@ultra.linux.cz>
+
+ The GNU C Library is free software; you can redistribute it and/or
+ modify it under the terms of the GNU Lesser General Public
+ License as published by the Free Software Foundation; either
+ version 2.1 of the License, or (at your option) any later version.
+
+ The GNU C Library is distributed in the hope that it will be useful,
+ but WITHOUT ANY WARRANTY; without even the implied warranty of
+ MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+ Lesser General Public License for more details.
+
+ You should have received a copy of the GNU Lesser General Public
+ License along with the GNU C Library; if not, see
+ <http://www.gnu.org/licenses/>. */
+
+/* For 0.1484375 + n/128.0, n=0..82 this table contains
+ first 64 bits of cosine, then at least 64 additional
+ bits and the same for sine.
+ 0.1484375+82.0/128.0 is the smallest number among above defined numbers
+ larger than pi/4.
+ Computed using MPFR.
+ */
+
+const long double __sincosl_table[] = {
+ /* x = 0.1484375 + 0/128. */
+ 0xf.d2f5320e1b7902100p-4L,
+ -0x6.4b225d06708635580p-68L,
+ 0x2.5dc50bc95711d0d80p-4L,
+ 0x1.787d108fd438cf5a0p-68L,
+ /* x = 0.1484375 + 1/128. */
+ 0xf.ce1a053e621438b00p-4L,
+ 0x6.d60c76e8c45bf0a80p-68L,
+ 0x2.7d66258bacd96a400p-4L,
+ -0x1.4cca4c9a3782a6bc0p-68L,
+ /* x = 0.1484375 + 2/128. */
+ 0xf.c8ffa01ba68074100p-4L,
+ 0x7.e05962b0d9fdf2000p-68L,
+ 0x2.9cfd49b8be4f66540p-4L,
+ -0x1.89354fe340fbd96c0p-68L,
+ /* x = 0.1484375 + 3/128. */
+ 0xf.c3a6170f767ac7300p-4L,
+ 0x5.d63d99a9d439e1d80p-68L,
+ 0x2.bc89f9f424de54840p-4L,
+ 0x1.de7ce03b2514952c0p-68L,
+ /* x = 0.1484375 + 4/128. */
+ 0xf.be0d7f7fef11e7100p-4L,
+ -0x5.5bc47540b095ba800p-68L,
+ 0x2.dc0bb80b49a97ffc0p-4L,
+ -0xc.b1722e07246208500p-72L,
+ /* x = 0.1484375 + 5/128. */
+ 0xf.b835efcf670dd2d00p-4L,
+ -0x1.90186db968115ec20p-68L,
+ 0x2.fb8205f75e56a2b40p-4L,
+ 0x1.6a1c4792f85625880p-68L,
+ /* x = 0.1484375 + 6/128. */
+ 0xf.b21f7f5c156696b00p-4L,
+ 0xa.c1fe28ac5fd766700p-76L,
+ 0x3.1aec65df552876f80p-4L,
+ 0x2.ece9a235671324700p-72L,
+ /* x = 0.1484375 + 7/128. */
+ 0xf.abca467fb3cb8f200p-4L,
+ -0x2.f960fe2715cc521c0p-68L,
+ 0x3.3a4a5a19d86246700p-4L,
+ 0x1.0f602c44df4fa5140p-68L,
+ /* x = 0.1484375 + 8/128. */
+ 0xf.a5365e8f1d3ca2800p-4L,
+ -0x4.1e24a289519b26800p-68L,
+ 0x3.599b652f40ec999c0p-4L,
+ 0x1.f12a0a4c8561de160p-68L,
+ /* x = 0.1484375 + 9/128. */
+ 0xf.9e63e1d9e8b6f6f00p-4L,
+ 0x2.e296bae5b5ed9c100p-68L,
+ 0x3.78df09db8c332ce00p-4L,
+ 0xd.2b53d865582e45200p-72L,
+ /* x = 0.1484375 + 10/128. */
+ 0xf.9752eba9fff6b9900p-4L,
+ -0x7.bd415254fab56cd00p-68L,
+ 0x3.9814cb10513453cc0p-4L,
+ -0x6.84de43e3595cc8500p-72L,
+ /* x = 0.1484375 + 11/128. */
+ 0xf.90039843324f9b900p-4L,
+ 0x4.0416c1984b6cbed00p-68L,
+ 0x3.b73c2bf6b4b9f6680p-4L,
+ 0xe.f9499c81f0d965100p-72L,
+ /* x = 0.1484375 + 12/128. */
+ 0xf.887604e2c39dbb200p-4L,
+ 0xe.4ec5825059a78a000p-72L,
+ 0x3.d654aff15cb457a00p-4L,
+ 0xf.ca854698aba330400p-72L,
+ /* x = 0.1484375 + 13/128. */
+ 0xf.80aa4fbef750ba800p-4L,
+ -0x7.c2cc346a06b075c00p-68L,
+ 0x3.f55dda9e62aed7500p-4L,
+ 0x1.3bd7b8e6a3d1635e0p-68L,
+ /* x = 0.1484375 + 14/128. */
+ 0xf.78a098069792dab00p-4L,
+ -0x4.3611bda6e483a5980p-68L,
+ 0x4.14572fd94556e6480p-4L,
+ -0xc.29dfd8ec7722b8400p-72L,
+ /* x = 0.1484375 + 15/128. */
+ 0xf.7058fde0788dfc800p-4L,
+ 0x5.b8fe88789e4f42500p-72L,
+ 0x4.334033bcd90d66080p-4L,
+ -0x3.0a0c93e2b47bbae40p-68L,
+ /* x = 0.1484375 + 16/128. */
+ 0xf.67d3a26af7d07aa00p-4L,
+ 0x4.bd6d42af8c0068000p-68L,
+ 0x4.52186aa5377ab2080p-4L,
+ 0x3.bf2524f52e3a06a80p-68L,
+ /* x = 0.1484375 + 17/128. */
+ 0xf.5f10a7bb77d3dfa00p-4L,
+ 0xc.1da8b578427832800p-72L,
+ 0x4.70df5931ae1d94600p-4L,
+ 0x7.6fe0dcff47fe31b80p-72L,
+ /* x = 0.1484375 + 18/128. */
+ 0xf.561030ddd7a789600p-4L,
+ 0xe.a9f4a32c652155500p-72L,
+ 0x4.8f948446abcd6b100p-4L,
+ -0x8.0334eff185e4d9100p-72L,
+ /* x = 0.1484375 + 19/128. */
+ 0xf.4cd261d3e6c15bb00p-4L,
+ 0x3.69c8758630d2ac000p-68L,
+ 0x4.ae37710fad27c8a80p-4L,
+ 0x2.9c4cf96c03519b9c0p-68L,
+ /* x = 0.1484375 + 20/128. */
+ 0xf.43575f94d4f6b2700p-4L,
+ 0x2.f5fb76b14d2a64ac0p-68L,
+ 0x4.ccc7a50127e1de100p-4L,
+ -0x3.494bf3cfd39ae0840p-68L,
+ /* x = 0.1484375 + 21/128. */
+ 0xf.399f500c9e9fd3800p-4L,
+ -0x5.166a8d9c254778900p-68L,
+ 0x4.eb44a5da74f600200p-4L,
+ 0x7.aaa090f0734e28880p-72L,
+ /* x = 0.1484375 + 22/128. */
+ 0xf.2faa5a1b74e82fd00p-4L,
+ 0x6.1fa05f9177380e900p-68L,
+ 0x5.09adf9a7b9a5a0f80p-4L,
+ -0x1.c75705c59f5e66be0p-68L,
+ /* x = 0.1484375 + 23/128. */
+ 0xf.2578a595224dd2e00p-4L,
+ 0x6.bfa2eb2f99cc67500p-68L,
+ 0x5.280326c3cf4818200p-4L,
+ 0x3.ba6bb08eac82c2080p-68L,
+ /* x = 0.1484375 + 24/128. */
+ 0xf.1b0a5b406b526d900p-4L,
+ -0x7.93aa0152372f23380p-68L,
+ 0x5.4643b3da29de9b380p-4L,
+ -0x2.8eaa110f0ccd04c00p-68L,
+ /* x = 0.1484375 + 25/128. */
+ 0xf.105fa4d66b607a600p-4L,
+ 0x7.d44e0427252044380p-68L,
+ 0x5.646f27e8bd65cbe00p-4L,
+ 0x3.a5d61ff0657229100p-68L,
+ /* x = 0.1484375 + 26/128. */
+ 0xf.0578ad01ede708000p-4L,
+ -0x5.c63f6239467b50100p-68L,
+ 0x5.82850a41e1dd46c80p-4L,
+ -0x9.fd15dbb3244403200p-76L,
+ /* x = 0.1484375 + 27/128. */
+ 0xe.fa559f5ec3aec3a00p-4L,
+ 0x4.eb03319278a2d4200p-68L,
+ 0x5.a084e28e35fda2780p-4L,
+ -0x9.202444aace28b3100p-72L,
+ /* x = 0.1484375 + 28/128. */
+ 0xe.eef6a879146af0c00p-4L,
+ -0x6.46a15d15f53f2c200p-72L,
+ 0x5.be6e38ce809554280p-4L,
+ 0x3.c14ee9da0d3648400p-68L,
+ /* x = 0.1484375 + 29/128. */
+ 0xe.e35bf5ccac8905300p-4L,
+ -0x3.26e2248cb2c5b81c0p-68L,
+ 0x5.dc40955d9084f4880p-4L,
+ 0x2.94675a2498de5d840p-68L,
+ /* x = 0.1484375 + 30/128. */
+ 0xe.d785b5c44741b4500p-4L,
+ -0x6.c3a943462cc75eb00p-68L,
+ 0x5.f9fb80f21b5364a00p-4L,
+ -0x3.bcdabf5af1dd3ad00p-68L,
+ /* x = 0.1484375 + 31/128. */
+ 0xe.cb7417b8d4ee3ff00p-4L,
+ -0x3.c8545bf8c55b70e00p-68L,
+ 0x6.179e84a09a5258a80p-4L,
+ -0x3.f164a0531fc1ada00p-68L,
+ /* x = 0.1484375 + 32/128. */
+ 0xe.bf274bf0bda4f6200p-4L,
+ 0x4.47e56a09362679900p-68L,
+ 0x6.352929dd264bd4480p-4L,
+ 0x2.02ea766325d8aa8c0p-68L,
+ /* x = 0.1484375 + 33/128. */
+ 0xe.b29f839f201fd1400p-4L,
+ -0x4.6c8697d86e9587100p-68L,
+ 0x6.529afa7d51b129600p-4L,
+ 0x3.1ec197c0a840a11c0p-68L,
+ /* x = 0.1484375 + 34/128. */
+ 0xe.a5dcf0e30cf03e700p-4L,
+ -0x6.8910f4e13d9aea080p-68L,
+ 0x6.6ff380ba014410a00p-4L,
+ -0x1.c65cdf4f5c05a02a0p-68L,
+ /* x = 0.1484375 + 35/128. */
+ 0xe.98dfc6c6be031e600p-4L,
+ 0xd.d3089cbdd18a75b00p-72L,
+ 0x6.8d324731433279700p-4L,
+ 0x3.bc712bcc4ccddc480p-68L,
+ /* x = 0.1484375 + 36/128. */
+ 0xe.8ba8393eca7821b00p-4L,
+ -0x5.a9c27cb6e49efee80p-68L,
+ 0x6.aa56d8e8249db4e80p-4L,
+ 0x3.60a761fe3f9e559c0p-68L,
+ /* x = 0.1484375 + 37/128. */
+ 0xe.7e367d2956cfb1700p-4L,
+ -0x4.955ee1abe632ffa80p-68L,
+ 0x6.c760c14c8585a5200p-4L,
+ -0x2.42cb99f5193ad5380p-68L,
+ /* x = 0.1484375 + 38/128. */
+ 0xe.708ac84d4172a3e00p-4L,
+ 0x2.737662213429e1400p-68L,
+ 0x6.e44f8c36eb10a1c80p-4L,
+ -0xa.d2f6c3ff0b2b84600p-72L,
+ /* x = 0.1484375 + 39/128. */
+ 0xe.62a551594b970a700p-4L,
+ 0x7.0b15d41d4c0e48400p-68L,
+ 0x7.0122c5ec5028c8d00p-4L,
+ -0xc.c540b02cbf333c800p-76L,
+ /* x = 0.1484375 + 40/128. */
+ 0xe.54864fe33e8575d00p-4L,
+ -0x5.40a42f1a30e4e5780p-68L,
+ 0x7.1dd9fb1ff46778500p-4L,
+ 0x3.acb970a9f6729c700p-68L,
+ /* x = 0.1484375 + 41/128. */
+ 0xe.462dfc670d421ab00p-4L,
+ 0x3.d1a15901228f146c0p-68L,
+ 0x7.3a74b8f52947b6800p-4L,
+ 0x1.baf6928eb3fb02180p-68L,
+ /* x = 0.1484375 + 42/128. */
+ 0xe.379c9045f29d51800p-4L,
+ -0x3.b7f755b683dfa84c0p-68L,
+ 0x7.56f28d011d9852880p-4L,
+ 0x2.44a75fc29c779bd80p-68L,
+ /* x = 0.1484375 + 43/128. */
+ 0xe.28d245c58baef7200p-4L,
+ 0x2.25e232abc003c4380p-68L,
+ 0x7.7353054ca72690d80p-4L,
+ -0x3.391e8e0266194c600p-68L,
+ /* x = 0.1484375 + 44/128. */
+ 0xe.19cf580eeec046b00p-4L,
+ -0x5.ebdd058b7f8131080p-68L,
+ 0x7.8f95b0560a9a3bd80p-4L,
+ -0x1.2084267e23c739ee0p-68L,
+ /* x = 0.1484375 + 45/128. */
+ 0xe.0a94032dbea7cee00p-4L,
+ -0x4.222625d0505267a80p-68L,
+ 0x7.abba1d12c17bfa200p-4L,
+ -0x2.6d0f26c09f2126680p-68L,
+ /* x = 0.1484375 + 46/128. */
+ 0xd.fb20840f3a9b36f00p-4L,
+ 0x7.ae2c515342890b600p-68L,
+ 0x7.c7bfdaf13e5ed1700p-4L,
+ 0x2.12f8a7525bfb113c0p-68L,
+ /* x = 0.1484375 + 47/128. */
+ 0xd.eb7518814a7a93200p-4L,
+ -0x4.433773ef632be3b00p-68L,
+ 0x7.e3a679daaf25c6780p-4L,
+ -0x1.abd434bfd72f69be0p-68L,
+ /* x = 0.1484375 + 48/128. */
+ 0xd.db91ff31879917300p-4L,
+ -0x4.2dbad2f5c7760ae80p-68L,
+ 0x7.ff6d8a34bd5e8fa80p-4L,
+ -0x2.b368b7d24aea62100p-68L,
+ /* x = 0.1484375 + 49/128. */
+ 0xd.cb7777ac420705100p-4L,
+ 0x6.8f31e3eb780ce9c80p-68L,
+ 0x8.1b149ce34caa5a500p-4L,
+ -0x1.9af072f602b295580p-68L,
+ /* x = 0.1484375 + 50/128. */
+ 0xd.bb25c25b8260c1500p-4L,
+ -0x9.1843671366e48f400p-72L,
+ 0x8.369b434a372da7f00p-4L,
+ -0x4.a3758e01c931e1f80p-68L,
+ /* x = 0.1484375 + 51/128. */
+ 0xd.aa9d2086082706400p-4L,
+ -0x2.1ae3f617aa166cd00p-72L,
+ 0x8.52010f4f080052100p-4L,
+ 0x3.78bd8dd614753d080p-68L,
+ /* x = 0.1484375 + 52/128. */
+ 0xd.99ddd44e44a43d500p-4L,
+ -0x2.b5c5c126adfbef900p-68L,
+ 0x8.6d45935ab396cb500p-4L,
+ -0x1.bde17dd211ab0caa0p-68L,
+ /* x = 0.1484375 + 53/128. */
+ 0xd.88e820b1526311e00p-4L,
+ -0x2.a9e1043f3e565ac80p-68L,
+ 0x8.8868625b4e1dbb200p-4L,
+ 0x3.13310133022527200p-68L,
+ /* x = 0.1484375 + 54/128. */
+ 0xd.77bc4985e93a60800p-4L,
+ -0x3.6279746f944394400p-68L,
+ 0x8.a3690fc5bfc11c000p-4L,
+ -0x6.aca1d8c657aed0b80p-68L,
+ /* x = 0.1484375 + 55/128. */
+ 0xd.665a937b4ef2b1f00p-4L,
+ 0x6.d51bad6d988a44180p-68L,
+ 0x8.be472f9776d809b00p-4L,
+ -0xd.477e8edbc29c29900p-72L,
+ /* x = 0.1484375 + 56/128. */
+ 0xd.54c3441844897fd00p-4L,
+ -0x7.07ac0f9aa0e459680p-68L,
+ 0x8.d902565817ee78400p-4L,
+ -0x6.431c32ed7f9fee680p-68L,
+ /* x = 0.1484375 + 57/128. */
+ 0xd.42f6a1b9f0168ce00p-4L,
+ -0xf.ce3d09c3726cfb200p-72L,
+ 0x8.f39a191b2ba612300p-4L,
+ -0x5.c05b0be2a5c002c00p-68L,
+ /* x = 0.1484375 + 58/128. */
+ 0xd.30f4f392c357ab000p-4L,
+ 0x6.61c5fa8a7d9b26600p-68L,
+ 0x9.0e0e0d81ca6787900p-4L,
+ 0x6.cc92c8ea8c2815c00p-68L,
+ /* x = 0.1484375 + 59/128. */
+ 0xd.1ebe81a95ee752e00p-4L,
+ 0x4.8a26bcd32d6e92300p-68L,
+ 0x9.285dc9bc45dd9ea00p-4L,
+ 0x3.d02457bcce59c4180p-68L,
+ /* x = 0.1484375 + 60/128. */
+ 0xd.0c5394d7722281900p-4L,
+ 0x5.e25736c0357470800p-68L,
+ 0x9.4288e48bd0335fc00p-4L,
+ 0x4.1c4cbd2920497a900p-68L,
+ /* x = 0.1484375 + 61/128. */
+ 0xc.f9b476c897c25c600p-4L,
+ -0x4.018af22c0cf715080p-68L,
+ 0x9.5c8ef544210ec0c00p-4L,
+ -0x6.e3b642d55f617ae80p-68L,
+ /* x = 0.1484375 + 62/128. */
+ 0xc.e6e171f92f2e27f00p-4L,
+ 0x3.2225327ec440ddb00p-68L,
+ 0x9.766f93cd18413a700p-4L,
+ -0x5.503e303903d754480p-68L,
+ /* x = 0.1484375 + 63/128. */
+ 0xc.d3dad1b5328a2e400p-4L,
+ 0x5.9f993f4f510881a00p-68L,
+ 0x9.902a58a45e27bed00p-4L,
+ 0x6.8412b426b675ed500p-68L,
+ /* x = 0.1484375 + 64/128. */
+ 0xc.c0a0e21709883a400p-4L,
+ -0xf.f6ee1ee5f811c4300p-76L,
+ 0x9.a9bedcdf01b38da00p-4L,
+ -0x6.c0c287df87e21d700p-68L,
+ /* x = 0.1484375 + 65/128. */
+ 0xc.ad33f00658fe5e800p-4L,
+ 0x2.04bbc0f3a66a0e6c0p-68L,
+ 0x9.c32cba2b14156ef00p-4L,
+ 0x5.256c4f857991ca680p-72L,
+ /* x = 0.1484375 + 66/128. */
+ 0xc.99944936cf48c8900p-4L,
+ 0x1.1ff93fe64b3ddb7a0p-68L,
+ 0x9.dc738ad14204e6900p-4L,
+ -0x6.53a7d2f07a7d9a700p-68L,
+ /* x = 0.1484375 + 67/128. */
+ 0xc.85c23c26ed7b6f000p-4L,
+ 0x1.4ef546c4792968220p-68L,
+ 0x9.f592e9b66a9cf9000p-4L,
+ 0x6.a3c7aa3c101998480p-68L,
+ /* x = 0.1484375 + 68/128. */
+ 0xc.71be181ecd6875d00p-4L,
+ -0x1.d25a9ea5fc335df80p-68L,
+ 0xa.0e8a725d33c828c00p-4L,
+ 0x1.1fa50fd9e9a15ffe0p-68L,
+ /* x = 0.1484375 + 69/128. */
+ 0xc.5d882d2ee48030c00p-4L,
+ 0x7.c07d28e981e348080p-68L,
+ 0xa.2759c0e79c3558200p-4L,
+ 0x5.27c32b55f5405c180p-68L,
+ /* x = 0.1484375 + 70/128. */
+ 0xc.4920cc2ec38fb8900p-4L,
+ 0x1.b38827db08884fc60p-68L,
+ 0xa.400072188acf49d00p-4L,
+ -0x2.94e8c7da1fc7cb900p-68L,
+ /* x = 0.1484375 + 71/128. */
+ 0xc.348846bbd36313400p-4L,
+ -0x7.001d401622ec7e600p-68L,
+ 0xa.587e23555bb080800p-4L,
+ 0x6.d02b9c662cdd29300p-68L,
+ /* x = 0.1484375 + 72/128. */
+ 0xc.1fbeef380e4ffdd00p-4L,
+ 0x5.a613ec8722f644000p-68L,
+ 0xa.70d272a76a8d4b700p-4L,
+ -0x2.5f136f8ed448b7480p-68L,
+ /* x = 0.1484375 + 73/128. */
+ 0xc.0ac518c8b6ae71100p-4L,
+ -0x4.5c85c1146f34ea500p-68L,
+ 0xa.88fcfebd9a8dd4800p-4L,
+ -0x1.d0c3891061dbc66e0p-68L,
+ /* x = 0.1484375 + 74/128. */
+ 0xb.f59b17550a4406800p-4L,
+ 0x7.5969296567cf3e380p-68L,
+ 0xa.a0fd66eddb9212300p-4L,
+ 0x2.c28520d3911b8a040p-68L,
+ /* x = 0.1484375 + 75/128. */
+ 0xb.e0413f84f2a771c00p-4L,
+ 0x6.14946a88cbf4da200p-68L,
+ 0xa.b8d34b36acd987200p-4L,
+ 0x1.0ed343ec65d7e3ae0p-68L,
+ /* x = 0.1484375 + 76/128. */
+ 0xb.cab7e6bfb2a14aa00p-4L,
+ -0x4.edd3a8b5c89413680p-68L,
+ 0xa.d07e4c409d08c5000p-4L,
+ -0x5.c56fa844f53db4780p-68L,
+ /* x = 0.1484375 + 77/128. */
+ 0xb.b4ff632a908f73f00p-4L,
+ -0x3.eae7c6346266c4b00p-68L,
+ 0xa.e7fe0b5fc786b2e00p-4L,
+ -0x6.991e2950ebf5b7780p-68L,
+ /* x = 0.1484375 + 78/128. */
+ 0xb.9f180ba77dd075100p-4L,
+ 0x6.28e135a9508299000p-68L,
+ 0xa.ff522a954f2ba1700p-4L,
+ -0x2.621023be91cc0a180p-68L,
+ /* x = 0.1484375 + 79/128. */
+ 0xb.890237d3bb3c28500p-4L,
+ -0x4.9eb5fac6fe9405f00p-68L,
+ 0xb.167a4c90d63c42400p-4L,
+ 0x4.cf5493b7cc23bd400p-68L,
+ /* x = 0.1484375 + 80/128. */
+ 0xb.72be40067aaf2c000p-4L,
+ 0x5.0dbdb7a14c3d7d500p-68L,
+ 0xb.2d7614b1f3aaa2500p-4L,
+ -0x2.0d291df5881e35c00p-68L,
+ /* x = 0.1484375 + 81/128. */
+ 0xb.5c4c7d4f7dae91600p-4L,
+ -0x5.3879330b4e5b67300p-68L,
+ 0xb.44452709a59752900p-4L,
+ 0x5.913765434a59d1100p-72L,
+ /* x = 0.1484375 + 82/128. */
+ 0xb.45ad4975b1294cb00p-4L,
+ -0x2.35b30bf1370dd5980p-68L,
+ 0xb.5ae7285bc10cf5100p-4L,
+ 0x5.753847e8f8b7a3100p-68L,
+};
diff --git a/sysdeps/x86_64/fpu/e_rem_pio2l.c b/sysdeps/x86_64/fpu/e_rem_pio2l.c
deleted file mode 100644
index 1347b04..0000000
--- a/sysdeps/x86_64/fpu/e_rem_pio2l.c
+++ /dev/null
@@ -1,3 +0,0 @@
-/* Empty. This file is only meant to avoid compiling the file with the
- same name in the libm-ieee754 directory. The code is not used since
- there is an assembler version for all users of this file. */
diff --git a/sysdeps/x86_64/fpu/k_cosl.c b/sysdeps/x86_64/fpu/k_cosl.c
deleted file mode 100644
index eea55a9..0000000
--- a/sysdeps/x86_64/fpu/k_cosl.c
+++ /dev/null
@@ -1 +0,0 @@
-/* Not needed. */
diff --git a/sysdeps/x86_64/fpu/k_sinl.c b/sysdeps/x86_64/fpu/k_sinl.c
deleted file mode 100644
index eea55a9..0000000
--- a/sysdeps/x86_64/fpu/k_sinl.c
+++ /dev/null
@@ -1 +0,0 @@
-/* Not needed. */
diff --git a/sysdeps/x86_64/fpu/s_cosl.S b/sysdeps/x86_64/fpu/s_cosl.S
deleted file mode 100644
index 6921cda..0000000
--- a/sysdeps/x86_64/fpu/s_cosl.S
+++ /dev/null
@@ -1,42 +0,0 @@
-/*
- * Written by J.T. Conklin <jtc@netbsd.org>.
- * Public domain.
- *
- * Adapted for `long double' by Ulrich Drepper <drepper@cygnus.com>.
- * Adapted for x86-64 by Andreas Jaeger <aj@suse.de>.
- * Fixed errno handling by Ulrich Drepper <drepper@redhat.com>.
- */
-
-#define __need_Emath
-#include <bits/errno.h>
-#include <machine/asm.h>
-
-ENTRY(__cosl)
- fldt 8(%rsp)
- fxam
- fstsw %ax
- movb $0x45, %dh
- andb %ah, %dh
- cmpb $0x05, %dh
- je 3f
-4: fcos
- fnstsw %ax
- testl $0x400,%eax
- jnz 1f
- ret
- .align ALIGNARG(4)
-1: fldpi
- fadd %st(0)
- fxch %st(1)
-2: fprem1
- fnstsw %ax
- testl $0x400,%eax
- jnz 2b
- fstp %st(1)
- fcos
- ret
-3: call __errno_location@PLT
- movl $EDOM, (%rax)
- jmp 4b
-END (__cosl)
-weak_alias (__cosl, cosl)
diff --git a/sysdeps/x86_64/fpu/s_sincosl.S b/sysdeps/x86_64/fpu/s_sincosl.S
deleted file mode 100644
index b394c04..0000000
--- a/sysdeps/x86_64/fpu/s_sincosl.S
+++ /dev/null
@@ -1,59 +0,0 @@
-/* Compute sine and cosine of argument.
- Copyright (C) 1997, 2000, 2001, 2005 Free Software Foundation, Inc.
- This file is part of the GNU C Library.
- Contributed by Ulrich Drepper <drepper@cygnus.com>, 1997.
-
- The GNU C Library is free software; you can redistribute it and/or
- modify it under the terms of the GNU Lesser General Public
- License as published by the Free Software Foundation; either
- version 2.1 of the License, or (at your option) any later version.
-
- The GNU C Library is distributed in the hope that it will be useful,
- but WITHOUT ANY WARRANTY; without even the implied warranty of
- MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
- Lesser General Public License for more details.
-
- You should have received a copy of the GNU Lesser General Public
- License along with the GNU C Library; if not, see
- <http://www.gnu.org/licenses/>. */
-
-#include <machine/asm.h>
-#include "bp-sym.h"
-#include "bp-asm.h"
-
-#define PARMS LINKAGE /* no space for saved regs */
-#define ANGLE PARMS
-#define SINP ANGLE+12
-#define COSP SINP+PTR_SIZE
-
- .text
-ENTRY (BP_SYM (__sincosl))
- ENTER
-
- fldt 8(%rsp)
- fsincos
- fnstsw %ax
- testl $0x400,%eax
- jnz 1f
- fstpt (%rsi)
- fstpt (%rdi)
-
- LEAVE
- retq
-
-1: fldpi
- fadd %st(0)
- fxch %st(1)
-2: fprem1
- fnstsw %ax
- testl $0x400,%eax
- jnz 2b
- fstp %st(1)
- fsincos
- fstpt (%rsi)
- fstpt (%rdi)
-
- LEAVE
- retq
-END (BP_SYM (__sincosl))
-weak_alias (BP_SYM (__sincosl), BP_SYM (sincosl))
diff --git a/sysdeps/x86_64/fpu/s_sinl.S b/sysdeps/x86_64/fpu/s_sinl.S
deleted file mode 100644
index 79fc4af..0000000
--- a/sysdeps/x86_64/fpu/s_sinl.S
+++ /dev/null
@@ -1,42 +0,0 @@
-/*
- * Written by J.T. Conklin <jtc@netbsd.org>.
- * Public domain.
- *
- * Adapted for `long double' by Ulrich Drepper <drepper@cygnus.com>.
- * Adapted for x86-64 by Andreas Jaeger <aj@suse.de>.
- * Fixed errno handling by Ulrich Drepper <drepper@redhat.com>.
- */
-
-#define __need_Emath
-#include <bits/errno.h>
-#include <machine/asm.h>
-
-ENTRY(__sinl)
- fldt 8(%rsp)
- fxam
- fstsw %ax
- movb $0x45, %dh
- andb %ah, %dh
- cmpb $0x05, %dh
- je 3f
-4: fsin
- fnstsw %ax
- testl $0x400,%eax
- jnz 1f
- ret
- .align ALIGNARG(4)
-1: fldpi
- fadd %st(0)
- fxch %st(1)
-2: fprem1
- fnstsw %ax
- testl $0x400,%eax
- jnz 2b
- fstp %st(1)
- fsin
- ret
-3: call __errno_location@PLT
- movl $EDOM, (%rax)
- jmp 4b
-END (__sinl)
-weak_alias (__sinl, sinl)
diff --git a/sysdeps/x86_64/fpu/libm-test-ulps b/sysdeps/x86_64/fpu/libm-test-ulps
index 33efe9d..88ae113 100644
--- a/sysdeps/x86_64/fpu/libm-test-ulps
+++ b/sysdeps/x86_64/fpu/libm-test-ulps
@@ -621,12 +621,19 @@ float: 1
ifloat: 1
ildouble: 1
ldouble: 1
+Test "cos_downward (10) == -0.8390715290764524522588639478240648345199":
+ildouble: 1
+ldouble: 1
Test "cos_downward (2) == -0.4161468365471423869975682295007621897660":
float: 1
ifloat: 1
+ildouble: 1
+ldouble: 1
Test "cos_downward (3) == -0.9899924966004454572715727947312613023937":
float: 1
ifloat: 1
+ildouble: 1
+ldouble: 1
Test "cos_downward (4) == -0.6536436208636119146391681830977503814241":
float: 1
ifloat: 1
@@ -688,6 +695,8 @@ ldouble: 1
Test "cos_towardzero (8) == -0.1455000338086135258688413818311946826093":
float: 1
ifloat: 1
+ildouble: 1
+ldouble: 1
# cos_upward
Test "cos_upward (10) == -0.8390715290764524522588639478240648345199":
@@ -701,6 +710,9 @@ ldouble: 1
Test "cos_upward (3) == -0.9899924966004454572715727947312613023937":
ildouble: 1
ldouble: 1
+Test "cos_upward (4) == -0.6536436208636119146391681830977503814241":
+ildouble: 1
+ldouble: 1
Test "cos_upward (5) == 0.2836621854632262644666391715135573083344":
ildouble: 1
ldouble: 1
@@ -712,9 +724,14 @@ ldouble: 1
Test "cos_upward (7) == 0.7539022543433046381411975217191820122183":
float: 1
ifloat: 1
+Test "cos_upward (8) == -0.1455000338086135258688413818311946826093":
+ildouble: 1
+ldouble: 1
Test "cos_upward (9) == -0.9111302618846769883682947111811653112463":
float: 2
ifloat: 2
+ildouble: 1
+ldouble: 1
# cosh_downward
Test "cosh_downward (22) == 1792456423.065795780980053377632656584997":
@@ -799,6 +816,8 @@ double: 1
float: 5
idouble: 1
ifloat: 5
+ildouble: 1
+ldouble: 1
Test "Imaginary part of: cpow (2 + 3 i, 4 + 0 i) == -119.0 - 120.0 i":
float: 2
ifloat: 2
@@ -1059,8 +1078,8 @@ double: 1
float: 1
idouble: 1
ifloat: 1
-ildouble: 1
-ldouble: 1
+ildouble: 2
+ldouble: 2
Test "j0 (0.75) == 0.864242275166648623555731103820923211":
float: 1
ifloat: 1
@@ -1077,8 +1096,8 @@ double: 1
float: 1
idouble: 1
ifloat: 1
-ildouble: 1
-ldouble: 1
+ildouble: 2
+ldouble: 2
Test "j0 (8.0) == 0.171650807137553906090869407851972001":
double: 2
float: 1
@@ -1106,8 +1125,8 @@ double: 1
float: 1
idouble: 1
ifloat: 1
-ildouble: 1
-ldouble: 1
+ildouble: 2
+ldouble: 2
Test "jn (0, 0.75) == 0.864242275166648623555731103820923211":
float: 1
ifloat: 1
@@ -1124,8 +1143,8 @@ double: 1
float: 1
idouble: 1
ifloat: 1
-ildouble: 1
-ldouble: 1
+ildouble: 2
+ldouble: 2
Test "jn (0, 8.0) == 0.171650807137553906090869407851972001":
double: 2
float: 1
@@ -1362,6 +1381,8 @@ ldouble: 1
Test "sin_tonearest (1) == 0.8414709848078965066525023216302989996226":
float: 1
ifloat: 1
+ildouble: 1
+ldouble: 1
Test "sin_tonearest (10) == -0.5440211108893698134047476618513772816836":
ildouble: 1
ldouble: 1
@@ -1387,6 +1408,8 @@ ldouble: 1
Test "sin_towardzero (4) == -0.7568024953079282513726390945118290941359":
float: 1
ifloat: 1
+ildouble: 1
+ldouble: 1
Test "sin_towardzero (5) == -0.9589242746631384688931544061559939733525":
float: 1
ifloat: 1
@@ -1409,14 +1432,27 @@ ldouble: 1
Test "sin_upward (1) == 0.8414709848078965066525023216302989996226":
float: 1
ifloat: 1
+ildouble: 1
+ldouble: 1
+Test "sin_upward (10) == -0.5440211108893698134047476618513772816836":
+ildouble: 1
+ldouble: 1
Test "sin_upward (2) == 0.9092974268256816953960198659117448427023":
float: 2
ifloat: 2
ildouble: 1
ldouble: 1
+Test "sin_upward (3) == 0.1411200080598672221007448028081102798469":
+ildouble: 1
+ldouble: 1
Test "sin_upward (4) == -0.7568024953079282513726390945118290941359":
float: 1
ifloat: 1
+ildouble: 2
+ldouble: 2
+Test "sin_upward (5) == -0.9589242746631384688931544061559939733525":
+ildouble: 1
+ldouble: 1
Test "sin_upward (6) == -0.2794154981989258728115554466118947596280":
ildouble: 1
ldouble: 1
@@ -2162,8 +2198,8 @@ double: 2
float: 2
idouble: 2
ifloat: 2
-ildouble: 1
-ldouble: 1
+ildouble: 2
+ldouble: 2
Function: "j1":
double: 1
@@ -2240,8 +2276,8 @@ ldouble: 1
Function: "sin_upward":
float: 2
ifloat: 2
-ildouble: 1
-ldouble: 1
+ildouble: 2
+ldouble: 2
Function: "sincos":
double: 1
--
Joseph S. Myers
joseph@codesourcery.com#include <stdio.h>
#include <mpfr.h>
int
main (void)
{
int j;
mpfr_t t, b, i, rs, rc, ts, tc, tsl, tcl;
mpfr_init2 (b, 64);
mpfr_init2 (i, 64);
mpfr_init2 (t, 64);
mpfr_set_str (b, "0.1484375", 0, MPFR_RNDN);
mpfr_set_str (i, "0x1p-7", 0, MPFR_RNDN);
mpfr_init2 (rs, 300);
mpfr_init2 (rc, 300);
mpfr_init2 (ts, 64);
mpfr_init2 (tc, 64);
mpfr_init2 (tsl, 64);
mpfr_init2 (tcl, 64);
for (j = 0; j <= 82; j++)
{
mpfr_mul_ui (t, i, j, MPFR_RNDN);
mpfr_add (t, t, b, MPFR_RNDN);
printf (" /* x = 0.1484375 + %d/128. */\n", j);
mpfr_cos (rc, t, MPFR_RNDN);
mpfr_sin (rs, t, MPFR_RNDN);
mpfr_set (tc, rc, MPFR_RNDN);
mpfr_set (ts, rs, MPFR_RNDN);
mpfr_sub (tcl, rc, tc, MPFR_RNDN);
mpfr_sub (tsl, rs, ts, MPFR_RNDN);
mpfr_printf (" %.17RaL,\n", tc);
mpfr_printf (" %.17RaL,\n", tcl);
mpfr_printf (" %.17RaL,\n", ts);
mpfr_printf (" %.17RaL,\n", tsl);
}
return 0;
}