View Javadoc
1   /*
2    * Licensed to the Hipparchus project under one or more
3    * contributor license agreements.  See the NOTICE file distributed with
4    * this work for additional information regarding copyright ownership.
5    * The Hipparchus project licenses this file to You under the Apache License, Version 2.0
6    * (the "License"); you may not use this file except in compliance with
7    * the License.  You may obtain a copy of the License at
8    *
9    *      https://www.apache.org/licenses/LICENSE-2.0
10   *
11   * Unless required by applicable law or agreed to in writing, software
12   * distributed under the License is distributed on an "AS IS" BASIS,
13   * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14   * See the License for the specific language governing permissions and
15   * limitations under the License.
16   */
17  package org.hipparchus.special.elliptic.carlson;
18  
19  import org.hipparchus.random.RandomGenerator;
20  import org.hipparchus.random.Well19937a;
21  import org.hipparchus.random.Well19937c;
22  import org.hipparchus.util.Binary64;
23  import org.hipparchus.util.FastMath;
24  import org.junit.Assert;
25  import org.junit.Test;
26  
27  public class CarlsonEllipticIntegralFieldTest {
28  
29      private Binary64 build(double d) {
30          return new Binary64(d);
31      }
32  
33      @Test
34      public void testNoConvergenceRf() {
35          Assert.assertTrue(CarlsonEllipticIntegral.rF(build(1), build(2), build(Double.NaN)).isNaN());
36      }
37  
38      @Test
39      public void testDlmfRf() {
40          Binary64 rf = CarlsonEllipticIntegral.rF(build(1), build(2), build(4));
41          Assert.assertEquals(0.6850858166, rf.getReal(), 1.0e-10);
42      }
43  
44      @Test
45      public void testCarlson1995rF() {
46  
47          Binary64 rf1 = CarlsonEllipticIntegral.rF(build(1), build(2), build(0));
48          Assert.assertEquals( 1.3110287771461, rf1.getReal(), 1.0e-13);
49  
50          Binary64 rf2 = CarlsonEllipticIntegral.rF(build(0.5), build(1), build(0));
51          Assert.assertEquals( 1.8540746773014, rf2.getReal(), 1.0e-13);
52  
53          Binary64 rf4 = CarlsonEllipticIntegral.rF(build(2), build(3), build(4));
54          Assert.assertEquals( 0.58408284167715, rf4.getReal(), 1.0e-13);
55  
56      }
57  
58      @Test
59      public void testCarlson1995ConsistencyRf() {
60          RandomGenerator random = new Well19937c(0x57f2689b3f4028b4l);
61          for (int i = 0; i < 10000; ++i) {
62              Binary64 x      = build(random.nextDouble() * 3);
63              Binary64 y      = build(random.nextDouble() * 3);
64              Binary64 lambda = build(random.nextDouble() * 3);
65              Binary64 mu     = x.multiply(y).divide(lambda);
66              Binary64 rfL    = CarlsonEllipticIntegral.rF(x.add(lambda), y.add(lambda), lambda);
67              Binary64 rfM    = CarlsonEllipticIntegral.rF(x.add(mu),     y.add(mu),     mu);
68              Binary64 rf0    = CarlsonEllipticIntegral.rF(x,             y,             Binary64.ZERO);
69              Assert.assertEquals(0.0, rfL.add(rfM).subtract(rf0).norm(), 2.0e-14);
70          }
71      }
72  
73      @Test
74      public void testNoConvergenceRc() {
75          Assert.assertTrue(CarlsonEllipticIntegral.rC(build(1), build(Double.NaN)).isNaN());
76      }
77  
78      @Test
79      public void testCarlson1995rC() {
80  
81          Binary64 rc1 = CarlsonEllipticIntegral.rC(build(0), build(0.25));
82          Assert.assertEquals(FastMath.PI, rc1.getReal(), 1.0e-15);
83  
84          Binary64 rc2 = CarlsonEllipticIntegral.rC(build(2.25), build(2));
85          Assert.assertEquals(FastMath.log(2), rc2.getReal(), 1.0e-15);
86  
87          Binary64 rc5 = CarlsonEllipticIntegral.rC(build(0.25), build(-2));
88          Assert.assertEquals(FastMath.log(2) / 3.0, rc5.getReal(), 1.0e-15);
89  
90      }
91  
92      @Test
93      public void testCarlson1995ConsistencyRc() {
94          RandomGenerator random = new Well19937c(0xf1170b6fc1a199cal);
95          for (int i = 0; i < 10000; ++i) {
96              Binary64 x      = build(random.nextDouble() * 3);
97              Binary64 lambda = build(random.nextDouble() * 3);
98              Binary64 mu     = x.multiply(x).divide(lambda);
99              Binary64 rcL    = CarlsonEllipticIntegral.rC(lambda,         x.add(lambda));
100             Binary64 rcM    = CarlsonEllipticIntegral.rC(mu,             x.add(mu));
101             Binary64 rc0    = CarlsonEllipticIntegral.rC(Binary64.ZERO, x);
102             Assert.assertEquals(0.0, rcL.add(rcM).subtract(rc0).norm(), 3.0e-14);
103         }
104     }
105 
106     @Test
107     public void testRfRc() {
108         RandomGenerator random = new Well19937a(0x7e8041334a8c20edl);
109         for (int i = 0; i < 10000; ++i) {
110             final Binary64 x = build(3 * random.nextDouble());
111             final Binary64 y = build(3 * random.nextDouble());
112             final Binary64 rf = CarlsonEllipticIntegral.rF(x, y, y);
113             final Binary64 rc = CarlsonEllipticIntegral.rC(x, y);
114             Assert.assertEquals(0.0, rf.subtract(rc).norm(), 4.0e-15);
115         }
116     }
117 
118     @Test
119     public void testNoConvergenceRj() {
120         Assert.assertTrue(CarlsonEllipticIntegral.rJ(build(1), build(1), build(1), build(Double.NaN)).isNaN());
121     }
122 
123     @Test
124     public void testCarlson1995rJ() {
125 
126         Binary64 rj01 = CarlsonEllipticIntegral.rJ(build(0), build(1), build(2), build(3));
127         Assert.assertEquals(0.77688623778582, rj01.getReal(), 1.0e-13);
128 
129         Binary64 rj02 = CarlsonEllipticIntegral.rJ(build(2), build(3), build(4), build(5));
130         Assert.assertEquals( 0.14297579667157, rj02.getReal(), 1.0e-13);
131 
132     }
133 
134     @Test
135     public void testCarlson1995ConsistencyRj() {
136         RandomGenerator random = new Well19937c(0x4af7bb722712e64el);
137         for (int i = 0; i < 10000; ++i) {
138             Binary64 x      = build(random.nextDouble() * 3);
139             Binary64 y      = build(random.nextDouble() * 3);
140             Binary64 p      = build(random.nextDouble() * 3);
141             Binary64 lambda = build(random.nextDouble() * 3);
142             Binary64 mu     = x.multiply(y).divide(lambda);
143             Binary64 a      = p.multiply(p).multiply(lambda.add(mu).add(x).add(y));
144             Binary64 b      = p.multiply(p.add(lambda)).multiply(p.add(mu));
145             Binary64 rjL    = CarlsonEllipticIntegral.rJ(x.add(lambda), y.add(lambda), lambda,         p.add(lambda));
146             Binary64 rjM    = CarlsonEllipticIntegral.rJ(x.add(mu),     y.add(mu),     mu,             p.add(mu));
147             Binary64 rj0    = CarlsonEllipticIntegral.rJ(x,             y,             Binary64.ZERO, p);
148             Binary64 rc     = CarlsonEllipticIntegral.rC(a, b);
149             Assert.assertEquals(0.0, rjL.add(rjM).subtract(rj0.subtract(rc.multiply(3))).norm(), 3.0e-13);
150         }
151     }
152 
153     @Test
154     public void testNoConvergenceRd() {
155         Assert.assertTrue(CarlsonEllipticIntegral.rD(build(1), build(1), build(Double.NaN)).isNaN());
156     }
157 
158     @Test
159     public void testCarlson1995rD() {
160 
161         Binary64 rd1 = CarlsonEllipticIntegral.rD(build(0), build(2), build(1));
162         Assert.assertEquals(1.7972103521034, rd1.getReal(), 1.0e-13);
163 
164         Binary64 rd2 = CarlsonEllipticIntegral.rD(build(2), build(3), build(4));
165         Assert.assertEquals( 0.16510527294261, rd2.getReal(), 1.0e-13);
166 
167     }
168 
169     @Test
170     public void testCarlson1995ConsistencyRd() {
171         RandomGenerator random = new Well19937c(0x17dea97eeb78206al);
172         for (int i = 0; i < 10000; ++i) {
173             Binary64 x      = build(random.nextDouble() * 3);
174             Binary64 y      = build(random.nextDouble() * 3);
175             Binary64 lambda = build(random.nextDouble() * 3);
176             Binary64 mu     = x.multiply(y).divide(lambda);
177             Binary64 rdL    = CarlsonEllipticIntegral.rD(lambda,         x.add(lambda), y.add(lambda));
178             Binary64 rdM    = CarlsonEllipticIntegral.rD(mu,             x.add(mu),     y.add(mu));
179             Binary64 rd0    = CarlsonEllipticIntegral.rD(Binary64.ZERO, x,             y);
180             Binary64 frac   = y.multiply(x.add(y).add(lambda).add(mu).sqrt()).reciprocal().multiply(3);
181             Assert.assertEquals(0.0, rdL.add(rdM).subtract(rd0.subtract(frac)).norm(), 9.0e-12);
182         }
183     }
184 
185     @Test
186     public void testRdNonSymmetry1() {
187         RandomGenerator random = new Well19937c(0x66db170b5ee1afc2l);
188         for (int i = 0; i < 10000; ++i) {
189             Binary64 x = build(random.nextDouble());
190             Binary64 y = build(random.nextDouble());
191             Binary64 z = build(random.nextDouble());
192             if (x.isZero() || y.isZero()) {
193                 continue;
194             }
195             // this is DLMF equation 19.21.7
196             Binary64 lhs = x.subtract(y).multiply(CarlsonEllipticIntegral.rD(y, z, x)).
197                             add(z.subtract(y).multiply(CarlsonEllipticIntegral.rD(x, y, z)));
198             Binary64 rhs = CarlsonEllipticIntegral.rF(x, y, z).subtract(y.divide(x.multiply(z)).sqrt()).multiply(3);
199             Assert.assertEquals(0.0, lhs.subtract(rhs).norm(), 1.0e-10);
200         }
201     }
202 
203     @Test
204     public void testRdNonSymmetry2() {
205         RandomGenerator random = new Well19937c(0x1a8994acc807438dl);
206         for (int i = 0; i < 10000; ++i) {
207             Binary64 x = build(random.nextDouble());
208             Binary64 y = build(random.nextDouble());
209             Binary64 z = build(random.nextDouble());
210             if (x.isZero() || y.isZero() || z.isZero()) {
211                 continue;
212             }
213             // this is DLMF equation 19.21.8
214             Binary64 lhs = CarlsonEllipticIntegral.rD(y, z, x).
215                             add(CarlsonEllipticIntegral.rD(z, x, y)).
216                             add(CarlsonEllipticIntegral.rD(x, y, z));
217             Binary64 rhs = x.multiply(y.multiply(z)).sqrt().reciprocal().multiply(3);
218             Assert.assertEquals(0.0, lhs.subtract(rhs).norm(), 2.0e-11);
219         }
220     }
221 
222     @Test
223     public void testCarlson1995rG() {
224 
225         Binary64 rg1 = CarlsonEllipticIntegral.rG(build(0), build(16), build(16));
226         Assert.assertEquals(FastMath.PI, rg1.getReal(), 1.0e-13);
227 
228         Binary64 rg2 = CarlsonEllipticIntegral.rG(build(2), build(3), build(4));
229         Assert.assertEquals(1.7255030280692, rg2.getReal(), 1.0e-13);
230 
231         Binary64 rg6 = CarlsonEllipticIntegral.rG(build(0), build(0.0796), build(4));
232         Assert.assertEquals( 1.0284758090288, rg6.getReal(), 1.0e-13);
233 
234     }
235 
236     @Test
237     public void testAlternateRG() {
238         RandomGenerator random = new Well19937c(0xa2946e4a55d133a6l);
239         for (int i = 0; i < 10000; ++i) {
240             Binary64 x = build(random.nextDouble() * 3);
241             Binary64 y = build(random.nextDouble() * 3);
242             Binary64 z = build(random.nextDouble() * 3);
243             Assert.assertEquals(0.0, CarlsonEllipticIntegral.rG(x, y, z).subtract(rgAlternateImplementation(x, y, z)).norm(), 2.0e-15);
244         }
245     }
246 
247     private Binary64 rgAlternateImplementation(final Binary64 x, final Binary64 y, final Binary64 z) {
248         // this implementation uses DLFM equation 19.21.11
249         return d(x, y, z).add(d(y, z, x)).add(d(z, x, y)).divide(6);
250     }
251 
252     private Binary64 d(final Binary64 u, final Binary64 v, final Binary64 w) {
253         return u.isZero() ? u : u.multiply(v.add(w)).multiply(new RdFieldDuplication<>(v, w, u).integral());
254     }
255 
256 }