darcy | 32db449 | 2009-01-26 19:49:26 -0800 | [diff] [blame] | 1 | /* |
| 2 | * Copyright 2003 Sun Microsystems, Inc. All Rights Reserved. |
| 3 | * DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER. |
| 4 | * |
| 5 | * This code is free software; you can redistribute it and/or modify it |
| 6 | * under the terms of the GNU General Public License version 2 only, as |
| 7 | * published by the Free Software Foundation. |
| 8 | * |
| 9 | * This code is distributed in the hope that it will be useful, but WITHOUT |
| 10 | * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or |
| 11 | * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License |
| 12 | * version 2 for more details (a copy is included in the LICENSE file that |
| 13 | * accompanied this code). |
| 14 | * |
| 15 | * You should have received a copy of the GNU General Public License version |
| 16 | * 2 along with this work; if not, write to the Free Software Foundation, |
| 17 | * Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA. |
| 18 | * |
| 19 | * Please contact Sun Microsystems, Inc., 4150 Network Circle, Santa Clara, |
| 20 | * CA 95054 USA or visit www.sun.com if you need additional information or |
| 21 | * have any questions. |
| 22 | */ |
| 23 | |
| 24 | /* |
| 25 | * @test |
| 26 | * @bug 4851638 4900189 4939441 |
| 27 | * @summary Tests for {Math, StrictMath}.expm1 |
| 28 | * @author Joseph D. Darcy |
| 29 | */ |
| 30 | |
| 31 | import sun.misc.DoubleConsts; |
| 32 | import sun.misc.FpUtils; |
| 33 | |
| 34 | /* |
| 35 | * The Taylor expansion of expxm1(x) = exp(x) -1 is |
| 36 | * |
| 37 | * 1 + x/1! + x^2/2! + x^3/3| + ... -1 = |
| 38 | * |
| 39 | * x + x^2/2! + x^3/3 + ... |
| 40 | * |
| 41 | * Therefore, for small values of x, expxm1 ~= x. |
| 42 | * |
| 43 | * For large values of x, expxm1(x) ~= exp(x) |
| 44 | * |
| 45 | * For large negative x, expxm1(x) ~= -1. |
| 46 | */ |
| 47 | |
| 48 | public class Expm1Tests { |
| 49 | |
| 50 | private Expm1Tests(){} |
| 51 | |
| 52 | static final double infinityD = Double.POSITIVE_INFINITY; |
| 53 | static final double NaNd = Double.NaN; |
| 54 | |
| 55 | static int testExpm1() { |
| 56 | int failures = 0; |
| 57 | |
| 58 | double [][] testCases = { |
| 59 | {Double.NaN, NaNd}, |
| 60 | {Double.longBitsToDouble(0x7FF0000000000001L), NaNd}, |
| 61 | {Double.longBitsToDouble(0xFFF0000000000001L), NaNd}, |
| 62 | {Double.longBitsToDouble(0x7FF8555555555555L), NaNd}, |
| 63 | {Double.longBitsToDouble(0xFFF8555555555555L), NaNd}, |
| 64 | {Double.longBitsToDouble(0x7FFFFFFFFFFFFFFFL), NaNd}, |
| 65 | {Double.longBitsToDouble(0xFFFFFFFFFFFFFFFFL), NaNd}, |
| 66 | {Double.longBitsToDouble(0x7FFDeadBeef00000L), NaNd}, |
| 67 | {Double.longBitsToDouble(0xFFFDeadBeef00000L), NaNd}, |
| 68 | {Double.longBitsToDouble(0x7FFCafeBabe00000L), NaNd}, |
| 69 | {Double.longBitsToDouble(0xFFFCafeBabe00000L), NaNd}, |
| 70 | {infinityD, infinityD}, |
| 71 | {-infinityD, -1.0}, |
| 72 | {-0.0, -0.0}, |
| 73 | {+0.0, +0.0}, |
| 74 | }; |
| 75 | |
| 76 | // Test special cases |
| 77 | for(int i = 0; i < testCases.length; i++) { |
| 78 | failures += testExpm1CaseWithUlpDiff(testCases[i][0], |
| 79 | testCases[i][1], 0, null); |
| 80 | } |
| 81 | |
| 82 | |
| 83 | // For |x| < 2^-54 expm1(x) ~= x |
| 84 | for(int i = DoubleConsts.MIN_SUB_EXPONENT; i <= -54; i++) { |
| 85 | double d = FpUtils.scalb(2, i); |
| 86 | failures += testExpm1Case(d, d); |
| 87 | failures += testExpm1Case(-d, -d); |
| 88 | } |
| 89 | |
| 90 | |
| 91 | // For values of y where exp(y) > 2^54, expm1(x) ~= exp(x). |
| 92 | // The least such y is ln(2^54) ~= 37.42994775023705; exp(x) |
| 93 | // overflows for x > ~= 709.8 |
| 94 | |
| 95 | // Use a 2-ulp error threshold to account for errors in the |
| 96 | // exp implementation; the increments of d in the loop will be |
| 97 | // exact. |
| 98 | for(double d = 37.5; d <= 709.5; d += 1.0) { |
| 99 | failures += testExpm1CaseWithUlpDiff(d, StrictMath.exp(d), 2, null); |
| 100 | } |
| 101 | |
| 102 | // For x > 710, expm1(x) should be infinity |
| 103 | for(int i = 10; i <= DoubleConsts.MAX_EXPONENT; i++) { |
| 104 | double d = FpUtils.scalb(2, i); |
| 105 | failures += testExpm1Case(d, infinityD); |
| 106 | } |
| 107 | |
| 108 | // By monotonicity, once the limit is reached, the |
| 109 | // implemenation should return the limit for all smaller |
| 110 | // values. |
| 111 | boolean reachedLimit [] = {false, false}; |
| 112 | |
| 113 | // Once exp(y) < 0.5 * ulp(1), expm1(y) ~= -1.0; |
| 114 | // The greatest such y is ln(2^-53) ~= -36.7368005696771. |
| 115 | for(double d = -36.75; d >= -127.75; d -= 1.0) { |
| 116 | failures += testExpm1CaseWithUlpDiff(d, -1.0, 1, |
| 117 | reachedLimit); |
| 118 | } |
| 119 | |
| 120 | for(int i = 7; i <= DoubleConsts.MAX_EXPONENT; i++) { |
| 121 | double d = -FpUtils.scalb(2, i); |
| 122 | failures += testExpm1CaseWithUlpDiff(d, -1.0, 1, reachedLimit); |
| 123 | } |
| 124 | |
| 125 | // Test for monotonicity failures near multiples of log(2). |
| 126 | // Test two numbers before and two numbers after each chosen |
| 127 | // value; i.e. |
| 128 | // |
| 129 | // pcNeighbors[] = |
| 130 | // {nextDown(nextDown(pc)), |
| 131 | // nextDown(pc), |
| 132 | // pc, |
| 133 | // nextUp(pc), |
| 134 | // nextUp(nextUp(pc))} |
| 135 | // |
| 136 | // and we test that expm1(pcNeighbors[i]) <= expm1(pcNeighbors[i+1]) |
| 137 | { |
| 138 | double pcNeighbors[] = new double[5]; |
| 139 | double pcNeighborsExpm1[] = new double[5]; |
| 140 | double pcNeighborsStrictExpm1[] = new double[5]; |
| 141 | |
| 142 | for(int i = -50; i <= 50; i++) { |
| 143 | double pc = StrictMath.log(2)*i; |
| 144 | |
| 145 | pcNeighbors[2] = pc; |
| 146 | pcNeighbors[1] = FpUtils.nextDown(pc); |
| 147 | pcNeighbors[0] = FpUtils.nextDown(pcNeighbors[1]); |
| 148 | pcNeighbors[3] = FpUtils.nextUp(pc); |
| 149 | pcNeighbors[4] = FpUtils.nextUp(pcNeighbors[3]); |
| 150 | |
| 151 | for(int j = 0; j < pcNeighbors.length; j++) { |
| 152 | pcNeighborsExpm1[j] = Math.expm1(pcNeighbors[j]); |
| 153 | pcNeighborsStrictExpm1[j] = StrictMath.expm1(pcNeighbors[j]); |
| 154 | } |
| 155 | |
| 156 | for(int j = 0; j < pcNeighborsExpm1.length-1; j++) { |
| 157 | if(pcNeighborsExpm1[j] > pcNeighborsExpm1[j+1] ) { |
| 158 | failures++; |
| 159 | System.err.println("Monotonicity failure for Math.expm1 on " + |
| 160 | pcNeighbors[j] + " and " + |
| 161 | pcNeighbors[j+1] + "\n\treturned " + |
| 162 | pcNeighborsExpm1[j] + " and " + |
| 163 | pcNeighborsExpm1[j+1] ); |
| 164 | } |
| 165 | |
| 166 | if(pcNeighborsStrictExpm1[j] > pcNeighborsStrictExpm1[j+1] ) { |
| 167 | failures++; |
| 168 | System.err.println("Monotonicity failure for StrictMath.expm1 on " + |
| 169 | pcNeighbors[j] + " and " + |
| 170 | pcNeighbors[j+1] + "\n\treturned " + |
| 171 | pcNeighborsStrictExpm1[j] + " and " + |
| 172 | pcNeighborsStrictExpm1[j+1] ); |
| 173 | } |
| 174 | |
| 175 | |
| 176 | } |
| 177 | |
| 178 | } |
| 179 | } |
| 180 | |
| 181 | return failures; |
| 182 | } |
| 183 | |
| 184 | public static int testExpm1Case(double input, |
| 185 | double expected) { |
| 186 | return testExpm1CaseWithUlpDiff(input, expected, 1, null); |
| 187 | } |
| 188 | |
| 189 | public static int testExpm1CaseWithUlpDiff(double input, |
| 190 | double expected, |
| 191 | double ulps, |
| 192 | boolean [] reachedLimit) { |
| 193 | int failures = 0; |
| 194 | double mathUlps = ulps, strictUlps = ulps; |
| 195 | double mathOutput; |
| 196 | double strictOutput; |
| 197 | |
| 198 | if (reachedLimit != null) { |
| 199 | if (reachedLimit[0]) |
| 200 | mathUlps = 0; |
| 201 | |
| 202 | if (reachedLimit[1]) |
| 203 | strictUlps = 0; |
| 204 | } |
| 205 | |
| 206 | failures += Tests.testUlpDiffWithLowerBound("Math.expm1(double)", |
| 207 | input, mathOutput=Math.expm1(input), |
| 208 | expected, mathUlps, -1.0); |
| 209 | failures += Tests.testUlpDiffWithLowerBound("StrictMath.expm1(double)", |
| 210 | input, strictOutput=StrictMath.expm1(input), |
| 211 | expected, strictUlps, -1.0); |
| 212 | if (reachedLimit != null) { |
| 213 | reachedLimit[0] |= (mathOutput == -1.0); |
| 214 | reachedLimit[1] |= (strictOutput == -1.0); |
| 215 | } |
| 216 | |
| 217 | return failures; |
| 218 | } |
| 219 | |
| 220 | public static void main(String argv[]) { |
| 221 | int failures = 0; |
| 222 | |
| 223 | failures += testExpm1(); |
| 224 | |
| 225 | if (failures > 0) { |
| 226 | System.err.println("Testing expm1 incurred " |
| 227 | + failures + " failures."); |
| 228 | throw new RuntimeException(); |
| 229 | } |
| 230 | } |
| 231 | } |