3

Using Cython, I try to do this:

cpdef myFun(double[:] array):
    cdef int[:] sortIndices = np.argsort(array, kind='mergesort')
    array = array[sortIndices]

The compiler complains:

Invalid index for memoryview specified, type int[:]

How do I index this memoryview using an integer array of some sort? Are only slices allowed? I could easily use 'array'-based indexing with the old NumPy array buffer support. (I just adapted my code to use memoryviews to see if it would improve performance, but it actually breaks...)

PDiracDelta
  • 2,348
  • 5
  • 21
  • 43

2 Answers2

4

@ead's suggestion of unrolling a loop yourself is a good one, however I'd be tempted to do this type of indexing on the underlying Numpy arrays, which you can access using the base attribute of a memoryview:

array = array.base[sortIndices]

or alternatively

array = np.asarray(array)[sortIndices]

This has the advantage of being quick to code and only needing minimal modifications from your working ndarray code. It has a few small disadvantages:

  • No Cython speedup since it's basically a Python object call - I hope this doesn't matter since Numpy indexing is generally pretty fast and the assumption would be that sortIndices is long enough to negate to Python object call.

  • The first version breaks if the underlying object isn't actually a Numpy array (thus the function is a little more constrained in what types it can take than the memoryview interface originally appears. You could work round that by using the second version, which should create a Numpy array wrapped round the memory of the memoryview.

DavidW
  • 29,336
  • 6
  • 55
  • 86
1

I'm afraid it is impossible to use a typed memory view as indices the way numpy can use int-arrays.

Cython's documentations states that

Memory views use Python slicing syntax in a similar way as NumPy.

"Similar" means it is the same for integers, slices (:), ellipsis (...) and None ( corresponds to numpy.newaxis('None')) but not for integer-arrays or boolean-arrays, which can be used as indices in numpy.

The Cython-code responsible for generating the access to typed memory views is generate_buffer_slice_code, and all you have to know is in the doc-string:

 """ 
 Slice a memoryviewslice. 
 indices     - list of index nodes. 
 If not a SliceNode, or NoneNode,  then it must be coercible to Py_ssize_t 
 ....
 """

So array is neither a SliceNode (i.e. :, ... or e.g. 0:33), nor None nor can it be coerced to Py_ssize_t and thus cannot be handled by Cython as it is.

My first thought was, that it will not be too hard to add the functionality to Cython's typed memory views. But this is probably not that easy to do in a consistent way: the result of the operation must be a new typed memory view (as we cannot change the array at hand in-place - the resulting array could have completely different dimensions) - but which type should be the underlying buffer? It is easier in the numpy-world, where everything is a numpy array but the underlying buffer of typed memory view could be numpy-array, an array.array, a c-array on the stack and so on.

Your best option right now is probably to roll out this reordering per hand (so you can explicitly choose the type of the underlying buffer) and replace the broken code through it or to fall back to numpy functionality for this operation as shown in @DavidWs answer.

ead
  • 32,758
  • 6
  • 90
  • 153