diff options
Diffstat (limited to 'src/main/java/org/apache/commons/math3/ode/sampling/StepNormalizer.java')
| -rw-r--r-- | src/main/java/org/apache/commons/math3/ode/sampling/StepNormalizer.java | 300 |
1 files changed, 300 insertions, 0 deletions
diff --git a/src/main/java/org/apache/commons/math3/ode/sampling/StepNormalizer.java b/src/main/java/org/apache/commons/math3/ode/sampling/StepNormalizer.java new file mode 100644 index 0000000..307c0d9 --- /dev/null +++ b/src/main/java/org/apache/commons/math3/ode/sampling/StepNormalizer.java @@ -0,0 +1,300 @@ +/* + * 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.math3.ode.sampling; + +import org.apache.commons.math3.exception.MaxCountExceededException; +import org.apache.commons.math3.util.FastMath; +import org.apache.commons.math3.util.Precision; + +/** + * This class wraps an object implementing {@link FixedStepHandler} + * into a {@link StepHandler}. + + * <p>This wrapper allows to use fixed step handlers with general + * integrators which cannot guaranty their integration steps will + * remain constant and therefore only accept general step + * handlers.</p> + * + * <p>The stepsize used is selected at construction time. The {@link + * FixedStepHandler#handleStep handleStep} method of the underlying + * {@link FixedStepHandler} object is called at normalized times. The + * normalized times can be influenced by the {@link StepNormalizerMode} and + * {@link StepNormalizerBounds}.</p> + * + * <p>There is no constraint on the integrator, it can use any time step + * it needs (time steps longer or shorter than the fixed time step and + * non-integer ratios are all allowed).</p> + * + * <p> + * <table border="1" align="center"> + * <tr BGCOLOR="#CCCCFF"><td colspan=6><font size="+2">Examples (step size = 0.5)</font></td></tr> + * <tr BGCOLOR="#EEEEFF"><font size="+1"><td>Start time</td><td>End time</td> + * <td>Direction</td><td>{@link StepNormalizerMode Mode}</td> + * <td>{@link StepNormalizerBounds Bounds}</td><td>Output</td></font></tr> + * <tr><td>0.3</td><td>3.1</td><td>forward</td><td>{@link StepNormalizerMode#INCREMENT INCREMENT}</td><td>{@link StepNormalizerBounds#NEITHER NEITHER}</td><td>0.8, 1.3, 1.8, 2.3, 2.8</td></tr> + * <tr><td>0.3</td><td>3.1</td><td>forward</td><td>{@link StepNormalizerMode#INCREMENT INCREMENT}</td><td>{@link StepNormalizerBounds#FIRST FIRST}</td><td>0.3, 0.8, 1.3, 1.8, 2.3, 2.8</td></tr> + * <tr><td>0.3</td><td>3.1</td><td>forward</td><td>{@link StepNormalizerMode#INCREMENT INCREMENT}</td><td>{@link StepNormalizerBounds#LAST LAST}</td><td>0.8, 1.3, 1.8, 2.3, 2.8, 3.1</td></tr> + * <tr><td>0.3</td><td>3.1</td><td>forward</td><td>{@link StepNormalizerMode#INCREMENT INCREMENT}</td><td>{@link StepNormalizerBounds#BOTH BOTH}</td><td>0.3, 0.8, 1.3, 1.8, 2.3, 2.8, 3.1</td></tr> + * <tr><td>0.3</td><td>3.1</td><td>forward</td><td>{@link StepNormalizerMode#MULTIPLES MULTIPLES}</td><td>{@link StepNormalizerBounds#NEITHER NEITHER}</td><td>0.5, 1.0, 1.5, 2.0, 2.5, 3.0</td></tr> + * <tr><td>0.3</td><td>3.1</td><td>forward</td><td>{@link StepNormalizerMode#MULTIPLES MULTIPLES}</td><td>{@link StepNormalizerBounds#FIRST FIRST}</td><td>0.3, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0</td></tr> + * <tr><td>0.3</td><td>3.1</td><td>forward</td><td>{@link StepNormalizerMode#MULTIPLES MULTIPLES}</td><td>{@link StepNormalizerBounds#LAST LAST}</td><td>0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.1</td></tr> + * <tr><td>0.3</td><td>3.1</td><td>forward</td><td>{@link StepNormalizerMode#MULTIPLES MULTIPLES}</td><td>{@link StepNormalizerBounds#BOTH BOTH}</td><td>0.3, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.1</td></tr> + * <tr><td>0.0</td><td>3.0</td><td>forward</td><td>{@link StepNormalizerMode#INCREMENT INCREMENT}</td><td>{@link StepNormalizerBounds#NEITHER NEITHER}</td><td>0.5, 1.0, 1.5, 2.0, 2.5, 3.0</td></tr> + * <tr><td>0.0</td><td>3.0</td><td>forward</td><td>{@link StepNormalizerMode#INCREMENT INCREMENT}</td><td>{@link StepNormalizerBounds#FIRST FIRST}</td><td>0.0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0</td></tr> + * <tr><td>0.0</td><td>3.0</td><td>forward</td><td>{@link StepNormalizerMode#INCREMENT INCREMENT}</td><td>{@link StepNormalizerBounds#LAST LAST}</td><td>0.5, 1.0, 1.5, 2.0, 2.5, 3.0</td></tr> + * <tr><td>0.0</td><td>3.0</td><td>forward</td><td>{@link StepNormalizerMode#INCREMENT INCREMENT}</td><td>{@link StepNormalizerBounds#BOTH BOTH}</td><td>0.0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0</td></tr> + * <tr><td>0.0</td><td>3.0</td><td>forward</td><td>{@link StepNormalizerMode#MULTIPLES MULTIPLES}</td><td>{@link StepNormalizerBounds#NEITHER NEITHER}</td><td>0.5, 1.0, 1.5, 2.0, 2.5, 3.0</td></tr> + * <tr><td>0.0</td><td>3.0</td><td>forward</td><td>{@link StepNormalizerMode#MULTIPLES MULTIPLES}</td><td>{@link StepNormalizerBounds#FIRST FIRST}</td><td>0.0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0</td></tr> + * <tr><td>0.0</td><td>3.0</td><td>forward</td><td>{@link StepNormalizerMode#MULTIPLES MULTIPLES}</td><td>{@link StepNormalizerBounds#LAST LAST}</td><td>0.5, 1.0, 1.5, 2.0, 2.5, 3.0</td></tr> + * <tr><td>0.0</td><td>3.0</td><td>forward</td><td>{@link StepNormalizerMode#MULTIPLES MULTIPLES}</td><td>{@link StepNormalizerBounds#BOTH BOTH}</td><td>0.0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0</td></tr> + * <tr><td>3.1</td><td>0.3</td><td>backward</td><td>{@link StepNormalizerMode#INCREMENT INCREMENT}</td><td>{@link StepNormalizerBounds#NEITHER NEITHER}</td><td>2.6, 2.1, 1.6, 1.1, 0.6</td></tr> + * <tr><td>3.1</td><td>0.3</td><td>backward</td><td>{@link StepNormalizerMode#INCREMENT INCREMENT}</td><td>{@link StepNormalizerBounds#FIRST FIRST}</td><td>3.1, 2.6, 2.1, 1.6, 1.1, 0.6</td></tr> + * <tr><td>3.1</td><td>0.3</td><td>backward</td><td>{@link StepNormalizerMode#INCREMENT INCREMENT}</td><td>{@link StepNormalizerBounds#LAST LAST}</td><td>2.6, 2.1, 1.6, 1.1, 0.6, 0.3</td></tr> + * <tr><td>3.1</td><td>0.3</td><td>backward</td><td>{@link StepNormalizerMode#INCREMENT INCREMENT}</td><td>{@link StepNormalizerBounds#BOTH BOTH}</td><td>3.1, 2.6, 2.1, 1.6, 1.1, 0.6, 0.3</td></tr> + * <tr><td>3.1</td><td>0.3</td><td>backward</td><td>{@link StepNormalizerMode#MULTIPLES MULTIPLES}</td><td>{@link StepNormalizerBounds#NEITHER NEITHER}</td><td>3.0, 2.5, 2.0, 1.5, 1.0, 0.5</td></tr> + * <tr><td>3.1</td><td>0.3</td><td>backward</td><td>{@link StepNormalizerMode#MULTIPLES MULTIPLES}</td><td>{@link StepNormalizerBounds#FIRST FIRST}</td><td>3.1, 3.0, 2.5, 2.0, 1.5, 1.0, 0.5</td></tr> + * <tr><td>3.1</td><td>0.3</td><td>backward</td><td>{@link StepNormalizerMode#MULTIPLES MULTIPLES}</td><td>{@link StepNormalizerBounds#LAST LAST}</td><td>3.0, 2.5, 2.0, 1.5, 1.0, 0.5, 0.3</td></tr> + * <tr><td>3.1</td><td>0.3</td><td>backward</td><td>{@link StepNormalizerMode#MULTIPLES MULTIPLES}</td><td>{@link StepNormalizerBounds#BOTH BOTH}</td><td>3.1, 3.0, 2.5, 2.0, 1.5, 1.0, 0.5, 0.3</td></tr> + * <tr><td>3.0</td><td>0.0</td><td>backward</td><td>{@link StepNormalizerMode#INCREMENT INCREMENT}</td><td>{@link StepNormalizerBounds#NEITHER NEITHER}</td><td>2.5, 2.0, 1.5, 1.0, 0.5, 0.0</td></tr> + * <tr><td>3.0</td><td>0.0</td><td>backward</td><td>{@link StepNormalizerMode#INCREMENT INCREMENT}</td><td>{@link StepNormalizerBounds#FIRST FIRST}</td><td>3.0, 2.5, 2.0, 1.5, 1.0, 0.5, 0.0</td></tr> + * <tr><td>3.0</td><td>0.0</td><td>backward</td><td>{@link StepNormalizerMode#INCREMENT INCREMENT}</td><td>{@link StepNormalizerBounds#LAST LAST}</td><td>2.5, 2.0, 1.5, 1.0, 0.5, 0.0</td></tr> + * <tr><td>3.0</td><td>0.0</td><td>backward</td><td>{@link StepNormalizerMode#INCREMENT INCREMENT}</td><td>{@link StepNormalizerBounds#BOTH BOTH}</td><td>3.0, 2.5, 2.0, 1.5, 1.0, 0.5, 0.0</td></tr> + * <tr><td>3.0</td><td>0.0</td><td>backward</td><td>{@link StepNormalizerMode#MULTIPLES MULTIPLES}</td><td>{@link StepNormalizerBounds#NEITHER NEITHER}</td><td>2.5, 2.0, 1.5, 1.0, 0.5, 0.0</td></tr> + * <tr><td>3.0</td><td>0.0</td><td>backward</td><td>{@link StepNormalizerMode#MULTIPLES MULTIPLES}</td><td>{@link StepNormalizerBounds#FIRST FIRST}</td><td>3.0, 2.5, 2.0, 1.5, 1.0, 0.5, 0.0</td></tr> + * <tr><td>3.0</td><td>0.0</td><td>backward</td><td>{@link StepNormalizerMode#MULTIPLES MULTIPLES}</td><td>{@link StepNormalizerBounds#LAST LAST}</td><td>2.5, 2.0, 1.5, 1.0, 0.5, 0.0</td></tr> + * <tr><td>3.0</td><td>0.0</td><td>backward</td><td>{@link StepNormalizerMode#MULTIPLES MULTIPLES}</td><td>{@link StepNormalizerBounds#BOTH BOTH}</td><td>3.0, 2.5, 2.0, 1.5, 1.0, 0.5, 0.0</td></tr> + * </table> + * </p> + * + * @see StepHandler + * @see FixedStepHandler + * @see StepNormalizerMode + * @see StepNormalizerBounds + * @since 1.2 + */ + +public class StepNormalizer implements StepHandler { + /** Fixed time step. */ + private double h; + + /** Underlying step handler. */ + private final FixedStepHandler handler; + + /** First step time. */ + private double firstTime; + + /** Last step time. */ + private double lastTime; + + /** Last state vector. */ + private double[] lastState; + + /** Last derivatives vector. */ + private double[] lastDerivatives; + + /** Integration direction indicator. */ + private boolean forward; + + /** The step normalizer bounds settings to use. */ + private final StepNormalizerBounds bounds; + + /** The step normalizer mode to use. */ + private final StepNormalizerMode mode; + + /** Simple constructor. Uses {@link StepNormalizerMode#INCREMENT INCREMENT} + * mode, and {@link StepNormalizerBounds#FIRST FIRST} bounds setting, for + * backwards compatibility. + * @param h fixed time step (sign is not used) + * @param handler fixed time step handler to wrap + */ + public StepNormalizer(final double h, final FixedStepHandler handler) { + this(h, handler, StepNormalizerMode.INCREMENT, + StepNormalizerBounds.FIRST); + } + + /** Simple constructor. Uses {@link StepNormalizerBounds#FIRST FIRST} + * bounds setting. + * @param h fixed time step (sign is not used) + * @param handler fixed time step handler to wrap + * @param mode step normalizer mode to use + * @since 3.0 + */ + public StepNormalizer(final double h, final FixedStepHandler handler, + final StepNormalizerMode mode) { + this(h, handler, mode, StepNormalizerBounds.FIRST); + } + + /** Simple constructor. Uses {@link StepNormalizerMode#INCREMENT INCREMENT} + * mode. + * @param h fixed time step (sign is not used) + * @param handler fixed time step handler to wrap + * @param bounds step normalizer bounds setting to use + * @since 3.0 + */ + public StepNormalizer(final double h, final FixedStepHandler handler, + final StepNormalizerBounds bounds) { + this(h, handler, StepNormalizerMode.INCREMENT, bounds); + } + + /** Simple constructor. + * @param h fixed time step (sign is not used) + * @param handler fixed time step handler to wrap + * @param mode step normalizer mode to use + * @param bounds step normalizer bounds setting to use + * @since 3.0 + */ + public StepNormalizer(final double h, final FixedStepHandler handler, + final StepNormalizerMode mode, + final StepNormalizerBounds bounds) { + this.h = FastMath.abs(h); + this.handler = handler; + this.mode = mode; + this.bounds = bounds; + firstTime = Double.NaN; + lastTime = Double.NaN; + lastState = null; + lastDerivatives = null; + forward = true; + } + + /** {@inheritDoc} */ + public void init(double t0, double[] y0, double t) { + + firstTime = Double.NaN; + lastTime = Double.NaN; + lastState = null; + lastDerivatives = null; + forward = true; + + // initialize the underlying handler + handler.init(t0, y0, t); + + } + + /** + * Handle the last accepted step + * @param interpolator interpolator for the last accepted step. For + * efficiency purposes, the various integrators reuse the same + * object on each call, so if the instance wants to keep it across + * all calls (for example to provide at the end of the integration a + * continuous model valid throughout the integration range), it + * should build a local copy using the clone method and store this + * copy. + * @param isLast true if the step is the last one + * @exception MaxCountExceededException if the interpolator throws one because + * the number of functions evaluations is exceeded + */ + public void handleStep(final StepInterpolator interpolator, final boolean isLast) + throws MaxCountExceededException { + // The first time, update the last state with the start information. + if (lastState == null) { + firstTime = interpolator.getPreviousTime(); + lastTime = interpolator.getPreviousTime(); + interpolator.setInterpolatedTime(lastTime); + lastState = interpolator.getInterpolatedState().clone(); + lastDerivatives = interpolator.getInterpolatedDerivatives().clone(); + + // Take the integration direction into account. + forward = interpolator.getCurrentTime() >= lastTime; + if (!forward) { + h = -h; + } + } + + // Calculate next normalized step time. + double nextTime = (mode == StepNormalizerMode.INCREMENT) ? + lastTime + h : + (FastMath.floor(lastTime / h) + 1) * h; + if (mode == StepNormalizerMode.MULTIPLES && + Precision.equals(nextTime, lastTime, 1)) { + nextTime += h; + } + + // Process normalized steps as long as they are in the current step. + boolean nextInStep = isNextInStep(nextTime, interpolator); + while (nextInStep) { + // Output the stored previous step. + doNormalizedStep(false); + + // Store the next step as last step. + storeStep(interpolator, nextTime); + + // Move on to the next step. + nextTime += h; + nextInStep = isNextInStep(nextTime, interpolator); + } + + if (isLast) { + // There will be no more steps. The stored one should be given to + // the handler. We may have to output one more step. Only the last + // one of those should be flagged as being the last. + boolean addLast = bounds.lastIncluded() && + lastTime != interpolator.getCurrentTime(); + doNormalizedStep(!addLast); + if (addLast) { + storeStep(interpolator, interpolator.getCurrentTime()); + doNormalizedStep(true); + } + } + } + + /** + * Returns a value indicating whether the next normalized time is in the + * current step. + * @param nextTime the next normalized time + * @param interpolator interpolator for the last accepted step, to use to + * get the end time of the current step + * @return value indicating whether the next normalized time is in the + * current step + */ + private boolean isNextInStep(double nextTime, + StepInterpolator interpolator) { + return forward ? + nextTime <= interpolator.getCurrentTime() : + nextTime >= interpolator.getCurrentTime(); + } + + /** + * Invokes the underlying step handler for the current normalized step. + * @param isLast true if the step is the last one + */ + private void doNormalizedStep(boolean isLast) { + if (!bounds.firstIncluded() && firstTime == lastTime) { + return; + } + handler.handleStep(lastTime, lastState, lastDerivatives, isLast); + } + + /** Stores the interpolated information for the given time in the current + * state. + * @param interpolator interpolator for the last accepted step, to use to + * get the interpolated information + * @param t the time for which to store the interpolated information + * @exception MaxCountExceededException if the interpolator throws one because + * the number of functions evaluations is exceeded + */ + private void storeStep(StepInterpolator interpolator, double t) + throws MaxCountExceededException { + lastTime = t; + interpolator.setInterpolatedTime(lastTime); + System.arraycopy(interpolator.getInterpolatedState(), 0, + lastState, 0, lastState.length); + System.arraycopy(interpolator.getInterpolatedDerivatives(), 0, + lastDerivatives, 0, lastDerivatives.length); + } +} |