如何扩展numpy

静态的和重复的是无聊的。充满活力的
随机性让人困惑。艺术介于两者之间。
--- 约翰·A·洛克
科学是一个微分方程。宗教是一个边界条件。
--- 艾伦图灵

编写扩展模块

虽然ndarray对象被设计为允许在python中快速计算，但它也被设计为通用的，并满足各种各样的计算需求。因此，如果绝对速度是必需的，那么就没有任何替代品可以替代一个精心设计的、专门针对应用程序和硬件的编译循环。这是NUMPY包含F2PY的原因之一，因此可以将简单的C/C++和（任意）FORTRAN代码直接链接到Python中的一种易于使用的机制。我们鼓励您使用和改进这种机制。本节的目的不是记录此工具，而是记录编写此工具所依赖的扩展模块的更基本步骤。

当扩展模块被写入、编译并安装到python路径（sys.path）的某个位置时，代码就可以像标准的python文件一样导入python。它将包含用C代码定义和编译的对象和方法。在python中实现这一点的基本步骤都有很好的文档记录，您可以在python本身的文档中找到更多信息，该文档可以在 www.python.org .

除了pythoncapi之外，NumPy还有一个完整而丰富的C-API，允许在C级进行复杂的操作。但是，对于大多数应用程序，通常只使用少数API调用。例如，如果您只需要提取一个指向内存的指针以及一些形状信息以传递给另一个计算例程，那么您将使用与尝试创建新的数组类型或为ndarray添加新的数据类型截然不同的调用。本章介绍最常用的API调用和宏。

所需子程序

只有一个函数必须在C代码中定义，以便Python将其用作扩展模块。必须调用函数init name其中name是python中模块的名称。必须声明此函数，以便它对例程之外的代码可见。除了添加所需的方法和常量之外，此子例程还必须包含如下调用 import_array() 和/或 import_ufunc() 取决于需要哪个C-API。一旦实际调用任何C-API子例程，忘记放置这些命令将显示为一个丑陋的分段错误（崩溃）。实际上，在一个文件中可以有多个in it name函数，在这种情况下，该文件将定义多个模块。但是，有一些技巧可以使其正确工作，这里没有介绍。

极小值 init{{name}} 方法如下：

PyMODINIT_FUNC
init{name}(void)
{
   (void)Py_InitModule({name}, mymethods);
   import_array();
}

mymethods必须是pymethoddef结构的数组（通常静态声明），其中包含方法名、实际的C函数、指示该方法是否使用关键字参数的变量以及docstring。这些将在下一节中解释。如果要向模块添加常量，则存储模块对象py-initmodule返回的值。向模块添加项的最一般方法是使用pymodule_getdict（module）获取模块字典。使用模块字典，您可以手动向模块添加您喜欢的任何内容。向模块添加对象的一个更简单的方法是使用另外三个不需要单独提取模块字典的python c-api调用中的一个。这些都记录在Python文档中，但为了方便起见，在这里重复：

int PyModule_AddObject(PyObject *module, char *name, PyObject *value)
int PyModule_AddIntConstant(PyObject *module, char *name, long value)
int PyModule_AddStringConstant(PyObject *module, char *name, char *value)

