8304028: Port fdlibm IEEEremainder to Java

Reviewed-by: bpb
This commit is contained in:
Joe Darcy 2023-03-31 19:48:03 +00:00
parent a565be4dc5
commit abfb900829
7 changed files with 716 additions and 129 deletions

View File

@ -3301,4 +3301,196 @@ class FdLibm {
return (jx >= 0)? z: -z;
}
}
static final class IEEEremainder {
private IEEEremainder() {throw new UnsupportedOperationException();}
static double compute(double x, double p) {
int hx, hp;
/*unsigned*/ int sx, lx, lp;
double p_half;
hx = __HI(x); // high word of x
lx = __LO(x); // low word of x
hp = __HI(p); // high word of p
lp = __LO(p); // low word of p
sx = hx & 0x8000_0000;
hp &= 0x7fff_ffff;
hx &= 0x7fff_ffff;
// purge off exception values
if ((hp | lp) == 0) {// p = 0
return (x*p)/(x*p);
}
if ((hx >= 0x7ff0_0000) || // not finite
((hp >= 0x7ff0_0000) && // p is NaN
(((hp - 0x7ff0_0000) | lp) != 0)))
return (x*p)/(x*p);
if (hp <= 0x7fdf_ffff) { // now x < 2p
x = __ieee754_fmod(x, p + p);
}
if (((hx - hp) | (lx - lp)) == 0) {
return 0.0*x;
}
x = Math.abs(x);
p = Math.abs(p);
if (hp < 0x0020_0000) {
if (x + x > p) {
x -= p;
if (x + x >= p) {
x -= p;
}
}
} else {
p_half = 0.5*p;
if (x > p_half) {
x -= p;
if (x >= p_half) {
x -= p;
}
}
}
return __HI(x, __HI(x)^sx);
}
private static double __ieee754_fmod(double x, double y) {
int n, hx, hy, hz, ix, iy, sx;
/*unsigned*/ int lx, ly, lz;
hx = __HI(x); // high word of x
lx = __LO(x); // low word of x
hy = __HI(y); // high word of y
ly = __LO(y); // low word of y
sx = hx & 0x8000_0000; // sign of x
hx ^= sx; // |x|
hy &= 0x7fff_ffff; // |y|
// purge off exception values
if ((hy | ly) == 0 || (hx >= 0x7ff0_0000)|| // y = 0, or x not finite
((hy | ((ly | -ly) >>> 31)) > 0x7ff0_0000)) // or y is NaN, unsigned shift
return (x*y)/(x*y);
if (hx <= hy) {
if ((hx < hy) || (Integer.compareUnsigned(lx, ly) < 0)) { // |x| < |y| return x
return x;
}
if (lx == ly) {
return signedZero(sx); // |x| = |y| return x*0
}
}
ix = ilogb(hx, lx);
iy = ilogb(hy, ly);
// set up {hx, lx}, {hy, ly} and align y to x
if (ix >= -1022)
hx = 0x0010_0000 | (0x000f_ffff & hx);
else { // subnormal x, shift x to normal
n = -1022 - ix;
if (n <= 31) {
hx = (hx << n) | (lx >>> (32 - n)); // unsigned shift
lx <<= n;
} else {
hx = lx << (n - 32);
lx = 0;
}
}
if (iy >= -1022)
hy = 0x0010_0000 | (0x000f_ffff & hy);
else { // subnormal y, shift y to normal
n = -1022 - iy;
if (n <= 31) {
hy = (hy << n)|(ly >>> (32 - n)); // unsigned shift
ly <<= n;
} else {
hy = ly << (n - 32);
ly = 0;
}
}
// fix point fmod
n = ix - iy;
while (n-- != 0) {
hz = hx - hy;
lz = lx - ly;
if (Integer.compareUnsigned(lx, ly) < 0) {
hz -= 1;
}
if (hz < 0){
hx = hx + hx +(lx >>> 31); // unsigned shift
lx = lx + lx;
} else {
if ((hz | lz) == 0) { // return sign(x)*0
return signedZero(sx);
}
hx = hz + hz + (lz >>> 31); // unsigned shift
lx = lz + lz;
}
}
hz = hx - hy;
lz = lx - ly;
if (Integer.compareUnsigned(lx, ly) < 0) {
hz -= 1;
}
if (hz >= 0) {
hx = hz;
lx = lz;
}
// convert back to floating value and restore the sign
if ((hx | lx) == 0) { // return sign(x)*0
return signedZero(sx);
}
while (hx < 0x0010_0000) { // normalize x
hx = hx + hx + (lx >>> 31); // unsigned shift
lx = lx + lx;
iy -= 1;
}
if (iy >= -1022) { // normalize output
hx = ((hx - 0x0010_0000) | ((iy + 1023) << 20));
x = __HI_LO(hx | sx, lx);
} else { // subnormal output
n = -1022 - iy;
if (n <= 20) {
lx = (lx >>> n)|(/*(unsigned)*/hx << (32 - n)); // unsigned shift
hx >>= n;
} else if (n <= 31) {
lx = (hx << (32 - n))|(lx >>> n); // unsigned shift
hx = sx;
} else {
lx = hx >>(n - 32);
hx = sx;
}
x = __HI_LO(hx | sx, lx);
x *= 1.0; // create necessary signal
}
return x; // exact output
}
/*
* Return a double zero with the same sign as the int argument.
*/
private static double signedZero(int sign) {
return +0.0 * ( (double)sign);
}
private static int ilogb(int hz, int lz) {
int iz, i;
if (hz < 0x0010_0000) { // subnormal z
if (hz == 0) {
for (iz = -1043, i = lz; i > 0; i <<= 1) {
iz -= 1;
}
} else {
for (iz = -1022, i = (hz << 11); i > 0; i <<= 1) {
iz -= 1;
}
}
} else {
iz = (hz >> 20) - 1023;
}
return iz;
}
}
}

