Dear Cythoneers,

I am simulating ordinary differential equation models with Pysundials
[1], a wrapper around the Sundials [2] suite of ODE solvers (written
in C). The ODE is specified by a Python function like:

def odefun(t, y, ydot, user_data):
   ydot[0] = ...
   ydot[1] = ...
   ...
   return 0

corresponding to the Sundials convention [3]:

typedef int (*CVRhsFn)( realtype t, N_Vector y, N_Vector ydot, void *user_data);

After an inspiring workshop with Dag Sverre Seljebotn [4], I see that
Cython could routinely speed up my simulation work by a factor of
hundreds. My benchmark used Numpy arrays for y and ydot, and then
cythonizing the ODE function was straightforward. However, CVODE (the
Sundials solver that I use) requires that y and ydot be of type
N_Vector [5]. Pysundials wraps that type, but at the cost of
inefficient Python indexing.

Somewhere inside pysundials.cvode.NVector is an array of double
screaming to get out [6]. However, Dag Sverre told me that passing it
to Cython may require converting between ctypes and Cython pointers.
Moreover, the solver makes lots of calls to odefun(), and these would
still be Python-based. Therefore, it is probably best to bypass
Pysundials altogether and call CVODE directly.

Typical use of CVODE [7] involves manual memory allocation and lots of
routine steps. Pysundials faithfully replicates this somewhat
unpythonic interface, which could make it fairly straightforward to
convert my code (which includes some higher-level wrappers) to call
the CVODE C library directly.

I typically use ODE functions autogenerated from a model repository
[8]. Although so far I've been generating Python code [9], with Cython
I'd use C directly [10].

So, to summarize:
I would be very grateful for any help in getting started with calling
the Sundials libraries from Cython. I'm running Python 2.5.2 with
Numpy 1.1.0 under Linux:
-bash-3.2$ dmesg | head -1
Linux version 2.6.18-92.1.13.el5 ([email protected]) (gcc
version 4.1.2 20071124 (Red Hat 4.1.2-42)) #1 SMP Wed Sep 24 19:32:05
EDT 2008

Subtasks include:
* Installing Sundials from [2] (Done; test programs run successfully.
Also, Pysundials uses the same Sundials installation.)
* (Optional: Experienced Cythoneers might perhaps find it informative
to see how Pysundials wraps things; in that case, install from [1]
(and watch out for the CPATH gotcha mentioned on that page))
* Telling Cython where to find cvode.h, nvector.h, sundials.h (see [11])
* Importing anything from Sundials (see [12]). Currently "python
setup.py build_ext --inplace" complains about not finding nvector.h
[13]
* Allocating an NVector with Sundials (see [3, 7]).
* Handle the NVector from Python: exchange data with a Numpy array,
and pass the NVector to a C or Cython function.
* Compiling a C function (for the ODE right-hand side) to work with
CVODE. I would modify [10] to fit the calling convention [3].

Given this, I think the other steps [7] will follow fairly easily.
Despite the lengthiness of this post, I suspect that interfacing with
Sundials would be quite straightforward if I knew which buttons to
push. The benefits would carry over to all my work with cellular ODE
models, and would hopefully have some reuse value for others. Thanks
in advance for any pointers (pardon the pun).

Best regards,
Jon Olav Vik

--
Jon Olav Vik
[email protected]
http://www.cigene.no/people/57-members/178-vik-jon-olav.html


[1] http://pysundials.sourceforge.net/
[2] https://computation.llnl.gov/casc/sundials/
[3] 
https://computation.llnl.gov/casc/sundials/documentation/cv_guide/node5.html#SECTION00561000000000000000
[4] http://wiki.cython.org/talks/notur2009
[5] https://computation.llnl.gov/casc/sundials/documentation/cv_guide/node7.html
[6]
In [25]: from pysundials import cvode
In [26]: y = cvode.NVector([3.14]) # state vector
In [38]: y.data.contents.content.contents.data[0]
Out[38]: 3.1400000000000001
In [39]: y.data.contents.content.contents.data[1]
Out[39]: 1.1297001364744253e-312
[7] 
https://computation.llnl.gov/casc/sundials/documentation/cv_guide/node5.html#SECTION00540000000000000000
[8] http://www.cellml.org/models
[9] https://tracker.physiomeproject.org/show_bug.cgi?id=1514#c28
[10] 
http://www.cellml.org/models/fitzhugh_1961_version05/generateRawSource?lang=C
[11] First attempt at setup.py:
"""Run with e.g. python setup.py build_ext --inplace"""
from distutils.core import setup
from distutils.extension import Extension
from Cython.Distutils import build_ext
setup(
   cmdclass = {'build_ext': build_ext},
   ext_modules = [
       Extension("hello",
           ["hello.pyx"],
           include_dirs=['/site/VERSIONS/sundials-2.3.0/include/cvode',
               '/site/VERSIONS/sundials-2.3.0/include/nvector',
               '/site/VERSIONS/sundials-2.3.0/include/sundials'],
           library_dirs=['/site/VERSIONS/sundials-2.3.0/lib'],
           libraries=['cvode']),
   ]
)
[12] First attempt at hello.pyx:
cdef extern from "cvode.h":
   int CVodeSetFdata(void *cvode_mem, void *f_data)
[13] Failed attempt at building:
-bash-3.2$ python setup.py build_ext --inplace
running build_ext
building 'hello' extension
gcc -pthread -fno-strict-aliasing -DNDEBUG -g -fwrapv -O3 -Wall
-Wstrict-prototy
              pes -fPIC -I/site/VERSIONS/sundials-2.3.0/include/cvode
-I/site/VERSIONS/sundial
                      s-2.3.0/include/nvector
-I/site/VERSIONS/sundials-2.3.0/include/sundials -I/site

/VERSIONS/compython-2.5/Linux/include/python2.5 -c hello.c -o
build/temp.linux-x
                86_64-2.5/hello.o
In file included from hello.c:135:
/site/VERSIONS/sundials-2.3.0/include/cvode/cvode.h:37:39: error:
sundials/sundi
            als_nvector.h: No such file or directory
In file included from hello.c:135:
/site/VERSIONS/sundials-2.3.0/include/cvode/cvode.h:177: error:
expected ')' bef
              ore 't'
...



-- 
Jon Olav Vik
[email protected]
http://www.cigene.no/people/57-members/178-vik-jon-olav.html
_______________________________________________
Cython-dev mailing list
[email protected]
http://codespeak.net/mailman/listinfo/cython-dev

Reply via email to