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

Zachary Pincus zachary.pincus@yale....
Tue Jun 24 21:50:20 CDT 2008

Hi all,

Thanks for the help!

> 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.

I'll try both of these -- I hadn't thought to use the solver's limit  
on iterations, which is a good idea...

> 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

My problem is similar, actually, but in 2D. I have an outline of a  
bacteria as a 2D polygon and want to lay down a coordinate system on  
it, with a midline running from pole to pole ("x-axis"), and well- 
defined lines normal to the midline ("y-axis"). The classic solution  
would be to use the medial axis transform, which I had been (using  
Martin Held's excellent if quirky VRONI code). However, the medial  
axis is only a partial midline, it may branch, and there's no  
guarantee that it is possible to find reasonable straight-line normals  
to the axis (you get a clearance radius at each point instead). So I  
formulated the problem as trying to find pairs of positions along the  
outline of the shape where the distance between the pairs is  
minimized, but adjacent points don't "clump" too much, which naturally  
led to this sort of spring-system picture. (This is a very brief  
explanation -- sorry if unclear.) It doesn't seem particularly elegant  
to solve it out like this, but I wanted to see if it worked at all...

Also, thanks for the references -- interesting indeed.

Rob and Lou:
> [try solving for dy/dt=0]

Good idea too! I'll try that one as well, and see which seems to work  

Thanks again, everyone,


More information about the SciPy-user mailing list