On Wed, 1 Mar 2023 05:28:34 GMT, Joe Darcy <da...@openjdk.org> wrote:
> Last and certainly not least in the port of FDLIBM to Java, the > transcendental methods for sin, cos, and tan. > > Some more tests are to be written in the StrictMath directory to verify that > the StrictMath algorihtm for sin/cos/tan is being used rather than a > different one. However, I wanted to get the rest of the change out for review > first. > > The sin/cos/tan methods are grouped together since they share the same > argument reduction logic. Argument reduction is the process of mapping an > argument of a function to an argument in a restricted range (and possibly > returning some function of the reduced argument). For sin, cos, and tan, > since they are fundamentally periodic with respect to a multiple of pi, > argument reduction is done to find the remainder of the original argument > with respect to pi/2. And for cos: $ diff -w Cos.c Cos.translit.java 1c1 < /* cos(x) --- > /** cos(x) 31,41c31,34 < < #include "fdlibm.h" < < #ifdef __STDC__ < double cos(double x) < #else < double cos(x) < double x; < #endif < { < double y[2],z=0.0; --- > static class Cos { > static double compute(double x) { > double[] y = new double[2]; > double z=0.0; 56c49 < n = __ieee754_rem_pio2(x,y); --- > n = RemPio2.__ieee754_rem_pio2(x,y); 58,60c51,53 < case 0: return __kernel_cos(y[0],y[1]); < case 1: return -__kernel_sin(y[0],y[1],1); < case 2: return -__kernel_cos(y[0],y[1]); --- > case 0: return Cos.__kernel_cos(y[0],y[1]); > case 1: return -Sin.__kernel_sin(y[0],y[1],1); > case 2: return -Cos.__kernel_cos(y[0],y[1]); 62c55 < return __kernel_sin(y[0],y[1],1); --- > return Sin.__kernel_sin(y[0],y[1],1); 66c59,60 < /* --- > > /** 100,107c94 < < #include "fdlibm.h" < < #ifdef __STDC__ < static const double < #else < static double < #endif --- > private static final double 116,123c103,104 < #ifdef __STDC__ < double __kernel_cos(double x, double y) < #else < double __kernel_cos(x, y) < double x,y; < #endif < { < double a,hz,z,r,qx; --- > static double __kernel_cos(double x, double y) { > double a,hz,z,r,qx = 0.0; 137,138c118,121 < __HI(qx) = ix-0x00200000; /* x/4 */ < __LO(qx) = 0; --- > //__HI(qx) = ix-0x00200000; /* x/4 */ > qx = __HI(qx, ix-0x00200000); > // __LO(qx) = 0; > qx = __LO(qx, 0); 144a128 > } $ diff -w Cos.translit.java Cos.fdlibm.java 31a32,33 > private Cos() {throw new UnsupportedOperationException();} > 37c39 < /* High word of x. */ --- > // High word of x. 40,48c42,48 < /* |x| ~< pi/4 */ < ix &= 0x7fffffff; < if(ix <= 0x3fe921fb) return __kernel_cos(x,z); < < /* cos(Inf or NaN) is NaN */ < else if (ix>=0x7ff00000) return x-x; < < /* argument reduction needed */ < else { --- > // |x| ~< pi/4 > ix &= 0x7fff_ffff; > if (ix <= 0x3fe9_21fb) { > return __kernel_cos(x, z); > } else if (ix >= 0x7ff0_0000) { // cos(Inf or NaN) is NaN > return x-x; > } else { // argument reduction needed 95,101c95,100 < one = 1.00000000000000000000e+00, /* 0x3FF00000, 0x00000000 */ < C1 = 4.16666666666666019037e-02, /* 0x3FA55555, 0x5555554C */ < C2 = -1.38888888888741095749e-03, /* 0xBF56C16C, 0x16C15177 */ < C3 = 2.48015872894767294178e-05, /* 0x3EFA01A0, 0x19CB1590 */ < C4 = -2.75573143513906633035e-07, /* 0xBE927E4F, 0x809C52AD */ < C5 = 2.08757232129817482790e-09, /* 0x3E21EE9E, 0xBDB4B1C4 */ < C6 = -1.13596475577881948265e-11; /* 0xBDA8FAE9, 0xBE8838D4 */ --- > C1 = 0x1.555555555554cp-5, // 4.16666666666666019037e-02 > C2 = -0x1.6c16c16c15177p-10, // -1.38888888888741095749e-03 > C3 = 0x1.a01a019cb159p-16, // 2.48015872894767294178e-05 > C4 = -0x1.27e4f809c52adp-22, // -2.75573143513906633035e-07 > C5 = 0x1.1ee9ebdb4b1c4p-29, // 2.08757232129817482790e-09 > C6 = -0x1.8fae9be8838d4p-37; // -1.13596475577881948265e-11 106,108c105,109 < ix = __HI(x)&0x7fffffff; /* ix = |x|'s high word*/ < if(ix<0x3e400000) { /* if x < 2**27 */ < if(((int)x)==0) return one; /* generate inexact */ --- > ix = __HI(x) & 0x7fff_ffff; // ix = |x|'s high word > if (ix < 0x3e40_0000) { // if x < 2**27 > if (((int)x) == 0) { // generate inexact > return 1.0; > } 112,115c113,116 < if(ix < 0x3FD33333) /* if |x| < 0.3 */ < return one - (0.5*z - (z*r - x*y)); < else { < if(ix > 0x3fe90000) { /* x > 0.78125 */ --- > if (ix < 0x3FD3_3333) { // if |x| < 0.3 > return 1.0 - (0.5*z - (z*r - x*y)); > } else { > if (ix > 0x3fe9_0000) { // x > 0.78125 118,121c119 < //__HI(qx) = ix-0x00200000; /* x/4 */ < qx = __HI(qx, ix-0x00200000); < // __LO(qx) = 0; < qx = __LO(qx, 0); --- > qx = __HI_LO(ix - 0x0020_0000, 0); 124c122 < a = one-qx; --- > a = 1.0 - qx; And tan: $ diff -w Tan.c Tan.translit.java 1c1 < /* tan(x) --- > /** tan(x) 30,40c30,33 < < #include "fdlibm.h" < < #ifdef __STDC__ < double tan(double x) < #else < double tan(x) < double x; < #endif < { < double y[2],z=0.0; --- > static class Tan { > static double compute(double x) { > double[] y= new double[2]; > double z=0.0; 55c48 < n = __ieee754_rem_pio2(x,y); --- > n = RemPio2.__ieee754_rem_pio2(x,y); 60c53,54 < /* __kernel_tan( x, y, k ) --- > > /** __kernel_tan( x, y, k ) 93,99c87 < < #include "fdlibm.h" < #ifdef __STDC__ < static const double < #else < static double < #endif --- > private static final double 119,125c107 < #ifdef __STDC__ < double __kernel_tan(double x, double y, int iy) < #else < double __kernel_tan(x, y, iy) < double x,y; int iy; < #endif < { --- > static double __kernel_tan(double x, double y, int iy) { 133c115 < return one / fabs(x); --- > return one / Math.abs(x); 141c123,124 < __LO(z) = 0; --- > // __LO(z) = 0; > z= __LO(z, 0); 144c127,128 < __LO(t) = 0; --- > //__LO(t) = 0; > t = __LO(t, 0); 179c163,164 < __LO(z) = 0; --- > // __LO(z) = 0; > z = __LO(z, 0); 182c167,168 < __LO(t) = 0; --- > // __LO(t) = 0; > t = __LO(t, 0); 186a173 > } $ diff -w Tan.translit.java Tan.fdlibm.java 30a31,32 > private Tan() {throw new UnsupportedOperationException();} > 36c38 < /* High word of x. */ --- > // High word of x. 39,47c41,47 < /* |x| ~< pi/4 */ < ix &= 0x7fffffff; < if(ix <= 0x3fe921fb) return __kernel_tan(x,z,1); < < /* tan(Inf or NaN) is NaN */ < else if (ix>=0x7ff00000) return x-x; /* NaN */ < < /* argument reduction needed */ < else { --- > // |x| ~< pi/4 > ix &= 0x7fff_ffff; > if (ix <= 0x3fe9_21fb) { > return __kernel_tan(x, z, 1); > } else if (ix >= 0x7ff0_0000) { // tan(Inf or NaN) is NaN > return x-x; // NaN > } else { // argument reduction needed 49,50c49 < return __kernel_tan(y[0],y[1],1-((n&1)<<1)); /* 1 -- n even < -1 -- n odd */ --- > return __kernel_tan(y[0], y[1], 1 - ((n & 1) << 1)); // 1 -- > n even; -1 -- n odd 88,90c87,88 < one = 1.00000000000000000000e+00, /* 0x3FF00000, 0x00000000 */ < pio4 = 7.85398163397448278999e-01, /* 0x3FE921FB, 0x54442D18 */ < pio4lo= 3.06161699786838301793e-17, /* 0x3C81A626, 0x33145C07 */ --- > pio4 = 0x1.921fb54442d18p-1, // 7.85398163397448278999e-01 > pio4lo= 0x1.1a62633145c07p-55, // 3.06161699786838301793e-17 92,104c90,102 < 3.33333333333334091986e-01, /* 0x3FD55555, 0x55555563 */ < 1.33333333333201242699e-01, /* 0x3FC11111, 0x1110FE7A */ < 5.39682539762260521377e-02, /* 0x3FABA1BA, 0x1BB341FE */ < 2.18694882948595424599e-02, /* 0x3F9664F4, 0x8406D637 */ < 8.86323982359930005737e-03, /* 0x3F8226E3, 0xE96E8493 */ < 3.59207910759131235356e-03, /* 0x3F6D6D22, 0xC9560328 */ < 1.45620945432529025516e-03, /* 0x3F57DBC8, 0xFEE08315 */ < 5.88041240820264096874e-04, /* 0x3F4344D8, 0xF2F26501 */ < 2.46463134818469906812e-04, /* 0x3F3026F7, 0x1A8D1068 */ < 7.81794442939557092300e-05, /* 0x3F147E88, 0xA03792A6 */ < 7.14072491382608190305e-05, /* 0x3F12B80F, 0x32F0A7E9 */ < -1.85586374855275456654e-05, /* 0xBEF375CB, 0xDB605373 */ < 2.59073051863633712884e-05, /* 0x3EFB2A70, 0x74BF7AD4 */ --- > 0x1.5555555555563p-2, // 3.33333333333334091986e-01 > 0x1.111111110fe7ap-3, // 1.33333333333201242699e-01 > 0x1.ba1ba1bb341fep-5, // 5.39682539762260521377e-02 > 0x1.664f48406d637p-6, // 2.18694882948595424599e-02 > 0x1.226e3e96e8493p-7, // 8.86323982359930005737e-03 > 0x1.d6d22c9560328p-9, // 3.59207910759131235356e-03 > 0x1.7dbc8fee08315p-10, // 1.45620945432529025516e-03 > 0x1.344d8f2f26501p-11, // 5.88041240820264096874e-04 > 0x1.026f71a8d1068p-12, // 2.46463134818469906812e-04 > 0x1.47e88a03792a6p-14, // 7.81794442939557092300e-05 > 0x1.2b80f32f0a7e9p-14, // 7.14072491382608190305e-05 > -0x1.375cbdb605373p-16, // -1.85586374855275456654e-05 > 0x1.b2a7074bf7ad4p-16, // 2.59073051863633712884e-05 110,117c108,115 < hx = __HI(x); /* high word of x */ < ix = hx&0x7fffffff; /* high word of |x| */ < if(ix<0x3e300000) { /* x < 2**-28 */ < if((int)x==0) { /* generate inexact */ < if (((ix | __LO(x)) | (iy + 1)) == 0) < return one / Math.abs(x); < else { < if (iy == 1) --- > hx = __HI(x); // high word of x > ix = hx&0x7fff_ffff; // high word of |x| > if (ix < 0x3e30_0000) { // x < 2**-28 > if ((int)x == 0) { // generate inexact > if (((ix | __LO(x)) | (iy + 1)) == 0) { > return 1.0 / Math.abs(x); > } else { > if (iy == 1) { 119c117 < else { /* compute -1 / (x+y) carefully */ --- > } else { // compute -1 / (x+y) carefully 123d120 < // __LO(z) = 0; 126,127c123 < t = a = -one / w; < //__LO(t) = 0; --- > t = a = -1.0 / w; 129c125 < s = one + t * z; --- > s = 1.0 + t * z; 135,136c131,135 < if(ix>=0x3FE59428) { /* |x|>=0.6744 */ < if(hx<0) {x = -x; y = -y;} --- > if (ix >= 0x3FE5_9428) { // |x| >= 0.6744 > if ( hx < 0) { > x = -x; > y = -y; > } 139c138,139 < x = z+w; y = 0.0; --- > x = z + w; > y = 0.0; 153c153 < if(ix>=0x3FE59428) { --- > if (ix >= 0x3FE5_9428) { 157,160c157,161 < if(iy==1) return w; < else { /* if allow error up to 2 ulp, < simply return -1.0/(x+r) here */ < /* compute -1.0/(x+r) accurately */ --- > if (iy == 1) { > return w; > } else { /* if were to allow error up to 2 ulp, > could simply return -1.0/(x + r) here */ > // compute -1.0/(x + r) accurately 163d163 < // __LO(z) = 0; 165,167c165,166 < v = r-(z - x); /* z+v = r+x */ < t = a = -1.0/w; /* a = -1.0/w */ < // __LO(t) = 0; --- > v = r - (z - x); // z + v = r + x > t = a = -1.0/w; // a = -1.0/w ------------- PR: https://git.openjdk.org/jdk/pull/12800