View File

@ -369,7 +369,9 @@ public final class StrictMath {
* @return the remainder when {@code f1} is divided by
* {@code f2}.
*/
public static native double IEEEremainder(double f1, double f2);
public static double IEEEremainder(double f1, double f2) {
return FdLibm.IEEEremainder.compute(f1, f2);
}
/**
* Returns the smallest (closest to negative infinity)

View File

@ -1,127 +0,0 @@
/*
* Copyright (c) 1994, 2016, Oracle and/or its affiliates. All rights reserved.
* DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER.
*
* This code is free software; you can redistribute it and/or modify it
* under the terms of the GNU General Public License version 2 only, as
* published by the Free Software Foundation. Oracle designates this
* particular file as subject to the "Classpath" exception as provided
* by Oracle in the LICENSE file that accompanied this code.
*
* This code 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 General Public License
* version 2 for more details (a copy is included in the LICENSE file that
* accompanied this code).
*
* You should have received a copy of the GNU General Public License version
* 2 along with this work; if not, write to the Free Software Foundation,
* Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
*
* Please contact Oracle, 500 Oracle Parkway, Redwood Shores, CA 94065 USA
* or visit www.oracle.com if you need additional information or have any
* questions.
*/
#include "jni.h"
#include "fdlibm.h"
#include "java_lang_StrictMath.h"
JNIEXPORT jdouble JNICALL
Java_java_lang_StrictMath_cos(JNIEnv *env, jclass unused, jdouble d)
{
return (jdouble) jcos((double)d);
}
JNIEXPORT jdouble JNICALL
Java_java_lang_StrictMath_sin(JNIEnv *env, jclass unused, jdouble d)
{
return (jdouble) jsin((double)d);
}
JNIEXPORT jdouble JNICALL
Java_java_lang_StrictMath_tan(JNIEnv *env, jclass unused, jdouble d)
{
return (jdouble) jtan((double)d);
}
JNIEXPORT jdouble JNICALL
Java_java_lang_StrictMath_asin(JNIEnv *env, jclass unused, jdouble d)
{
return (jdouble) jasin((double)d);
}
JNIEXPORT jdouble JNICALL
Java_java_lang_StrictMath_acos(JNIEnv *env, jclass unused, jdouble d)
{
return (jdouble) jacos((double)d);
}
JNIEXPORT jdouble JNICALL
Java_java_lang_StrictMath_atan(JNIEnv *env, jclass unused, jdouble d)
{
return (jdouble) jatan((double)d);
}
JNIEXPORT jdouble JNICALL
Java_java_lang_StrictMath_log(JNIEnv *env, jclass unused, jdouble d)
{
return (jdouble) jlog((double)d);
}
JNIEXPORT jdouble JNICALL
Java_java_lang_StrictMath_log10(JNIEnv *env, jclass unused, jdouble d)
{
return (jdouble) jlog10((double)d);
}
JNIEXPORT jdouble JNICALL
Java_java_lang_StrictMath_sqrt(JNIEnv *env, jclass unused, jdouble d)
{
return (jdouble) jsqrt((double)d);
}
JNIEXPORT jdouble JNICALL
Java_java_lang_StrictMath_atan2(JNIEnv *env, jclass unused, jdouble d1, jdouble d2)
{
return (jdouble) jatan2((double)d1, (double)d2);
}
JNIEXPORT jdouble JNICALL
Java_java_lang_StrictMath_IEEEremainder(JNIEnv *env, jclass unused,
jdouble dividend,
jdouble divisor)
{
return (jdouble) jremainder(dividend, divisor);
}
JNIEXPORT jdouble JNICALL
Java_java_lang_StrictMath_cosh(JNIEnv *env, jclass unused, jdouble d)
{
return (jdouble) jcosh((double)d);
}
JNIEXPORT jdouble JNICALL
Java_java_lang_StrictMath_sinh(JNIEnv *env, jclass unused, jdouble d)
{
return (jdouble) jsinh((double)d);
}
JNIEXPORT jdouble JNICALL
Java_java_lang_StrictMath_tanh(JNIEnv *env, jclass unused, jdouble d)
{
return (jdouble) jtanh((double)d);
}
JNIEXPORT jdouble JNICALL
Java_java_lang_StrictMath_log1p(JNIEnv *env, jclass unused, jdouble d)
{
return (jdouble) jlog1p((double)d);
}
JNIEXPORT jdouble JNICALL
Java_java_lang_StrictMath_expm1(JNIEnv *env, jclass unused, jdouble d)
{
return (jdouble) jexpm1((double)d);
}

View File

@ -0,0 +1,209 @@
/*
* Copyright (c) 2004, 2023, Oracle and/or its affiliates. All rights reserved.
* DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER.
*
* This code is free software; you can redistribute it and/or modify it
* under the terms of the GNU General Public License version 2 only, as
* published by the Free Software Foundation.
*
* This code 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 General Public License
* version 2 for more details (a copy is included in the LICENSE file that
* accompanied this code).
*
* You should have received a copy of the GNU General Public License version
* 2 along with this work; if not, write to the Free Software Foundation,
* Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
*
* Please contact Oracle, 500 Oracle Parkway, Redwood Shores, CA 94065 USA
* or visit www.oracle.com if you need additional information or have any
* questions.
*/
/*
* @test
* @bug 8304028
* @summary Tests for {Math, StrictMath}.IEEEremainder
*/
public class IeeeRemainderTests {
private IeeeRemainderTests(){}
public static void main(String... args) {
int failures = 0;
failures += testIeeeRemainderSpecials();
failures += testIeeeRemainderZeroResult();
failures += testIeeeRemainderOneResult();
failures += testIeeeRemainderRounding();
if (failures > 0) {
System.err.println("Testing IEEEremainder incurred "
+ failures + " failures.");
throw new RuntimeException();
}
}
private static final double NaNd = Double.NaN;
private static final double MIN_VALUE = Double.MIN_VALUE;
private static final double MIN_NORM = Double.MIN_NORMAL;
private static final double MAX_VALUE = Double.MAX_VALUE;
private static final double InfinityD = Double.POSITIVE_INFINITY;
/**
* Special cases from the spec interspersed with test cases.
*/
private static int testIeeeRemainderSpecials() {
int failures = 0;
/*
* If either argument is NaN, or the first argument is
* infinite, or the second argument is positive zero or
* negative zero, then the result is NaN.
*
*/
for(double nan : Tests.NaNs) {
failures += testIEEEremainderCase(nan, 1.0, NaNd);
failures += testIEEEremainderCase(1.0, nan, NaNd);
}
double [][] nanResultCases = {
{ InfinityD, InfinityD},
{-InfinityD, InfinityD},
{ InfinityD, 1.0},
{-InfinityD, 1.0},
{ InfinityD, NaNd},
{-InfinityD, NaNd},
{ InfinityD, 0.0},
{-InfinityD, 0.0},
{ InfinityD, -0.0},
{-InfinityD, -0.0},
};
for(double[] testCase : nanResultCases) {
failures += testIEEEremainderCase(testCase[0], testCase[1], NaNd);
}
/*
* If the first argument is finite and the second argument is
* infinite, then the result is the same as the first
* argument.
*
*/
double [] specialCases = {
+0.0,
+MIN_VALUE,
+MIN_NORM,
+MAX_VALUE,
-0.0,
-MIN_VALUE,
-MIN_NORM,
-MAX_VALUE,
};
double [] infinities = {
+InfinityD,
-InfinityD
};
for (double specialCase : specialCases) {
for (double infinity: infinities) {
failures += testIEEEremainderCase(specialCase, infinity, specialCase);
}
}
return failures;
}
private static int testIeeeRemainderZeroResult() {
int failures = 0;
double [] testCases = {
+MIN_VALUE,
+MIN_NORM,
+MAX_VALUE*0.5,
-MIN_VALUE,
-MIN_NORM,
-MAX_VALUE*0.5,
};
for (double testCase : testCases) {
/*
* "If the remainder is zero, its sign is the same as the sign of the first argument."
*/
failures += testIEEEremainderCase(testCase*2.0, +testCase, Math.copySign(0.0, testCase));
failures += testIEEEremainderCase(testCase*2.0, -testCase, Math.copySign(0.0, testCase));
}
return failures;
}
/*
* Construct test cases where the remainder is one.
*/
private static int testIeeeRemainderOneResult() {
int failures = 0;
double [][] testCases = {
{4.0, 3.0},
{10_001.0, 5000.0},
{15_001.0, 5000.0},
{10_000.0, 9999.0},
{0x1.0p52 + 1.0, 0x1.0p52},
{0x1.fffffffffffffp52, 0x1.ffffffffffffep52},
};
for (var testCase : testCases) {
failures += testIEEEremainderCase(testCase[0], testCase[1], 1.0);
}
return failures;
}
/*
* Test cases that differ in rounding between % and IEEEremainder.
*/
private static int testIeeeRemainderRounding() {
int failures = 0;
double [][] testCases = {
{3.0, 2.0, -1.0},
{3.0, -2.0, -1.0},
};
for (var testCase : testCases) {
failures += testIEEEremainderCase(testCase[0], testCase[1], testCase[2]);
}
return failures;
}
/*
* For exact cases, built-in % remainder and IEEE remainder should
* be the same since the rounding mode in the implicit divide
* doesn't come into play.
*/
private static double remainder(double a, double b) {
return a % b;
}
private static int testIEEEremainderCase(double input1, double input2, double expected) {
int failures = 0;
failures += Tests.test("StrictMath.IEEEremainder", input1, input2, StrictMath::IEEEremainder, expected);
failures += Tests.test("Math.IEEEremainder", input1, input2, Math::IEEEremainder, expected);
return failures;
}
}

View File

@ -23,7 +23,7 @@
/*
* @test
* @bug 8301833 8302026 8301444 8302028 8302040 8302027
* @bug 8301833 8302026 8301444 8302028 8302040 8302027 8304028
* @build Tests
* @build FdlibmTranslit
* @build ExhaustingTests
@ -135,6 +135,7 @@ public class ExhaustingTests {
BinaryTestCase[] testCases = {
new BinaryTestCase("hypot", FdlibmTranslit::hypot, StrictMath::hypot, 20, 20),
new BinaryTestCase("atan2", FdlibmTranslit::atan2, StrictMath::atan2, 20, 20),
new BinaryTestCase("IEEEremainder", FdlibmTranslit::IEEEremainder, StrictMath::IEEEremainder, 20, 20),
};
for (var testCase : testCases) {

View File

@ -142,6 +142,10 @@ public class FdlibmTranslit {
return Tanh.compute(x);
}
public static double IEEEremainder(double f1, double f2) {
return IEEEremainder.compute(f1, f2);
}
// -----------------------------------------------------------------------------------------
/** sin(x)
@ -2583,4 +2587,171 @@ public class FdlibmTranslit {
return (jx>=0)? z: -z;
}
}
private static final class IEEEremainder {
private static final double zero = 0.0;
private static double one = 1.0;
private static double[] Zero = {0.0, -0.0,};
static double compute(double x, double p) {
int hx,hp;
/*unsigned*/ int sx,lx,lp;
double p_half;
hx = __HI(x); /* high word of x */
lx = __LO(x); /* low word of x */
hp = __HI(p); /* high word of p */
lp = __LO(p); /* low word of p */
sx = hx&0x80000000;
hp &= 0x7fffffff;
hx &= 0x7fffffff;
/* purge off exception values */
if((hp|lp)==0) return (x*p)/(x*p); /* p = 0 */
if((hx>=0x7ff00000)|| /* x not finite */
((hp>=0x7ff00000)&& /* p is NaN */
(((hp-0x7ff00000)|lp)!=0)))
return (x*p)/(x*p);
if (hp<=0x7fdfffff) x = __ieee754_fmod(x,p+p); /* now x < 2p */
if (((hx-hp)|(lx-lp))==0) return zero*x;
x = Math.abs(x);
p = Math.abs(p);
if (hp<0x00200000) {
if(x+x>p) {
x-=p;
if(x+x>=p) x -= p;
}
} else {
p_half = 0.5*p;
if(x>p_half) {
x-=p;
if(x>=p_half) x -= p;
}
}
// __HI(x) ^= sx;
x = __HI(x, __HI(x)^sx);
return x;
}
private static double __ieee754_fmod(double x, double y) {
int n,hx,hy,hz,ix,iy,sx,i;
/*unsigned*/ int lx,ly,lz;
hx = __HI(x); /* high word of x */
lx = __LO(x); /* low word of x */
hy = __HI(y); /* high word of y */
ly = __LO(y); /* low word of y */
sx = hx&0x80000000; /* sign of x */
hx ^=sx; /* |x| */
hy &= 0x7fffffff; /* |y| */
/* purge off exception values */
if((hy|ly)==0||(hx>=0x7ff00000)|| /* y=0,or x not finite */
((hy|((ly|-ly)>>>31))>0x7ff00000)) /* or y is NaN */ // unsigned shift
return (x*y)/(x*y);
if(hx<=hy) {
// if((hx<hy)||(lx<ly)) return x; /* |x|<|y| return x */
if((hx<hy)||(Integer.compareUnsigned(lx,ly) < 0)) return x; /* |x|<|y| return x */
if(lx==ly)
return Zero[/*(unsigned)*/sx>>>31]; /* |x|=|y| return x*0*/ // unsigned shift
}
/* determine ix = ilogb(x) */
if(hx<0x00100000) { /* subnormal x */
if(hx==0) {
for (ix = -1043, i=lx; i>0; i<<=1) ix -=1;
} else {
for (ix = -1022,i=(hx<<11); i>0; i<<=1) ix -=1;
}
} else ix = (hx>>20)-1023;
/* determine iy = ilogb(y) */
if(hy<0x00100000) { /* subnormal y */
if(hy==0) {
for (iy = -1043, i=ly; i>0; i<<=1) iy -=1;
} else {
for (iy = -1022,i=(hy<<11); i>0; i<<=1) iy -=1;
}
} else iy = (hy>>20)-1023;
/* set up {hx,lx}, {hy,ly} and align y to x */
if(ix >= -1022)
hx = 0x00100000|(0x000fffff&hx);
else { /* subnormal x, shift x to normal */
n = -1022-ix;
if(n<=31) {
hx = (hx<<n)|(lx >>> (32-n)); // unsigned shift
lx <<= n;
} else {
hx = lx<<(n-32);
lx = 0;
}
}
if(iy >= -1022)
hy = 0x00100000|(0x000fffff&hy);
else { /* subnormal y, shift y to normal */
n = -1022-iy;
if(n<=31) {
hy = (hy<<n)|(ly >>> (32-n)); // unsigned shift
ly <<= n;
} else {
hy = ly<<(n-32);
ly = 0;
}
}
/* fix point fmod */
n = ix - iy;
while(n-- != 0) {
hz=hx-hy;lz=lx-ly;
// if(lx<ly) hz -= 1;
if(Integer.compareUnsigned(lx, ly) < 0) hz -= 1;
if(hz<0){hx = hx+hx+(lx >>> 31); lx = lx+lx;} // unsigned shift
else {
if((hz|lz)==0) /* return sign(x)*0 */
return Zero[/*(unsigned)*/sx>>>31]; // unsigned shift
hx = hz+hz+(lz >>> 31); // unsigned shift
lx = lz+lz;
}
}
hz=hx-hy;lz=lx-ly;
// if(lx<ly) hz -= 1;
if(Integer.compareUnsigned(lx, ly) < 0) hz -= 1;
if(hz>=0) {hx=hz;lx=lz;}
/* convert back to floating value and restore the sign */
if((hx|lx)==0) /* return sign(x)*0 */
return Zero[/*(unsigned)*/sx >>> 31]; // unsigned shift
while(hx<0x00100000) { /* normalize x */
hx = hx+hx+(lx >>> 31); lx = lx+lx; // unsigned shift
iy -= 1;
}
if(iy>= -1022) { /* normalize output */
hx = ((hx-0x00100000)|((iy+1023)<<20));
// __HI(x) = hx|sx;
x = __HI(x, hx|sx);
// __LO(x) = lx;
x = __LO(x, lx);
} else { /* subnormal output */
n = -1022 - iy;
if(n<=20) {
lx = (lx >>> n)|(/*(unsigned)*/hx<<(32-n)); // unsigned shift
hx >>= n;
} else if (n<=31) {
lx = (hx<<(32-n))|(lx >>> n); // unsigned shift
hx = sx;
} else {
lx = hx>>(n-32); hx = sx;
}
// __HI(x) = hx|sx;
x = __HI(x, hx|sx);
// __LO(x) = lx;
x = __LO(x, lx);
x *= one; /* create necessary signal */
}
return x; /* exact output */
}
}
}

