[SciPy-user] integrate ODE to steady-state?
Nathan Bell
wnbell at gmail.com
Tue Jun 24 20:51:57 EDT 2008
On Tue, Jun 24, 2008 at 6:45 PM, Zachary Pincus <zachary.pincus at yale.edu> wrote:
>
> Any thoughts on what the easiest way to do this with scipy.integrate
> is? Ideally, I'd just want the solver to take as large steps as
> possible until things converge, and so I don't really care about the
> "time" values. One option would be to use odeint and just tell it to
> integrate to a distant time-point when I'm sure things will be in
> equilibrium, but that seems dorky, wasteful, and potentially incorrect.
>
Regarding numerical methods, you should read about the "Forward Euler"
integrator, which is very simple and intuitive. For applications
where accuracy is unimportant and a reasonable time step is known a
priori, Forward Euler is sufficient.
Also, explicit ODE integrators (like Forward Euler and Runge-Kutta)
are generally cheap, so the dominant cost of the integration is
usually computing the time derivatives for each variable.
You are right that you'll want to take as-large-as-possible timesteps
when integrating through fictitious time. However, even though the
timescale is arbitrary, the allowable timestep is limited by both the
accuracy of the trajectory (probably not a concern in your case) and
the stability of the system. In short, there's no free lunch :)
I don't know what interface scipy's ode solvers support, but I would
suggest setting t_final to be a large value and instructing the solver
to use at most N iterations. I would use an "adaptive" method like
RK45 (or whatever scipy supports) which automatically finds a timestep
that achieves numerical stability and the desired accuracy.
Alternatively, you could use your second proposal, and terminate
relaxation when the functional is sufficiently small.
BTW, are you using a Lennard-Jones like potential to distribute points
? I did something similar to distribute points evenly over a mesh
surface using a simple integrator:
http://graphics.cs.uiuc.edu/~wnbell/publications/2005-05-SCA-Granular/BeYiMu2005.pdf
An interesting alternative is the CVT:
http://www.math.psu.edu/qdu/Res/Pic/gallery3.html
--
Nathan Bell wnbell at gmail.com
http://graphics.cs.uiuc.edu/~wnbell/
More information about the SciPy-User
mailing list