9

I have some C code that has the following declaration:

int myfunc(int m, int n, const double **a, double **b, double *c);

So a is a constant 2D array, b is a 2D array, and c is a 1D array, all dynamically allocated. b and c do not need to be anything specifically before they are passed to myfunc, and should be understood as output information. For the purposes of this question, I'm not allowed to change the declaration of myfunc.

Question 1: How do I convert a given numpy array a_np into an array a with the format required by this C function, so that I can call this C function in Cython with a?

Question 2: Are the declarations for b and c below correct, or do they need to be in some other format for the C function to understand them as a 2D and 1D array (respectively)?

My attempt:

myfile.pxd

cdef extern from "myfile.h":
    int myfunc(int p, int q, const double **a, double **b, double *c)

mytest.pyx

cimport cython
cimport myfile
import numpy as np
cimport numpy as np

p = 3
q = 4
cdef:
    double** a = np.random.random([p,q])
    double** b
    double* c

myfile.myfunc(p, q, a, b, c)

Then in iPython I run

import pyximport; pyximport.install()
import mytest

The line with the definition of a gives me the error message Cannot convert Python object to 'double **'. I don't get any error messages regarding b or c, but since I'm unable to run the C function at this time, I'm not sure the declarations of b and c are written correctly (that is, in a way that will enable the C function to output a 2D and a 1D array, respectively).

Other attempts: I've also tried following the solution here, but this doesn't work with the double-asterisk type of arrays I have in the myfunc declaration. The solution here does not apply to my task because I can't change the declaration of myfunc.

Improve this question
3
  • When you say "dynamically allocated", you mean outside myfunc? And since you're trying pass numpy arrays into myfunc, that is irrelevant and you just need to convert those numpy arrays into the suitable argument format (double and single pointers to double), correct? Commented Nov 23, 2016 at 1:51
  • @Evert First of all, let me warn you that I'm not very knowledgable about C. I'm just trying to use myfunc to compute the arrays b and c and I don't need them to be dynamically allocated or anything special. I only called them "dynamically allocated" because that's the format I thought the double and single pointers to double necessitated. In short, yes, you are correct. Commented Nov 23, 2016 at 9:57
  • Using double** doesn't match well with numpy. See stackoverflow.com/questions/27681814/… for some discussion. Commented Nov 23, 2016 at 17:50

2 Answers 2

Reset to default
13

Create a helper array in cython

To get a double** from a numpy array, you can create a helper-array of pointers in your *.pyx file. Further more, you have to make sure that the numpy array has the correct memory layout. (It might involve creating a copy)

Fortran order

If your C-function expects fortran order (all x-coordinates in one list, all y coordinates in another list, all z-coordinates in a third list, if your array a corresponds to a list of points in 3D space)

N,M = a.shape
# Make sure the array a has the correct memory layout (here F-order)
cdef np.ndarray[double, ndim=2, mode="fortran"] a_cython =
                         np.asarray(a, dtype = float, order="F")
#Create our helper array
cdef double** point_to_a = <double **>malloc(M * sizeof(double*))
if not point_to_a: raise MemoryError
try:
    #Fillup the array with pointers
    for i in range(M): 
        point_to_a[i] = &a_cython[0, i]
    # Call the C function that expects a double**
    myfunc(... &point_to_a[0], ...)
finally:
    free(point_to_a)

C-order

If your C-function expects C-order ([x1,y1,z1] is the first list, [x2,y2,z2] the second list for a list of 3D points):

N,M = a.shape
# Make sure the array a has the correct memory layout (here C-order)
cdef np.ndarray[double, ndim=2, mode="c"] a_cython =
                         np.asarray(a, dtype = float, order="C")
#Create our helper array
cdef double** point_to_a = <double **>malloc(N * sizeof(double*))
if not point_to_a: raise MemoryError
try:
    for i in range(N): 
        point_to_a[i] = &a_cython[i, 0]
    # Call the C function that expects a double**
    myfunc(... &point_to_a[0], ...)
finally:
    free(point_to_a)
Improve this answer
Sign up to request clarification or add additional context in comments.

5 Comments

You have no idea how long I've been searching for this answer. Thank you for posting it. But where is a_cython defined?
cdef np.ndarray[double, ndim=2, mode="fortran"] a_cython = np.asarray(a, dtype = float, order="F")
Thanks I must have not scrolled over far enough on the SO app.
Very helpful answer, thanks +1! Note that the Cython docs recommends using PyMem_Malloc() and PyMem_Free() instead of malloc() and free().
Furthermore, using the more generic typed memoryviews (e.g. cdef [:,::1] a_view = np.ascontiguousarray(a) instead of cdef np.ndarray[.....] a_view = ...) comes with improved readability and other advantages. See also this tutorial or this (duplicate) post.
0

Reply 1: You can pass NumPy array via Cython to C using the location of the start of the array (see code below).

Reply 2: Your declarations seem correct but I don't use this approach of explicit memory management. You can use NumPy to declare cdef-ed arrays.

Use

cdef double[:,::1] a = np.random.random([p, q])
cdef double[:,::1] b = np.empty([p, q])
cdef double[::1] b = np.empty(q)

Then pass &a[0], the location of the start of the array, to your C function. The ::1 is to ensure contiguousness.

A good reference for this is Jake Vanderplas' blog: https://jakevdp.github.io/blog/2012/08/08/memoryview-benchmarks/

Finally, typically one creates functions in Cython and calls them in Python, so your Python code would be:

import pyximport; pyximport.install()
import mytest
mytest.mywrappedfunc()

where mywrappedfunc is a Python (def and not cdef) function defined in the module that can do the array declaration show above.

Improve this answer

3 Comments

Thanks for your reply, but this doesn't work... I get the following errors: Cannot take address of memoryview slice for the array a, Cannot assign type 'double[:, ::1]' to 'double **' for b, and Cannot assign type 'double[::1]' to 'double *' for c.
Yes, I replaced my cdef statements with the cdef statements you provided.
Hi, I understand my confusion. I use cythonized Fortran code, for which it works. You are using in C a 2D array that is a 1D array of pointers to 1D arrays, whereas in Fortran you have a direct pointer to the start of the array as the argument. Options: 1. Set at compile time all (except the first) dimensions of your C array. 2. Use a 1D array with manual "sub-indexing" of the data in C. 3. If you cannot change the C code, you must build the array as an "array of pointers to 1D arrays" at the Cython level and pass this instead. The syntax &c[0] should work, I just checked it.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.