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

Ryan Krauss ryanlists at gmail.com
Fri Aug 3 15:41:00 EDT 2012


Fortran, so fast, yet so painful.  Once I got it working, it was 94
times faster than my pure Python version.

Thanks to Jim and Pauli for helping me find my error.  Ironically, I
was thinking like a C programmer.  Just because a Fortran function
declares its return value data type doesn't mean all calling functions
or subroutines will know the data of the function when they call it.

I am still open to Cython suggestions.  I don't want to bring more F77
code into the world.....

On Fri, Aug 3, 2012 at 2:16 PM, Ryan Krauss <ryanlists at gmail.com> wrote:
> Thanks for the suggestions.
>
>> - Do not use small numpy arrays inside loops.
>>   Use C constructs instead.
>
> This is where I ran into trouble with my knowledge of C.  I have
> several 3x1 arrays that I need to pass into the dxdt function,
> multiply by scalars, and add together.  I don't know how to do that
> cleanly in C.  For example:
> x_out = x_prev + 1.0/6*(g1 + 2*g2 + 2*g3 + g4)
> where x_prev, g1, g2, g3, and g4 are all 3x1.
>
> A little googling lead me to valarray's, but I don't know if that is
> the best approach or how to use them within Cython.
>
> How would you do basic math on small arrays in pure C?
>
>
>
> On Fri, Aug 3, 2012 at 1:56 PM, Pauli Virtanen <pav at iki.fi> wrote:
>> 03.08.2012 19:02, Ryan Krauss kirjoitti:
>> [clip]
>>> Can anyone either make my Cython or Fortran approaches work or point
>>> me in a different direction?
>>
>> Regarding Cython: run
>>
>>         cython -a runge_kuta.pyx
>>
>> and check the created HTML file. Slow points are highlighted with yellow.
>>
>> Regarding this case:
>>
>> - `cdef`, not `def` for the dxdt_* function
>>
>> - from libc.math import exp
>>
>> - Do not use small numpy arrays inside loops.
>>   Use C constructs instead.
>>
>> - Use @cython.cdivision(True), @cython.boundscheck(False)
>>
>>
>>
>> PS. Runge-Kutta
>>
>> --
>> Pauli Virtanen
>>
>> _______________________________________________
>> SciPy-User mailing list
>> SciPy-User at scipy.org
>> http://mail.scipy.org/mailman/listinfo/scipy-user



More information about the SciPy-User mailing list