[newlib-cygwin] Fix error in float trig. function range reduction

[newlib-cygwin] Fix error in float trig. function range reduction

[Corinna Vinschen]

https://sourceware.org/git/gitweb.cgi?p=newlib-cygwin.git;h=d4bcecb3e91f73091195adbfdcb0d4287cdf1753

commit d4bcecb3e91f73091195adbfdcb0d4287cdf1753
Author: Fabian Schriever <fabian.schriever@gtd-gmbh.de>
Date:   Tue Mar 3 14:49:06 2020 +0100

    Fix error in float trig. function range reduction
    
    The single-precision trigonometric functions show rather high errors in
    specific ranges starting at about 30000 radians. For example the sinf
    procedure produces an error of 7626.55 ULP with the input
    5.195880078125e+04 (0x474AF6CD) (compared with MPFR in 128bit
    precision). For the test we used 100k values evenly spaced in the range
    of [30k, 70k]. The issues are periodic at higher ranges.
    
    This error was introduced when the double precision range reduction was
    first converted to float. The shift by 8 bits always returns 0 as iq is
    never higher than 255.
    
    The fix reduces the error of the example above to 0.45 ULP, highest
    error within the test set fell to 1.31 ULP, which is not perfect, but
    still a significant improvement. Testing other previously erroneous
    ranges no longer show particularly large accuracy errors.

Diff:
---
 newlib/libm/math/kf_rem_pio2.c | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/newlib/libm/math/kf_rem_pio2.c b/newlib/libm/math/kf_rem_pio2.c
index 1573ca9..3ca7cc2 100644
--- a/newlib/libm/math/kf_rem_pio2.c
+++ b/newlib/libm/math/kf_rem_pio2.c
@@ -101,7 +101,7 @@ recompute:
 	    iq[jz-1] -= i<<(8-q0);
 	    ih = iq[jz-1]>>(7-q0);
 	} 
-	else if(q0==0) ih = iq[jz-1]>>8;
+	else if(q0==0) ih = iq[jz-1]>>7;
 	else if(z>=(float)0.5) ih=2;
 
 	if(ih>0) {	/* q > 0.5 */