1   /*
2    * Licensed to the Apache Software Foundation (ASF) 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 ASF 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    *      http://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  
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   * This class is used to handle steps for the test problems
27   * integrated during the junit tests for the ODE integrators.
28   */
29  public class TestProblemHandler
30    implements StepHandler {
31  
32    /** Associated problem. */
33    private TestProblemAbstract problem;
34  
35    /** Maximal errors encountered during the integration. */
36    private double maxValueError;
37    private double maxTimeError;
38  
39    /** Error at the end of the integration. */
40    private double lastError;
41  
42    /** Time at the end of integration. */
43    private double lastTime;
44  
45    /** ODE solver used. */
46    private ODEIntegrator integrator;
47  
48    /** Expected start for step. */
49    private double expectedStepStart;
50  
51    /**
52     * Simple constructor.
53     * @param problem problem for which steps should be handled
54     * @param integrator ODE solver used
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          // multistep integrators do not handle the first steps themselves
80          // so we have to make sure the integrator we look at has really started its work
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      // store the error at the last step
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     // walk through the step
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       // update the errors
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    * Get the maximal value error encountered during integration.
120    * @return maximal value error
121    */
122   public double getMaximalValueError() {
123     return maxValueError;
124   }
125 
126   /**
127    * Get the maximal time error encountered during integration.
128    * @return maximal time error
129    */
130   public double getMaximalTimeError() {
131     return maxTimeError;
132   }
133 
134   /**
135    * Get the error at the end of the integration.
136    * @return error at the end of the integration
137    */
138   public double getLastError() {
139     return lastError;
140   }
141 
142   /**
143    * Get the time at the end of the integration.
144    * @return time at the end of the integration.
145    */
146   public double getLastTime() {
147     return lastTime;
148   }
149 
150 }