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 org.apache.commons.math.ode.DerivativeException;
21 import org.apache.commons.math.ode.ODEIntegrator;
22 import org.apache.commons.math.ode.sampling.StepHandler;
23 import org.apache.commons.math.ode.sampling.StepInterpolator;
24
25
26
27
28
29 public class TestProblemHandler
30 implements StepHandler {
31
32
33 private TestProblemAbstract problem;
34
35
36 private double maxValueError;
37 private double maxTimeError;
38
39
40 private double lastError;
41
42
43 private double lastTime;
44
45
46 private ODEIntegrator integrator;
47
48
49 private double expectedStepStart;
50
51
52
53
54
55
56 public TestProblemHandler(TestProblemAbstract problem, ODEIntegrator integrator) {
57 this.problem = problem;
58 this.integrator = integrator;
59 reset();
60 }
61
62 public boolean requiresDenseOutput() {
63 return true;
64 }
65
66 public void reset() {
67 maxValueError = 0;
68 maxTimeError = 0;
69 lastError = 0;
70 expectedStepStart = Double.NaN;
71 }
72
73 public void handleStep(StepInterpolator interpolator,
74 boolean isLast)
75 throws DerivativeException {
76
77 double start = integrator.getCurrentStepStart();
78 if (Math.abs((start - problem.getInitialTime()) / integrator.getCurrentSignedStepsize()) > 0.001) {
79
80
81 if (!Double.isNaN(expectedStepStart)) {
82 maxTimeError = Math.max(maxTimeError, Math.abs(start - expectedStepStart));
83 }
84 expectedStepStart = start + integrator.getCurrentSignedStepsize();
85 }
86
87 double pT = interpolator.getPreviousTime();
88 double cT = interpolator.getCurrentTime();
89 double[] errorScale = problem.getErrorScale();
90
91
92 if (isLast) {
93 double[] interpolatedY = interpolator.getInterpolatedState();
94 double[] theoreticalY = problem.computeTheoreticalState(cT);
95 for (int i = 0; i < interpolatedY.length; ++i) {
96 double error = Math.abs(interpolatedY[i] - theoreticalY[i]);
97 lastError = Math.max(error, lastError);
98 }
99 lastTime = cT;
100 }
101
102
103 for (int k = 0; k <= 20; ++k) {
104
105 double time = pT + (k * (cT - pT)) / 20;
106 interpolator.setInterpolatedTime(time);
107 double[] interpolatedY = interpolator.getInterpolatedState();
108 double[] theoreticalY = problem.computeTheoreticalState(interpolator.getInterpolatedTime());
109
110
111 for (int i = 0; i < interpolatedY.length; ++i) {
112 double error = errorScale[i] * Math.abs(interpolatedY[i] - theoreticalY[i]);
113 maxValueError = Math.max(error, maxValueError);
114 }
115 }
116 }
117
118
119
120
121
122 public double getMaximalValueError() {
123 return maxValueError;
124 }
125
126
127
128
129
130 public double getMaximalTimeError() {
131 return maxTimeError;
132 }
133
134
135
136
137
138 public double getLastError() {
139 return lastError;
140 }
141
142
143
144
145
146 public double getLastTime() {
147 return lastTime;
148 }
149
150 }