1 /*
2  * Copyright (C) 2015 The Android Open Source Project
3  *
4  * Licensed under the Apache License, Version 2.0 (the "License");
5  * you may not use this file except in compliance with the License.
6  * You may obtain a copy of the License at
7  *
8  *   http://www.apache.org/licenses/LICENSE-2.0
9  *
10  * Unless required by applicable law or agreed to in writing, software
11  * distributed under the License is distributed on an "AS IS" BASIS,
12  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13  * See the License for the specific language governing permissions and
14  * limitations under the License.
15  */
16 
17 // THIS TYPICALLY TAKES > 4 MINUTES TO RUN!
18 // It should generate no output during that time.
19 
20 package com.hp.creals;
21 
22 import junit.framework.AssertionFailedError;
23 import junit.framework.TestCase;
24 
25 import java.math.BigInteger;
26 import java.util.Random;
27 
28 public class SlowCRTest extends TestCase {
check(boolean x, String s)29     private static void check(boolean x, String s) {
30         if (!x) throw new AssertionFailedError(s);
31     }
32     final static int TEST_PREC = -200; // 200 bits to the right of
33                                        // binary point.
34     final static int NRANDOM = 100;    // Number of random values to
35                                        // test.  Bigger ==> slower
checkEq(CR x, CR y, String s)36     private static void checkEq(CR x, CR y, String s) {
37         check(x.compareTo(y, TEST_PREC) == 0, s);
38     }
checkApprEq(double x, double y, String s)39     private static void checkApprEq(double x, double y, String s) {
40         check(Math.abs(x - y) <= 0.000001, s);
41     }
42     final static BigInteger MASK =
43             BigInteger.ONE.shiftLeft(-TEST_PREC).subtract(BigInteger.ONE);
isApprInt(CR x)44     private static boolean isApprInt(CR x) {
45         BigInteger appr = x.get_appr(TEST_PREC);
46         return appr.and(MASK).signum() == 0;
47     }
48 
49     final static CR ZERO = CR.valueOf(0);
50     final static CR ONE = CR.valueOf(1);
51     final static CR TWO = CR.valueOf(2);
52     final static CR BIG = CR.valueOf(200).exp();
53     final static CR SMALL = BIG.inverse();
54     final static CR HALF_PI = CR.PI.divide(CR.valueOf(2));
55 
56     final static UnaryCRFunction ATAN = UnaryCRFunction.atanFunction;
57     final static UnaryCRFunction TAN = UnaryCRFunction.tanFunction;
58     final static UnaryCRFunction COSINE = UnaryCRFunction.sinFunction
59                                              .monotoneDerivative(ZERO, CR.PI);
60     final static UnaryCRFunction ARCSINE =
61                  UnaryCRFunction.sinFunction.inverseMonotone(HALF_PI.negate(),
62                                                              HALF_PI);
63 
64     // Perform some consistency checks on trig functions at x
65     // We assume that x is within floating point range.
checkTrig(CR x)66     private static void checkTrig(CR x) {
67         double xAsDouble = x.doubleValue();
68         if (Math.abs(xAsDouble) < 1000000.0) {
69             checkApprEq(x.sin().doubleValue(), Math.sin(xAsDouble),
70                           "sin float compare:" + xAsDouble);
71             checkApprEq(x.cos().doubleValue(), Math.cos(xAsDouble),
72                           "cos float compare:" + xAsDouble);
73             checkApprEq(TAN.execute(x).doubleValue(), Math.tan(xAsDouble),
74                           "tan float compare:" + xAsDouble);
75             checkApprEq(ATAN.execute(x).doubleValue(), Math.atan(xAsDouble),
76                           "atan float compare:" + xAsDouble);
77         }
78         if (Math.abs(xAsDouble) < 1.0) {
79             checkApprEq(x.asin().doubleValue(), Math.asin(xAsDouble),
80                           "asin float compare:" + xAsDouble);
81             checkApprEq(x.acos().doubleValue(), Math.acos(xAsDouble),
82                           "acos float compare:" + xAsDouble);
83             checkEq(ARCSINE.execute(x), x.asin(),
84                           "inverse(sin) compare:" + xAsDouble);
85         }
86         if (xAsDouble < 3.1415926535 && xAsDouble > 0.0) {
87             checkApprEq(COSINE.execute(x).doubleValue(), Math.cos(xAsDouble),
88                           "deriv(sin) float compare:"  + xAsDouble);
89             checkEq(COSINE.execute(x), x.cos(),
90                           "deriv(sin) float compare:"  + xAsDouble);
91         }
92         // Check that sin(x+v) = sin(x)cos(v) + cos(x)sin(v)
93         // for a couple of different values of v.
94         for (int i = 1; i <= 5; ++i) {
95             CR v = CR.valueOf(i);
96             checkEq(
97                 x.add(v).sin(),
98                 x.sin().multiply(v.cos()).add(x.cos().multiply(v.sin())),
99                 "Angle sum formula failed for " + xAsDouble + " + " + i);
100         }
101         checkEq(x.cos().multiply(x.cos()).add(x.sin().multiply(x.sin())),
102                     CR.valueOf(1),
103                     "sin(x)^2 + cos(x)^2 != 1:" + xAsDouble);
104         // Check that inverses are consistent
105         checkEq(x, TAN.execute(ATAN.execute(x)),
106                       "tan(atan(" + xAsDouble + ")");
107         CR xcos = x.cos();
108         CR tmp = xcos.acos();
109         // Result or its inverse should differ from x by an
110         // exact multiple of pi.
111         check(isApprInt(tmp.subtract(x).divide(CR.PI))
112               || isApprInt(tmp.add(x).divide(CR.PI)),
113               "acos(cos):" + xAsDouble);
114         CR xsin = x.sin();
115         tmp = ARCSINE.execute(xsin);
116         CR tmp2 = xsin.asin();
117         checkEq(tmp, tmp2, "Asin(sin) computations differ:" + xAsDouble);
118         // Result or its inverse should differ from x by an
119         // exact multiple of pi.
120         check(isApprInt(tmp.subtract(x).divide(CR.PI))
121               || isApprInt(tmp.add(x).divide(CR.PI)),
122               "acos(cos):" + xAsDouble);
123     }
124 
checkExpLn(CR x)125     private static void checkExpLn(CR x) {
126         double xAsDouble = x.doubleValue();
127         if (Math.abs(xAsDouble) < 10.0) {
128             checkApprEq(x.exp().doubleValue(), Math.exp(xAsDouble),
129                           "exp float compare:" + xAsDouble);
130         }
131         if (Math.abs(xAsDouble) <= 1000.0) {
132             checkEq(x, x.exp().ln(), "ln(exp) failed:" + xAsDouble);
133             checkEq(x.multiply(CR.valueOf(2)).exp(),
134 
135                         x.exp().multiply(x.exp()),
136                         "exp^2 failed:" + xAsDouble);
137         }
138         if (xAsDouble > 0.000000001) {
139             checkApprEq(x.ln().doubleValue(), Math.log(xAsDouble),
140                           "exp float compare:" + xAsDouble);
141             checkEq(x, x.ln().exp(), "exp(ln) failed:" + xAsDouble);
142             checkEq(x.ln().divide(CR.valueOf(2)), x.sqrt().ln(),
143                         "ln(sqrt) failed:" + xAsDouble);
144             // Check that ln(xv) = ln(x) + ln(v) for various v
145             for (int i = 1; i <= 5; ++i) {
146                 CR v = CR.valueOf(i);
147                 checkEq(
148                     x.ln().add(v.ln()),
149                     x.multiply(v).ln(),
150                     "ln(product) formula failed for:" + xAsDouble + "," + i);
151             }
152         }
153     }
154 
checkBasic(CR x)155     private static void checkBasic(CR x) {
156         checkEq(x.abs().sqrt().multiply(x.abs().sqrt()), x.abs(),
157                     "sqrt*sqrt:" + x.doubleValue());
158         if (!x.get_appr(TEST_PREC).equals(BigInteger.ZERO)) {
159             checkEq(x.inverse().inverse(), x,
160                         "inverse(inverse):" + x.doubleValue());
161         }
162     }
163 
testSlowTrig()164     public void testSlowTrig() {
165         checkEq(ZERO.acos(), CR.PI.divide(TWO), "acos(0)");
166         checkEq(ONE.acos(), ZERO, "acos(1)");
167         checkEq(ONE.negate().acos(), CR.PI, "acos(-1)");
168         checkEq(ZERO.asin(), ZERO, "asin(0)");
169         checkEq(ONE.asin(), CR.PI.divide(TWO), "asin(1)");
170         checkEq(ONE.negate().asin(), CR.PI.divide(TWO).negate(), "asin(-1)");
171         checkTrig(ZERO);
172         CR BIG = CR.valueOf(200).exp();
173         checkTrig(BIG);
174         checkTrig(BIG.negate());
175         checkTrig(SMALL);
176         checkTrig(SMALL.negate());
177         checkTrig(CR.PI);
178         checkTrig(CR.PI.subtract(SMALL));
179         checkTrig(CR.PI.add(SMALL));
180         checkTrig(CR.PI.negate());
181         checkTrig(CR.PI.negate().subtract(SMALL));
182         checkTrig(CR.PI.negate().add(SMALL));
183         Random r = new Random();  // Random seed!
184         for (int i = 0; i < NRANDOM; ++i) {
185             double d = Math.exp(2.0 * r.nextDouble() - 1.0);
186             if (r.nextBoolean()) d = -d;
187             final CR x = CR.valueOf(d);
188             checkTrig(x);
189         }
190         // And a few big ones
191         for (int i = 0; i < 10; ++i) {
192             double d = Math.exp(200.0 * r.nextDouble());
193             if (r.nextBoolean()) d = -d;
194             final CR x = CR.valueOf(d);
195             checkTrig(x);
196         }
197     }
198 
testSlowExpLn()199     public void testSlowExpLn() {
200         checkEq(CR.valueOf(1).ln(), CR.valueOf(0), "ln(1) != 0");
201         checkExpLn(CR.valueOf(0));
202         CR BIG = CR.valueOf(200).exp();
203         checkExpLn(BIG);
204         checkExpLn(BIG.negate());
205         checkExpLn(SMALL);
206         checkExpLn(SMALL.negate());
207         checkExpLn(CR.PI);
208         checkExpLn(ONE);
209         checkExpLn(ONE.subtract(SMALL));
210         checkExpLn(ONE.negate().subtract(SMALL));
211         Random r = new Random();  // Random seed!
212         for (int i = 0; i < NRANDOM; ++i) {
213             double d = Math.exp(10.0 * r.nextDouble() - 1.0);
214             if (r.nextBoolean()) d = -d;
215             final CR x = CR.valueOf(d);
216             checkExpLn(x);
217         }
218         // And a few big ones
219         for (int i = 0; i < 10; ++i) {
220             double d = Math.exp(200.0 * r.nextDouble());
221             if (r.nextBoolean()) d = -d;
222             final CR x = CR.valueOf(d);
223             checkExpLn(x);
224         }
225     }
226 
testSlowBasic()227     public void testSlowBasic() {
228         checkEq(ZERO.sqrt(), ZERO, "sqrt(0)");
229         checkEq(ZERO.abs(), ZERO, "abs(0)");
230         Random r = new Random();  // Random seed!
231         for (int i = 0; i < NRANDOM; ++i) {
232             double d = Math.exp(10.0 * r.nextDouble() - 1.0);
233             if (r.nextBoolean()) d = -d;
234             final CR x = CR.valueOf(d);
235             checkBasic(x);
236         }
237         // And a few very big ones, but within IEEE double range
238         for (int i = 0; i < 10; ++i) {
239             double d = Math.exp(600.0 * r.nextDouble());
240             if (r.nextBoolean()) d = -d;
241             final CR x = CR.valueOf(d);
242             checkBasic(x);
243         }
244     }
245 }
246