diff --git a/src/main/java/org/apache/commons/math/ode/AbstractIntegrator.java b/src/main/java/org/apache/commons/math/ode/AbstractIntegrator.java index 60b2214d0..538311190 100644 --- a/src/main/java/org/apache/commons/math/ode/AbstractIntegrator.java +++ b/src/main/java/org/apache/commons/math/ode/AbstractIntegrator.java @@ -28,6 +28,8 @@ import java.util.TreeSet; import org.apache.commons.math.ConvergenceException; import org.apache.commons.math.MaxEvaluationsExceededException; +import org.apache.commons.math.analysis.solvers.BrentSolver; +import org.apache.commons.math.analysis.solvers.UnivariateRealSolver; import org.apache.commons.math.exception.MathUserException; import org.apache.commons.math.exception.util.LocalizedFormats; import org.apache.commons.math.ode.events.EventException; @@ -123,7 +125,18 @@ public abstract class AbstractIntegrator implements FirstOrderIntegrator { final double maxCheckInterval, final double convergence, final int maxIterationCount) { - eventsStates.add(new EventState(handler, maxCheckInterval, convergence, maxIterationCount)); + addEventHandler(handler, maxCheckInterval, convergence, + maxIterationCount, new BrentSolver(convergence)); + } + + /** {@inheritDoc} */ + public void addEventHandler(final EventHandler handler, + final double maxCheckInterval, + final double convergence, + final int maxIterationCount, + final UnivariateRealSolver solver) { + eventsStates.add(new EventState(handler, maxCheckInterval, convergence, + maxIterationCount, solver)); } /** {@inheritDoc} */ diff --git a/src/main/java/org/apache/commons/math/ode/ODEIntegrator.java b/src/main/java/org/apache/commons/math/ode/ODEIntegrator.java index bf62fd93d..e7a99f9a5 100644 --- a/src/main/java/org/apache/commons/math/ode/ODEIntegrator.java +++ b/src/main/java/org/apache/commons/math/ode/ODEIntegrator.java @@ -19,6 +19,8 @@ package org.apache.commons.math.ode; import java.util.Collection; +import org.apache.commons.math.analysis.solvers.UnivariateRealSolver; +import org.apache.commons.math.analysis.solvers.BrentSolver; import org.apache.commons.math.ode.events.EventHandler; import org.apache.commons.math.ode.sampling.StepHandler; @@ -62,7 +64,9 @@ public interface ODEIntegrator { */ void clearStepHandlers(); - /** Add an event handler to the integrator. + /** Add an event handler to the integrator. Uses a {@link BrentSolver} + * with an absolute accuracy equal to the given convergence threshold, + * as root-finding algorithm to detect the state events. * @param handler event handler * @param maxCheckInterval maximal time interval between switching * function checks (this interval prevents missing sign changes in @@ -76,6 +80,23 @@ public interface ODEIntegrator { void addEventHandler(EventHandler handler, double maxCheckInterval, double convergence, int maxIterationCount); + /** Add an event handler to the integrator. + * @param handler event handler + * @param maxCheckInterval maximal time interval between switching + * function checks (this interval prevents missing sign changes in + * case the integration steps becomes very large) + * @param convergence convergence threshold in the event time search + * @param maxIterationCount upper limit of the iteration count in + * the event time search + * @param solver The root-finding algorithm to use to detect the state + * events. + * @see #getEventHandlers() + * @see #clearEventHandlers() + */ + void addEventHandler(EventHandler handler, double maxCheckInterval, + double convergence, int maxIterationCount, + UnivariateRealSolver solver); + /** Get all the event handlers that have been added to the integrator. * @return an unmodifiable collection of the added events handlers * @see #addEventHandler(EventHandler, double, double, int) 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 2c97ba876..4cc74d789 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 @@ -19,7 +19,7 @@ package org.apache.commons.math.ode.events; import org.apache.commons.math.ConvergenceException; import org.apache.commons.math.analysis.UnivariateRealFunction; -import org.apache.commons.math.analysis.solvers.BrentSolver; +import org.apache.commons.math.analysis.solvers.UnivariateRealSolver; import org.apache.commons.math.exception.MathInternalError; import org.apache.commons.math.exception.MathUserException; import org.apache.commons.math.ode.sampling.StepInterpolator; @@ -81,6 +81,9 @@ public class EventState { /** Next action indicator. */ private int nextAction; + /** Root-finding algorithm to use to detect state events. */ + private final UnivariateRealSolver solver; + /** Simple constructor. * @param handler event handler * @param maxCheckInterval maximal time interval between switching @@ -89,13 +92,16 @@ public class EventState { * @param convergence convergence threshold in the event time search * @param maxIterationCount upper limit of the iteration count in * the event time search + * @param solver Root-finding algorithm to use to detect state events */ public EventState(final EventHandler handler, final double maxCheckInterval, - final double convergence, final int maxIterationCount) { + final double convergence, final int maxIterationCount, + final UnivariateRealSolver solver) { this.handler = handler; this.maxCheckInterval = maxCheckInterval; this.convergence = FastMath.abs(convergence); this.maxIterationCount = maxIterationCount; + this.solver = solver; // some dummy values ... t0 = Double.NaN; @@ -235,7 +241,6 @@ public class EventState { } } }; - final BrentSolver solver = new BrentSolver(convergence); if (ga * gb >= 0) { // this is a corner case: diff --git a/src/main/java/org/apache/commons/math/ode/nonstiff/GraggBulirschStoerIntegrator.java b/src/main/java/org/apache/commons/math/ode/nonstiff/GraggBulirschStoerIntegrator.java index 67dcfc03f..c56531684 100644 --- a/src/main/java/org/apache/commons/math/ode/nonstiff/GraggBulirschStoerIntegrator.java +++ b/src/main/java/org/apache/commons/math/ode/nonstiff/GraggBulirschStoerIntegrator.java @@ -17,6 +17,7 @@ package org.apache.commons.math.ode.nonstiff; +import org.apache.commons.math.analysis.solvers.UnivariateRealSolver; import org.apache.commons.math.exception.MathUserException; import org.apache.commons.math.ode.FirstOrderDifferentialEquations; import org.apache.commons.math.ode.IntegratorException; @@ -347,8 +348,10 @@ public class GraggBulirschStoerIntegrator extends AdaptiveStepsizeIntegrator { public void addEventHandler(final EventHandler function, final double maxCheckInterval, final double convergence, - final int maxIterationCount) { - super.addEventHandler(function, maxCheckInterval, convergence, maxIterationCount); + final int maxIterationCount, + final UnivariateRealSolver solver) { + super.addEventHandler(function, maxCheckInterval, convergence, + maxIterationCount, solver); // reinitialize the arrays initializeArrays(); diff --git a/src/site/xdoc/changes.xml b/src/site/xdoc/changes.xml index 19ecb450c..2a283b3a5 100644 --- a/src/site/xdoc/changes.xml +++ b/src/site/xdoc/changes.xml @@ -52,6 +52,9 @@ The type attribute can be add,update,fix,remove. If the output is not quite correct, check for invisible trailing spaces! --> + + Added a way to specify a custom root solver to find events in ODE integration. + Fixed error in javadoc describing Integer distribution inverse cumulative probability API. 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 57b0fa165..18f16cad3 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 @@ -19,6 +19,7 @@ package org.apache.commons.math.ode.events; import org.apache.commons.math.ConvergenceException; +import org.apache.commons.math.analysis.solvers.BrentSolver; import org.apache.commons.math.exception.MathUserException; import org.apache.commons.math.ode.sampling.AbstractStepInterpolator; import org.apache.commons.math.ode.sampling.DummyStepInterpolator; @@ -47,7 +48,9 @@ public class EventStateTest { }; final double tolerance = 0.1; - EventState es = new EventState(closeEventsGenerator, 1.5 * gap, tolerance, 10); + EventState es = new EventState(closeEventsGenerator, 1.5 * gap, + tolerance, 10, + new BrentSolver(tolerance)); AbstractStepInterpolator interpolator = new DummyStepInterpolator(new double[0], new double[0], true);