[SciPy-user] integrate ODE to steady-state?
Tue Jun 24 19:51:57 CDT 2008
On Tue, Jun 24, 2008 at 6:45 PM, Zachary Pincus <email@example.com> 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:
An interesting alternative is the CVT:
Nathan Bell firstname.lastname@example.org
More information about the SciPy-user