[SciPy-User] help speeding up a Runge-Kuta algorithm (cython, f2py, ...)

[Ryan Krauss]

I upgraded to cython 0.16 and made a bit more progress.

I don't know if this is headed in the right direction or not, but
based on the memview.pyx example I changed
double _x[::1]
to
np.float64_t[::1] _x
and did the same thing with
cdef double out[::1] = np.zeros(3)

I seem to be closer to compiling successfully, but now have this error:
Error compiling Cython file:
------------------------------------------------------------
...

import numpy as np
cimport numpy as np
from libc.math cimport exp, fabs

cdef inline void dxdt_runge_kuta(double *x "const double *",
                               ^
------------------------------------------------------------

runge_kuta_v2.pyx:8:32: Function argument cannot have C name specification
(carrot points to the last a in runge_kuta

Thanks again,

Ryan


On Mon, Aug 6, 2012 at 9:02 AM, Ryan Krauss <ryanlists at gmail.com> wrote:
> So, I get the same error when I try to compile Stula's memview.pyx
> example.  I think I have too old of a version of cython:
>
> Cython version 0.15.1
>
> Let me look into that...
>
> On Mon, Aug 6, 2012 at 8:51 AM, Ryan Krauss <ryanlists at gmail.com> wrote:
>> Thanks to Sturla for helping me get this working in Cython.
>>
>> I am trying to compile the code to compare it against fortran for
>> speed.  I have run into two bugs so far (I mentioned that my C skills
>> are weak).
>>
>> The first has to do with the "const trick":
>> Error compiling Cython file:
>> ------------------------------------------------------------
>> ...
>> cdef inline void dxdt_runge_kuta(double *x "const double *",
>>                                  double voltage "const double",
>>                                  double *dxdt):
>>     #cdef double J = 0.0011767297528720126 "const double"
>>     cdef double J = 0.0011767297528720126
>>     cdef double alpha0 = 4.1396263800000002 "const double"
>>                                                                  ^
>> ------------------------------------------------------------
>>
>> runge_kuta_v2.pyx:12:44: Syntax error in C variable declaration
>>
>> I don't know what the problem is here, so for now I just got rid of
>> all the "const double" statements. (In case the formatting doesn't
>> come through, the little error carrot ^ points to the space between
>> the last number and the quote.
>>
>> After getting rid of all the "const double" expressions (just to see
>> if everything else would compile), I got this:
>> Error compiling Cython file:
>> ------------------------------------------------------------
>> ...
>>     dxdt[0] = vel
>>     dxdt[1] = accel
>>     dxdt[2] = dzdt
>>
>>
>> def runge_kuta_one_step(double _x[::1], Py_ssize_t factor, double volts,
>>                                                   ^
>> ------------------------------------------------------------
>>
>> runge_kuta_v2.pyx:31:34: Expected an identifier or literal
>>
>> The carrot points to the first square bracket.
>>
>> Thanks,
>>
>> Ryan
>>
>>
>> On Sat, Aug 4, 2012 at 6:28 PM, Sturla Molden <sturla at molden.no> wrote:
>>> Not tested and debugged, but to me it looks like something like this might
>>> be what you want.
>>>
>>> Sturla
>>>
>>>
>>> Den 03.08.2012 19:02, skrev Ryan Krauss:
>>>
>>> I need help speeding up some code I wrote to perform a Runge-Kuta
>>> integration.  I need to do the integration as part of a real-time
>>> control algorithm, so it needs to be fairly fast.
>>> scipy.integrate.odeint does too much error checking to be fast enough.
>>>  My pure Python version was just a little too slow, so I tried coding
>>> it up in Cython.  I have only used Cython once before, so I don't know
>>> if I did it correctly (the .pyx file is attached).
>>>
>>> The code runs just fine, but there is almost no speed up.  I think the
>>> core issue is that my dxdt_runge_kuta function gets called about 4000
>>> times per second, so most of my overhead is in the function calls (I
>>> think).  I am running my real-time control algorithm at 500 Hz and I
>>> need at least 2 Runge-Kuta integration steps per real-time steps for
>>> numeric stability.  And the Runge-Kuta algorithm needs to evaluate the
>>> derivative 4 times per times step.  So, 500 Hz * 2 * 4 = 4000 calls
>>> per second.
>>>
>>> I also tried coding this up in fortran and using f2py, but I am
>>> getting a type mismatch error I don't understand.  I have a function
>>> that declares its return values as double precision:
>>>
>>> double precision function dzdt(x,voltage)
>>>
>>> and I declare the variable I want to store the returned value in to
>>> also be double precision:
>>>
>>> double precision F,z,vel,accel,zdot1,zdot2,zdot3,zdot4
>>>
>>> zdot1 = dzdt(x_prev,volts)
>>>
>>> but some how it is not happy.
>>>
>>>
>>> My C skills are pretty weak (the longer I use Python, the more C I
>>> forget, and I didn't know that much to start with).  I started looking
>>> into Boost as well as using f2py on C code, but I got stuck.
>>>
>>>
>>> Can anyone either make my Cython or Fortran approaches work or point
>>> me in a different direction?
>>>
>>> Thanks,
>>>
>>> Ryan
>>>
>>>
>>>
>>> _______________________________________________
>>> SciPy-User mailing list
>>> SciPy-User at scipy.org
>>> http://mail.scipy.org/mailman/listinfo/scipy-user
>>>
>>>
>>>
>>> _______________________________________________
>>> SciPy-User mailing list
>>> SciPy-User at scipy.org
>>> http://mail.scipy.org/mailman/listinfo/scipy-user
>>>