org.apache.commons.math.ode.DormandPrince54StepInterpolator Maven / Gradle / Ivy
/*
* Licensed to the Apache Software Foundation (ASF) under one or more
* contributor license agreements. See the NOTICE file distributed with
* this work for additional information regarding copyright ownership.
* The ASF licenses this file to You under the Apache License, Version 2.0
* (the "License"); you may not use this file except in compliance with
* the License. You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*/
package org.apache.commons.math.ode;
/**
* This class represents an interpolator over the last step during an
* ODE integration for the 5(4) Dormand-Prince integrator.
*
* @see DormandPrince54Integrator
*
* @version $Revision: 620312 $ $Date: 2008-02-10 12:28:59 -0700 (Sun, 10 Feb 2008) $
* @since 1.2
*/
class DormandPrince54StepInterpolator
extends RungeKuttaStepInterpolator {
/** Simple constructor.
* This constructor builds an instance that is not usable yet, the
* {@link #reinitialize} method should be called before using the
* instance in order to initialize the internal arrays. This
* constructor is used only in order to delay the initialization in
* some cases. The {@link EmbeddedRungeKuttaIntegrator} uses the
* prototyping design pattern to create the step interpolators by
* cloning an uninitialized model and latter initializing the copy.
*/
public DormandPrince54StepInterpolator() {
super();
v1 = null;
v2 = null;
v3 = null;
v4 = null;
vectorsInitialized = false;
}
/** Copy constructor.
* @param interpolator interpolator to copy from. The copy is a deep
* copy: its arrays are separated from the original arrays of the
* instance
*/
public DormandPrince54StepInterpolator(DormandPrince54StepInterpolator interpolator) {
super(interpolator);
if (interpolator.v1 == null) {
v1 = null;
v2 = null;
v3 = null;
v4 = null;
vectorsInitialized = false;
} else {
v1 = (double[]) interpolator.v1.clone();
v2 = (double[]) interpolator.v2.clone();
v3 = (double[]) interpolator.v3.clone();
v4 = (double[]) interpolator.v4.clone();
vectorsInitialized = interpolator.vectorsInitialized;
}
}
/** Really copy the finalized instance.
* @return a copy of the finalized instance
*/
protected StepInterpolator doCopy() {
return new DormandPrince54StepInterpolator(this);
}
/** Reinitialize the instance
* @param equations set of differential equations being integrated
* @param y reference to the integrator array holding the state at
* the end of the step
* @param yDotK reference to the integrator array holding all the
* intermediate slopes
* @param forward integration direction indicator
*/
public void reinitialize(FirstOrderDifferentialEquations equations,
double[] y, double[][] yDotK, boolean forward) {
super.reinitialize(equations, y, yDotK, forward);
v1 = null;
v2 = null;
v3 = null;
v4 = null;
vectorsInitialized = false;
}
/** Store the current step time.
* @param t current time
*/
public void storeTime(double t) {
super.storeTime(t);
vectorsInitialized = false;
}
/** Compute the state at the interpolated time.
* @param theta normalized interpolation abscissa within the step
* (theta is zero at the previous time step and one at the current time step)
* @param oneMinusThetaH time gap between the interpolated time and
* the current time
* @throws DerivativeException this exception is propagated to the caller if the
* underlying user function triggers one
*/
protected void computeInterpolatedState(double theta,
double oneMinusThetaH)
throws DerivativeException {
if (! vectorsInitialized) {
if (v1 == null) {
v1 = new double[interpolatedState.length];
v2 = new double[interpolatedState.length];
v3 = new double[interpolatedState.length];
v4 = new double[interpolatedState.length];
}
// no step finalization is needed for this interpolator
// we need to compute the interpolation vectors for this time step
for (int i = 0; i < interpolatedState.length; ++i) {
v1[i] = h * (a70 * yDotK[0][i] + a72 * yDotK[2][i] + a73 * yDotK[3][i] +
a74 * yDotK[4][i] + a75 * yDotK[5][i]);
v2[i] = h * yDotK[0][i] - v1[i];
v3[i] = v1[i] - v2[i] - h * yDotK[6][i];
v4[i] = h * (d0 * yDotK[0][i] + d2 * yDotK[2][i] + d3 * yDotK[3][i] +
d4 * yDotK[4][i] + d5 * yDotK[5][i] + d6 * yDotK[6][i]);
}
vectorsInitialized = true;
}
// interpolate
double eta = oneMinusThetaH / h;
for (int i = 0; i < interpolatedState.length; ++i) {
interpolatedState[i] = currentState[i] -
eta * (v1[i] - theta * (v2[i] + theta * (v3[i] + eta * v4[i])));
}
}
/** First vector for interpolation. */
private double[] v1;
/** Second vector for interpolation. */
private double[] v2;
/** Third vector for interpolation. */
private double[] v3;
/** Fourth vector for interpolation. */
private double[] v4;
/** Initialization indicator for the interpolation vectors. */
private boolean vectorsInitialized;
/** Last row of the Butcher-array internal weights, element 0. */
private static final double a70 = 35.0 / 384.0;
// element 1 is zero, so it is neither stored nor used
/** Last row of the Butcher-array internal weights, element 2. */
private static final double a72 = 500.0 / 1113.0;
/** Last row of the Butcher-array internal weights, element 3. */
private static final double a73 = 125.0 / 192.0;
/** Last row of the Butcher-array internal weights, element 4. */
private static final double a74 = -2187.0 / 6784.0;
/** Last row of the Butcher-array internal weights, element 5. */
private static final double a75 = 11.0 / 84.0;
/** Shampine (1986) Dense output, element 0. */
private static final double d0 = -12715105075.0 / 11282082432.0;
// element 1 is zero, so it is neither stored nor used
/** Shampine (1986) Dense output, element 2. */
private static final double d2 = 87487479700.0 / 32700410799.0;
/** Shampine (1986) Dense output, element 3. */
private static final double d3 = -10690763975.0 / 1880347072.0;
/** Shampine (1986) Dense output, element 4. */
private static final double d4 = 701980252875.0 / 199316789632.0;
/** Shampine (1986) Dense output, element 5. */
private static final double d5 = -1453857185.0 / 822651844.0;
/** Shampine (1986) Dense output, element 6. */
private static final double d6 = 69997945.0 / 29380423.0;
/** Serializable version identifier */
private static final long serialVersionUID = 4104157279605906956L;
}
© 2015 - 2025 Weber Informatics LLC | Privacy Policy