Fixes for bugs in sinpi/cospi/tanpi

patch to fix half-cycle trigonometric functions

Paul Zimmermann, a MPFR developer, contacted me about large errors in
the half-cycle trigonometric functions. I'

Fixes for bugs in sinpi/cospi/tanpi

patch to fix half-cycle trigonometric functions

Paul Zimmermann, a MPFR developer, contacted me about large errors in
the half-cycle trigonometric functions. I've have investigated these
issues and developed the attached patch. The float, double, and ld80
(long double) changes have been tested.

Caveat emptor: The ld128 changes have not been compiled. The ld128
changes have not been tested. I do not have access to a system that
uses ld128 floating point.

Here is an itemized list of changes:

* lib/msun/src/math_private.h:
. Add fast floor macros to compute the integer part of |x| for
0 <= |x| 01xp(N-1), where N is the precision of the type of x.
These macros are used in the half-cycle trigonometric functions
(e.g., sinpi(x)).
. The FFLOOR80 macros is used with the Intel 80-bit extended double
functions. This macors corrects an off-by-one error, which led to
enormous error for |x| > 0x1p32.

* lib/msun/src/s_cospif.c:
* lib/msun/src/s_cospi.c:
* lib/msun/ld80/s_cospil.c:
. Update Copyright years.
. Use FFLOOR*() macro to get integer part of |x|.
. Correct handle the range 0x1p(N-1) <= |x| < 0x1pN. Here, one needs
to determine if the integral value of |x| is even or odd to choose
+1 or -1. If |x| >= 0x1pN, always return +1.

* lib/msun/src/s_sinpif.c:
* lib/msun/src/s_sinpi.c:
* lib/msun/ld80/s_sinpil.c:
. Update Copyright years.
. Use FFLOOR*() macro to get integer part of |x|.

* lib/msun/src/s_tanpif.c:
* lib/msun/src/s_tanpi.c:
* lib/msun/ld80/s_tanpil.c:
. Update Copyright years.
. For +-0.5, return +-inf. Previously, tanpi[fl]() returned an NaN.
. Use FFLOOR*() to get integer part of |x|. Need to determine if the
integer part is even or odd. This is used to set +-0 for |x|
integral
and +-inf for (n+1/2).
. For 0x1p(N-1) <= |x| < 0x1pN need to determine if x is an even or
odd
integer to select +0 or -0. For |x| >= 0x1pN, it is always an even
integer, select 0.
. Note, tanpi[fl](x) is an odd function, so one needs to consider
tanpi[fl](-|x|) = - tanpi[fl](|x|).

* lib/msun/ld128/s_cospil.c:
* lib/msun/ld128/s_sinpil.c:
* lib/msun/ld128/s_tanpil.c:
. Update Copyright years.
. These routines use an FFLOOR128 macros, which likely should be
replaced by a bit twiddling algorithm.
. The same considerations above are applied to 0x1p112 <= |x| <
0x1p113,
and |x| >= 0x1p113 cases.
. Note, even and odd determination used fmodl(x,2.), which is likely
slow.

PR: 272742
MFC after: 1 week

show more ...

sinpi[fl] etc: Fix the ld128 implementations

Mark Murray graciously provided access to an aarch64 system
to test the ld128 implementations. This patch address
* Misuses of copysignl() in sinpil() a

sinpi[fl] etc: Fix the ld128 implementations

Mark Murray graciously provided access to an aarch64 system
to test the ld128 implementations. This patch address
* Misuses of copysignl() in sinpil() and tanpil().
* Redo the splitting of argument 'x' into an integer part and
remainder. The remainder must satify 0 <= r < 1.
* Update the reduction of the integer part to something that can
easily be seen as even or odd, e.g., sin(pi*x) = (-1)^n*sin(pi*r)
with n <= 2^112 and we an reduce n by subtracting integer powers
of 2.
* In s_cospil.c, fix typos where 'x' is used where 'ax', the
remainder, is required.
* In tanpil(), fix the use of an uninitialized variable, ax = fabsl(ax),
ax should be x in fabsl().

One item of note, in the limited tested on aarch64, the max ULP
for sinpil() and cospil() were less than 1.1 ULP, which is higher
that the desired max ULP less than 1. This was traced to the
kernel for cosl() in the fundamental interval [0,pi/4].
The coefficients in the minmax polynomial likely need refinement.

PR: 218514
MFC after: 1 week

show more ...

[LIBM] implementations of sinpi[fl], cospi[fl], and tanpi[fl]

Both IEEE-754 2008 and ISO/IEC TS 18661-4 define the half-cycle
trignometric functions cospi, sinpi, and tanpi. The attached
patch impl

[LIBM] implementations of sinpi[fl], cospi[fl], and tanpi[fl]

Both IEEE-754 2008 and ISO/IEC TS 18661-4 define the half-cycle
trignometric functions cospi, sinpi, and tanpi. The attached
patch implements cospi[fl], sinpi[fl], and tanpi[fl]. Limited
testing on the cospi and sinpi reveal a max ULP less than 0.89;
while tanpi is more problematic with a max ULP less than 2.01
in the interval [0,0.5]. The algorithms used in these functions
are documented in {ks}_cospi.c, {ks}_sinpi.c, and s_tanpi.c.

Note. I no longer have access to a system with ld128 and
adequate support to compile and test the ld128 implementations
of these functions. Given the almost complete lack of input from
others on improvements to libm, I doubt that anyone cares. If
someone does care, the ld128 files contain a number of FIXME comments,
and in particular, while the polynomial coefficients are given
I did not update the polynomial algorithms to properly use the
coefficients.

PR: 218514
MFC after: 2 weeks

show more ...