public class GraggBulirschStoerIntegrator extends AdaptiveStepsizeIntegrator
The Gragg-Bulirsch-Stoer algorithm is one of the most efficient ones currently available for smooth problems. It uses Richardson extrapolation to estimate what would be the solution if the step size could be decreased down to zero.
This method changes both the step size and the order during
integration, in order to minimize computation cost. It is
particularly well suited when a very high precision is needed. The
limit where this method becomes more efficient than high-order
embedded Runge-Kutta methods like Dormand-Prince 8(5,3) depends on the problem. Results given in the
Hairer, Norsett and Wanner book show for example that this limit
occurs for accuracy around 1e-6 when integrating Saltzam-Lorenz
equations (the authors note this problem is extremely sensitive
to the errors in the first integration steps), and around 1e-11
for a two dimensional celestial mechanics problems with seven
bodies (pleiades problem, involving quasi-collisions for which
automatic step size control is essential).
This implementation is basically a reimplementation in Java of the odex fortran code by E. Hairer and G. Wanner. The redistribution policy for this code is available here, for convenience, it is reproduced below.
| Copyright (c) 2004, Ernst Hairer |
Redistribution and use in source and binary forms, with or
without modification, are permitted provided that the following
conditions are met:
|
| THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. |
mainSetDimension, scalAbsoluteTolerance, scalRelativeTolerance, vecAbsoluteTolerance, vecRelativeTolerance| Constructor and Description |
|---|
GraggBulirschStoerIntegrator(double minStep,
double maxStep,
double[] vecAbsoluteTolerance,
double[] vecRelativeTolerance)
Simple constructor.
|
GraggBulirschStoerIntegrator(double minStep,
double maxStep,
double scalAbsoluteTolerance,
double scalRelativeTolerance)
Simple constructor.
|
| Modifier and Type | Method and Description |
|---|---|
ODEStateAndDerivative |
integrate(ExpandableODE equations,
ODEState initialState,
double finalTime)
Integrate the differential equations up to the given time.
|
void |
setControlFactors(double control1,
double control2,
double control3,
double control4)
Set the step size control factors.
|
void |
setInterpolationControl(boolean useInterpolationErrorForControl,
int mudifControlParameter)
Set the interpolation order control parameter.
|
void |
setOrderControl(int maximalOrder,
double control1,
double control2)
Set the order control parameters.
|
void |
setStabilityCheck(boolean performStabilityCheck,
int maxNumIter,
int maxNumChecks,
double stepsizeReductionFactor)
Set the stability check controls.
|
filterStep, getMaxStep, getMinStep, initializeStep, resetInternalState, sanityChecks, setInitialStepSize, setStepSizeControl, setStepSizeControlacceptStep, addEventHandler, addEventHandler, addStepHandler, clearEventHandlers, clearStepHandlers, computeDerivatives, getCurrentSignedStepsize, getCurrentStepStart, getEquations, getEvaluations, getEvaluationsCounter, getEventHandlers, getMaxEvaluations, getName, getStepHandlers, getStepSize, getStepStart, initIntegration, isLastStep, resetOccurred, setIsLastStep, setMaxEvaluations, setStateInitialized, setStepSize, setStepStartclone, equals, finalize, getClass, hashCode, notify, notifyAll, toString, wait, wait, waitintegrate, integratepublic GraggBulirschStoerIntegrator(double minStep,
double maxStep,
double scalAbsoluteTolerance,
double scalRelativeTolerance)
minStep - minimal step (sign is irrelevant, regardless of
integration direction, forward or backward), the last step can
be smaller than thismaxStep - maximal step (sign is irrelevant, regardless of
integration direction, forward or backward), the last step can
be smaller than thisscalAbsoluteTolerance - allowed absolute errorscalRelativeTolerance - allowed relative errorpublic GraggBulirschStoerIntegrator(double minStep,
double maxStep,
double[] vecAbsoluteTolerance,
double[] vecRelativeTolerance)
minStep - minimal step (must be positive even for backward
integration), the last step can be smaller than thismaxStep - maximal step (must be positive even for backward
integration)vecAbsoluteTolerance - allowed absolute errorvecRelativeTolerance - allowed relative errorpublic void setStabilityCheck(boolean performStabilityCheck,
int maxNumIter,
int maxNumChecks,
double stepsizeReductionFactor)
The stability check is performed on the first few iterations of the extrapolation scheme. If this test fails, the step is rejected and the stepsize is reduced.
By default, the test is performed, at most during two iterations at each step, and at most once for each of these iterations. The default stepsize reduction factor is 0.5.
performStabilityCheck - if true, stability check will be performed,
if false, the check will be skippedmaxNumIter - maximal number of iterations for which checks are
performed (the number of iterations is reset to default if negative
or null)maxNumChecks - maximal number of checks for each iteration
(the number of checks is reset to default if negative or null)stepsizeReductionFactor - stepsize reduction factor in case of
failure (the factor is reset to default if lower than 0.0001 or
greater than 0.9999)public void setControlFactors(double control1,
double control2,
double control3,
double control4)
The new step size hNew is computed from the old one h by:
hNew = h * stepControl2 / (err/stepControl1)^(1/(2k + 1))where err is the scaled error and k the iteration number of the extrapolation scheme (counting from 0). The default values are 0.65 for stepControl1 and 0.94 for stepControl2.
The step size is subject to the restriction:
stepControl3^(1/(2k + 1))/stepControl4 <= hNew/h <= 1/stepControl3^(1/(2k + 1))The default values are 0.02 for stepControl3 and 4.0 for stepControl4.
control1 - first stepsize control factor (the factor is
reset to default if lower than 0.0001 or greater than 0.9999)control2 - second stepsize control factor (the factor
is reset to default if lower than 0.0001 or greater than 0.9999)control3 - third stepsize control factor (the factor is
reset to default if lower than 0.0001 or greater than 0.9999)control4 - fourth stepsize control factor (the factor
is reset to default if lower than 1.0001 or greater than 999.9)public void setOrderControl(int maximalOrder,
double control1,
double control2)
The Gragg-Bulirsch-Stoer method changes both the step size and the order during integration, in order to minimize computation cost. Each extrapolation step increases the order by 2, so the maximal order that will be used is always even, it is twice the maximal number of columns in the extrapolation table.
order is decreased if w(k - 1) <= w(k) * orderControl1 order is increased if w(k) <= w(k - 1) * orderControl2
where w is the table of work per unit step for each order (number of function calls divided by the step length), and k is the current order.
The default maximal order after construction is 18 (i.e. the maximal number of columns is 9). The default values are 0.8 for orderControl1 and 0.9 for orderControl2.
maximalOrder - maximal order in the extrapolation table (the
maximal order is reset to default if order <= 6 or odd)control1 - first order control factor (the factor is
reset to default if lower than 0.0001 or greater than 0.9999)control2 - second order control factor (the factor
is reset to default if lower than 0.0001 or greater than 0.9999)public void setInterpolationControl(boolean useInterpolationErrorForControl,
int mudifControlParameter)
useInterpolationErrorForControl - if true, interpolation error is used
for stepsize controlmudifControlParameter - interpolation order control parameter (the parameter
is reset to default if <= 0 or >= 7)public ODEStateAndDerivative integrate(ExpandableODE equations, ODEState initialState, double finalTime) throws MathIllegalArgumentException, MathIllegalStateException
This method solves an Initial Value Problem (IVP).
Since this method stores some internal state variables made
available in its public interface during integration (ODEIntegrator.getCurrentSignedStepsize()), it is not thread-safe.
equations - differential equations to integrateinitialState - initial state (time, primary and secondary state vectors)finalTime - target time for the integration
(can be set to a value smaller than t0 for backward integration)finalTime if
integration reached its target, but may be different if some ODEEventHandler stops it at some point.MathIllegalArgumentException - if integration step is too smallMathIllegalStateException - if the number of functions evaluations is exceededCopyright © 2016–2020 Hipparchus.org. All rights reserved.