1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18 package org.apache.commons.math.ode;
19
20 import junit.framework.*;
21
22 import org.apache.commons.math.ode.ClassicalRungeKuttaIntegrator;
23 import org.apache.commons.math.ode.DerivativeException;
24 import org.apache.commons.math.ode.FirstOrderIntegrator;
25 import org.apache.commons.math.ode.IntegratorException;
26 import org.apache.commons.math.ode.StepHandler;
27 import org.apache.commons.math.ode.StepInterpolator;
28 import org.apache.commons.math.ode.SwitchingFunction;
29
30 public class ClassicalRungeKuttaIntegratorTest
31 extends TestCase {
32
33 public ClassicalRungeKuttaIntegratorTest(String name) {
34 super(name);
35 }
36
37 public void testSanityChecks() {
38 try {
39 TestProblem1 pb = new TestProblem1();
40 new ClassicalRungeKuttaIntegrator(0.01).integrate(pb,
41 0.0, new double[pb.getDimension()+10],
42 1.0, new double[pb.getDimension()]);
43 fail("an exception should have been thrown");
44 } catch(DerivativeException de) {
45 fail("wrong exception caught");
46 } catch(IntegratorException ie) {
47 }
48 try {
49 TestProblem1 pb = new TestProblem1();
50 new ClassicalRungeKuttaIntegrator(0.01).integrate(pb,
51 0.0, new double[pb.getDimension()],
52 1.0, new double[pb.getDimension()+10]);
53 fail("an exception should have been thrown");
54 } catch(DerivativeException de) {
55 fail("wrong exception caught");
56 } catch(IntegratorException ie) {
57 }
58 try {
59 TestProblem1 pb = new TestProblem1();
60 new ClassicalRungeKuttaIntegrator(0.01).integrate(pb,
61 0.0, new double[pb.getDimension()],
62 0.0, new double[pb.getDimension()]);
63 fail("an exception should have been thrown");
64 } catch(DerivativeException de) {
65 fail("wrong exception caught");
66 } catch(IntegratorException ie) {
67 }
68 }
69
70 public void testDecreasingSteps()
71 throws DerivativeException, IntegratorException {
72
73 TestProblemAbstract[] problems = TestProblemFactory.getProblems();
74 for (int k = 0; k < problems.length; ++k) {
75
76 double previousError = Double.NaN;
77 for (int i = 4; i < 10; ++i) {
78
79 TestProblemAbstract pb = (TestProblemAbstract) problems[k].clone();
80 double step = (pb.getFinalTime() - pb.getInitialTime())
81 * Math.pow(2.0, -i);
82
83 FirstOrderIntegrator integ = new ClassicalRungeKuttaIntegrator(step);
84 TestProblemHandler handler = new TestProblemHandler(pb, integ);
85 integ.setStepHandler(handler);
86 SwitchingFunction[] functions = pb.getSwitchingFunctions();
87 for (int l = 0; l < functions.length; ++l) {
88 integ.addSwitchingFunction(functions[l],
89 Double.POSITIVE_INFINITY, 1.0e-6 * step, 1000);
90 }
91 integ.integrate(pb, pb.getInitialTime(), pb.getInitialState(),
92 pb.getFinalTime(), new double[pb.getDimension()]);
93
94 double error = handler.getMaximalValueError();
95 if (i > 4) {
96 assertTrue(error < Math.abs(previousError));
97 }
98 previousError = error;
99 assertEquals(0, handler.getMaximalTimeError(), 1.0e-12);
100 }
101
102 }
103
104 }
105
106 public void testSmallStep()
107 throws DerivativeException, IntegratorException {
108
109 TestProblem1 pb = new TestProblem1();
110 double step = (pb.getFinalTime() - pb.getInitialTime()) * 0.001;
111
112 FirstOrderIntegrator integ = new ClassicalRungeKuttaIntegrator(step);
113 TestProblemHandler handler = new TestProblemHandler(pb, integ);
114 integ.setStepHandler(handler);
115 integ.integrate(pb, pb.getInitialTime(), pb.getInitialState(),
116 pb.getFinalTime(), new double[pb.getDimension()]);
117
118 assertTrue(handler.getLastError() < 2.0e-13);
119 assertTrue(handler.getMaximalValueError() < 4.0e-12);
120 assertEquals(0, handler.getMaximalTimeError(), 1.0e-12);
121 assertEquals("classical Runge-Kutta", integ.getName());
122 }
123
124 public void testBigStep()
125 throws DerivativeException, IntegratorException {
126
127 TestProblem1 pb = new TestProblem1();
128 double step = (pb.getFinalTime() - pb.getInitialTime()) * 0.2;
129
130 FirstOrderIntegrator integ = new ClassicalRungeKuttaIntegrator(step);
131 TestProblemHandler handler = new TestProblemHandler(pb, integ);
132 integ.setStepHandler(handler);
133 integ.integrate(pb, pb.getInitialTime(), pb.getInitialState(),
134 pb.getFinalTime(), new double[pb.getDimension()]);
135
136 assertTrue(handler.getLastError() > 0.0004);
137 assertTrue(handler.getMaximalValueError() > 0.005);
138 assertEquals(0, handler.getMaximalTimeError(), 1.0e-12);
139
140 }
141
142 public void testKepler()
143 throws DerivativeException, IntegratorException {
144
145 final TestProblem3 pb = new TestProblem3(0.9);
146 double step = (pb.getFinalTime() - pb.getInitialTime()) * 0.0003;
147
148 FirstOrderIntegrator integ = new ClassicalRungeKuttaIntegrator(step);
149 integ.setStepHandler(new KeplerHandler(pb));
150 integ.integrate(pb,
151 pb.getInitialTime(), pb.getInitialState(),
152 pb.getFinalTime(), new double[pb.getDimension()]);
153 }
154
155 private static class KeplerHandler implements StepHandler {
156 public KeplerHandler(TestProblem3 pb) {
157 this.pb = pb;
158 reset();
159 }
160 public boolean requiresDenseOutput() {
161 return false;
162 }
163 public void reset() {
164 maxError = 0;
165 }
166 public void handleStep(StepInterpolator interpolator,
167 boolean isLast) {
168
169 double[] interpolatedY = interpolator.getInterpolatedState ();
170 double[] theoreticalY = pb.computeTheoreticalState(interpolator.getCurrentTime());
171 double dx = interpolatedY[0] - theoreticalY[0];
172 double dy = interpolatedY[1] - theoreticalY[1];
173 double error = dx * dx + dy * dy;
174 if (error > maxError) {
175 maxError = error;
176 }
177 if (isLast) {
178
179
180
181 assertTrue(maxError > 0.005);
182 }
183 }
184 private double maxError = 0;
185 private TestProblem3 pb;
186 }
187
188 public static Test suite() {
189 return new TestSuite(ClassicalRungeKuttaIntegratorTest.class);
190 }
191
192 }