Przechodząc C++ wektora NumPy przez Cython bez kopiowania i dbanie o zarządzaniu pamięcią automatycznie

Question

Radzenie sobie z przetwarzaniem dużych matryc (NxM z 1K < = N < = 20K & 10K < = M < = 200K), często muszą przechodzić Numpy macierze od C++ do Cython, aby wykonać zadanie i działa zgodnie z oczekiwaniami & bez kopiowania.

Jednakże, zdarza się, że potrzebne do zainicjowania przebiegu wyprzedzającego i matryca, C++ i przekazać je do NumPy (Python 3.6). Załóżmy, że macierze są linearyzowane (więc rozmiar jest N * M i jest to matryca 1D - col/row major nie ma tutaj znaczenia). Po informacji tutaj: exposing C-computed arrays in Python without data copies & modyfikując go pod kątem zgodności z C++, jestem w stanie przekazać tablicę C++.

Problem to, jeśli chcę użyć wektora std zamiast inicjować tablicę, otrzymam błąd segmentacji. Na przykład, biorąc pod uwagę następujące pliki:

fast.h

#include <iostream> 
#include <vector> 

using std::cout; using std::endl; using std::vector; 
int* doit(int length); 

fast.cpp

#include "fast.h" 
int* doit(int length) { 
    // Something really heavy 
    cout << "C++: doing it fast " << endl; 

    vector<int> WhyNot; 

    // Heavy stuff - like reading a big file and preprocessing it 
    for(int i=0; i<length; ++i) 
     WhyNot.push_back(i); // heavy stuff 

    cout << "C++: did it really fast" << endl; 
    return &WhyNot[0]; // or WhyNot.data() 
} 

faster.pyx

cimport numpy as np 
import numpy as np 
from libc.stdlib cimport free 
from cpython cimport PyObject, Py_INCREF 

np.import_array() 

cdef extern from "fast.h": 
    int* doit(int length) 

cdef class ArrayWrapper: 
    cdef void* data_ptr 
    cdef int size 

    cdef set_data(self, int size, void* data_ptr): 
     self.data_ptr = data_ptr 
     self.size = size 

    def __array__(self): 
     print ("Cython: __array__ called") 
     cdef np.npy_intp shape[1] 
     shape[0] = <np.npy_intp> self.size 
     ndarray = np.PyArray_SimpleNewFromData(1, shape, 
               np.NPY_INT, self.data_ptr) 
     print ("Cython: __array__ done") 
     return ndarray 

    def __dealloc__(self): 
     print("Cython: __dealloc__ called") 
     free(<void*>self.data_ptr) 
     print("Cython: __dealloc__ done") 


def faster(length): 
    print("Cython: calling C++ function to do it") 
    cdef int *array = doit(length) 
    print("Cython: back from C++") 
    cdef np.ndarray ndarray 
    array_wrapper = ArrayWrapper() 
    array_wrapper.set_data(length, <void*> array) 
    print("Ctyhon: array wrapper set") 
    ndarray = np.array(array_wrapper, copy=False) 
    ndarray.base = <PyObject*> array_wrapper 
    Py_INCREF(array_wrapper) 
    print("Cython: all done - returning") 
    return ndarray 

setup.py

from distutils.core import setup 
from distutils.extension import Extension 
from Cython.Distutils import build_ext 
import numpy 

ext_modules = [Extension(
    "faster", 
    ["faster.pyx", "fast.cpp"], 
    language='c++', 
    extra_compile_args=["-std=c++11"], 
    extra_link_args=["-std=c++11"] 
)] 

setup(
    cmdclass = {'build_ext': build_ext}, 
    ext_modules = ext_modules, 
    include_dirs=[numpy.get_include()] 
) 

Jeśli budować, to z

python setup.py build_ext --inplace 

i uruchomić interpreter Pythona 3.6, jeśli wprowadź następujące chcesz dostać seg winy po kilku próbach.

>>> from faster import faster 
>>> a = faster(1000000) 
Cython: calling C++ function to do it 
C++: doing it fast 
C++: did it really fast 
Cython: back from C++ 
Ctyhon: array wrapper set 
Cython: __array__ called 
Cython: __array__ done 
Cython: all done - returning 
>>> a = faster(1000000) 
Cython: calling C++ function to do it 
C++: doing it fast 
C++: did it really fast 
Cython: back from C++ 
Ctyhon: array wrapper set 
Cython: __array__ called 
Cython: __array__ done 
Cython: all done - returning 
Cython: __dealloc__ called 
Segmentation fault (core dumped) 

