[SciPy-dev] Help with f2py

[Pearu Peterson]

On Wed, March 5, 2008 7:51 am, Russel Howe wrote:
> I am attempting to add support for lsodi from odepack to
> scipy.integrate.ode.  I have almost succeeded, but to be complete I
> would like to be able to recover values of ydot during integration.  The
> documentation for lsodi recommends a call to intdy to return these
> values.  How do I write a pyf file for this function?
>
>   From the lsodi documentation:
>
> c   call intdy (t, k, rwork(21), nyh, dky, iflag)
> c
> c the input parameters are..
> c
> c t         = value of independent variable where answers are desired
> c             (normally the same as the t last returned by lsodi).
> c             for valid results, t must lie between tcur - hu and tcur.
> c             (see optional outputs for tcur and hu.)
> c k         = integer order of the derivative desired.  k must satisfy
> c             0 .le. k .le. nqcur, where nqcur is the current order
> c             (see optional outputs).  the capability corresponding
> c             to k = 0, i.e. computing y(t), is already provided
> c             by lsodi directly.  since nqcur .ge. 1, the first
> c             derivative dy/dt is always available with intdy.
> c rwork(21) = the base address of the history array yh.
> c nyh       = column length of yh, equal to the initial value of neq.
>
>
> t, k, and nyh are simple but how do I pass the base address of an array?

For example, use rwork[21-1:] where rwork is the array passed on to lsodi
rwork argument. Make sure that no copies are made to rwork when calling
lsodi (use contiguous array with proper type).

>        Further digging shows the the yh array has dimensions (neq, nqcur)
> both of which are accessible at call time but unknown at compile time.

Then in the wrapper of intdy use `*` for yh dimensions.

> I use this information to construct an apropriate 2-d array to pass to
> intdy, but I don't see how to pass the full dimensions to f2py.

What do you mean?

HTH,
Pearu