[SciPy-user] integrate ODE to steady-state?

Rob Clewley rob.clewley@gmail....
Tue Jun 24 19:55:42 CDT 2008


Zach,

It somewhat depends on your equations, but in principle you don't need
to actually integrate the ODEs at all. An equilibrium corresponds to
no motion, i.e. when the right-hand sides of all the ODEs are
zero-valued. So for a system dx/dt = f(x), where x is a vector, you
just need to solve f(x) = 0. This may have multiple solutions, of
course, if f is nonlinear. There's no scipy code that will do all this
for you but fsolve will do in most cases. It's best to have a good
initial guess as a starting condition for the search. So, integration
can help you find that. If f is strongly nonlinear there could be many
equilibria and your task will be to identify which initial conditions
lead to which equilibria. This is not a trivial task - you may need an
exhaustive search of initial conditions (e.g. based on sampling your
phase space of x) to get started. I have some naive code that
pre-samples the phase space and then uses fsolve, and works OK on some
nonlinear problems. It's the find_fixedpoints function in PyDSTool,
but it's extremely easy to remove the PyDSTool dependence :)

If this approach turns out to be too numerically problematic, you'll
just have to go back to integrating for long times...

Rob

On Tue, Jun 24, 2008 at 7:45 PM, Zachary Pincus <zachary.pincus@yale.edu> wrote:
> Hi all,
>
> So, after a brief bout of the stupids (thanks Robert), I have
> formulated my optimization problem as a physical system governed by an
> ODE, and I wish to learn the equilibrium configuration of the system.
>
> 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.
>
> Alternately, I could use the ode class and keep asking it to integrate
> small time-steps until the RMS change drops below a threshold. There,
> still, I'd need to choose a reasonable time-step, and also the inner
> loop would be in python instead of fortran.
>
> Any recommendations? (Or a I again being daft? I never really took a
> class in numerical methods, so sorry for dim-bulb questions!)
>
> Zach
> _______________________________________________
> SciPy-user mailing list
> SciPy-user@scipy.org
> http://projects.scipy.org/mailman/listinfo/scipy-user
>


More information about the SciPy-user mailing list