[SciPy-user] Optimization & Parallelization of, integrate.odeint

[Anne Archibald]

On 08/05/07, Lars Friedrich <lfriedri at imtek.de> wrote:

> Anne Archibald wrote:
> > If you write it efficiently, with minimal loops, you will not lose a
> > great deal of performance by using python.
> Ok. I have the problem, that my ODEs are usually very simple like
>
> x'(t) = input(t) - k(x) * x(t)
>
> where k(x) is a nonlinear function of x and x is a scalar(!). So I don't
> have the chance to avoid loops, because there are none.

Well, the answer depends on exactly why the simple approach is too
slow: are you solving the equation many times? are you following the
equation at many times? does the solution have many sharp features the
solver is working hard to track? is the problem stiff (forcing the
solver to work as if there were sharp features even if there aren't -
scipy's solvers should work fine here but check their full output to
see if it's happening)? is the nonlinear function k(x) expensive to
evaluate?

There are ways to speed it up in most of the above cases, but they are
very different.

> I understand. This means, that the solver will to the unavoidable
> "looping" over the system function for me. But if the system function is
> a basic python function, this will cause some overhead, that makes the
> simulation slow, or did I get something wrong here?

Yes. It's not looping per se that's the problem, it's that execution
of python code is fairly slow. And note, "fairly slow" just means
"much slower than C code accomplishing the same task". It may be fast
enough for your application.

> Ok. What do you recommend? (Fortran is no option for me ;-) Numexpr is
> in the sandbox, so I started to look at weave. But it seems to me that I
> need some work to start using weave. (setting up the compiler?,
> understanding how weave works...) This is no problem, if I know, that I
> am on the right way. So can you recommend the combination weave/odeint
> or is there a better way? Do you think this will speed up simulations
> with *scalar* or *low-dimensional* states?

Again, it depends what the problem is. Using odeint will certainly
incur *some* overhead for every call of the RHS function, simply
because the call has to go through python. But weave will allow you to
speed up the RHS itself. There's also a pure python optimizer (pyrex?
pyro? not sure what it's called) that's supposed to work fairly well.

In short: the answer depends totally on why the obvious approach is too slow.

Anne