所有这三个功能都需要 模块 对象（py-initmodule的返回值）。这个 name 是标记模块中的值的字符串。根据调用的函数， 价值 参数不是常规对象 (PyModule_AddObject 窃取对它的引用）、整数常量或字符串常量。

定义函数

传递给py in it module函数的第二个参数是一个结构，它使在模块中定义函数变得容易。在上面给出的示例中，mymethods结构将在文件的前面（通常在init name子例程之前）定义为：

static PyMethodDef mymethods[] = {
    { nokeywordfunc,nokeyword_cfunc,
      METH_VARARGS,
      Doc string},
    { keywordfunc, keyword_cfunc,
      METH_VARARGS|METH_KEYWORDS,
      Doc string},
    {NULL, NULL, 0, NULL} /* Sentinel */
}

mymethods数组中的每个条目都是 PyMethodDef 结构包含1）python名称，2）实现函数的C函数，3）指示此函数是否接受关键字的标志，以及4）函数的docstring。通过向该表中添加更多条目，可以为单个模块定义任意数量的函数。最后一个条目必须全部为空，如图所示，才能充当哨兵。python查找此条目以了解模块的所有函数都已定义。

完成扩展模块必须做的最后一件事是实际编写执行所需函数的代码。有两种函数：不接受关键字参数的函数和接受关键字参数的函数。

不带关键字参数的函数

不接受关键字参数的函数应编写为：

static PyObject*
nokeyword_cfunc (PyObject *dummy, PyObject *args)
{
    /* convert Python arguments */
    /* do function */
    /* return something */
}

此上下文中未使用伪参数，可以安全地忽略它。这个 args 参数包含作为元组传入函数的所有参数。此时您可以做任何您想做的事情，但通常管理输入参数的最简单方法是调用 PyArg_ParseTuple （args，格式化字符串，地址到变量…）或 PyArg_UnpackTuple （元组，“名称”，最小值，最大值，…）。pythoncapi参考手册第5.5节（解析参数和构建值）中包含了关于如何使用第一个函数的详细说明。您应该特别注意使用转换函数在Python对象和C对象之间传递的“O&”格式。所有其他格式函数都可以（大部分）看作是这个一般规则的特例。NumPy C-API中定义了几个可能有用的转换器函数。特别是 PyArray_DescrConverter 函数对于支持任意数据类型规范非常有用。此函数将任何有效的数据类型python对象转换为 PyArray_Descr* 对象。记住传递应该填写的C变量的地址。

有很多关于如何使用的例子 PyArg_ParseTuple 整个numpy源代码。标准用法如下：

PyObject *input;
PyArray_Descr *dtype;
if (!PyArg_ParseTuple(args, "OO&", &input,
                      PyArray_DescrConverter,
                      &dtype)) return NULL;

重要的是要记住 借来 使用“o”格式字符串时对对象的引用。但是，转换器功能通常需要某种形式的内存处理。在本例中，如果转换成功， D型 将保存对 PyArray_Descr* 对象，而 输入 将持有借来的参考资料。因此，如果此转换与另一个转换（例如整数）混合，并且数据类型转换成功，但整数转换失败，则需要在返回之前释放对数据类型对象的引用计数。一个典型的方法是设置 D型 到 NULL 打电话之前 PyArg_ParseTuple 然后使用 Py_XDECREF 在 D型 返回前。

在处理输入参数之后，将编写实际执行工作的代码（可能根据需要调用其他函数）。C函数的最后一步是返回一些东西。如果遇到错误，则 NULL 应返回（确保实际设置了错误）。如果不应返回任何内容，则递增 Py_None 然后把它还给我。如果应该返回单个对象，则返回该对象（确保您首先拥有对它的引用）。如果应该返回多个对象，那么需要返回一个元组。这个 Py_BuildValue （format_string，c_variables…）函数使从c变量构建python对象的元组变得容易。请特别注意格式字符串中的“n”和“o”之间的区别，否则很容易造成内存泄漏。“o”格式字符串增加 PyObject* 它对应的C变量，而“n”格式字符串窃取对相应的 PyObject* C变量。如果已经为对象创建了一个引用，并且只想将该引用赋给元组，那么应该使用“n”。如果只有对对象的借用引用，并且需要创建一个对象来为元组提供引用，那么应该使用“o”。

带关键字参数的函数

这些函数与没有关键字参数的函数非常相似。唯一的区别是函数签名是：

static PyObject*
keyword_cfunc (PyObject *dummy, PyObject *args, PyObject *kwds)
{
...
}

kwds参数包含一个python字典，其键是关键字参数的名称，其值是相应的关键字参数值。这本词典可以按你认为合适的方式处理。然而，处理它的最简单方法是替换 PyArg_ParseTuple （args，format_string，addresses…）函数，调用 PyArg_ParseTupleAndKeywords （args、kwds、format_-string、char * KWLIST[]，地址…）。此函数的kwlist参数是 NULL -终止了提供所需关键字参数的字符串数组。格式为“字符串”的每个条目应有一个字符串。如果传入无效的关键字参数，使用此函数将引发类型错误。

有关此函数的更多帮助，请参见Python文档中“扩展和嵌入”教程的第1.8节（扩展函数的关键字参数）。

参考计数

编写扩展模块时最大的困难是引用计数。这是一个重要的原因，流行的f2py，编织，赛通，ctypes等。如果处理错误的引用计数，则可能会从内存泄漏到分段错误。我所知道的正确处理参考计数的唯一策略是血液、汗液和眼泪。首先，将每个python变量都有一个引用计数强制输入您的大脑。然后，您可以准确地了解每个函数对对象的引用计数所做的操作，以便在需要时正确地使用decref和incref。引用计数可以真正测试您对编程过程的耐心和勤奋。尽管描述很严酷，但大多数引用计数的情况都非常简单，最常见的困难是，由于某些错误，在从例程早期退出之前，不在对象上使用减量法。第二，是传递给将要窃取引用的函数或宏的对象上没有引用的常见错误（ e.g. PyTuple_SET_ITEM 以及大多数功能 PyArray_Descr 对象）。

通常，当变量被创建或是某个函数的返回值时，您会得到一个新的变量引用（但是，也有一些突出的异常，比如从元组或字典中获取一个项）。当你拥有证明人时，你有责任确保 Py_DECREF （var）在不再需要变量时调用（并且没有其他函数“窃取”其引用）。另外，如果您要将一个python对象传递给一个将“窃取”引用的函数，那么您需要确保您拥有它（或使用 Py_INCREF 以获取您自己的参考）。您还将遇到借用参考资料的概念。借用引用的函数不会更改对象的引用计数，也不希望“保留”引用。它只是暂时使用对象。当你使用 PyArg_ParseTuple 或 PyArg_UnpackTuple 您将收到对元组中对象的借用引用，并且不应在函数内更改它们的引用计数。通过练习，你可以学会正确地进行参考计数，但一开始可能会令人沮丧。

引用计数错误的一个常见来源是 Py_BuildValue 功能。注意“n”格式字符和“o”格式字符之间的区别。如果在子例程（如输出数组）中创建一个新对象，并将其传递回一个返回值的元组，则最有可能在 Py_BuildValue . “o”字符将使引用计数增加一。这将给调用者留下两个新数组的引用计数。当变量被删除并且引用计数减1时，仍然会有额外的引用计数，并且数组将永远不会被释放。您将有一个引用计数引起的内存泄漏。使用“n”字符将避免这种情况，因为它将返回具有单个引用计数的对象（在元组内）。

处理数组对象

numpy的大多数扩展模块都需要访问ndarray对象（或其子类之一）的内存。要做到这一点，最简单的方法不需要你对numpy的内部了解太多。方法是

  1. 确保处理的是一个行为良好的数组（按机器字节顺序和单个段对齐），该数组的类型和维数正确。

    1. 通过使用 PyArray_FromAny 或者一个基于它的宏。

    2. 通过使用 PyArray_NewFromDescr 或者基于它的更简单的宏或函数。

  2. 获取数组的形状和指向其实际数据的指针。

  3. 将数据和形状信息传递给实际执行计算的子例程或代码的其他部分。

  4. 如果您正在编写算法，那么我建议您使用数组中包含的跨步信息来访问数组的元素（即 PyArray_GetPtr 宏使这变得无痛）。然后，您可以放宽您的需求，以避免强制执行单个段数组和可能导致的数据复制。

下面的小节将介绍这些子主题中的每一个。

转换任意序列对象

从任何可以转换为数组的python对象获取数组的主要例程是 PyArray_FromAny . 这个函数非常灵活，有许多输入参数。几个宏使基本函数的使用更加容易。 PyArray_FROM_OTF 可以说是这些宏中最常用的。它允许您将任意python对象转换为特定内置数据类型的数组。（ e.g. float），同时指定一组特定的需求（ e.g. 连续、对齐和可写）。语法是

PyArray_FROM_OTF

从任何python对象返回一个ndarray， obj ，可以转换为数组。返回数组中的维度数由对象确定。返回数组的所需数据类型在中提供 铅字 应该是枚举类型之一。这个 要求 对于返回的数组，可以是标准数组标志的任意组合。下面将更详细地解释这些论点中的每一个。成功时将收到对数组的新引用。失败论 NULL 返回并设置异常。

obj

对象可以是任何可转换为ndarray的python对象。如果对象已经是满足需求的ndarray（的子类），则返回新的引用。否则，将构造一个新的数组。内容 obj 除非使用数组接口，否则将复制到新数组，这样就不必复制数据。可以转换为数组的对象包括：1）任何嵌套的序列对象，2）任何暴露数组接口的对象，3）具有 __array__ 方法（应该返回一个ndarray）和4）任何标量对象（成为零维数组）。将通过以其他方式符合要求的ndarray的子类。如果要确保基类的ndarray，请使用 NPY_ARRAY_ENSUREARRAY 在需求标志中。只有在必要时才复印。如果你想保证一份，那就交给我 NPY_ARRAY_ENSURECOPY 到需求标志。

铅字

枚举类型之一或 NPY_NOTYPE 如果应该从对象本身确定数据类型。基于C的名称可以使用：

或者，可以使用平台上支持的位宽度名称。例如：

只有在不丢失精度的情况下才能将对象转换为所需的类型。否则 NULL 将返回并引发错误。使用 NPY_ARRAY_FORCECAST 在“需求”标志中重写此行为。

要求

一个ndarray的内存模型允许在每个维度中任意前进到数组的下一个元素。但是，通常需要与需要C-连续或Fortran连续内存布局的代码进行接口。此外，如果尝试将指针取消引用到数组数据中，则ndarray可能会错位（元素的地址不是元素大小的整数倍），这可能会导致程序崩溃（或至少工作得更慢）。这两个问题都可以通过将python对象转换为一个更“行为良好”的数组来解决，以满足您的特定使用。

需求标志允许指定什么类型的数组是可接受的。如果传入的对象不满足此要求，则会制作一个副本，以便thre返回的对象满足这些要求。这些ndarray可以使用非常通用的内存指针。此标志允许指定返回的数组对象的所需属性。所有标志都在详细的API章节中解释。最常用的标志是 NPY_ARRAY_IN_ARRAY ， NPY_OUT_ARRAY 和 NPY_ARRAY_INOUT_ARRAY ：

NPY_ARRAY_IN_ARRAY

此标志对于必须按C-连续顺序排列并对齐的数组很有用。这些类型的数组通常是某些算法的输入数组。

NPY_ARRAY_OUT_ARRAY

此标志用于指定以C-连续顺序排列、对齐并可以写入的数组。这样的数组通常作为输出返回（尽管通常这样的输出数组是从零开始创建的）。

NPY_ARRAY_INOUT_ARRAY

此标志用于指定将用于输入和输出的数组。 PyArray_ResolveWritebackIfCopy 必须在之前调用 Py_DECREF 在接口例程的末尾，将临时数据写回传入的原始数组。使用 NPY_ARRAY_WRITEBACKIFCOPY 或 NPY_ARRAY_UPDATEIFCOPY 标志要求输入对象已经是一个数组（因为其他对象不能以这种方式自动更新）。如果发生错误，请使用 PyArray_DiscardWritebackIfCopy （obj）在设置了这些标志的数组上。这会将基础基数组设置为可写，而不会导致内容被复制回原始数组。

其他可以作为附加要求的有用标志包括：

NPY_ARRAY_FORCECAST

强制转换为所需的类型，即使它不能在不丢失信息的情况下完成。

NPY_ARRAY_ENSURECOPY

确保生成的数组是原始数组的副本。

NPY_ARRAY_ENSUREARRAY

确保生成的对象是实际的ndarray，而不是子类。

注解

数组是否进行字节交换取决于数组的数据类型。本机字节顺序数组总是由 PyArray_FROM_OTF 所以没有必要 NPY_ARRAY_NOTSWAPPED 需求参数中的标志。也没有办法从这个例程中获得字节交换数组。

创建全新的日历

通常，必须从扩展模块代码中创建新的数组。可能需要一个输出数组，而您不希望调用者必须提供它。可能只需要一个临时数组来保存中间计算。无论需要什么，都有简单的方法来获取所需数据类型的ndarray对象。最常用的功能是 PyArray_NewFromDescr . 所有的数组创建函数都要经过这段大量重复使用的代码。由于它的灵活性，使用起来可能有些混乱。因此，存在更简单的形式，更易于使用。这些表格是 PyArray_SimpleNew 一系列函数，通过为常见用例提供默认值来简化接口。

获取ndarray内存并访问ndarray的元素

如果obj是ndarray (PyArrayObject* ，则该空区域指向该数据数组的数据区域。 * pointer PyArray_DATA (obj) or the char* 指针 PyArray_BYTES （OBJ）。请记住（通常）此数据区域可能不会根据数据类型进行对齐，它可能表示字节交换数据，和/或它可能不可写。如果数据区域是按本机字节顺序对齐的，那么如何获取数组的特定元素只由npy-intp变量数组决定。 PyArray_STRIDES （OBJ）。特别是，这个整数的C数组显示了 字节 必须添加到当前元素指针才能到达每个维度中的下一个元素。对于小于4维的数组 PyArray_GETPTR{{k}} （obj，…）宏，其中k是整数1、2、3或4，使使用数组步幅更容易。参数….在数组中表示非负整数索引。例如，假设 E 是一个三维数组。指向元素的（void*）指针 E[i,j,k] 获得作为 PyArray_GETPTR3 （E，I，J，K）。

如前所述，C样式的连续数组和Fortran样式的连续数组具有特定的跨线模式。两个数组标志 (NPY_ARRAY_C_CONTIGUOUS 和 NPY_ARRAY_F_CONTIGUOUS ）指示特定数组的删除模式是否与C样式continuous或Fortran样式continuous匹配。无论划线模式是否与标准C或Fortran模式匹配，都可以使用 PyArray_IS_C_CONTIGUOUS （OBJ）和 PyArray_ISFORTRAN （obj）分别。大多数第三方库需要连续的数组。但是，通常支持通用跨步并不困难。我鼓励您尽可能在自己的代码中使用删除信息，并保留包装第三方代码的单个段要求。使用与数据阵列一起提供的划线信息而不需要连续的划线，可以减少必须进行的复制。

例子

下面的示例演示如何编写一个接受两个输入参数（将转换为数组）和一个输出参数（必须是数组）的包装器。函数返回none并更新输出数组。注意更新了对numpy v1.14及更高版本的writebackifcopy语义的使用

static PyObject *
example_wrapper(PyObject *dummy, PyObject *args)
{
    PyObject *arg1=NULL, *arg2=NULL, *out=NULL;
    PyObject *arr1=NULL, *arr2=NULL, *oarr=NULL;

    if (!PyArg_ParseTuple(args, "OOO!", &arg1, &arg2,
        &PyArray_Type, &out)) return NULL;

    arr1 = PyArray_FROM_OTF(arg1, NPY_DOUBLE, NPY_ARRAY_IN_ARRAY);
    if (arr1 == NULL) return NULL;
    arr2 = PyArray_FROM_OTF(arg2, NPY_DOUBLE, NPY_ARRAY_IN_ARRAY);
    if (arr2 == NULL) goto fail;
#if NPY_API_VERSION >= 0x0000000c
    oarr = PyArray_FROM_OTF(out, NPY_DOUBLE, NPY_ARRAY_INOUT_ARRAY2);
#else
    oarr = PyArray_FROM_OTF(out, NPY_DOUBLE, NPY_ARRAY_INOUT_ARRAY);
#endif
    if (oarr == NULL) goto fail;

    /* code that makes use of arguments */
    /* You will probably need at least
       nd = PyArray_NDIM(<..>)    -- number of dimensions
       dims = PyArray_DIMS(<..>)  -- npy_intp array of length nd
                                     showing length in each dim.
       dptr = (double *)PyArray_DATA(<..>) -- pointer to data.

       If an error occurs goto fail.
     */

    Py_DECREF(arr1);
    Py_DECREF(arr2);
#if NPY_API_VERSION >= 0x0000000c
    PyArray_ResolveWritebackIfCopy(oarr);
#endif
    Py_DECREF(oarr);
    Py_INCREF(Py_None);
    return Py_None;

 fail:
    Py_XDECREF(arr1);
    Py_XDECREF(arr2);
#if NPY_API_VERSION >= 0x0000000c
    PyArray_DiscardWritebackIfCopy(oarr);
#endif
    Py_XDECREF(oarr);
    return NULL;
}
使用numpy c-api 使用python作为胶水