How to process a 2d array efficiently

503 views
Skip to first unread message

Raphael C

unread,
Dec 23, 2016, 8:48:52 AM12/23/16
to cython-users
I have some code I am trying to make fast using cython. It takes in a 2d numpy array of floating point numbers in a matrix called M. When I run cython -a on it I see that one crucial line is still being interpreted in python. The offending line is

v[i] -= 2*d[j]*M[j][i]

in the code below.

from __future__ import division
import numpy as np
cimport numpy as np
DTYPE_int = np.int
ctypedef np.int_t DTYPE_int_t
DTYPE_float = np.float64
ctypedef np.float64_t DTYPE_float_t
def permfunc(np.ndarray [DTYPE_float_t, ndim =2] M):
    cdef int n = M.shape[0]
    cdef np.ndarray[DTYPE_float_t, ndim =1] d = np.ones(n, dtype=DTYPE_float)
    cdef int j =  0
    cdef int s = 1
    cdef np.ndarray [DTYPE_int_t, ndim =1] f = np.arange(n, dtype=DTYPE_int)
    cdef np.ndarray [DTYPE_float_t, ndim =1] v = M.sum(axis=0)
    cdef DTYPE_float_t p = 1
    cdef int i
    cdef double prod
    for i in range(n):
        p *= v[i]
    while (j < n-1):
        for i in range(n):
            v[i] -= 2*d[j]*M[j][i]
        d[j] = -d[j]
        s = -s
        prod = 1
        for i in range(n):
            prod *= v[i]
        p += s*prod
        f[0] = 0
        f[j] = f[j+1]
        f[j+1] = j+1
        j = f[0]
    return p/2**(n-1) 


What am I doing wrong?

Raphael

Gregor Thalhammer

unread,
Dec 23, 2016, 8:59:41 AM12/23/16
to cython...@googlegroups.com
Am 23.12.2016 um 11:51 schrieb Raphael C <drr...@gmail.com>:

I have some code I am trying to make fast using cython. It takes in a 2d numpy array of floating point numbers in a matrix called M. When I run cython -a on it I see that one crucial line is still being interpreted in python. The offending line is

v[i] -= 2*d[j]*M[j][i]

try

v[i] -= 2*d[j]*M[j,i]

Does this improve speed?
Gregor

--

---
You received this message because you are subscribed to the Google Groups "cython-users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to cython-users...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

Raphael C

unread,
Dec 23, 2016, 1:03:33 PM12/23/16
to cython...@googlegroups.com
On 23 December 2016 at 13:59, Gregor Thalhammer
<gregor.t...@gmail.com> wrote:
>
> Am 23.12.2016 um 11:51 schrieb Raphael C <drr...@gmail.com>:
>
> I have some code I am trying to make fast using cython. It takes in a 2d
> numpy array of floating point numbers in a matrix called M. When I run
> cython -a on it I see that one crucial line is still being interpreted in
> python. The offending line is
>
>> v[i] -= 2*d[j]*M[j][i]
>
>
> try
>
> v[i] -= 2*d[j]*M[j,i]
>
> Does this improve speed?


Amazingly that seems to be the answer!

Thank you

Raphael

Chris Barker

unread,
Dec 27, 2016, 3:05:30 PM12/27/16
to cython-users
On Fri, Dec 23, 2016 at 6:00 AM, Raphael C <drr...@gmail.com> wrote:
On 23 December 2016 at 13:59, Gregor Thalhammer
>> v[i] -= 2*d[j]*M[j][i]

vs:
 
> v[i] -= 2*d[j]*M[j,i]
>
> Does this improve speed?


Amazingly that seems to be the answer!

good catch Gregor!

but it's not amazing :-)

M[j][i] 

is the C way to index "multi-dimensional" arrays, but in P/Cython what that does is index M as a one-d array, returning a 1-d "row" as a ndarray, and then indexing that. crating that intermediate ndarray is going to slow things down.

M[i,j]

does the 2-d indexing directly -- Cython "Understands" n-d arrays so can do that.

YOU may be able to eak out a tiny bit more speed by:

@cython.boundscheck(False)

using an unsigned type for your indexes.

maybe flattening the 2-D array and doing a single loop through it -- faster that pointer arithmetic for each access -- though I haven't looked to see if that's at all easy to do with your algorithm.

-CHB



--

Christopher Barker, Ph.D.
Oceanographer

Emergency Response Division
NOAA/NOS/OR&R            (206) 526-6959   voice
7600 Sand Point Way NE   (206) 526-6329   fax
Seattle, WA  98115       (206) 526-6317   main reception

Chris....@noaa.gov
Reply all
Reply to author
Forward
0 new messages