View File

@ -0,0 +1,139 @@
/*
* Copyright (c) 2003, 2023, Oracle and/or its affiliates. All rights reserved.
* DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER.
*
* This code is free software; you can redistribute it and/or modify it
* under the terms of the GNU General Public License version 2 only, as
* published by the Free Software Foundation.
*
* This code 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 General Public License
* version 2 for more details (a copy is included in the LICENSE file that
* accompanied this code).
*
* You should have received a copy of the GNU General Public License version
* 2 along with this work; if not, write to the Free Software Foundation,
* Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.
*
* Please contact Oracle, 500 Oracle Parkway, Redwood Shores, CA 94065 USA
* or visit www.oracle.com if you need additional information or have any
* questions.
*/
/*
* @test
* @bug 8304028
* @key randomness
* @summary Tests for StrictMath.IEEEremainder
* @library /test/lib
* @build jdk.test.lib.RandomFactory
* @build Tests
* @build FdlibmTranslit
* @build IeeeRemainderTests
* @run main IeeeRemainderTests
*/
import jdk.test.lib.RandomFactory;
/**
* The tests in ../Math/IeeeRemainderTests.java test properties that
* should hold for any IEEEremainder implementation, including the
* FDLIBM-based one required for StrictMath.IEEEremainder. Therefore,
* the test cases in ../Math/IEEEremainderTests.java are run against
* both the Math and StrictMath versions of IEEEremainder. The role
* of this test is to verify that the FDLIBM IEEEremainder algorithm
* is being used by running golden file tests on values that may vary
* from one conforming IEEEremainder implementation to another.
*/
public class IeeeRemainderTests {
private IeeeRemainderTests(){}
public static void main(String... args) {
int failures = 0;
failures += testAgainstTranslit();
if (failures > 0) {
System.err.println("Testing IEEEremainder incurred "
+ failures + " failures.");
throw new RuntimeException();
}
}
// Initialize shared random number generator
private static java.util.Random random = RandomFactory.getRandom();
/**
* Test StrictMath.IEEEremainder against transliteration port of IEEEremainder.
*/
private static int testAgainstTranslit() {
int failures = 0;
// The exact special cases for infinity, NaN, zero,
// etc. inputs are checked in the Math tests.
// Test exotic NaN bit patterns
double[][] exoticNaNs = {
{Double.longBitsToDouble(0x7FF0_0000_0000_0001L), 0.0},
{0.0, Double.longBitsToDouble(0x7FF0_0000_0000_0001L)},
{Double.longBitsToDouble(0xFFF_00000_0000_0001L), 0.0},
{0.0, Double.longBitsToDouble(0xFFF0_0000_0000_0001L)},
{Double.longBitsToDouble(0x7FF_00000_7FFF_FFFFL), 0.0},
{0.0, Double.longBitsToDouble(0x7FF0_7FFF_0000_FFFFL)},
{Double.longBitsToDouble(0xFFF_00000_7FFF_FFFFL), 0.0},
{0.0, Double.longBitsToDouble(0xFFF0_7FFF_0000_FFFFL)},
};
for (double[] exoticNaN: exoticNaNs) {
failures += testIEEEremainderCase(exoticNaN[0], exoticNaN[1],
FdlibmTranslit.IEEEremainder(exoticNaN[0], exoticNaN[1]));
}
// Probe near decision points in the FDLIBM algorithm.
double[][] decisionPoints = {
{0x1.fffffp1022, 100.0},
{0x1.fffffp1022, 0x1.fffffp1022},
{2.0*0x1.0p-1022, 0x1.0p-1022},
{2.0*0x1.0p-1022, 0x1.0p-1023},
};
for (var decisionPoint : decisionPoints) {
double x = decisionPoint[0];
double p = decisionPoint[1];
double increment_x = Math.ulp(x);
double increment_p = Math.ulp(p);
x = x - 64*increment_x;
p = p - 64*increment_p;
for (int i = 0; i < 128; i++, x += increment_x) {
for (int j = 0; j < 126; j++, p += increment_p) {
failures += testIEEEremainderCase( x, p, FdlibmTranslit.IEEEremainder( x, p));
failures += testIEEEremainderCase(-x, p, FdlibmTranslit.IEEEremainder(-x, p));
failures += testIEEEremainderCase( x, -p, FdlibmTranslit.IEEEremainder( x, -p));
failures += testIEEEremainderCase(-x, -p, FdlibmTranslit.IEEEremainder(-x, -p));
}
}
}
// Check random values
for (int k = 0; k < 200; k++ ) {
double x = random.nextDouble();
double p = random.nextDouble();
failures += testIEEEremainderCase(x, p, FdlibmTranslit.IEEEremainder(x, p));
}
return failures;
}
private static int testIEEEremainderCase(double input1, double input2, double expected) {
int failures = 0;
failures += Tests.test("StrictMath.IEEEremainder(double)", input1, input2,
StrictMath::IEEEremainder, expected);
return failures;
}
}