diff --git a/src/main/java/org/apache/commons/math/ode/events/CombinedEventsManager.java b/src/main/java/org/apache/commons/math/ode/events/CombinedEventsManager.java
index 0831e4ab8..08d873862 100644
--- a/src/main/java/org/apache/commons/math/ode/events/CombinedEventsManager.java
+++ b/src/main/java/org/apache/commons/math/ode/events/CombinedEventsManager.java
@@ -135,11 +135,8 @@ public class CombinedEventsManager {
if (! initialized) {
// initialize the events states
- final double t0 = interpolator.getPreviousTime();
- interpolator.setInterpolatedTime(t0);
- final double [] y = interpolator.getInterpolatedState();
for (EventState state : states) {
- state.reinitializeBegin(t0, y);
+ state.reinitializeBegin(interpolator);
}
initialized = true;
@@ -170,6 +167,10 @@ public class CombinedEventsManager {
return first != null;
} catch (EventException se) {
+ final Throwable cause = se.getCause();
+ if ((cause != null) && (cause instanceof DerivativeException)) {
+ throw (DerivativeException) cause;
+ }
throw new IntegratorException(se);
} catch (ConvergenceException ce) {
throw new IntegratorException(ce);
diff --git a/src/main/java/org/apache/commons/math/ode/events/EventState.java b/src/main/java/org/apache/commons/math/ode/events/EventState.java
index 173265084..e8c10f285 100644
--- a/src/main/java/org/apache/commons/math/ode/events/EventState.java
+++ b/src/main/java/org/apache/commons/math/ode/events/EventState.java
@@ -140,18 +140,49 @@ public class EventState {
}
/** Reinitialize the beginning of the step.
- * @param tStart value of the independent time variable at the
- * beginning of the step
- * @param yStart array containing the current value of the state vector
- * at the beginning of the step
+ * @param interpolator valid for the current step
* @exception EventException if the event handler
* value cannot be evaluated at the beginning of the step
*/
- public void reinitializeBegin(final double tStart, final double[] yStart)
+ public void reinitializeBegin(final StepInterpolator interpolator)
throws EventException {
- t0 = tStart;
- g0 = handler.g(tStart, yStart);
- g0Positive = g0 >= 0;
+ try {
+ // excerpt from MATH-421 issue:
+ // If an ODE solver is setup with an EventHandler that return STOP
+ // when the even is triggered, the integrator stops (which is exactly
+ // the expected behavior). If however the user want to restart the
+ // solver from the final state reached at the event with the same
+ // configuration (expecting the event to be triggered again at a
+ // later time), then the integrator may fail to start. It can get stuck
+ // at the previous event.
+
+ // The use case for the bug MATH-421 is fairly general, so events occurring
+ // less than epsilon after the solver start in the first step should be ignored,
+ // where epsilon is the convergence threshold of the event. The sign of the g
+ // function should be evaluated after this initial ignore zone, not exactly at
+ // beginning (if there are no event at the very beginning g(t0) and g(t0+epsilon)
+ // have the same sign, so this does not hurt ; if there is an event at the very
+ // beginning, g(t0) and g(t0+epsilon) have opposite signs and we want to start
+ // with the second one. Of course, the sign of epsilon depend on the integration
+ // direction (forward or backward). This explains what is done below.
+
+ final double ignoreZone = interpolator.isForward() ? getConvergence() : -getConvergence();
+ t0 = interpolator.getPreviousTime() + ignoreZone;
+ interpolator.setInterpolatedTime(t0);
+ g0 = handler.g(t0, interpolator.getInterpolatedState());
+ if (g0 == 0) {
+ // extremely rare case: there is a zero EXACTLY at end of ignore zone
+ // we will use the opposite of sign at step beginning to force ignoring this zero
+ final double tStart = interpolator.getPreviousTime();
+ interpolator.setInterpolatedTime(tStart);
+ g0Positive = handler.g(tStart, interpolator.getInterpolatedState()) <= 0;
+ } else {
+ g0Positive = g0 >= 0;
+ }
+
+ } catch (DerivativeException de) {
+ throw new EventException(de);
+ }
}
/** Evaluate the impact of the proposed step on the event handler.
diff --git a/src/site/xdoc/changes.xml b/src/site/xdoc/changes.xml
index b7a560856..909b85091 100644
--- a/src/site/xdoc/changes.xml
+++ b/src/site/xdoc/changes.xml
@@ -52,6 +52,9 @@ The
When an event is triggered, the event time, current state and an indicator diff --git a/src/test/java/org/apache/commons/math/ode/events/EventStateTest.java b/src/test/java/org/apache/commons/math/ode/events/EventStateTest.java index e4baf3bd8..dff6d6b5b 100644 --- a/src/test/java/org/apache/commons/math/ode/events/EventStateTest.java +++ b/src/test/java/org/apache/commons/math/ode/events/EventStateTest.java @@ -49,11 +49,16 @@ public class EventStateTest { final double tolerance = 0.1; EventState es = new EventState(closeEventsGenerator, 1.5 * gap, tolerance, 10); - double t0 = r1 - 0.5 * gap; - es.reinitializeBegin(t0, new double[0]); AbstractStepInterpolator interpolator = new DummyStepInterpolator(new double[0], new double[0], true); - interpolator.storeTime(t0); + interpolator.storeTime(r1 - 2.5 * gap); + interpolator.shift(); + interpolator.storeTime(r1 - 1.5 * gap); + es.reinitializeBegin(interpolator); + + interpolator.shift(); + interpolator.storeTime(r1 - 0.5 * gap); + Assert.assertFalse(es.evaluateStep(interpolator)); interpolator.shift(); interpolator.storeTime(0.5 * (r1 + r2));