我有一些生成和操作矩陣數(shù)組的 C++ 代碼Eigen。最后我想在 python 中使用這些矩陣，并認(rèn)為這可能是pybind11.基本上我想要在 python 中返回的是兩個(gè)嵌套列表/numpy 數(shù)組 mat_a(I, 4, 4)和mat_b(J, K, 4, 4). 因?yàn)槲冶仨氃?C++ 中做很多線性代數(shù)的東西，所以我想使用 Eigen，我使用的數(shù)據(jù)結(jié)構(gòu)是 std::array<std::array<Eigen::Matrix4f, 2>, 3>>> mat_b // for J=3, K=2. 現(xiàn)在的問題是如何有效地將它傳遞給python？此外，我想對(duì)多個(gè)輸入x = [x_0, x_1, ..., x_N] 執(zhí)行這些計(jì)算，結(jié)果mat_a(N, I, 4, 4)超出mat_b(N, J, K, 4, 4)預(yù)期。每個(gè)計(jì)算都是獨(dú)立的，但我認(rèn)為用 C++x_i重寫這個(gè)循環(huán)可能會(huì)更快。x_i另一方面，如果我們?cè)?C++ 中只有固定大小的數(shù)組，任務(wù)會(huì)變得更容易，這個(gè)循環(huán)也可以轉(zhuǎn)移到 python。這是我的問題的一些虛擬代碼（I=5，J=3，K=2）：// example.cpp#include <pybind11/pybind11.h>#include <pybind11/eigen.h>#include <pybind11/stl.h>#include <pybind11/functional.h>#include <pybind11/stl_bind.h>#include <array>#include <vector>#include <Eigen/Dense>Eigen::Matrix4f get_dummy(){ Eigen::Matrix4f mat_a; mat_a << 1, 2, 3, 4, 5, 6, 7, 8, 9, 8, 7, 6, 5, 4, 3, 2; return mat_a;}std::pair< std::vector<std::array<Eigen::Matrix4f, 5> >, std::vector<std::array<std::array<Eigen::Matrix4f, 2>, 3> > > get_matrices(std::vector<float> & x){ std::vector<std::array<Eigen::Matrix4f, 5> > mat_a(x.size()); std::vector< std::array< std::array< Eigen::Matrix4f, 2>, 3> > mat_b(x.size()); // for (u_int i=0; i< x.size(); i++) // do_stuff(x[i], mat_a[i], mat_b[i]); mat_a[0][0] = get_dummy(); return std::make_pair(mat_a, mat_b); }PYBIND11_MODULE(example, m) { m.def("get_dummy", &get_dummy, pybind11::return_value_policy::reference_internal); m.def("get_matrices", &get_matrices, pybind11::return_value_policy::reference_internal);}我通過以下方式編譯代碼：c++ -O3 -Wall -shared -std=c++14 -fPIC `python3 -m pybind11 --includes` example.cpp -o example`python3-config --extension-suffix`
2 回答
寶慕林4294392
TA貢獻(xiàn)2021條經(jīng)驗(yàn) 獲得超8個(gè)贊
您最好的選擇可能是在 python 端創(chuàng)建數(shù)據(jù)，以便對(duì)其進(jìn)行重新計(jì)數(shù)和垃圾收集。
test.py
import example
import numpy as np
array = np.zeros((3, 2, 4, 4), 'f4')
example.do_math(array, 3, 2)
print(array[0, 0])
例子.cpp
#define PY_SSIZE_T_CLEAN
#include <Python.h>
#include <Eigen/Dense>
Eigen::Matrix4f get_dummy() {
Eigen::Matrix4f mat_a;
mat_a << 1, 2, 3, 4,
5, 6, 7, 8,
9, 8, 7, 6,
5, 4, 3, 2;
return mat_a;
}
PyObject * example_meth_do_math(PyObject * self, PyObject * args, PyObject * kwargs) {
static char * keywords[] = {"array", "rows", "cols", NULL};
PyObject * array;
int rows, cols;
if (!PyArg_ParseTupleAndKeywords(args, kwargs, "Oii", keywords, &array, &rows, &cols)) {
return NULL;
}
Py_buffer view = {};
if (PyObject_GetBuffer(array, &view, PyBUF_SIMPLE)) {
return NULL;
}
Eigen::Matrix4f * ptr = (Eigen::Matrix4f *)view.buf;
for (int i = 0; i < rows; ++i) {
for (int j = 0; j < cols; ++j) {
ptr[i * cols + j] = get_dummy();
}
}
PyBuffer_Release(&view);
Py_RETURN_NONE;
}
PyMethodDef module_methods[] = {
{"do_math", (PyCFunction)example_meth_do_math, METH_VARARGS | METH_KEYWORDS, NULL},
{},
};
PyModuleDef module_def = {PyModuleDef_HEAD_INIT, "example", NULL, -1, module_methods};
extern "C" PyObject * PyInit_example() {
PyObject * module = PyModule_Create(&module_def);
return module;
}
setup.py
from setuptools import Extension, setup
ext = Extension(
name='example',
sources=['./example.cpp'],
extra_compile_args=['-fpermissive'],
include_dirs=['.'], # add the path of Eigen
library_dirs=[],
libraries=[],
)
setup(
name='example',
version='0.1.0',
ext_modules=[ext],
)
從這里添加第二個(gè)參數(shù)并使用兩個(gè)數(shù)組進(jìn)行計(jì)算應(yīng)該是微不足道的。
您可以使用python setup.py develop.
如果你想分發(fā)它，你可以創(chuàng)建一個(gè) wheel 文件python setup.py bdist_wheel。
我曾經(jīng)numpy創(chuàng)建數(shù)據(jù)，這確保了數(shù)據(jù)的底層內(nèi)存是 C 連續(xù)的。
這個(gè)例子很簡單，它使用一個(gè) Matrix4f 指針來迭代一個(gè) 3x2 矩陣數(shù)組。隨意將 轉(zhuǎn)換ptr為Eigen::Array<Eigen::Matrix4f>, 3, 2>。您不能將其強(qiáng)制轉(zhuǎn)換為 an std::vector，因?yàn)?an 的內(nèi)部數(shù)據(jù)std::vector包含指針。
請(qǐng)注意，std::vector<std::array<...>>內(nèi)存中沒有單個(gè)連續(xù)數(shù)組。改用Eigen::Array。
編輯：
這是一個(gè)使用Eigen Array Map:
PyObject * example_meth_do_math(PyObject * self, PyObject * args, PyObject * kwargs) {
static char * keywords[] = {"array", NULL};
PyObject * array;
if (!PyArg_ParseTupleAndKeywords(args, kwargs, "O", keywords, &array)) {
return NULL;
}
Py_buffer view = {};
if (PyObject_GetBuffer(array, &view, PyBUF_SIMPLE)) {
return NULL;
}
Eigen::Map<Eigen::Array<Eigen::Matrix4f, 2, 3>> array_map((Eigen::Matrix4f *)view.buf, 2, 3);
for (int i = 0; i < 2; ++i) {
for (int j = 0; j < 3; ++j) {
array_map(i, j) = get_dummy();
}
}
PyBuffer_Release(&view);
Py_RETURN_NONE;
}
至尊寶的傳說
TA貢獻(xiàn)1789條經(jīng)驗(yàn) 獲得超10個(gè)贊
線性代數(shù)不會(huì)那么流暢（在那里很難擊敗 Eigen），但會(huì)類似于您在 numpy 中所做的（np.dot(A,B)例如。
如果您想堅(jiān)持使用 Eigen，請(qǐng)注意使用 STL 有一些技術(shù)細(xì)節(jié)。由于您std::array不再能夠包含固定數(shù)量的矩陣，因此當(dāng)您移動(dòng)到std::vector您會(huì)遇到對(duì)齊問題（誠然，我不完全理解）。很快就會(huì)為您提供 xtensor 的有效實(shí)現(xiàn)。
添加回答
舉報(bào)
0/150
提交
取消