[SciPy-user] unexpected behavior in odeint

[Travis Oliphant]

>
> I have a question about the intended behavior of the
> scipy.integrate.odeint function.
>
> The documented process for constructing an ODE solution using the
> odeint() function is to pass a sequence (NumPy array) of times t at
> which a solution y(t) is desired.  Looking into the underlying lsoda
> subroutine, and the __odepack.h wrapper, however, I see that this
> multiple time-step interface is provided in the Python wrapper; that
> is, lsoda would appear to be constructed to execute a single time step
> (over and over again), and the multiple time stepping is a convenience
> provided in the Python overlayer.
>

This is true....

> The problem I am having is that lsoda does not take a time step
> exactly as requested, but instead takes a step to ensure that the
> current time is at least as large as the requested time value.  As a
> result, in the process of feeding an array of requested time points,
> lsoda can "outrun" the driver and produce multiple steps with the same
> value of the time tcur.  For example, I may request a solution at times
> t=arange(0., 10., 0.1).  odeint guarantees that each step will be at
> least as large as the requested time, but in some cases it would
> appear that the current time actually already exceeds the requested
> time.  In a particular example of mine, I find:
>


> t = t_requested[65:68] = array([ 6.5,         6.6,         6.7])
> tcur = t_actual[65:68] = array([ 6.56972519,  6.70778359,  6.70778359])
>
> What we see is that on step 66, odeint has incremented itself to a
> time which exceeds the requested value for step 67, so that it does
> not step to a larger time during step 67.  odeint does, however,
> return a different function value at steps 66 and 67, so that the
> result function y(t) is multivalued for some t (which is not especially
> appealing to me).

The problem that is happening is that I did not copy the array t before
passing it to the _odepack.odeint C-interface.  This interface hands the
array t to the lsoda routine which causes your t array to be overwritten
with possible overshot values for t0.

My understanding, though is that lsoda uses interpolation to give the
output value at the requested time.

So, I think your y output vector corresponds to the correct values but
your t vector may have been overwritten.  What does t_original compared to
the ouput look like?

>
> Am I misusing this multiple time-step feature of odeint?  Is the
> behavior that I have describe intended?  Should I be exercising caution
> in the construction of time-sequences to pass to odeint?  An alternate
> strategy would be to take the current time returned by lsoda, and
> set the next requested time based on that (but in such a case, one
> might get fewer time slices than originally requested).

I had not intended for the t vector to be overwritten (i don't think it is
a problem all the time -- but you've uncovered a case).

Is it possible for you to send an example fprime function, intial
conditions and time-points so I could try to replicate this behavior.

Thanks,

-Travis