ptolemy.domains.continuous.kernel.solver
Class ExplicitRK45Solver

java.lang.Object
  ptolemy.domains.continuous.kernel.ContinuousODESolver
      ptolemy.domains.continuous.kernel.solver.ExplicitRK45Solver

public class ExplicitRK45Solver
extends ContinuousODESolver

This class implements a fourth-order Runge-Kutta ODE solving method. The algorithm was introduced in "A Variable Order Runge-Kutta Method for Initial Value Problems with Rapidly Varying Right-Hand Sides" by J. R. Cash and Alan H. Karp, ACM Transactions on Mathematical Software, vol 16, pp. 201-222, 1990. For completeness, a brief explanation of the algorithm is explained below.

For an ODE of the form:

 dx(t)/dt = f(x(t), t), x(0) = x0
 

it does the following:

 K0 = f(x(n), tn);
 K1 = f(x(n) + 0.2*K0*h, tn + 0.2*h);
 K2 = f(x(n) + (3.0/40*K0 + 9.0/40*K1)*h, tn + 0.3*h);
 K3 = f(x(n) + (0.3*K0 - 0.9*K1 + 1.2*K2)*h, tn + 0.6*h);
 K4 = f(x(n) + (-11/54*K0 + 5.0/2*K1 -70/27*K2 + 35/27*K3)*h, tn + 1.0*h);
 K5 = f(x(n) + (1631/55296*K0 + 175/512*K1 + 575/13824*K2 + 3544275/110592*K3 + 253/4096*K4)*h, tn + 7/8*h);
 x(n+1) = x(n)+(37/378*K0 + 250/621*K2 + 125.0/594*K3 + 512.0/1771*K5)*h;
 

, and error control:

 LTE = [(37.0/378 - 2825.0/27648)*K0 + (250.0/621 - 18575.0/48384)*K2 +
 (125.0/594 - 13525.0/55296)*K3 + (0.0 - 277.0/14336)*K4 +
 (512.0/1771 - 0.25)*K5]*h.
 

If the LTE is less than the error tolerance, then this step size h is considered successful, and the next integration step size h' is predicted as:

 h' = h * Math.pow((ErrorTolerance/LTE), 1.0/5.0)
 

This is a fourth order method, but uses a fifth order procedure to estimate the local truncation error.

It takes 6 steps for this solver to resolve a state with an integration step size.

Since:

Ptolemy II 6.0

Version:

$Id: ExplicitRK45Solver.java 57040 2010-01-27 20:52:32Z cxh $

Author:

Haiyang Zheng, Edward A. Lee

Accepted Rating:

Green (hyzheng)

Proposed Rating:

Green (hyzheng)

Field Summary

private static double[][]	_B B coefficients
private static double[]	_E E coefficients
private static int	_ERROR_INDEX The index of the error stored in the auxiliary variables.
private static int	_ORDER The order of the algorithm.
private int	_roundCount The round counter.
protected static double[]	_TIME_INCREMENTS The ratio of time increments within one integration step.

Fields inherited from class ptolemy.domains.continuous.kernel.ContinuousODESolver

_director

Constructor Summary

ExplicitRK45Solver()

Method Summary

protected int	_getRound() Return the current round.
protected double	_getRoundTimeIncrement() Get the current round factor.
protected boolean	_isStepFinished() Return true if the current integration step is finished.
protected void	_reset() Reset the solver, indicating to it that we are starting an integration step.
protected void	_setRound(int round) Set the round for the next integration step.
int	getIntegratorAuxVariableCount() Return the number of time increments plus one (to store the truncation error).
void	integratorIntegrate(ContinuousIntegrator integrator) Fire the given integrator.
boolean	integratorIsAccurate(ContinuousIntegrator integrator) Return true if the integration is accurate for the given integrator.
double	integratorSuggestedStepSize(ContinuousIntegrator integrator) Predict the next step size for the integrators executed under this solver.

Methods inherited from class ptolemy.domains.continuous.kernel.ContinuousODESolver

_debug, _isDebugging, _makeSolverOf

Methods inherited from class java.lang.Object

clone, equals, finalize, getClass, hashCode, notify, notifyAll, toString, wait, wait, wait

Field Detail

_TIME_INCREMENTS

protected static final double[] _TIME_INCREMENTS

The ratio of time increments within one integration step.

_B

private static final double[][] _B

B coefficients

_E

private static final double[] _E

E coefficients

_ERROR_INDEX

private static final int _ERROR_INDEX

The index of the error stored in the auxiliary variables.

_ORDER

private static final int _ORDER

The order of the algorithm.

See Also:

Constant Field Values

_roundCount

private int _roundCount

The round counter.

Constructor Detail

ExplicitRK45Solver

public ExplicitRK45Solver()

Method Detail

getIntegratorAuxVariableCount

public final int getIntegratorAuxVariableCount()

Return the number of time increments plus one (to store the truncation error).

Specified by:

getIntegratorAuxVariableCount in class ContinuousODESolver

Returns:

The number of time increments plus one.

integratorIntegrate

public void integratorIntegrate(ContinuousIntegrator integrator)
                         throws IllegalActionException

Fire the given integrator. This method performs the ODE solving algorithm described in the class comment.

Specified by:

integratorIntegrate in class ContinuousODESolver

Parameters:

integrator - The integrator of that calls this method.

Throws:

IllegalActionException - If there is no director, or can not read input, or can not send output.

integratorIsAccurate

public boolean integratorIsAccurate(ContinuousIntegrator integrator)

Return true if the integration is accurate for the given integrator. This estimates the local truncation error for that integrator and compare it with the error tolerance.

Specified by:

integratorIsAccurate in class ContinuousODESolver

Parameters:

integrator - The integrator of that calls this method.

Returns:

True if the integration is successful.

integratorSuggestedStepSize

public double integratorSuggestedStepSize(ContinuousIntegrator integrator)

Predict the next step size for the integrators executed under this solver. This uses the algorithm in the class comments to predict the next step size based on the current estimation of the local truncation error.

Specified by:

integratorSuggestedStepSize in class ContinuousODESolver

Parameters:

integrator - The integrator that calls this method.

Returns:

The next step size suggested by the given integrator.

_getRound

protected int _getRound()

Return the current round.

Specified by:

_getRound in class ContinuousODESolver

Returns:

The current round.

_getRoundTimeIncrement

protected final double _getRoundTimeIncrement()

Get the current round factor. If the step is finished, then return 1.0.

Specified by:

_getRoundTimeIncrement in class ContinuousODESolver

Returns:

The current round factor.

See Also:

ContinuousODESolver._isStepFinished()

_isStepFinished

protected final boolean _isStepFinished()

Return true if the current integration step is finished. This method will return true if _incrementRound() has been called 6 or more times since _reset().

Specified by:

_isStepFinished in class ContinuousODESolver

Returns:

Return true if the solver has finished an integration step.

See Also:

_reset()

_reset

protected final void _reset()

Reset the solver, indicating to it that we are starting an integration step. This method resets the round counter.

Specified by:

_reset in class ContinuousODESolver

_setRound

protected void _setRound(int round)

Set the round for the next integration step.

Specified by:

_setRound in class ContinuousODESolver

Parameters:

round - The round for the next integration step.

See Also:

ContinuousODESolver.getIntegratorAuxVariableCount()