1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18 package org.apache.commons.math4.legacy.ode.nonstiff;
19
20 import org.apache.commons.math4.legacy.exception.DimensionMismatchException;
21 import org.apache.commons.math4.legacy.exception.MaxCountExceededException;
22 import org.apache.commons.math4.legacy.exception.NoBracketingException;
23 import org.apache.commons.math4.legacy.exception.NumberIsTooSmallException;
24 import org.apache.commons.math4.legacy.ode.FirstOrderIntegrator;
25 import org.apache.commons.math4.legacy.ode.TestProblem1;
26 import org.apache.commons.math4.legacy.ode.TestProblem3;
27 import org.apache.commons.math4.legacy.ode.TestProblem4;
28 import org.apache.commons.math4.legacy.ode.TestProblem5;
29 import org.apache.commons.math4.legacy.ode.TestProblemAbstract;
30 import org.apache.commons.math4.legacy.ode.TestProblemHandler;
31 import org.apache.commons.math4.legacy.ode.events.EventHandler;
32 import org.apache.commons.math4.legacy.ode.sampling.StepHandler;
33 import org.apache.commons.math4.legacy.ode.sampling.StepInterpolator;
34 import org.apache.commons.math4.core.jdkmath.JdkMath;
35 import org.junit.Assert;
36 import org.junit.Test;
37
38
39 public class DormandPrince54IntegratorTest {
40
41 @Test(expected=DimensionMismatchException.class)
42 public void testDimensionCheck()
43 throws DimensionMismatchException, NumberIsTooSmallException,
44 MaxCountExceededException, NoBracketingException {
45 TestProblem1 pb = new TestProblem1();
46 DormandPrince54Integrator integrator = new DormandPrince54Integrator(0.0, 1.0,
47 1.0e-10, 1.0e-10);
48 integrator.integrate(pb,
49 0.0, new double[pb.getDimension()+10],
50 1.0, new double[pb.getDimension()+10]);
51 }
52
53 @Test(expected=NumberIsTooSmallException.class)
54 public void testMinStep()
55 throws DimensionMismatchException, NumberIsTooSmallException,
56 MaxCountExceededException, NoBracketingException {
57
58 TestProblem1 pb = new TestProblem1();
59 double minStep = 0.1 * (pb.getFinalTime() - pb.getInitialTime());
60 double maxStep = pb.getFinalTime() - pb.getInitialTime();
61 double[] vecAbsoluteTolerance = { 1.0e-15, 1.0e-16 };
62 double[] vecRelativeTolerance = { 1.0e-15, 1.0e-16 };
63
64 FirstOrderIntegrator integ = new DormandPrince54Integrator(minStep, maxStep,
65 vecAbsoluteTolerance,
66 vecRelativeTolerance);
67 TestProblemHandler handler = new TestProblemHandler(pb, integ);
68 integ.addStepHandler(handler);
69 integ.integrate(pb,
70 pb.getInitialTime(), pb.getInitialState(),
71 pb.getFinalTime(), new double[pb.getDimension()]);
72 Assert.fail("an exception should have been thrown");
73 }
74
75 @Test
76 public void testSmallLastStep()
77 throws DimensionMismatchException, NumberIsTooSmallException,
78 MaxCountExceededException, NoBracketingException {
79
80 TestProblemAbstract pb = new TestProblem5();
81 double minStep = 1.25;
82 double maxStep = JdkMath.abs(pb.getFinalTime() - pb.getInitialTime());
83 double scalAbsoluteTolerance = 6.0e-4;
84 double scalRelativeTolerance = 6.0e-4;
85
86 AdaptiveStepsizeIntegrator integ =
87 new DormandPrince54Integrator(minStep, maxStep,
88 scalAbsoluteTolerance,
89 scalRelativeTolerance);
90
91 DP54SmallLastHandler handler = new DP54SmallLastHandler(minStep);
92 integ.addStepHandler(handler);
93 integ.setInitialStepSize(1.7);
94 integ.integrate(pb,
95 pb.getInitialTime(), pb.getInitialState(),
96 pb.getFinalTime(), new double[pb.getDimension()]);
97 Assert.assertTrue(handler.wasLastSeen());
98 Assert.assertEquals("Dormand-Prince 5(4)", integ.getName());
99 }
100
101 @Test
102 public void testBackward()
103 throws DimensionMismatchException, NumberIsTooSmallException,
104 MaxCountExceededException, NoBracketingException {
105
106 TestProblem5 pb = new TestProblem5();
107 double minStep = 0;
108 double maxStep = pb.getFinalTime() - pb.getInitialTime();
109 double scalAbsoluteTolerance = 1.0e-8;
110 double scalRelativeTolerance = 0.01 * scalAbsoluteTolerance;
111
112 FirstOrderIntegrator integ = new DormandPrince54Integrator(minStep, maxStep,
113 scalAbsoluteTolerance,
114 scalRelativeTolerance);
115 TestProblemHandler handler = new TestProblemHandler(pb, integ);
116 integ.addStepHandler(handler);
117 integ.integrate(pb, pb.getInitialTime(), pb.getInitialState(),
118 pb.getFinalTime(), new double[pb.getDimension()]);
119
120 Assert.assertTrue(handler.getLastError() < 2.0e-7);
121 Assert.assertTrue(handler.getMaximalValueError() < 2.0e-7);
122 Assert.assertEquals(0, handler.getMaximalTimeError(), 1.0e-12);
123 Assert.assertEquals("Dormand-Prince 5(4)", integ.getName());
124 }
125
126 private static final class DP54SmallLastHandler implements StepHandler {
127
128 DP54SmallLastHandler(double minStep) {
129 lastSeen = false;
130 this.minStep = minStep;
131 }
132
133 @Override
134 public void init(double t0, double[] y0, double t) {
135 }
136
137 @Override
138 public void handleStep(StepInterpolator interpolator, boolean isLast) {
139 if (isLast) {
140 lastSeen = true;
141 double h = interpolator.getCurrentTime() - interpolator.getPreviousTime();
142 Assert.assertTrue(JdkMath.abs(h) < minStep);
143 }
144 }
145
146 public boolean wasLastSeen() {
147 return lastSeen;
148 }
149
150 private boolean lastSeen;
151 private double minStep;
152 }
153
154 @Test
155 public void testIncreasingTolerance()
156 throws DimensionMismatchException, NumberIsTooSmallException,
157 MaxCountExceededException, NoBracketingException {
158
159 int previousCalls = Integer.MAX_VALUE;
160 for (int i = -12; i < -2; ++i) {
161 TestProblem1 pb = new TestProblem1();
162 double minStep = 0;
163 double maxStep = pb.getFinalTime() - pb.getInitialTime();
164 double scalAbsoluteTolerance = JdkMath.pow(10.0, i);
165 double scalRelativeTolerance = 0.01 * scalAbsoluteTolerance;
166
167 EmbeddedRungeKuttaIntegrator integ =
168 new DormandPrince54Integrator(minStep, maxStep,
169 scalAbsoluteTolerance, scalRelativeTolerance);
170 TestProblemHandler handler = new TestProblemHandler(pb, integ);
171 integ.setSafety(0.8);
172 integ.setMaxGrowth(5.0);
173 integ.setMinReduction(0.3);
174 integ.addStepHandler(handler);
175 integ.integrate(pb,
176 pb.getInitialTime(), pb.getInitialState(),
177 pb.getFinalTime(), new double[pb.getDimension()]);
178 Assert.assertEquals(0.8, integ.getSafety(), 1.0e-12);
179 Assert.assertEquals(5.0, integ.getMaxGrowth(), 1.0e-12);
180 Assert.assertEquals(0.3, integ.getMinReduction(), 1.0e-12);
181
182
183
184
185 Assert.assertTrue(handler.getMaximalValueError() < (0.7 * scalAbsoluteTolerance));
186 Assert.assertEquals(0, handler.getMaximalTimeError(), 1.0e-12);
187
188 int calls = pb.getCalls();
189 Assert.assertEquals(integ.getEvaluations(), calls);
190 Assert.assertTrue(calls <= previousCalls);
191 previousCalls = calls;
192 }
193 }
194
195 @Test
196 public void testEvents()
197 throws DimensionMismatchException, NumberIsTooSmallException,
198 MaxCountExceededException, NoBracketingException {
199
200 TestProblem4 pb = new TestProblem4();
201 double minStep = 0;
202 double maxStep = pb.getFinalTime() - pb.getInitialTime();
203 double scalAbsoluteTolerance = 1.0e-8;
204 double scalRelativeTolerance = 0.01 * scalAbsoluteTolerance;
205
206 FirstOrderIntegrator integ = new DormandPrince54Integrator(minStep, maxStep,
207 scalAbsoluteTolerance,
208 scalRelativeTolerance);
209 TestProblemHandler handler = new TestProblemHandler(pb, integ);
210 integ.addStepHandler(handler);
211 EventHandler[] functions = pb.getEventsHandlers();
212 double convergence = 1.0e-8 * maxStep;
213 for (int l = 0; l < functions.length; ++l) {
214 integ.addEventHandler(functions[l],
215 Double.POSITIVE_INFINITY, convergence, 1000);
216 }
217 Assert.assertEquals(functions.length, integ.getEventHandlers().size());
218 integ.integrate(pb,
219 pb.getInitialTime(), pb.getInitialState(),
220 pb.getFinalTime(), new double[pb.getDimension()]);
221
222 Assert.assertTrue(handler.getMaximalValueError() < 5.0e-6);
223 Assert.assertEquals(0, handler.getMaximalTimeError(), convergence);
224 Assert.assertEquals(12.0, handler.getLastTime(), convergence);
225 integ.clearEventHandlers();
226 Assert.assertEquals(0, integ.getEventHandlers().size());
227 }
228
229 @Test
230 public void testKepler()
231 throws DimensionMismatchException, NumberIsTooSmallException,
232 MaxCountExceededException, NoBracketingException {
233
234 final TestProblem3 pb = new TestProblem3(0.9);
235 double minStep = 0;
236 double maxStep = pb.getFinalTime() - pb.getInitialTime();
237 double scalAbsoluteTolerance = 1.0e-8;
238 double scalRelativeTolerance = scalAbsoluteTolerance;
239
240 FirstOrderIntegrator integ = new DormandPrince54Integrator(minStep, maxStep,
241 scalAbsoluteTolerance,
242 scalRelativeTolerance);
243 integ.addStepHandler(new KeplerHandler(pb));
244 integ.integrate(pb,
245 pb.getInitialTime(), pb.getInitialState(),
246 pb.getFinalTime(), new double[pb.getDimension()]);
247
248 Assert.assertEquals(integ.getEvaluations(), pb.getCalls());
249 Assert.assertTrue(pb.getCalls() < 2800);
250 }
251
252 @Test
253 public void testVariableSteps()
254 throws DimensionMismatchException, NumberIsTooSmallException,
255 MaxCountExceededException, NoBracketingException {
256
257 final TestProblem3 pb = new TestProblem3(0.9);
258 double minStep = 0;
259 double maxStep = pb.getFinalTime() - pb.getInitialTime();
260 double scalAbsoluteTolerance = 1.0e-8;
261 double scalRelativeTolerance = scalAbsoluteTolerance;
262
263 FirstOrderIntegrator integ = new DormandPrince54Integrator(minStep, maxStep,
264 scalAbsoluteTolerance,
265 scalRelativeTolerance);
266 integ.addStepHandler(new VariableHandler());
267 double stopTime = integ.integrate(pb, pb.getInitialTime(), pb.getInitialState(),
268 pb.getFinalTime(), new double[pb.getDimension()]);
269 Assert.assertEquals(pb.getFinalTime(), stopTime, 1.0e-10);
270 }
271
272 private static final class KeplerHandler implements StepHandler {
273 KeplerHandler(TestProblem3 pb) {
274 this.pb = pb;
275 }
276 @Override
277 public void init(double t0, double[] y0, double t) {
278 nbSteps = 0;
279 maxError = 0;
280 }
281 @Override
282 public void handleStep(StepInterpolator interpolator, boolean isLast)
283 throws MaxCountExceededException {
284
285 ++nbSteps;
286 for (int a = 1; a < 10; ++a) {
287
288 double prev = interpolator.getPreviousTime();
289 double curr = interpolator.getCurrentTime();
290 double interp = ((10 - a) * prev + a * curr) / 10;
291 interpolator.setInterpolatedTime(interp);
292
293 double[] interpolatedY = interpolator.getInterpolatedState ();
294 double[] theoreticalY = pb.computeTheoreticalState(interpolator.getInterpolatedTime());
295 double dx = interpolatedY[0] - theoreticalY[0];
296 double dy = interpolatedY[1] - theoreticalY[1];
297 double error = dx * dx + dy * dy;
298 if (error > maxError) {
299 maxError = error;
300 }
301 }
302 if (isLast) {
303 Assert.assertTrue(maxError < 7.0e-10);
304 Assert.assertTrue(nbSteps < 400);
305 }
306 }
307 private int nbSteps;
308 private double maxError;
309 private TestProblem3 pb;
310 }
311
312 private static final class VariableHandler implements StepHandler {
313 VariableHandler() {
314 firstTime = true;
315 minStep = 0;
316 maxStep = 0;
317 }
318 @Override
319 public void init(double t0, double[] y0, double t) {
320 firstTime = true;
321 minStep = 0;
322 maxStep = 0;
323 }
324 @Override
325 public void handleStep(StepInterpolator interpolator,
326 boolean isLast) {
327
328 double step = JdkMath.abs(interpolator.getCurrentTime()
329 - interpolator.getPreviousTime());
330 if (firstTime) {
331 minStep = JdkMath.abs(step);
332 maxStep = minStep;
333 firstTime = false;
334 } else {
335 if (step < minStep) {
336 minStep = step;
337 }
338 if (step > maxStep) {
339 maxStep = step;
340 }
341 }
342
343 if (isLast) {
344 Assert.assertTrue(minStep < (1.0 / 450.0));
345 Assert.assertTrue(maxStep > (1.0 / 4.2));
346 }
347 }
348 private boolean firstTime;
349 private double minStep;
350 private double maxStep;
351 }
352 }