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.nonstiff;
19
20 import java.io.IOException;
21 import java.io.ObjectInput;
22 import java.io.ObjectOutput;
23
24 import org.apache.commons.math.ode.AbstractIntegrator;
25 import org.apache.commons.math.ode.sampling.AbstractStepInterpolator;
26
27 /** This class represents an interpolator over the last step during an
28 * ODE integration for Runge-Kutta and embedded Runge-Kutta integrators.
29 *
30 * @see RungeKuttaIntegrator
31 * @see EmbeddedRungeKuttaIntegrator
32 *
33 * @version $Revision: 783103 $ $Date: 2009-06-09 15:33:19 -0400 (Tue, 09 Jun 2009) $
34 * @since 1.2
35 */
36
37 abstract class RungeKuttaStepInterpolator
38 extends AbstractStepInterpolator {
39
40 /** Simple constructor.
41 * This constructor builds an instance that is not usable yet, the
42 * {@link #reinitialize} method should be called before using the
43 * instance in order to initialize the internal arrays. This
44 * constructor is used only in order to delay the initialization in
45 * some cases. The {@link RungeKuttaIntegrator} and {@link
46 * EmbeddedRungeKuttaIntegrator} classes use the prototyping design
47 * pattern to create the step interpolators by cloning an
48 * uninitialized model and latter initializing the copy.
49 */
50 protected RungeKuttaStepInterpolator() {
51 super();
52 yDotK = null;
53 integrator = null;
54 }
55
56 /** Copy constructor.
57
58 * <p>The copied interpolator should have been finalized before the
59 * copy, otherwise the copy will not be able to perform correctly any
60 * interpolation and will throw a {@link NullPointerException}
61 * later. Since we don't want this constructor to throw the
62 * exceptions finalization may involve and since we don't want this
63 * method to modify the state of the copied interpolator,
64 * finalization is <strong>not</strong> done automatically, it
65 * remains under user control.</p>
66
67 * <p>The copy is a deep copy: its arrays are separated from the
68 * original arrays of the instance.</p>
69
70 * @param interpolator interpolator to copy from.
71
72 */
73 public RungeKuttaStepInterpolator(final RungeKuttaStepInterpolator interpolator) {
74
75 super(interpolator);
76
77 if (interpolator.currentState != null) {
78 final int dimension = currentState.length;
79
80 yDotK = new double[interpolator.yDotK.length][];
81 for (int k = 0; k < interpolator.yDotK.length; ++k) {
82 yDotK[k] = new double[dimension];
83 System.arraycopy(interpolator.yDotK[k], 0,
84 yDotK[k], 0, dimension);
85 }
86
87 } else {
88 yDotK = null;
89 }
90
91 // we cannot keep any reference to the equations in the copy
92 // the interpolator should have been finalized before
93 integrator = null;
94
95 }
96
97 /** Reinitialize the instance
98 * <p>Some Runge-Kutta integrators need fewer functions evaluations
99 * than their counterpart step interpolators. So the interpolator
100 * should perform the last evaluations they need by themselves. The
101 * {@link RungeKuttaIntegrator RungeKuttaIntegrator} and {@link
102 * EmbeddedRungeKuttaIntegrator EmbeddedRungeKuttaIntegrator}
103 * abstract classes call this method in order to let the step
104 * interpolator perform the evaluations it needs. These evaluations
105 * will be performed during the call to <code>doFinalize</code> if
106 * any, i.e. only if the step handler either calls the {@link
107 * AbstractStepInterpolator#finalizeStep finalizeStep} method or the
108 * {@link AbstractStepInterpolator#getInterpolatedState
109 * getInterpolatedState} method (for an interpolator which needs a
110 * finalization) or if it clones the step interpolator.</p>
111 * @param integrator integrator being used
112 * @param y reference to the integrator array holding the state at
113 * the end of the step
114 * @param yDotK reference to the integrator array holding all the
115 * intermediate slopes
116 * @param forward integration direction indicator
117 */
118 public void reinitialize(final AbstractIntegrator integrator,
119 final double[] y, final double[][] yDotK, final boolean forward) {
120 reinitialize(y, forward);
121 this.yDotK = yDotK;
122 this.integrator = integrator;
123 }
124
125 /** {@inheritDoc} */
126 @Override
127 public void writeExternal(final ObjectOutput out)
128 throws IOException {
129
130 // save the state of the base class
131 writeBaseExternal(out);
132
133 // save the local attributes
134 final int n = (currentState == null) ? -1 : currentState.length;
135 final int kMax = (yDotK == null) ? -1 : yDotK.length;
136 out.writeInt(kMax);
137 for (int k = 0; k < kMax; ++k) {
138 for (int i = 0; i < n; ++i) {
139 out.writeDouble(yDotK[k][i]);
140 }
141 }
142
143 // we do not save any reference to the equations
144
145 }
146
147 /** {@inheritDoc} */
148 @Override
149 public void readExternal(final ObjectInput in)
150 throws IOException {
151
152 // read the base class
153 final double t = readBaseExternal(in);
154
155 // read the local attributes
156 final int n = (currentState == null) ? -1 : currentState.length;
157 final int kMax = in.readInt();
158 yDotK = (kMax < 0) ? null : new double[kMax][];
159 for (int k = 0; k < kMax; ++k) {
160 yDotK[k] = (n < 0) ? null : new double[n];
161 for (int i = 0; i < n; ++i) {
162 yDotK[k][i] = in.readDouble();
163 }
164 }
165
166 integrator = null;
167
168 if (currentState != null) {
169 // we can now set the interpolated time and state
170 setInterpolatedTime(t);
171 } else {
172 interpolatedTime = t;
173 }
174
175 }
176
177 /** Slopes at the intermediate points */
178 protected double[][] yDotK;
179
180 /** Reference to the integrator. */
181 protected AbstractIntegrator integrator;
182
183 }