Apache
/
commons-math
mirror of https://github.com/apache/commons-math.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

MATH-905 

Avoid overflow on the whole range covered by the equivalent functions in
the standard "Math" class.


git-svn-id: https://svn.apache.org/repos/asf/commons/proper/math/trunk@1413600 13f79535-47bb-0310-9956-ffa450edef68
This commit is contained in:
Gilles Sadowski 2012-11-26 13:25:27 +00:00
parent 9e75d6faa9
commit 46ed0af1c2
2 changed files with 79 additions and 13 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

48
src/main/java/org/apache/commons/math3/util/FastMath.java
View File

@ -78,6 +78,8 @@ import java.io.PrintStream;
* @since 2.2
*/
public class FastMath {
/** StrictMath.log(Double.MAX_VALUE): {@value} */
private static final double LOG_MAX_VALUE = StrictMath.log(Double.MAX_VALUE);
/** Archimede's constant PI, ratio of circle circumference to diameter. */
public static final double PI = 105414357.0 / 33554432.0 + 1.984187159361080883e-9;
@ -389,15 +391,25 @@ public class FastMath {
// for numbers with magnitude 20 or so,
// exp(-z) can be ignored in comparison with exp(z)
if (x > 20.0) {
return exp(x)/2.0;
if (x > 20) {
if (x >= LOG_MAX_VALUE) {
// Avoid overflow (MATH-905).
final double t = exp(0.5 * x);
return (0.5 * t) * t;
} else {
return 0.5 * exp(x);
}
} else if (x < -20) {
if (x <= -LOG_MAX_VALUE) {
// Avoid overflow (MATH-905).
final double t = exp(-0.5 * x);
return (0.5 * t) * t;
} else {
return 0.5 * exp(-x);
}
}
if (x < -20) {
return exp(-x)/2.0;
}
double hiPrec[] = new double[2];
final double hiPrec[] = new double[2];
if (x < 0.0) {
x = -x;
}
@ -449,12 +461,22 @@ public class FastMath {
// for values of z larger than about 20,
// exp(-z) can be ignored in comparison with exp(z)
if (x > 20.0) {
return exp(x)/2.0;
}
if (x < -20) {
return -exp(-x)/2.0;
if (x > 20) {
if (x >= LOG_MAX_VALUE) {
// Avoid overflow (MATH-905).
final double t = exp(0.5 * x);
return (0.5 * t) * t;
} else {
return 0.5 * exp(x);
}
} else if (x < -20) {
if (x <= -LOG_MAX_VALUE) {
// Avoid overflow (MATH-905).
final double t = exp(-0.5 * x);
return (-0.5 * t) * t;
} else {
return -0.5 * exp(-x);
}
}
if (x == 0) {

44
src/test/java/org/apache/commons/math3/util/FastMathTest.java
View File

@ -157,6 +157,50 @@ public class FastMathTest {
}
@Test
public void testMath905LargePositive() {
final double start = StrictMath.log(Double.MAX_VALUE);
final double endT = StrictMath.sqrt(2) * StrictMath.sqrt(Double.MAX_VALUE);
final double end = 2 * StrictMath.log(endT);
double maxErr = 0;
for (double x = start; x < end; x += 1e-3) {
final double tst = FastMath.cosh(x);
final double ref = Math.cosh(x);
maxErr = FastMath.max(maxErr, FastMath.abs(ref - tst) / FastMath.ulp(ref));
}
Assert.assertEquals(0, maxErr, 3);
for (double x = start; x < end; x += 1e-3) {
final double tst = FastMath.sinh(x);
final double ref = Math.sinh(x);
maxErr = FastMath.max(maxErr, FastMath.abs(ref - tst) / FastMath.ulp(ref));
}
Assert.assertEquals(0, maxErr, 3);
}
@Test
public void testMath905LargeNegative() {
final double start = -StrictMath.log(Double.MAX_VALUE);
final double endT = StrictMath.sqrt(2) * StrictMath.sqrt(Double.MAX_VALUE);
final double end = -2 * StrictMath.log(endT);
double maxErr = 0;
for (double x = start; x > end; x -= 1e-3) {
final double tst = FastMath.cosh(x);
final double ref = Math.cosh(x);
maxErr = FastMath.max(maxErr, FastMath.abs(ref - tst) / FastMath.ulp(ref));
}
Assert.assertEquals(0, maxErr, 3);
for (double x = start; x > end; x -= 1e-3) {
final double tst = FastMath.sinh(x);
final double ref = Math.sinh(x);
maxErr = FastMath.max(maxErr, FastMath.abs(ref - tst) / FastMath.ulp(ref));
}
Assert.assertEquals(0, maxErr, 3);
}
@Test
public void testHyperbolicInverses() {
double maxErr = 0;