Kilka rzeczy do uwaga:

• Jeśli używasz tablicę zamiast wektora (w fast.cpp) to będzie działać jak czar!
• Jeśli zadzwonisz pod numer faster(1000000) i umieścisz wynik w innym niż variable a, zadziała to.

Jeśli wpiszesz mniejsza liczba jak faster(10) chcesz uzyskać więcej szczegółowych informacji, takich jak:

Cython: calling C++ function to do it 
C++: doing it fast 
C++: did it really fast 
Cython: back from C++ 
Ctyhon: array wrapper set 
Cython: __array__ called 
Cython: __array__ done 
Cython: all done - returning 
Cython: __dealloc__ called <--- Perhaps this happened too early or late? 
*** Error in 'python': double free or corruption (fasttop): 0x0000000001365570 *** 
======= Backtrace: ========= 
More info here .... 

To naprawdę zastanawiające, że dlaczego tak się nie dzieje z tablicami? Nieważne co!

Bardzo często używam wektorów i chciałbym móc z nich korzystać w tych scenariuszach.

Answer 1

myślę @ odpowiedź FlorianWeimer zapewnia przyzwoite rozwiązanie (przydzielić vector i przekazać, że w swojej funkcji C++), ale powinno być możliwe t o zwróć wektor z doit i unikaj kopii za pomocą konstruktora ruchu.

from libcpp.vector cimport vector 

cdef extern from "<utility>" namespace "std" nogil: 
    T move[T](T) # don't worry that this doesn't quite match the c++ signature 

cdef extern from "fast.h": 
    vector[int] doit(int length) 

# define ArrayWrapper as holding in a vector 
cdef class ArrayWrapper: 
    cdef vector[int] vec 
    cdef Py_ssize_t shape[1] 
    cdef Py_ssize_t strides[1] 

    # constructor and destructor are fairly unimportant now since 
    # vec will be destroyed automatically. 

    cdef set_data(self, vector[int]& data): 
     self.vec = move(data) 

    # now implement the buffer protocol for the class 
    # which makes it generally useful to anything that expects an array 
    def __getbuffer__(self, Py_buffer *buffer, int flags): 
     # relevant documentation http://cython.readthedocs.io/en/latest/src/userguide/buffer.html#a-matrix-class 
     cdef Py_ssize_t itemsize = sizeof(self.vec[0]) 

     self.shape[0] = self.vec.size() 
     self.strides[0] = sizeof(int) 
     buffer.buf = <char *>&(self.vec[0]) 
     buffer.format = 'i' 
     buffer.internal = NULL 
     buffer.itemsize = itemsize 
     buffer.len = self.v.size() * itemsize # product(shape) * itemsize 
     buffer.ndim = 1 
     buffer.obj = self 
     buffer.readonly = 0 
     buffer.shape = self.shape 
     buffer.strides = self.strides 
     buffer.suboffsets = NULL 

Powinieneś być w stanie użyć go jako:

cdef vector[int] array = doit(length) 
cdef ArrayWrapper w 
w.set_data(array) # "array" itself is invalid from here on 
numpy_array = np.asarray(w) 

Answer 2

Po powrocie z obiektu doit obiekt WhyNot wychodzi poza zakres, a elementy tablicy są zwolnione.Oznacza to, że &WhyNot[0] nie jest już poprawnym wskaźnikiem. Obiekt WhyNot należy przechowywać w innym miejscu, prawdopodobnie w miejscu udostępnionym przez wywołującego.

Jednym ze sposobów, aby to zrobić, to podzielić doit na trzy funkcje, które przydziela doit_allocate wektor i zwraca wskaźnik do niego, doit jak wcześniej (ale z argumentem, który otrzymuje wskaźnik do zdefiniowanej przez wektor , and doit_free` który dealokujący wektor

coś takiego:.

vector<int> * 
doit_allocate() 
{ 
    return new vector<int>; 
} 

int * 
doit(vector<int> *WhyNot, int length) 
{ 
    // Something really heavy 
    cout << "C++: doing it fast " << endl; 

    // Heavy stuff - like reading a big file and preprocessing it 
    for(int i=0; i<length; ++i) 
     WhyNot->push_back(i); // heavy stuff 

    cout << "C++: did it really fast" << endl; 
    return WhyNot->front(); 
} 

void 
doit_free(vector<int> *WhyNot) 
{ 
    delete WhyNot; 
}