将三维numpy数组传递给C

fjnneemd 于 2023-06-28 发布在其他

关注(0)|答案(4)|浏览(160)

为了提高速度，我正在为我的Python程序编写一个C扩展，并且在尝试传入一个三维numpy数组时遇到了一些非常奇怪的行为。它适用于2维数组，但我确信我在试图让它适用于3维的指针上搞砸了一些事情。但奇怪的是。如果我只是传入一个3-D数组，它会崩溃并显示一个总线错误。如果（在Python中）我首先将变量创建为2D数组，然后用3D数组覆盖它，它可以完美地工作。如果变量首先是空数组，然后是3D数组，则会崩溃并出现Seg Fault。怎么可能发生这种事？
还有，谁能帮我弄一个3D数组？或者我应该给予并传入一个2D数组并自己重塑它？
下面是我的C代码：

static PyObject* func(PyObject* self, PyObject* args) {
  PyObject *list2_obj;
  PyObject *list3_obj;
  if (!PyArg_ParseTuple(args, "OO", &list2_obj, &list3_obj))
    return NULL;

  double **list2;
  double ***list3;

  //Create C arrays from numpy objects:
  int typenum = NPY_DOUBLE;
  PyArray_Descr *descr;
  descr = PyArray_DescrFromType(typenum);
  npy_intp dims[3];
  if (PyArray_AsCArray(&list2_obj, (void **)&list2, dims, 2, descr) < 0 || PyArray_AsCArray(&list3_obj, (void ***)&list3, dims, 3, descr) < 0) {
    PyErr_SetString(PyExc_TypeError, "error converting to c array");
    return NULL;
  }
  printf("2D: %f, 3D: %f.\n", list2[3][1], list3[1][0][2]);
}

下面是调用上述函数的Python代码：

import cmod, numpy
l2 = numpy.array([[1.0,2.0,3.0], [4.0,5.0,6.0], [7.0,8.0,9.0], [3.0, 5.0, 0.0]])

l3 = numpy.array([[2,7, 1], [6, 3, 9], [1, 10, 13], [4, 2, 6]])  # Line A
l3 = numpy.array([])                                             # Line B

l3 = numpy.array([[[2,7, 1, 11], [6, 3, 9, 12]],
                 [[1, 10, 13, 15], [4, 2, 6, 2]]])

cmod.func(l2, l3)

因此，如果我同时注解掉A行和B行，它将崩溃，并出现Bus错误。如果行A在那里，但行B被注解掉，它可以正确运行，没有错误。如果行B在那里，但行A被注解掉，它会打印正确的数字，但随后出现Seg错误。最后，如果两条线都存在，它也会打印正确的数字，然后是Seg故障。这到底是怎么回事

**EDIT：**Ok.哇。所以我在Python中使用了int，但在C中称之为double。这在1D和2D阵列上工作得很好。但不是3D。所以我把Python中l3的定义改成了floats，现在一切都很好用了（Thank you very much Bi Rico）。

但现在，更奇怪的行为与线A和B！现在，如果两行都被注解掉，程序就可以工作了。如果B行存在，但A被注解掉，则它可以工作，如果两者都未注解，则同上。但是如果A行存在，而B被注解掉，我又会得到那个奇怪的总线错误。我真的很想在将来避免这些，所以有人知道为什么Python变量的声明会有这种影响吗？

**编辑2：**好吧，尽管这些错误很疯狂，但它们都是由于我传入的三维numpy数组。如果我只传入1-D或2-D数组，它的行为与预期的一样，并且对其他Python变量的操作没有任何作用。这让我相信问题出在Python的引用计数上。在C代码中，引用计数的减少比3-D数组的减少要多，当该函数返回时，Python试图清理对象，并试图删除一个NULL指针。这只是我的猜测，我已经尝试了Py_INCREF();一切我能想到的都无济于事。我想我将使用一个2D数组并在C中重塑它。

来源：https://stackoverflow.com/questions/15576136/passing-3-dimensional-numpy-array-to-c

4条答案

按热度按时间

aiqt4smr1#

我通常使用PyArray_GETPTR直接访问numpy数组元素，而不是转换为c样式的数组（请参阅https：//numpy.org/doc/stable/reference/c-api/array.html#data-access）。
例如，要访问double类型的三维numpy数组的元素，请使用double elem=*((double *)PyArray_GETPTR3(list3_obj,i,j,k))。
对于您的应用程序，您可以使用PyArray_NDIM检测每个数组的正确维数，然后使用适当版本的PyArray_GETPTR访问元素。

