summaryrefslogtreecommitdiff
path: root/src/main/java/org/apache/commons/math/ode/sampling
diff options
context:
space:
mode:
Diffstat (limited to 'src/main/java/org/apache/commons/math/ode/sampling')
-rw-r--r--src/main/java/org/apache/commons/math/ode/sampling/AbstractStepInterpolator.java519
-rw-r--r--src/main/java/org/apache/commons/math/ode/sampling/DummyStepHandler.java101
-rw-r--r--src/main/java/org/apache/commons/math/ode/sampling/DummyStepInterpolator.java150
-rw-r--r--src/main/java/org/apache/commons/math/ode/sampling/FixedStepHandler.java62
-rw-r--r--src/main/java/org/apache/commons/math/ode/sampling/NordsieckStepInterpolator.java291
-rw-r--r--src/main/java/org/apache/commons/math/ode/sampling/StepHandler.java78
-rw-r--r--src/main/java/org/apache/commons/math/ode/sampling/StepInterpolator.java132
-rw-r--r--src/main/java/org/apache/commons/math/ode/sampling/StepNormalizer.java161
-rw-r--r--src/main/java/org/apache/commons/math/ode/sampling/package.html60
9 files changed, 1554 insertions, 0 deletions
diff --git a/src/main/java/org/apache/commons/math/ode/sampling/AbstractStepInterpolator.java b/src/main/java/org/apache/commons/math/ode/sampling/AbstractStepInterpolator.java
new file mode 100644
index 0000000..5cb2979
--- /dev/null
+++ b/src/main/java/org/apache/commons/math/ode/sampling/AbstractStepInterpolator.java
@@ -0,0 +1,519 @@
+/*
+ * 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.sampling;
+
+import java.io.IOException;
+import java.io.ObjectInput;
+import java.io.ObjectOutput;
+
+import org.apache.commons.math.ode.DerivativeException;
+
+/** This abstract class represents an interpolator over the last step
+ * during an ODE integration.
+ *
+ * <p>The various ODE integrators provide objects extending this class
+ * to the step handlers. The handlers can use these objects to
+ * retrieve the state vector at intermediate times between the
+ * previous and the current grid points (dense output).</p>
+ *
+ * @see org.apache.commons.math.ode.FirstOrderIntegrator
+ * @see org.apache.commons.math.ode.SecondOrderIntegrator
+ * @see StepHandler
+ *
+ * @version $Revision: 1073158 $ $Date: 2011-02-21 22:46:52 +0100 (lun. 21 févr. 2011) $
+ * @since 1.2
+ *
+ */
+
+public abstract class AbstractStepInterpolator
+ implements StepInterpolator {
+
+ /** current time step */
+ protected double h;
+
+ /** current state */
+ protected double[] currentState;
+
+ /** interpolated time */
+ protected double interpolatedTime;
+
+ /** interpolated state */
+ protected double[] interpolatedState;
+
+ /** interpolated derivatives */
+ protected double[] interpolatedDerivatives;
+
+ /** global previous time */
+ private double globalPreviousTime;
+
+ /** global current time */
+ private double globalCurrentTime;
+
+ /** soft previous time */
+ private double softPreviousTime;
+
+ /** soft current time */
+ private double softCurrentTime;
+
+ /** indicate if the step has been finalized or not. */
+ private boolean finalized;
+
+ /** integration direction. */
+ private boolean forward;
+
+ /** indicator for dirty state. */
+ private boolean dirtyState;
+
+
+ /** 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. As an example, the {@link
+ * org.apache.commons.math.ode.nonstiff.EmbeddedRungeKuttaIntegrator}
+ * class uses the prototyping design pattern to create the step
+ * interpolators by cloning an uninitialized model and latter
+ * initializing the copy.
+ */
+ protected AbstractStepInterpolator() {
+ globalPreviousTime = Double.NaN;
+ globalCurrentTime = Double.NaN;
+ softPreviousTime = Double.NaN;
+ softCurrentTime = Double.NaN;
+ h = Double.NaN;
+ interpolatedTime = Double.NaN;
+ currentState = null;
+ interpolatedState = null;
+ interpolatedDerivatives = null;
+ finalized = false;
+ this.forward = true;
+ this.dirtyState = true;
+ }
+
+ /** Simple constructor.
+ * @param y reference to the integrator array holding the state at
+ * the end of the step
+ * @param forward integration direction indicator
+ */
+ protected AbstractStepInterpolator(final double[] y, final boolean forward) {
+
+ globalPreviousTime = Double.NaN;
+ globalCurrentTime = Double.NaN;
+ softPreviousTime = Double.NaN;
+ softCurrentTime = Double.NaN;
+ h = Double.NaN;
+ interpolatedTime = Double.NaN;
+
+ currentState = y;
+ interpolatedState = new double[y.length];
+ interpolatedDerivatives = new double[y.length];
+
+ finalized = false;
+ this.forward = forward;
+ this.dirtyState = true;
+
+ }
+
+ /** Copy constructor.
+
+ * <p>The copied interpolator should have been finalized before the
+ * copy, otherwise the copy will not be able to perform correctly
+ * any derivative computation and will throw a {@link
+ * NullPointerException} later. Since we don't want this constructor
+ * to throw the exceptions finalization may involve and since we
+ * don't want this method to modify the state of the copied
+ * interpolator, finalization is <strong>not</strong> done
+ * automatically, it remains under user control.</p>
+ *
+ * <p>The copy is a deep copy: its arrays are separated from the
+ * original arrays of the instance.</p>
+ *
+ * @param interpolator interpolator to copy from.
+ *
+ */
+ protected AbstractStepInterpolator(final AbstractStepInterpolator interpolator) {
+
+ globalPreviousTime = interpolator.globalPreviousTime;
+ globalCurrentTime = interpolator.globalCurrentTime;
+ softPreviousTime = interpolator.softPreviousTime;
+ softCurrentTime = interpolator.softCurrentTime;
+ h = interpolator.h;
+ interpolatedTime = interpolator.interpolatedTime;
+
+ if (interpolator.currentState != null) {
+ currentState = interpolator.currentState.clone();
+ interpolatedState = interpolator.interpolatedState.clone();
+ interpolatedDerivatives = interpolator.interpolatedDerivatives.clone();
+ } else {
+ currentState = null;
+ interpolatedState = null;
+ interpolatedDerivatives = null;
+ }
+
+ finalized = interpolator.finalized;
+ forward = interpolator.forward;
+ dirtyState = interpolator.dirtyState;
+
+ }
+
+ /** Reinitialize the instance
+ * @param y reference to the integrator array holding the state at
+ * the end of the step
+ * @param isForward integration direction indicator
+ */
+ protected void reinitialize(final double[] y, final boolean isForward) {
+
+ globalPreviousTime = Double.NaN;
+ globalCurrentTime = Double.NaN;
+ softPreviousTime = Double.NaN;
+ softCurrentTime = Double.NaN;
+ h = Double.NaN;
+ interpolatedTime = Double.NaN;
+
+ currentState = y;
+ interpolatedState = new double[y.length];
+ interpolatedDerivatives = new double[y.length];
+
+ finalized = false;
+ this.forward = isForward;
+ this.dirtyState = true;
+
+ }
+
+ /** {@inheritDoc} */
+ public StepInterpolator copy() throws DerivativeException {
+
+ // finalize the step before performing copy
+ finalizeStep();
+
+ // create the new independent instance
+ return doCopy();
+
+ }
+
+ /** Really copy the finalized instance.
+ * <p>This method is called by {@link #copy()} after the
+ * step has been finalized. It must perform a deep copy
+ * to have an new instance completely independent for the
+ * original instance.
+ * @return a copy of the finalized instance
+ */
+ protected abstract StepInterpolator doCopy();
+
+ /** Shift one step forward.
+ * Copy the current time into the previous time, hence preparing the
+ * interpolator for future calls to {@link #storeTime storeTime}
+ */
+ public void shift() {
+ globalPreviousTime = globalCurrentTime;
+ softPreviousTime = globalPreviousTime;
+ softCurrentTime = globalCurrentTime;
+ }
+
+ /** Store the current step time.
+ * @param t current time
+ */
+ public void storeTime(final double t) {
+
+ globalCurrentTime = t;
+ softCurrentTime = globalCurrentTime;
+ h = globalCurrentTime - globalPreviousTime;
+ setInterpolatedTime(t);
+
+ // the step is not finalized anymore
+ finalized = false;
+
+ }
+
+ /** Restrict step range to a limited part of the global step.
+ * <p>
+ * This method can be used to restrict a step and make it appear
+ * as if the original step was smaller. Calling this method
+ * <em>only</em> changes the value returned by {@link #getPreviousTime()},
+ * it does not change any other property
+ * </p>
+ * @param softPreviousTime start of the restricted step
+ * @since 2.2
+ */
+ public void setSoftPreviousTime(final double softPreviousTime) {
+ this.softPreviousTime = softPreviousTime;
+ }
+
+ /** Restrict step range to a limited part of the global step.
+ * <p>
+ * This method can be used to restrict a step and make it appear
+ * as if the original step was smaller. Calling this method
+ * <em>only</em> changes the value returned by {@link #getCurrentTime()},
+ * it does not change any other property
+ * </p>
+ * @param softCurrentTime end of the restricted step
+ * @since 2.2
+ */
+ public void setSoftCurrentTime(final double softCurrentTime) {
+ this.softCurrentTime = softCurrentTime;
+ }
+
+ /**
+ * Get the previous global grid point time.
+ * @return previous global grid point time
+ * @since 2.2
+ */
+ public double getGlobalPreviousTime() {
+ return globalPreviousTime;
+ }
+
+ /**
+ * Get the current global grid point time.
+ * @return current global grid point time
+ * @since 2.2
+ */
+ public double getGlobalCurrentTime() {
+ return globalCurrentTime;
+ }
+
+ /**
+ * Get the previous soft grid point time.
+ * @return previous soft grid point time
+ * @see #setSoftPreviousTime(double)
+ */
+ public double getPreviousTime() {
+ return softPreviousTime;
+ }
+
+ /**
+ * Get the current soft grid point time.
+ * @return current soft grid point time
+ * @see #setSoftCurrentTime(double)
+ */
+ public double getCurrentTime() {
+ return softCurrentTime;
+ }
+
+ /** {@inheritDoc} */
+ public double getInterpolatedTime() {
+ return interpolatedTime;
+ }
+
+ /** {@inheritDoc} */
+ public void setInterpolatedTime(final double time) {
+ interpolatedTime = time;
+ dirtyState = true;
+ }
+
+ /** {@inheritDoc} */
+ public boolean isForward() {
+ return forward;
+ }
+
+ /** Compute the state and derivatives at the interpolated time.
+ * This is the main processing method that should be implemented by
+ * the derived classes to perform the interpolation.
+ * @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 abstract void computeInterpolatedStateAndDerivatives(double theta,
+ double oneMinusThetaH)
+ throws DerivativeException;
+
+ /** {@inheritDoc} */
+ public double[] getInterpolatedState() throws DerivativeException {
+
+ // lazy evaluation of the state
+ if (dirtyState) {
+ final double oneMinusThetaH = globalCurrentTime - interpolatedTime;
+ final double theta = (h == 0) ? 0 : (h - oneMinusThetaH) / h;
+ computeInterpolatedStateAndDerivatives(theta, oneMinusThetaH);
+ dirtyState = false;
+ }
+
+ return interpolatedState;
+
+ }
+
+ /** {@inheritDoc} */
+ public double[] getInterpolatedDerivatives() throws DerivativeException {
+
+ // lazy evaluation of the state
+ if (dirtyState) {
+ final double oneMinusThetaH = globalCurrentTime - interpolatedTime;
+ final double theta = (h == 0) ? 0 : (h - oneMinusThetaH) / h;
+ computeInterpolatedStateAndDerivatives(theta, oneMinusThetaH);
+ dirtyState = false;
+ }
+
+ return interpolatedDerivatives;
+
+ }
+
+ /**
+ * Finalize the step.
+ *
+ * <p>Some embedded Runge-Kutta integrators need fewer functions
+ * evaluations than their counterpart step interpolators. These
+ * interpolators should perform the last evaluations they need by
+ * themselves only if they need them. This method triggers these
+ * extra evaluations. It can be called directly by the user step
+ * handler and it is called automatically if {@link
+ * #setInterpolatedTime} is called.</p>
+ *
+ * <p>Once this method has been called, <strong>no</strong> other
+ * evaluation will be performed on this step. If there is a need to
+ * have some side effects between the step handler and the
+ * differential equations (for example update some data in the
+ * equations once the step has been done), it is advised to call
+ * this method explicitly from the step handler before these side
+ * effects are set up. If the step handler induces no side effect,
+ * then this method can safely be ignored, it will be called
+ * transparently as needed.</p>
+ *
+ * <p><strong>Warning</strong>: since the step interpolator provided
+ * to the step handler as a parameter of the {@link
+ * StepHandler#handleStep handleStep} is valid only for the duration
+ * of the {@link StepHandler#handleStep handleStep} call, one cannot
+ * simply store a reference and reuse it later. One should first
+ * finalize the instance, then copy this finalized instance into a
+ * new object that can be kept.</p>
+ *
+ * <p>This method calls the protected <code>doFinalize</code> method
+ * if it has never been called during this step and set a flag
+ * indicating that it has been called once. It is the <code>
+ * doFinalize</code> method which should perform the evaluations.
+ * This wrapping prevents from calling <code>doFinalize</code> several
+ * times and hence evaluating the differential equations too often.
+ * Therefore, subclasses are not allowed not reimplement it, they
+ * should rather reimplement <code>doFinalize</code>.</p>
+ *
+ * @throws DerivativeException this exception is propagated to the
+ * caller if the underlying user function triggers one
+ */
+ public final void finalizeStep()
+ throws DerivativeException {
+ if (! finalized) {
+ doFinalize();
+ finalized = true;
+ }
+ }
+
+ /**
+ * Really finalize the step.
+ * The default implementation of this method does nothing.
+ * @throws DerivativeException this exception is propagated to the
+ * caller if the underlying user function triggers one
+ */
+ protected void doFinalize()
+ throws DerivativeException {
+ }
+
+ /** {@inheritDoc} */
+ public abstract void writeExternal(ObjectOutput out)
+ throws IOException;
+
+ /** {@inheritDoc} */
+ public abstract void readExternal(ObjectInput in)
+ throws IOException, ClassNotFoundException;
+
+ /** Save the base state of the instance.
+ * This method performs step finalization if it has not been done
+ * before.
+ * @param out stream where to save the state
+ * @exception IOException in case of write error
+ */
+ protected void writeBaseExternal(final ObjectOutput out)
+ throws IOException {
+
+ if (currentState == null) {
+ out.writeInt(-1);
+ } else {
+ out.writeInt(currentState.length);
+ }
+ out.writeDouble(globalPreviousTime);
+ out.writeDouble(globalCurrentTime);
+ out.writeDouble(softPreviousTime);
+ out.writeDouble(softCurrentTime);
+ out.writeDouble(h);
+ out.writeBoolean(forward);
+
+ if (currentState != null) {
+ for (int i = 0; i < currentState.length; ++i) {
+ out.writeDouble(currentState[i]);
+ }
+ }
+
+ out.writeDouble(interpolatedTime);
+
+ // we do not store the interpolated state,
+ // it will be recomputed as needed after reading
+
+ // finalize the step (and don't bother saving the now true flag)
+ try {
+ finalizeStep();
+ } catch (DerivativeException e) {
+ IOException ioe = new IOException(e.getLocalizedMessage());
+ ioe.initCause(e);
+ throw ioe;
+ }
+
+ }
+
+ /** Read the base state of the instance.
+ * This method does <strong>neither</strong> set the interpolated
+ * time nor state. It is up to the derived class to reset it
+ * properly calling the {@link #setInterpolatedTime} method later,
+ * once all rest of the object state has been set up properly.
+ * @param in stream where to read the state from
+ * @return interpolated time be set later by the caller
+ * @exception IOException in case of read error
+ */
+ protected double readBaseExternal(final ObjectInput in)
+ throws IOException {
+
+ final int dimension = in.readInt();
+ globalPreviousTime = in.readDouble();
+ globalCurrentTime = in.readDouble();
+ softPreviousTime = in.readDouble();
+ softCurrentTime = in.readDouble();
+ h = in.readDouble();
+ forward = in.readBoolean();
+ dirtyState = true;
+
+ if (dimension < 0) {
+ currentState = null;
+ } else {
+ currentState = new double[dimension];
+ for (int i = 0; i < currentState.length; ++i) {
+ currentState[i] = in.readDouble();
+ }
+ }
+
+ // we do NOT handle the interpolated time and state here
+ interpolatedTime = Double.NaN;
+ interpolatedState = (dimension < 0) ? null : new double[dimension];
+ interpolatedDerivatives = (dimension < 0) ? null : new double[dimension];
+
+ finalized = true;
+
+ return in.readDouble();
+
+ }
+
+}
diff --git a/src/main/java/org/apache/commons/math/ode/sampling/DummyStepHandler.java b/src/main/java/org/apache/commons/math/ode/sampling/DummyStepHandler.java
new file mode 100644
index 0000000..cc759b9
--- /dev/null
+++ b/src/main/java/org/apache/commons/math/ode/sampling/DummyStepHandler.java
@@ -0,0 +1,101 @@
+/*
+ * 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.sampling;
+
+/**
+ * This class is a step handler that does nothing.
+
+ * <p>This class is provided as a convenience for users who are only
+ * interested in the final state of an integration and not in the
+ * intermediate steps. Its handleStep method does nothing.</p>
+ *
+ * <p>Since this class has no internal state, it is implemented using
+ * the Singleton design pattern. This means that only one instance is
+ * ever created, which can be retrieved using the getInstance
+ * method. This explains why there is no public constructor.</p>
+ *
+ * @see StepHandler
+ * @version $Revision: 811827 $ $Date: 2009-09-06 17:32:50 +0200 (dim. 06 sept. 2009) $
+ * @since 1.2
+ */
+
+public class DummyStepHandler implements StepHandler {
+
+ /** Private constructor.
+ * The constructor is private to prevent users from creating
+ * instances (Singleton design-pattern).
+ */
+ private DummyStepHandler() {
+ }
+
+ /** Get the only instance.
+ * @return the only instance
+ */
+ public static DummyStepHandler getInstance() {
+ return LazyHolder.INSTANCE;
+ }
+
+ /** Determines whether this handler needs dense output.
+ * Since this handler does nothing, it does not require dense output.
+ * @return always false
+ */
+ public boolean requiresDenseOutput() {
+ return false;
+ }
+
+ /** Reset the step handler.
+ * Initialize the internal data as required before the first step is
+ * handled.
+ */
+ public void reset() {
+ }
+
+ /**
+ * Handle the last accepted step.
+ * This method does nothing in this class.
+ * @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
+ */
+ public void handleStep(final StepInterpolator interpolator, final boolean isLast) {
+ }
+
+ // CHECKSTYLE: stop HideUtilityClassConstructor
+ /** Holder for the instance.
+ * <p>We use here the Initialization On Demand Holder Idiom.</p>
+ */
+ private static class LazyHolder {
+ /** Cached field instance. */
+ private static final DummyStepHandler INSTANCE = new DummyStepHandler();
+ }
+ // CHECKSTYLE: resume HideUtilityClassConstructor
+
+ /** Handle deserialization of the singleton.
+ * @return the singleton instance
+ */
+ private Object readResolve() {
+ // return the singleton instance
+ return LazyHolder.INSTANCE;
+ }
+
+}
diff --git a/src/main/java/org/apache/commons/math/ode/sampling/DummyStepInterpolator.java b/src/main/java/org/apache/commons/math/ode/sampling/DummyStepInterpolator.java
new file mode 100644
index 0000000..af28c13
--- /dev/null
+++ b/src/main/java/org/apache/commons/math/ode/sampling/DummyStepInterpolator.java
@@ -0,0 +1,150 @@
+/*
+ * 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.sampling;
+
+import java.io.IOException;
+import java.io.ObjectInput;
+import java.io.ObjectOutput;
+
+/** This class is a step interpolator that does nothing.
+ *
+ * <p>This class is used when the {@link StepHandler "step handler"}
+ * set up by the user does not need step interpolation. It does not
+ * recompute the state when {@link AbstractStepInterpolator#setInterpolatedTime
+ * setInterpolatedTime} is called. This implies the interpolated state
+ * is always the state at the end of the current step.</p>
+ *
+ * @see StepHandler
+ *
+ * @version $Revision: 1037327 $ $Date: 2010-11-20 21:57:37 +0100 (sam. 20 nov. 2010) $
+ * @since 1.2
+ */
+
+public class DummyStepInterpolator
+ extends AbstractStepInterpolator {
+
+ /** Serializable version identifier. */
+ private static final long serialVersionUID = 1708010296707839488L;
+
+ /** Current derivative. */
+ private double[] currentDerivative;
+
+ /** Simple constructor.
+ * This constructor builds an instance that is not usable yet, the
+ * <code>AbstractStepInterpolator.reinitialize</code> protected 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. As an example, the {@link
+ * org.apache.commons.math.ode.nonstiff.EmbeddedRungeKuttaIntegrator} uses
+ * the prototyping design pattern to create the step interpolators by
+ * cloning an uninitialized model and latter initializing the copy.
+ */
+ public DummyStepInterpolator() {
+ super();
+ currentDerivative = null;
+ }
+
+ /** Simple constructor.
+ * @param y reference to the integrator array holding the state at
+ * the end of the step
+ * @param yDot reference to the integrator array holding the state
+ * derivative at some arbitrary point within the step
+ * @param forward integration direction indicator
+ */
+ public DummyStepInterpolator(final double[] y, final double[] yDot, final boolean forward) {
+ super(y, forward);
+ currentDerivative = yDot;
+ }
+
+ /** 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 DummyStepInterpolator(final DummyStepInterpolator interpolator) {
+ super(interpolator);
+ currentDerivative = interpolator.currentDerivative.clone();
+ }
+
+ /** Really copy the finalized instance.
+ * @return a copy of the finalized instance
+ */
+ @Override
+ protected StepInterpolator doCopy() {
+ return new DummyStepInterpolator(this);
+ }
+
+ /** Compute the state at the interpolated time.
+ * In this class, this method does nothing: the interpolated state
+ * is always the state at the end of the current step.
+ * @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
+ */
+ @Override
+ protected void computeInterpolatedStateAndDerivatives(final double theta, final double oneMinusThetaH) {
+ System.arraycopy(currentState, 0, interpolatedState, 0, currentState.length);
+ System.arraycopy(currentDerivative, 0, interpolatedDerivatives, 0, currentDerivative.length);
+ }
+
+ /** Write the instance to an output channel.
+ * @param out output channel
+ * @exception IOException if the instance cannot be written
+ */
+ @Override
+ public void writeExternal(final ObjectOutput out)
+ throws IOException {
+
+ // save the state of the base class
+ writeBaseExternal(out);
+
+ if (currentDerivative != null) {
+ for (int i = 0; i < currentDerivative.length; ++i) {
+ out.writeDouble(currentDerivative[i]);
+ }
+ }
+
+ }
+
+ /** Read the instance from an input channel.
+ * @param in input channel
+ * @exception IOException if the instance cannot be read
+ */
+ @Override
+ public void readExternal(final ObjectInput in)
+ throws IOException {
+
+ // read the base class
+ final double t = readBaseExternal(in);
+
+ if (currentState == null) {
+ currentDerivative = null;
+ } else {
+ currentDerivative = new double[currentState.length];
+ for (int i = 0; i < currentDerivative.length; ++i) {
+ currentDerivative[i] = in.readDouble();
+ }
+ }
+
+ // we can now set the interpolated time and state
+ setInterpolatedTime(t);
+
+ }
+
+}
diff --git a/src/main/java/org/apache/commons/math/ode/sampling/FixedStepHandler.java b/src/main/java/org/apache/commons/math/ode/sampling/FixedStepHandler.java
new file mode 100644
index 0000000..190c999
--- /dev/null
+++ b/src/main/java/org/apache/commons/math/ode/sampling/FixedStepHandler.java
@@ -0,0 +1,62 @@
+/*
+ * 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.sampling;
+
+import org.apache.commons.math.ode.DerivativeException;
+
+/**
+ * This interface represents a handler that should be called after
+ * each successful fixed step.
+
+ * <p>This interface should be implemented by anyone who is interested
+ * in getting the solution of an ordinary differential equation at
+ * fixed time steps. Objects implementing this interface should be
+ * wrapped within an instance of {@link StepNormalizer} that itself
+ * is used as the general {@link StepHandler} by the integrator. The
+ * {@link StepNormalizer} object is called according to the integrator
+ * internal algorithms and it calls objects implementing this
+ * interface as necessary at fixed time steps.</p>
+ *
+ * @see StepHandler
+ * @see StepNormalizer
+ * @version $Revision: 1073158 $ $Date: 2011-02-21 22:46:52 +0100 (lun. 21 févr. 2011) $
+ * @since 1.2
+ */
+
+public interface FixedStepHandler {
+
+ /**
+ * Handle the last accepted step
+ * @param t time of the current step
+ * @param y state vector at t. For efficiency purposes, the {@link
+ * StepNormalizer} class reuses the same array 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 complete array of all
+ * steps), it should build a local copy store this copy.
+ * @param yDot derivatives of the state vector state vector at t.
+ * For efficiency purposes, the {@link StepNormalizer} class reuses
+ * the same array 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 complete array of all
+ * steps), it should build a local copy store this copy.
+ * @param isLast true if the step is the last one
+ * @throws DerivativeException if some error condition is encountered
+ */
+ void handleStep(double t, double[] y, double[] yDot, boolean isLast) throws DerivativeException;
+
+}
diff --git a/src/main/java/org/apache/commons/math/ode/sampling/NordsieckStepInterpolator.java b/src/main/java/org/apache/commons/math/ode/sampling/NordsieckStepInterpolator.java
new file mode 100644
index 0000000..2090276
--- /dev/null
+++ b/src/main/java/org/apache/commons/math/ode/sampling/NordsieckStepInterpolator.java
@@ -0,0 +1,291 @@
+/*
+ * 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.sampling;
+
+import java.io.IOException;
+import java.io.ObjectInput;
+import java.io.ObjectOutput;
+import java.util.Arrays;
+
+import org.apache.commons.math.ode.DerivativeException;
+import org.apache.commons.math.linear.Array2DRowRealMatrix;
+import org.apache.commons.math.util.FastMath;
+
+/**
+ * This class implements an interpolator for integrators using Nordsieck representation.
+ *
+ * <p>This interpolator computes dense output around the current point.
+ * The interpolation equation is based on Taylor series formulas.
+ *
+ * @see org.apache.commons.math.ode.nonstiff.AdamsBashforthIntegrator
+ * @see org.apache.commons.math.ode.nonstiff.AdamsMoultonIntegrator
+ * @version $Revision: 1073158 $ $Date: 2011-02-21 22:46:52 +0100 (lun. 21 févr. 2011) $
+ * @since 2.0
+ */
+
+public class NordsieckStepInterpolator extends AbstractStepInterpolator {
+
+ /** Serializable version identifier */
+ private static final long serialVersionUID = -7179861704951334960L;
+
+ /** State variation. */
+ protected double[] stateVariation;
+
+ /** Step size used in the first scaled derivative and Nordsieck vector. */
+ private double scalingH;
+
+ /** Reference time for all arrays.
+ * <p>Sometimes, the reference time is the same as previousTime,
+ * sometimes it is the same as currentTime, so we use a separate
+ * field to avoid any confusion.
+ * </p>
+ */
+ private double referenceTime;
+
+ /** First scaled derivative. */
+ private double[] scaled;
+
+ /** Nordsieck vector. */
+ private Array2DRowRealMatrix nordsieck;
+
+ /** Simple constructor.
+ * This constructor builds an instance that is not usable yet, the
+ * {@link AbstractStepInterpolator#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.
+ */
+ public NordsieckStepInterpolator() {
+ }
+
+ /** 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 NordsieckStepInterpolator(final NordsieckStepInterpolator interpolator) {
+ super(interpolator);
+ scalingH = interpolator.scalingH;
+ referenceTime = interpolator.referenceTime;
+ if (interpolator.scaled != null) {
+ scaled = interpolator.scaled.clone();
+ }
+ if (interpolator.nordsieck != null) {
+ nordsieck = new Array2DRowRealMatrix(interpolator.nordsieck.getDataRef(), true);
+ }
+ if (interpolator.stateVariation != null) {
+ stateVariation = interpolator.stateVariation.clone();
+ }
+ }
+
+ /** {@inheritDoc} */
+ @Override
+ protected StepInterpolator doCopy() {
+ return new NordsieckStepInterpolator(this);
+ }
+
+ /** Reinitialize the instance.
+ * <p>Beware that all arrays <em>must</em> be references to integrator
+ * arrays, in order to ensure proper update without copy.</p>
+ * @param y reference to the integrator array holding the state at
+ * the end of the step
+ * @param forward integration direction indicator
+ */
+ @Override
+ public void reinitialize(final double[] y, final boolean forward) {
+ super.reinitialize(y, forward);
+ stateVariation = new double[y.length];
+ }
+
+ /** Reinitialize the instance.
+ * <p>Beware that all arrays <em>must</em> be references to integrator
+ * arrays, in order to ensure proper update without copy.</p>
+ * @param time time at which all arrays are defined
+ * @param stepSize step size used in the scaled and nordsieck arrays
+ * @param scaledDerivative reference to the integrator array holding the first
+ * scaled derivative
+ * @param nordsieckVector reference to the integrator matrix holding the
+ * nordsieck vector
+ */
+ public void reinitialize(final double time, final double stepSize,
+ final double[] scaledDerivative,
+ final Array2DRowRealMatrix nordsieckVector) {
+ this.referenceTime = time;
+ this.scalingH = stepSize;
+ this.scaled = scaledDerivative;
+ this.nordsieck = nordsieckVector;
+
+ // make sure the state and derivatives will depend on the new arrays
+ setInterpolatedTime(getInterpolatedTime());
+
+ }
+
+ /** Rescale the instance.
+ * <p>Since the scaled and Nordiseck arrays are shared with the caller,
+ * this method has the side effect of rescaling this arrays in the caller too.</p>
+ * @param stepSize new step size to use in the scaled and nordsieck arrays
+ */
+ public void rescale(final double stepSize) {
+
+ final double ratio = stepSize / scalingH;
+ for (int i = 0; i < scaled.length; ++i) {
+ scaled[i] *= ratio;
+ }
+
+ final double[][] nData = nordsieck.getDataRef();
+ double power = ratio;
+ for (int i = 0; i < nData.length; ++i) {
+ power *= ratio;
+ final double[] nDataI = nData[i];
+ for (int j = 0; j < nDataI.length; ++j) {
+ nDataI[j] *= power;
+ }
+ }
+
+ scalingH = stepSize;
+
+ }
+
+ /**
+ * Get the state vector variation from current to interpolated state.
+ * <p>This method is aimed at computing y(t<sub>interpolation</sub>)
+ * -y(t<sub>current</sub>) accurately by avoiding the cancellation errors
+ * that would occur if the subtraction were performed explicitly.</p>
+ * <p>The returned vector is a reference to a reused array, so
+ * it should not be modified and it should be copied if it needs
+ * to be preserved across several calls.</p>
+ * @return state vector at time {@link #getInterpolatedTime}
+ * @see #getInterpolatedDerivatives()
+ * @throws DerivativeException if this call induces an automatic
+ * step finalization that throws one
+ */
+ public double[] getInterpolatedStateVariation()
+ throws DerivativeException {
+ // compute and ignore interpolated state
+ // to make sure state variation is computed as a side effect
+ getInterpolatedState();
+ return stateVariation;
+ }
+
+ /** {@inheritDoc} */
+ @Override
+ protected void computeInterpolatedStateAndDerivatives(final double theta, final double oneMinusThetaH) {
+
+ final double x = interpolatedTime - referenceTime;
+ final double normalizedAbscissa = x / scalingH;
+
+ Arrays.fill(stateVariation, 0.0);
+ Arrays.fill(interpolatedDerivatives, 0.0);
+
+ // apply Taylor formula from high order to low order,
+ // for the sake of numerical accuracy
+ final double[][] nData = nordsieck.getDataRef();
+ for (int i = nData.length - 1; i >= 0; --i) {
+ final int order = i + 2;
+ final double[] nDataI = nData[i];
+ final double power = FastMath.pow(normalizedAbscissa, order);
+ for (int j = 0; j < nDataI.length; ++j) {
+ final double d = nDataI[j] * power;
+ stateVariation[j] += d;
+ interpolatedDerivatives[j] += order * d;
+ }
+ }
+
+ for (int j = 0; j < currentState.length; ++j) {
+ stateVariation[j] += scaled[j] * normalizedAbscissa;
+ interpolatedState[j] = currentState[j] + stateVariation[j];
+ interpolatedDerivatives[j] =
+ (interpolatedDerivatives[j] + scaled[j] * normalizedAbscissa) / x;
+ }
+
+ }
+
+ /** {@inheritDoc} */
+ @Override
+ public void writeExternal(final ObjectOutput out)
+ throws IOException {
+
+ // save the state of the base class
+ writeBaseExternal(out);
+
+ // save the local attributes
+ out.writeDouble(scalingH);
+ out.writeDouble(referenceTime);
+
+ final int n = (currentState == null) ? -1 : currentState.length;
+ if (scaled == null) {
+ out.writeBoolean(false);
+ } else {
+ out.writeBoolean(true);
+ for (int j = 0; j < n; ++j) {
+ out.writeDouble(scaled[j]);
+ }
+ }
+
+ if (nordsieck == null) {
+ out.writeBoolean(false);
+ } else {
+ out.writeBoolean(true);
+ out.writeObject(nordsieck);
+ }
+
+ // we don't save state variation, it will be recomputed
+
+ }
+
+ /** {@inheritDoc} */
+ @Override
+ public void readExternal(final ObjectInput in)
+ throws IOException, ClassNotFoundException {
+
+ // read the base class
+ final double t = readBaseExternal(in);
+
+ // read the local attributes
+ scalingH = in.readDouble();
+ referenceTime = in.readDouble();
+
+ final int n = (currentState == null) ? -1 : currentState.length;
+ final boolean hasScaled = in.readBoolean();
+ if (hasScaled) {
+ scaled = new double[n];
+ for (int j = 0; j < n; ++j) {
+ scaled[j] = in.readDouble();
+ }
+ } else {
+ scaled = null;
+ }
+
+ final boolean hasNordsieck = in.readBoolean();
+ if (hasNordsieck) {
+ nordsieck = (Array2DRowRealMatrix) in.readObject();
+ } else {
+ nordsieck = null;
+ }
+
+ if (hasScaled && hasNordsieck) {
+ // we can now set the interpolated time and state
+ stateVariation = new double[n];
+ setInterpolatedTime(t);
+ } else {
+ stateVariation = null;
+ }
+
+ }
+
+}
diff --git a/src/main/java/org/apache/commons/math/ode/sampling/StepHandler.java b/src/main/java/org/apache/commons/math/ode/sampling/StepHandler.java
new file mode 100644
index 0000000..53c2fbe
--- /dev/null
+++ b/src/main/java/org/apache/commons/math/ode/sampling/StepHandler.java
@@ -0,0 +1,78 @@
+/*
+ * 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.sampling;
+
+import org.apache.commons.math.ode.DerivativeException;
+
+/**
+ * This interface represents a handler that should be called after
+ * each successful step.
+ *
+ * <p>The ODE integrators compute the evolution of the state vector at
+ * some grid points that depend on their own internal algorithm. Once
+ * they have found a new grid point (possibly after having computed
+ * several evaluation of the derivative at intermediate points), they
+ * provide it to objects implementing this interface. These objects
+ * typically either ignore the intermediate steps and wait for the
+ * last one, store the points in an ephemeris, or forward them to
+ * specialized processing or output methods.</p>
+ *
+ * @see org.apache.commons.math.ode.FirstOrderIntegrator
+ * @see org.apache.commons.math.ode.SecondOrderIntegrator
+ * @see StepInterpolator
+ * @version $Revision: 1073158 $ $Date: 2011-02-21 22:46:52 +0100 (lun. 21 févr. 2011) $
+ * @since 1.2
+ */
+
+public interface StepHandler {
+
+ /** Determines whether this handler needs dense output.
+ * <p>This method allows the integrator to avoid performing extra
+ * computation if the handler does not need dense output. If this
+ * method returns false, the integrator will call the {@link
+ * #handleStep} method with a {@link DummyStepInterpolator} rather
+ * than a custom interpolator.</p>
+ * @return true if the handler needs dense output
+ */
+ boolean requiresDenseOutput();
+
+ /** Reset the step handler.
+ * Initialize the internal data as required before the first step is
+ * handled.
+ */
+ void reset();
+
+ /**
+ * 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, as the
+ * {@link org.apache.commons.math.ode.ContinuousOutputModel
+ * ContinuousOutputModel} class does), it should build a local copy
+ * using the clone method of the interpolator and store this copy.
+ * Keeping only a reference to the interpolator and reusing it will
+ * result in unpredictable behavior (potentially crashing the application).
+ * @param isLast true if the step is the last one
+ * @exception DerivativeException if user code called from step interpolator
+ * finalization triggers one
+ */
+ void handleStep(StepInterpolator interpolator, boolean isLast) throws DerivativeException;
+
+}
diff --git a/src/main/java/org/apache/commons/math/ode/sampling/StepInterpolator.java b/src/main/java/org/apache/commons/math/ode/sampling/StepInterpolator.java
new file mode 100644
index 0000000..ceb5a86
--- /dev/null
+++ b/src/main/java/org/apache/commons/math/ode/sampling/StepInterpolator.java
@@ -0,0 +1,132 @@
+/*
+ * 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.sampling;
+
+import java.io.Externalizable;
+
+import org.apache.commons.math.ode.DerivativeException;
+
+/** This interface represents an interpolator over the last step
+ * during an ODE integration.
+ *
+ * <p>The various ODE integrators provide objects implementing this
+ * interface to the step handlers. These objects are often custom
+ * objects tightly bound to the integrator internal algorithms. The
+ * handlers can use these objects to retrieve the state vector at
+ * intermediate times between the previous and the current grid points
+ * (this feature is often called dense output).</p>
+ * <p>One important thing to note is that the step handlers may be so
+ * tightly bound to the integrators that they often share some internal
+ * state arrays. This imply that one should <em>never</em> use a direct
+ * reference to a step interpolator outside of the step handler, either
+ * for future use or for use in another thread. If such a need arise, the
+ * step interpolator <em>must</em> be copied using the dedicated
+ * {@link #copy()} method.
+ * </p>
+ *
+ * @see org.apache.commons.math.ode.FirstOrderIntegrator
+ * @see org.apache.commons.math.ode.SecondOrderIntegrator
+ * @see StepHandler
+ * @version $Revision: 1073158 $ $Date: 2011-02-21 22:46:52 +0100 (lun. 21 févr. 2011) $
+ * @since 1.2
+ */
+
+public interface StepInterpolator extends Externalizable {
+
+ /**
+ * Get the previous grid point time.
+ * @return previous grid point time
+ */
+ double getPreviousTime();
+
+ /**
+ * Get the current grid point time.
+ * @return current grid point time
+ */
+ double getCurrentTime();
+
+ /**
+ * Get the time of the interpolated point.
+ * If {@link #setInterpolatedTime} has not been called, it returns
+ * the current grid point time.
+ * @return interpolation point time
+ */
+ double getInterpolatedTime();
+
+ /**
+ * Set the time of the interpolated point.
+ * <p>Setting the time outside of the current step is now allowed, but
+ * should be used with care since the accuracy of the interpolator will
+ * probably be very poor far from this step. This allowance has been
+ * added to simplify implementation of search algorithms near the
+ * step endpoints.</p>
+ * <p>Setting the time changes the instance internal state. If a
+ * specific state must be preserved, a copy of the instance must be
+ * created using {@link #copy()}.</p>
+ * @param time time of the interpolated point
+ */
+ void setInterpolatedTime(double time);
+
+ /**
+ * Get the state vector of the interpolated point.
+ * <p>The returned vector is a reference to a reused array, so
+ * it should not be modified and it should be copied if it needs
+ * to be preserved across several calls.</p>
+ * @return state vector at time {@link #getInterpolatedTime}
+ * @see #getInterpolatedDerivatives()
+ * @exception DerivativeException if user code called from step interpolator
+ * finalization triggers one
+ */
+ double[] getInterpolatedState() throws DerivativeException;
+
+ /**
+ * Get the derivatives of the state vector of the interpolated point.
+ * <p>The returned vector is a reference to a reused array, so
+ * it should not be modified and it should be copied if it needs
+ * to be preserved across several calls.</p>
+ * @return derivatives of the state vector at time {@link #getInterpolatedTime}
+ * @see #getInterpolatedState()
+ * @exception DerivativeException if user code called from step interpolator
+ * finalization triggers one
+ * @since 2.0
+ */
+ double[] getInterpolatedDerivatives() throws DerivativeException;
+
+ /** Check if the natural integration direction is forward.
+ * <p>This method provides the integration direction as specified by
+ * the integrator itself, it avoid some nasty problems in
+ * degenerated cases like null steps due to cancellation at step
+ * initialization, step control or discrete events
+ * triggering.</p>
+ * @return true if the integration variable (time) increases during
+ * integration
+ */
+ boolean isForward();
+
+ /** Copy the instance.
+ * <p>The copied instance is guaranteed to be independent from the
+ * original one. Both can be used with different settings for
+ * interpolated time without any side effect.</p>
+ * @return a deep copy of the instance, which can be used independently.
+ * @exception DerivativeException if user code called from step interpolator
+ * finalization triggers one
+ * @see #setInterpolatedTime(double)
+ */
+ StepInterpolator copy() throws DerivativeException;
+
+}
diff --git a/src/main/java/org/apache/commons/math/ode/sampling/StepNormalizer.java b/src/main/java/org/apache/commons/math/ode/sampling/StepNormalizer.java
new file mode 100644
index 0000000..73bcc23
--- /dev/null
+++ b/src/main/java/org/apache/commons/math/ode/sampling/StepNormalizer.java
@@ -0,0 +1,161 @@
+/*
+ * 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.sampling;
+
+import org.apache.commons.math.ode.DerivativeException;
+import org.apache.commons.math.util.FastMath;
+
+/**
+ * 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 the beginning time of
+ * the integration t0 and also at times t0+h, t0+2h, ... If the
+ * integration range is an integer multiple of the stepsize, then the
+ * last point handled will be the endpoint of the integration tend, if
+ * not, the last point will belong to the interval [tend - h ;
+ * tend].</p>
+ *
+ * <p>There is no constraint on the integrator, it can use any
+ * timestep it needs (time steps longer or shorter than the fixed time
+ * step and non-integer ratios are all allowed).</p>
+ *
+ * @see StepHandler
+ * @see FixedStepHandler
+ * @version $Revision: 1073158 $ $Date: 2011-02-21 22:46:52 +0100 (lun. 21 févr. 2011) $
+ * @since 1.2
+ */
+
+public class StepNormalizer implements StepHandler {
+
+ /** Fixed time step. */
+ private double h;
+
+ /** Underlying step handler. */
+ private final FixedStepHandler handler;
+
+ /** Last step time. */
+ private double lastTime;
+
+ /** Last State vector. */
+ private double[] lastState;
+
+ /** Last Derivatives vector. */
+ private double[] lastDerivatives;
+
+ /** Integration direction indicator. */
+ private boolean forward;
+
+ /** Simple constructor.
+ * @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 = FastMath.abs(h);
+ this.handler = handler;
+ reset();
+ }
+
+ /** Determines whether this handler needs dense output.
+ * This handler needs dense output in order to provide data at
+ * regularly spaced steps regardless of the steps the integrator
+ * uses, so this method always returns true.
+ * @return always true
+ */
+ public boolean requiresDenseOutput() {
+ return true;
+ }
+
+ /** Reset the step handler.
+ * Initialize the internal data as required before the first step is
+ * handled.
+ */
+ public void reset() {
+ lastTime = Double.NaN;
+ lastState = null;
+ lastDerivatives = null;
+ forward = true;
+ }
+
+ /**
+ * 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
+ * @throws DerivativeException this exception is propagated to the
+ * caller if the underlying user function triggers one
+ */
+ public void handleStep(final StepInterpolator interpolator, final boolean isLast)
+ throws DerivativeException {
+
+ if (lastState == null) {
+
+ 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;
+ }
+
+ }
+
+ double nextTime = lastTime + h;
+ boolean nextInStep = forward ^ (nextTime > interpolator.getCurrentTime());
+ while (nextInStep) {
+
+ // output the stored previous step
+ handler.handleStep(lastTime, lastState, lastDerivatives, false);
+
+ // store the next step
+ lastTime = nextTime;
+ interpolator.setInterpolatedTime(lastTime);
+ System.arraycopy(interpolator.getInterpolatedState(), 0,
+ lastState, 0, lastState.length);
+ System.arraycopy(interpolator.getInterpolatedDerivatives(), 0,
+ lastDerivatives, 0, lastDerivatives.length);
+
+ nextTime += h;
+ nextInStep = forward ^ (nextTime > interpolator.getCurrentTime());
+
+ }
+
+ if (isLast) {
+ // there will be no more steps,
+ // the stored one should be flagged as being the last
+ handler.handleStep(lastTime, lastState, lastDerivatives, true);
+ }
+
+ }
+
+}
diff --git a/src/main/java/org/apache/commons/math/ode/sampling/package.html b/src/main/java/org/apache/commons/math/ode/sampling/package.html
new file mode 100644
index 0000000..46de4ef
--- /dev/null
+++ b/src/main/java/org/apache/commons/math/ode/sampling/package.html
@@ -0,0 +1,60 @@
+<html>
+<!--
+ 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.
+ -->
+ <!-- $Revision: 613620 $ -->
+<body>
+<p>
+This package provides classes to handle sampling steps during
+Ordinary Differential Equations integration.
+</p>
+
+<p>
+In addition to computing the evolution of the state vector at some grid points, all
+ODE integrators also build up interpolation models of this evolution <em>inside</em> the
+last computed step. If users are interested in these interpolators, they can register a
+{@link org.apache.commons.math.ode.sampling.StepHandler StepHandler} instance using the
+{@link org.apache.commons.math.ode.FirstOrderIntegrator#addStepHandler addStepHandler}
+method which is supported by all integrators. The integrator will call this instance
+at the end of each accepted step and provide it the interpolator. The user can do
+whatever he wants with this interpolator, which computes both the state and its
+time-derivative. A typical use of step handler is to provide some output to monitor
+the integration process.
+</p>
+
+<p>
+In a sense, this is a kind of Inversion Of Control: rather than having the master
+application driving the slave integrator by providing the target end value for
+the free variable, we get a master integrator scheduling the free variable
+evolution and calling the slave application callbacks that were registered at
+configuration time.
+</p>
+
+<p>
+Since some integrators may use variable step size, the generic {@link
+org.apache.commons.math.ode.sampling.StepHandler StepHandler} interface can be called
+either at regular or irregular rate. This interface allows to navigate to any location
+within the last computed step, thanks to the provided {@link
+org.apache.commons.math.ode.sampling.StepInterpolator StepInterpolator} object.
+If regular output is desired (for example in order to write an ephemeris file), then
+the simpler {@link org.apache.commons.math.ode.sampling.FixedStepHandler FixedStepHandler}
+interface can be used. Objects implementing this interface should be wrapped within a
+{@link org.apache.commons.math.ode.sampling.StepNormalizer StepNormalizer} instance
+in order to be registered to the integrator.
+</p>
+
+</body>
+</html>