[SciPy-user] stopping criterion for integrate.odeint?

Fernando Perez Fernando.Perez at colorado.edu
Wed Mar 23 18:34:53 CST 2005


Ryan Gutenkunst wrote:

> Consider the integration of dy/dt where y can be a vector. We want to 
> know exactly when the condition c(y) = 0 is satisfied. We integrate 
> dy/dt from 0 to some time T, requesting a reasonably fine grid of times. 
> Then we go back and calculate c(y, t) for each time point in our 
> trajectory. If two adjacent time points (t1 and t2) have different signs 
> for c(y, t), we know our condition was fulfilled sometime in between 
> those times.
> 
> Now we can integrate backward from t2 to find the exact time the 
> condition was fulfilled. To do so, we make a change of variables from t 
> to c. Namely, we integrate the equations defined by:
>    [dy/dc = dy/dt / dc/dt, dt/dc = 1 / dc/dt]
> from c = c(t2) to c = 0. (Note that this is one more equation that we 
> were integrating before.) The initial conditions are [y(t2), 0]. This 
> integration terminates with the values [y(tc), t2 - tc] where tc is the 
> time when c crossed zero. We can then insert these values into our 
> trajectory, and terminate the integration process if we want.

In case anyone wants a full reference:

  - M. Henon, Physica 5D, 412 (1982).

This is a standard trick in the construction of Poincare maps in dynamical 
systems.  The above paper by Henon is the place where I first read about it, 
though it is something simple enough that it has probably been 'reinvented' 
numerous times and in many different contexts.

Best,

f



More information about the SciPy-user mailing list