赞(0）回复(0）举报 2023-06-28

rlcwz9us2#

我已经在评论中提到了这一点，但我希望稍微冲洗一下有助于使它更清楚。
当你在C中使用numpy数组时，最好明确数组的类型。具体地说，看起来你把指针声明为double ***list3，但是在你的python代码中创建l3的方式，你会得到一个dtype为npy_intp的数组（我认为）。你可以通过在创建数组时显式地使用dtype来解决这个问题。

import cmod, numpy
l2 = numpy.array([[1.0,2.0,3.0],
                  [4.0,5.0,6.0],
                  [7.0,8.0,9.0],
                  [3.0, 5.0, 0.0]], dtype="double")

l3 = numpy.array([[[2,7, 1, 11], [6, 3, 9, 12]],
                  [[1, 10, 13, 15], [4, 2, 6, 2]]], dtype="double")

cmod.func(l2, l3)

另外要注意的是，由于python的工作方式，“行A”和“行B”几乎不可能对C代码产生任何影响。我知道这似乎与你的经验相冲突，但我很确定这一点。
我对此不太确定，但根据我使用C语言的经验，总线错误和段错误不是确定性的。它们取决于内存分配、对齐和地址。在某些情况下，代码似乎运行了10次，第11次运行时失败，即使没有任何变化。
您是否考虑过使用cython？我知道这不是每个人的选择，但如果它是一个选项，你可以得到几乎C级的加速使用typed memoryviews。

赞(0）回复(0）举报 2023-06-28

zpjtge223#

根据http://docs.scipy.org/doc/numpy/reference/c-api.array.html?highlight=pyarray_ascarray#PyArray_AsCArray:
注意对于二维和三维数组，C样式数组的模拟并不完整。例如，模拟的指针数组不能传递给需要特定的、静态定义的2-d和3-d数组的子程序。要传递给需要这类输入的函数，必须静态定义所需的数组并复制数据。
我认为这意味着PyArray_AsCArray返回一个内存块，其中的数据按C顺序排列。但是，要访问该数据，需要更多信息（请参见http://www.phy225.dept.shef.ac.uk/mediawiki/index.php/Arrays,_dynamic_array_allocation）。这可以通过提前知道维度，声明数组，然后以正确的顺序复制数据来实现。然而，我怀疑更一般的情况更有用：你不知道尺寸，直到他们返回。我认为下面的代码将创建必要的C指针框架，以允许数据被寻址。

static PyObject* func(PyObject* self, PyObject* args) {
    PyObject *list2_obj;
    PyObject *list3_obj;
    if (!PyArg_ParseTuple(args, "OO", &list2_obj, &list3_obj)) return NULL;

    double **list2;
    double ***list3;

    // For the final version
    double **final_array2;
    double **final_array2;

    // For loops
    int i,j;

    //Create C arrays from numpy objects:
    int typenum = NPY_DOUBLE;
    PyArray_Descr *descr;
    descr = PyArray_DescrFromType(typenum);

    // One per array coming back ...
    npy_intp dims2[2];
    npy_intp dims3[3];

    if (PyArray_AsCArray(&list2_obj, (void **)&list2, dims2, 2, descr) < 0 || PyArray_AsCArray(&list3_obj, (void ***)&list3, dims3, 3, descr) < 0) {
        PyErr_SetString(PyExc_TypeError, "error converting to c array");
        return NULL;
    }

    // Create the pointer arrays needed to access the data

    // 2D array
    final_array2 = calloc(dim2[0], sizeof(double *));
    for (i=0; i<dim[0]; i++) final_array2[i] = list2 + dim2[1]*sizeof(double);

    // 2D array
    final_array3    = calloc(dim3[0], sizeof(double **));
    final_array3[0] = calloc(dim3[0]*dim3[1], sizeof(double *));
    for (i=0; i<dim[0]; i++) {
         final_array3[i] = list2 + dim3[1]*sizeof(double *);
         for (j=0; j<dim[1]; j++) {
             final_array[i][j] = final_array[i] + dim3[2]*sizeof(double);
         }
    }

    printf("2D: %f, 3D: %f.\n", final_array2[3][1], final_array3[1][0][2]);
    // Do stuff with the arrays

    // When ready to complete, free the array access stuff
    free(final_array2);

    free(final_array3[0]);
    free(final_array3);

    // I would guess you also need to free the stuff allocated by PyArray_AsCArray, if so:
    free(list2);
    free(list3);
}

我找不到npy_intp的定义，上面假设它与int相同。如果不是，则需要在执行代码之前将dim2和dim3转换为int数组。

赞(0）回复(0）举报 2023-06-28

jvlzgdj94#

在numpy C-API中有一个bug，现在应该修复：
https://github.com/numpy/numpy/pull/5314

赞(0）回复(0）举报 2023-06-28