使用 NumPy 数组的 C 扩展¶
| 日期 | 2010-06-20（最后修改），2006-12-08（创建） |
|---|
我编写了几个处理 NumPy 数组的 C 扩展。它们很简单，但似乎运行良好。它们将向您展示如何将 Python 变量和 NumPy 数组传递给您的 C 代码。一旦您学会了如何做，它就非常简单。我怀疑它们足以满足大多数数值代码的需求。我将其写成草稿，并提供了代码和文档文件。我发现对于我的数值需求，我实际上只需要传递有限的一组东西（整数、浮点数、字符串和 NumPy 数组）。如果这是您的类别，此代码可能会对您有所帮助。
我已经测试了这些例程，到目前为止一切顺利，但我不能保证任何事情。我对这方面还比较陌生。如果您发现任何错误，请在 SciPy 邮件列表中发布消息。
下面是包含代码和文档的 tar 包的链接。
我最近更新了一些信息并添加了更多示例。下面提供的文档是原始文档，仍然有用。下面的链接包含最新的文档和源代码。
-- Lou Pecora
NumPy 和 Python 的 C 扩展¶
作者：Lou Pecora - 2006-12-07（草稿版本 0.1）
概述¶
简介 - 一些背景¶
在我使用 Python 的过程中，我遇到了一个典型的问题：我需要加快代码的特定部分。我不是 Python 专家，也不是任何类型的编码/计算机专家。我使用 Python 进行数值计算，并且大量使用 Numeric/NumPy。几乎每本 Python 书籍或教程都会告诉你，当你需要一个快速运行的例程时，就构建 Python 的 C 扩展。C 扩展是可以用 C 代码编写的，可以编译并链接到一个共享库，该库可以像任何 Python 模块一样导入，你可以像调用 Python 函数一样调用指定的 C 例程。
听起来不错，但我有一些保留意见。它看起来并不简单（在某种程度上确实如此）。所以我寻找其他解决方案。我找到了它们。它们是诸如 SWIG、Pyrex、ctypes、Psyco 和 Weave 之类的方案。当我尝试这些方案时，我经常能使给出的简单示例正常工作（但并非全部）。但是，当我尝试将它们应用于 NumPy 时，我遇到了障碍。然后就会涉及到类型映射或其他混合结构。我并不是在贬低这些方法，但我始终无法弄清楚它们，也无法在我的代码上运行，尽管我阅读了大量的在线教程，并从各种 Python 支持组和邮件列表中获得了有益的建议。
所以我决定看看是否可以自己编写 C 扩展。我得到了帮助，形式是 2000 年左右关于使用 Numeric 的一些简单的 C 扩展示例，这些示例来自康奈尔大学的 Tom Loredo。这些示例一直保存在我的硬盘驱动器中，直到 5 年后，出于绝望，我将它们拿出来，并使用他的示例，我能够快速地将几个 C 扩展组合在一起，这些扩展（至少对我来说）处理了我想要加速的所有情况（到目前为止）。这些情况主要涉及传递 Python 整数、浮点数（=C 双精度浮点数）、字符串以及 NumPy 1D 和 2D 浮点数和整数数组。我很少需要将其他任何东西传递给 C 例程来进行计算。如果你和我处境相同，那么我整理的这个包可能会对你有所帮助。一旦你开始使用，它就会变得相当容易。
请注意，尽管我深受 Tom Loredo 的帮助，但他对我的代码或说明中的任何错误不承担任何责任。同样，此代码仅供研究使用。它只通过我的开发和使用进行了测试。它没有保证，也不提供任何担保。不要在可能造成生命、肢体、财产或金钱损失或任何你或其他人珍视的东西的威胁的地方使用此代码。
我在一台使用系统 OS X 10.4（基本 BSD Unix）、Python 2.4、NumPy 0.9x 和 GNU 编译器和链接器 gcc 的 Macintosh G4 笔记本电脑上开发了这些 C 扩展及其 Python 包装器。我认为我在这里告诉你的大部分内容很容易移植到 Linux 和其他 Unix 系统，而不仅仅是 Mac。我不确定 Windows。我希望我的低级方法也能让 Windows 用户轻松使用。
扩展的代码（C 和 Python）可能看起来很多，但它非常重复。一旦你掌握了某个扩展函数的主要方案，你就会看到它在所有其他函数中一遍又一遍地重复出现，只是为了处理不同的参数或向调用例程返回不同的对象而略有变化。不要被代码吓倒。好消息是，对于许多数值用途，扩展将遵循相同的格式，因此你可以快速重用你已经为新项目编写的代码。专注于一个扩展函数并详细了解它（事实上，我将在下面这样做）。一旦你理解了它，其他例程就几乎是显而易见的。对于与包一起提供的几个实用程序函数也是如此。它们帮助你创建、测试和操作数据，并且它们也具有很多重复性。实用程序函数也非常短小简单，所以不用担心。
NumPy 扩展的一般方案¶
这将在下面详细介绍，但首先我想让你了解每个扩展是如何组织的。
在 C 源文件中的 C 扩展函数之前必须完成三件事。
你必须包含 Python 和 NumPy 头文件。
每个扩展必须在文件开头定义的结构中命名。这是一个用于从 Python 函数访问扩展的名称。
接下来，进行一组初始化调用以设置 Python 和 NumPy 调用和接口。对于所有涉及 NumPy 和 Python 的扩展，它将是相同的，除非你添加扩展以访问 NumPy 数组之外的其他 Python 包或类。我在这里不会介绍任何这些内容（因为我不知道）。因此，init 调用可以复制到每个扩展文件。
每个 C 扩展将具有以下形式。
参数将始终相同：（
PyObject *self，PyObject *args） - 如果你不知道这些到底是什么，不要担心。它们是指向一般 Python 对象的指针，它们是由你将从 NumPy 和 Python 本身使用的头文件自动提供的。你只需要知道这些。args 通过一个解析它们的函数调用进行处理，并将它们分配给 C 定义的对象。
接下来，该解析的结果可能会由一个实用程序例程进行检查，该例程进入表示对象的结构，并确保数据是正确的类型（浮点数或整数、一维或二维数组等）。虽然我包含了一些这些 C 级别的检查，但你会发现我认为最好在用于包装 C 扩展的 Python 函数中完成它们。它们在 Python 中也更容易完成。在我的调用 Python 包装器中有很多数据检查。通常，这不会导致太多开销，因为你不会在某个循环中数百万次调用这些扩展，而是将它们用作通往 C 或 C++ 例程的入口，以执行长时间、复杂、重复、专门的计算。
在一些可能的数据检查之后，C 数据类型被初始化为指向 NumPy 数组的数据部分，并借助实用程序函数。
提取下一维度的信息，以便您了解列数、行数、向量维度等。
现在您可以使用 C 数组来操作 NumPy 数组中的数据。C 数组和来自上述解析的 C 数据指向原始的 Python/NumPy 数据，因此您所做的更改会在扩展返回 Python 后影响数组值。您可以将数组传递给执行计算等的其它 C 函数。请记住，您正在操作原始的 NumPy 矩阵和向量。
在计算完成后，您必须释放为 NumPy 数组的 C 数据构造分配的任何内存。这再次由实用程序函数完成。只有在您分配了内存来处理数组（例如在矩阵例程中）时才需要执行此步骤，但如果您没有分配内存（例如在向量例程中），则不需要执行此步骤。
最后，您通过返回 Python 值或 NumPy 数组来返回到 Python 调用函数。我有一些 C 扩展，其中展示了这两种方法的示例。
Python 包装器函数¶
最好通过调用一个然后调用扩展的 Python 函数来调用 C 扩展。这被称为 Python 包装器函数。它使您的代码更具 Python 风格（例如，您可以轻松地使用关键字），并且正如我上面指出的那样，它允许您在将函数参数和数据传递给 C 扩展和其他 C 函数进行大型计算之前，轻松地检查它们是否正确。这似乎是一个不必要的额外步骤，但它值得这样做。
代码¶
在本节中，我将参考源文件 C_arraytest.h、C_arraytest.c、C_arraytest.py 和 C_arraytest.mak 中的代码。您应该将这些文件放在手边（可能打印出来），以便您可以按照下面对代码的解释进行操作。
C 代码 - 带有实用程序的一个详细示例¶
首先，我将使用来自 C_arraytest.h、C_arraytest.c 的代码示例，用于名为 matsq 的例程。此函数以 (NumPy) 矩阵 A、整数 i 和 (Python) 浮点数 y 作为输入，并输出一个返回 (NumPy) 矩阵 B，其每个分量等于输入矩阵分量的平方乘以整数乘以浮点数。数学上
$B_{ij} = i y (A_{ij})^2$
调用 matsq 例程的 Python 代码为 A=matsq(B,i,y)。以下是相关代码
头文件，C_arraytest.h
/* Header to test of C modules for arrays for Python: C_test.c */ /* ==== Prototypes =================================== */ // .... Python callable Vector functions .................. static PyObject *vecfcn1(PyObject *self, PyObject *args); static PyObject *vecsq(PyObject *self, PyObject *args); /* .... C vector utility functions ..................*/ PyArrayObject *pyvector(PyObject *objin); double *pyvector_to_Carrayptrs(PyArrayObject *arrayin); int not_doublevector(PyArrayObject *vec); /* .... Python callable Matrix functions ..................*/ static PyObject *rowx2(PyObject *self, PyObject *args); static PyObject *rowx2_v2(PyObject *self, PyObject *args); static PyObject *matsq(PyObject *self, PyObject *args); static PyObject *contigmat(PyObject *self, PyObject *args); /* .... C matrix utility functions ..................*/ PyArrayObject *pymatrix(PyObject *objin); double **pymatrix_to_Carrayptrs(PyArrayObject *arrayin); double **ptrvector(long n); void free_Carrayptrs(double **v); int not_doublematrix(PyArrayObject *mat); /* .... Python callable integer 2D array functions ..................*/ static PyObject *intfcn1(PyObject *self, PyObject *args); /* .... C 2D int array utility functions ..................*/ PyArrayObject *pyint2Darray(PyObject *objin); int **pyint2Darray_to_Carrayptrs(PyArrayObject *arrayin); int **ptrintvector(long n); void free_Cint2Darrayptrs(int **v); int not_int2Darray(PyArrayObject *mat);
源文件，C_arraytest.c
/* A file to test imorting C modules for handling arrays to Python */ #include "Python.h" #include "arrayobject.h" #include "C_arraytest.h" #include/* #### Globals #################################### */ /* ==== Set up the methods table ====================== */ static PyMethodDef _C_arraytestMethods[] = { {"vecfcn1", vecfcn1, METH_VARARGS}, {"vecsq", vecsq, METH_VARARGS}, {"rowx2", rowx2, METH_VARARGS}, {"rowx2_v2", rowx2_v2, METH_VARARGS}, {"matsq", matsq, METH_VARARGS}, {"contigmat", contigmat, METH_VARARGS}, {"intfcn1", intfcn1, METH_VARARGS}, {NULL, NULL} /* Sentinel - marks the end of this structure */ }; /* ==== Initialize the C_test functions ====================== */ // Module name must be _C_arraytest in compile and linked void init_C_arraytest() { (void) Py_InitModule("_C_arraytest", _C_arraytestMethods); import_array(); // Must be present for NumPy. Called first after above line. } /* #### Vector Extensions ############################## */ /* ==== vector function - manipulate vector in place ====================== Multiply the input by 2 x dfac and put in output Interface: vecfcn1(vec1, vec2, str1, d1) vec1, vec2 are NumPy vectors, str1 is Python string, d1 is Python float (double) Returns integer 1 if successful */ static PyObject *vecfcn1(PyObject *self, PyObject *args) { PyArrayObject *vecin, *vecout; // The python objects to be extracted from the args double *cin, *cout; // The C vectors to be created to point to the // python vectors, cin and cout point to the row // of vecin and vecout, respectively int i,j,n; const char *str; double dfac; /* Parse tuples separately since args will differ between C fcns */ if (!PyArg_ParseTuple(args, "O!O!sd", &PyArray_Type, &vecin, &PyArray_Type, &vecout, &str, &dfac)) return NULL; if (NULL == vecin) return NULL; if (NULL == vecout) return NULL; // Print out input string printf("Input string: %s\n", str); /* Check that objects are 'double' type and vectors Not needed if python wrapper function checks before call to this routine */ if (not_doublevector(vecin)) return NULL; if (not_doublevector(vecout)) return NULL; /* Change contiguous arrays into C * arrays */ cin=pyvector_to_Carrayptrs(vecin); cout=pyvector_to_Carrayptrs(vecout); /* Get vector dimension. */ n=vecin->dimensions[0]; /* Operate on the vectors */ for ( i=0; i dimensions[0]; /* Make a new double vector of same dimension */ vecout=(PyArrayObject *) PyArray_FromDims(1,dims,NPY_DOUBLE); /* Change contiguous arrays into C *arrays */ cin=pyvector_to_Carrayptrs(vecin); cout=pyvector_to_Carrayptrs(vecout); /* Do the calculation. */ for ( i=0; i dimensions[0]; return (double *) arrayin->data; /* pointer to arrayin data as double */ } /* ==== Check that PyArrayObject is a double (Float) type and a vector ============== return 1 if an error and raise exception */ int not_doublevector(PyArrayObject *vec) { if (vec->descr->type_num != NPY_DOUBLE || vec->nd != 1) { PyErr_SetString(PyExc_ValueError, "In not_doublevector: array must be of type Float and 1 dimensional (n)."); return 1; } return 0; } /* #### Matrix Extensions ############################## */ /* ==== Row x 2 function - manipulate matrix in place ====================== Multiply the 2nd row of the input by 2 and put in output interface: rowx2(mat1, mat2) mat1 and mat2 are NumPy matrices Returns integer 1 if successful */ static PyObject *rowx2(PyObject *self, PyObject *args) { PyArrayObject *matin, *matout; // The python objects to be extracted from the args double **cin, **cout; // The C matrices to be created to point to the // python matrices, cin and cout point to the rows // of matin and matout, respectively int i,j,n,m; /* Parse tuples separately since args will differ between C fcns */ if (!PyArg_ParseTuple(args, "O!O!", &PyArray_Type, &matin, &PyArray_Type, &matout)) return NULL; if (NULL == matin) return NULL; if (NULL == matout) return NULL; /* Check that objects are 'double' type and matrices Not needed if python wrapper function checks before call to this routine */ if (not_doublematrix(matin)) return NULL; if (not_doublematrix(matout)) return NULL; /* Change contiguous arrays into C ** arrays (Memory is Allocated!) */ cin=pymatrix_to_Carrayptrs(matin); cout=pymatrix_to_Carrayptrs(matout); /* Get matrix dimensions. */ n=matin->dimensions[0]; m=matin->dimensions[1]; /* Operate on the matrices */ for ( i=0; i dimensions[0]; m=matin->dimensions[1]; /* Operate on the matrices */ for ( i=0; i dimensions[0]; m=dims[1]=matin->dimensions[1]; /* Make a new double matrix of same dims */ matout=(PyArrayObject *) PyArray_FromDims(2,dims,NPY_DOUBLE); /* Change contiguous arrays into C ** arrays (Memory is Allocated!) */ cin=pymatrix_to_Carrayptrs(matin); cout=pymatrix_to_Carrayptrs(matout); /* Do the calculation. */ for ( i=0; i dimensions[0]; m=dims[1]=matin->dimensions[1]; ncomps=n*m; /* Make a new double matrix of same dims */ matout=(PyArrayObject *) PyArray_FromDims(2,dims,NPY_DOUBLE); /* Change contiguous arrays into C * arrays pointers to PyArrayObject data */ cin=pyvector_to_Carrayptrs(matin); cout=pyvector_to_Carrayptrs(matout); /* Do the calculation. */ printf("In contigmat, cout (as contiguous memory) =\n"); for ( i=0; i dimensions[0]; m=arrayin->dimensions[1]; c=ptrvector(n); a=(double *) arrayin->data; /* pointer to arrayin data as double */ for ( i=0; i descr->type_num != NPY_DOUBLE || mat->nd != 2) { PyErr_SetString(PyExc_ValueError, "In not_doublematrix: array must be of type Float and 2 dimensional (n x m)."); return 1; } return 0; } /* #### Integer 2D Array Extensions ############################## */ /* ==== Integer function - manipulate integer 2D array in place ====================== Replace >=0 integer with 1 and < 0 integer with 0 and put in output interface: intfcn1(int1, afloat) int1 is a NumPy integer 2D array, afloat is a Python float Returns integer 1 if successful */ static PyObject *intfcn1(PyObject *self, PyObject *args) { PyArrayObject *intin, *intout; // The python objects to be extracted from the args int **cin, **cout; // The C integer 2D arrays to be created to point to the // python integer 2D arrays, cin and cout point to the rows // of intin and intout, respectively int i,j,n,m, dims[2]; double afloat; /* Parse tuples separately since args will differ between C fcns */ if (!PyArg_ParseTuple(args, "O!d", &PyArray_Type, &intin, &afloat)) return NULL; if (NULL == intin) return NULL; printf("In intfcn1, the input Python float = %e, a C double\n",afloat); /* Check that object input is int type and a 2D array Not needed if python wrapper function checks before call to this routine */ if (not_int2Darray(intin)) return NULL; /* Get the dimensions of the input */ n=dims[0]=intin->dimensions[0]; m=dims[1]=intin->dimensions[1]; /* Make a new int array of same dims */ intout=(PyArrayObject *) PyArray_FromDims(2,dims,NPY_LONG); /* Change contiguous arrays into C ** arrays (Memory is Allocated!) */ cin=pyint2Darray_to_Carrayptrs(intin); cout=pyint2Darray_to_Carrayptrs(intout); /* Do the calculation. */ for ( i=0; i = 0) { cout[i][j]= 1; } else { cout[i][j]= 0; } } } printf("In intfcn1, the output array is,\n\n"); for ( i=0; i dimensions[0]; m=arrayin->dimensions[1]; c=ptrintvector(n); a=(int *) arrayin->data; /* pointer to arrayin data as int */ for ( i=0; i descr->type_num != NPY_LONG || mat->nd != 2) { PyErr_SetString(PyExc_ValueError, "In not_int2Darray: array must be of type int and 2 dimensional (n x m)."); return 1; } return 0; } // EOF
现在，让我们看看源代码的较小部分。
头文件¶
您必须包含以下头文件，其中 Python.h **始终** 是第一个包含的头文件。
#include "Python.h"
#include "arrayobject.h"
我还包含头文件 C_arraytest.h，其中包含 matsq 函数的原型。
static PyObject *matsq(PyObject *self, PyObject *args);
函数声明前面的 static 关键字使此函数对您的扩展模块私有。链接器将无法看到它。这样，您可以对所有扩展模块使用相同的直观函数名（例如 sum、check、trace），而不会在链接时发生名称冲突。函数的类型是 PyObject *，因为它将始终返回到 Python 调用函数，因此您可以（实际上必须）返回一个 Python 对象。参数始终相同，
PyObject *self and PyObject *args
第一个 self 从未使用，但由于 Python 传递参数的方式而必不可少。第二个 args 是指向 Python 元组的指针，该元组包含函数的所有参数（B、i、x）。
方法定义¶
这将设置一个函数名称表，该表将成为您的 Python 代码与您的 C 扩展之间的接口。C 扩展模块的名称将为 _C_arraytest（注意前导下划线）。正确获取名称在每次使用时都很重要，因为对代码中使用模块名称有严格的要求。名称首先出现在方法定义表中，作为表名称的第一部分
static PyMethodDef _C_arraytestMethods[] = {
...,
{"matsq", matsq, METH_VARARGS},
...,
{NULL, NULL} /* Sentinel - marks the end of this structure */
};
我使用省略号 (...) 来忽略与该函数无关的其他代码。METH_VARARGS 参数告诉编译器您将以通常的方式传递参数，而无需使用关键字，如上面的示例 A=matsq(B,i,x) 中所示。有一些方法可以使用 Python 关键字，但我还没有尝试过。该表应始终以 {NULL, NULL} 结尾，这只是一个“标记”，用于标记表的结尾。
初始化¶
这些函数告诉 Python 解释器在加载模块时调用什么。请注意，模块的名称 (_C_arraytest) 必须直接位于初始化结构名称中的 init 之后。
void init_C_arraytest() {
(void) Py_InitModule("_C_arraytest", _C_arraytestMethods);
import_array(); // Must be present for NumPy. Called first after above line.
}
顺序很重要，您必须首先调用这两个初始化函数。
matsqfunction 代码¶
现在，这是您将从 Python 代码中调用的实际函数。我将把它拆开并解释每个部分。
函数名称和类型
static PyObject *matsq(PyObject *self, PyObject *args)
您可以看到它们与 C_arraytest.h 中的原型匹配。
局部变量
PyArrayObject *matin, *matout;
double **cin, **cout, dfactor;
int i,j,n,m, dims[2], ifactor;
PyArrayObjects 是在 NumPy 头文件中定义的结构，它们将被分配指向实际输入和输出 NumPy 数组 (A 和 B) 的指针。C 数组 cin 和 cout 是 C 指针，它们将指向 (最终) NumPy 数组中的实际数据，并允许您操作它。变量 dfactor 将是 Python 浮点数 y，ifactor 将是 Python 整数 i，变量 i、j、n 和 m 将是循环变量 (i 和 j) 以及 A 和 B 中的矩阵维度 (n=行数，m=列数)。数组 dims 将用于从 NumPy 数组访问 n 和 m。所有这些都在下面发生。首先，我们必须从 args 元组中提取输入变量 (A、i、y)。这是通过调用完成的，
/* Parse tuples separately since args will differ between C fcns */
if (!PyArg_ParseTuple(args, "O!id",
&PyArray_Type, &matin, &ifactor, &dfactor)) return NULL;
PyArg_ParseTuple 函数接受 args 元组，并使用下一个出现的格式字符串 ("O!id") 将元组的每个成员分配给一个 C 变量。请注意，您必须通过引用传递所有 C 变量。即使 C 变量是指向字符串的指针，这也是正确的 (参见 vecfcn1 例程中的代码)。格式字符串告诉解析函数使用什么类型的变量。Python 的常用变量都有字母名称 (例如，s 代表字符串，i 代表整数，d 代表 (double - Python 浮点数))。您可以在 Guido 的教程 (https://docs.pythonlang.cn/ext/ext.html) 中找到这些变量以及更多变量的列表。对于像 NumPy 数组一样不在标准 Python 中的数据类型，您使用 O! 符号，它告诉解析器查找类型结构 (在本例中为 NumPy 结构 PyArray_Type) 来帮助它转换将被分配给指向 NumPy 数组结构的局部变量 (matin) 的元组成员。请注意，这些也是通过引用传递的。必须保持顺序并与您想要的 Python 函数的调用接口匹配。格式字符串定义了接口，如果您没有从 Python 调用函数，因此参数数量与格式字符串中的数量匹配，您将收到错误。这很好，因为它将指出问题所在。
如果这不起作用，我们将返回 NULL，这将导致 Python 异常。
if (NULL == matin) return NULL;
接下来，我们检查输入矩阵是否确实是 NumPy 类型 double 的矩阵。此测试也在此 C 扩展的 Python 包装器中进行。最好在那里进行，但我在这里包含此测试以向您展示您可以在 C 扩展中进行测试，并且可以“进入”NumPy 结构以提取其参数。实用程序函数 not_doublematrix 稍后解释。
/* Check that object input is 'double' type and a matrix
Not needed if python wrapper function checks before call to this routine */
if (not_doublematrix(matin)) return NULL;
以下是如何进入 NumPy 结构以获取矩阵 matin 的维度并将它们分配给上面提到的局部变量的示例。
/* Get the dimensions of the input */
n=dims[0]=matin->dimensions[0];
m=dims[1]=matin->dimensions[1];
现在，我们使用这些矩阵参数在我们的 C 扩展中生成一个新的 NumPy 矩阵 matout（我们的输出）。PyArray_FromDims(2,dims,NPY_DOUBLE) 是 NumPy 提供的实用程序函数（不是我），它的参数告诉 NumPy NumPy 对象的秩（2）、每个维度的尺寸（dims）和数据类型（NPY_DOUBLE）。其他 C 扩展中还有创建不同 NumPy 数组的其他示例。
/* Make a new double matrix of same dims */
matout=(PyArrayObject *) PyArray_FromDims(2,dims,NPY_DOUBLE);
为了实际进行我们的计算，我们需要 C 结构来处理我们的数据，因此我们生成两个 C 二维数组（cin 和 cout），它们将分别指向 matin 和 matout 中的数据。请注意，这里分配了内存，因为我们需要创建一个指向 C doubles 的指针数组，以便我们可以像使用两个索引的常规 C 矩阵一样寻址 cin 和 cout。此内存必须在此 C 扩展结束时释放。这种内存分配并不总是必要的。请参阅 C 扩展中用于 NumPy 向量操作和将 NumPy 矩阵视为连续数组（如 NumPy 中那样）的例程（例程 contigmat）。
/* Change contiguous arrays into C ** arrays (Memory is Allocated!) */
cin=pymatrix_to_Carrayptrs(matin);
cout=pymatrix_to_Carrayptrs(matout);
最后，我们到达可以操作矩阵并进行计算的点。以下是执行原始方程运算 $B_{ij}= i y (A_{ij})^2$ 的部分。请注意，我们正在直接操作传递到此扩展的原始 NumPy 数组 A 和 B 中的数据。因此，您在此对 cin 或 cout 的组件所做的任何操作都将对原始矩阵进行，并且当您返回 Python 代码时，这些操作将出现在那里。
/* Do the calculation. */
for ( i=0; i<n; i++) {
for ( j=0; j<m; j++) {
cout[i][j]= ifactor*dfactor*cin[i][j]*cin[i][j];
}
}
我们准备返回 Python 调用例程，但首先我们要释放为 cin 和 cout 分配的内存。
/* Free memory, close file and return */
free_Carrayptrs(cin);
free_Carrayptrs(cout);
现在我们返回计算结果。
return PyArray_Return(matout);
如果您查看其他 C 扩展，您会发现您还可以使用另一个 Python 提供的函数 Py_BuildValue("i", 1) 返回常规 Python 变量（如整数），其中字符串 "i" 告诉函数数据类型，第二个参数（这里为 1）是要返回的数据值。如果您决定不返回任何内容，您 '''必须''' 返回 Python 关键字 None，如下所示
Py_INCREF(Py_None);
return Py_None;
Py_INCREF 函数增加了对 None 的引用次数（记住 Python 在没有更多对数据的引用时会自动收集分配的内存）。您必须在 C 扩展中小心这一点。有关更多信息，请参阅 Guido 的教程。
实用程序函数¶
以下是矩阵实用函数的简要描述。它们基本上是不言自明的。向量和整数数组实用函数非常相似。
第一个实用函数没有在任何 C 扩展中使用，但我将其包含在内，因为一位热心人士将其与一些代码一起发送，并且它确实展示了如何将 Python 对象转换为 NumPy 数组。我还没有尝试过。使用风险自负。
PyArrayObject *pymatrix(PyObject *objin) {
return (PyArrayObject *) PyArray_ContiguousFromObject(objin,
NPY_DOUBLE, 2,2);
}
下一个函数创建用于指向 NumPy 矩阵行的 C 数组。这会分配指向 NumPy 数据的指针数组。NumPy 数据是连续的，并且使用步幅 (m) 来访问每一行。此函数调用 ptrvector(n)，它执行实际的内存分配。请记住在使用完此函数后释放内存。
double **pymatrix_to_Carrayptrs(PyArrayObject *arrayin) {
double **c, *a;
int i,n,m;
n=arrayin->dimensions[0];
m=arrayin->dimensions[1];
c=ptrvector(n);
a=(double *) arrayin->data; /* pointer to arrayin data as double */
for ( i=0; i<n; i++) {
c[i]=a+i*m;
}
return c;
}
这是为 C 指针数组分配内存的地方。它是一个非常标准的数组内存分配器。
double **ptrvector(long n) {
double **v;
v=(double **)malloc((size_t) (n*sizeof(double)));
if (!v) {
printf("In **ptrvector. Allocation of memory for double array failed.");
exit(0);
}
return v;
}
这是释放内存的例程。
void free_Carrayptrs(double **v) {
free((char*) v);
}
注意：有一个标准的 C-API 用于将 Python 对象转换为 C 样式的指针数组，称为 PyArray_AsCArray
这是一个实用函数，用于检查由解析生成的物体是否为 NumPy 矩阵。您可以看到它如何进入 NumPy 对象结构。
int not_doublematrix(PyArrayObject *mat) {
if (mat->descr->type_num != NPY_DOUBLE || mat->nd != 2) {
PyErr_SetString(PyExc_ValueError,
"In not_doublematrix: array must be of type Float and 2 dimensional (n x m).");
return 1; }
return 0;
}C 代码 – 其他变体¶
正如我在引言中提到的，这些函数是重复的。所有其他函数都遵循非常相似的模式。它们在方法结构中有一行，它们具有相同的参数，它们解析参数，它们可能会在解析后检查 C 结构，它们设置用于操作的指向输入数据的 C 变量，它们执行实际的计算，它们释放内存（如果需要）并且它们返回一些内容给 Python（要么是 None，要么是 Python 对象）。我将只提一下上述矩阵 C 扩展 matsq 中代码的一些差异。
vecfcn1
解析函数的格式字符串指定来自 Python 的变量是一个字符串 (s)。
if (!PyArg_ParseTuple(args, "O!O!sd", &PyArray_Type, &vecin,
&PyArray_Type, &vecout, &str, &dfac)) return NULL;
在本地 C 数组的指针分配中没有分配内存，因为它们已经是向量。
cin=pyvector_to_Carrayptrs(vecin);
cout=pyvector_to_Carrayptrs(vecout);
返回值为 int = 1，如果成功。这将作为 Python int 返回。
return Py_BuildValue("i", 1);
rowx2
在这个例程中，我们利用参数元组列表中按引用传递的特性，通过对输出进行就地修改来“回传”输出。将此与在 matsq 中从 C 扩展返回数组进行比较。两种方法都能完成任务。
contigmat
这里，矩阵数据被视为一个长向量（就像将行首尾相连地堆叠起来一样）。如果您在 C++ 中有数组类，这些类将数据存储为一个长向量，然后使用步长来访问它，就像一个数组（二维、三维或其他），那么这将很有用。即使 matin 和 matout 是“矩阵”，我们也将其视为向量，并使用向量工具来获取我们的 C 指针 cin 和 cout。
/* Change contiguous arrays into C * arrays pointers to PyArrayObject data */
cin=pyvector_to_Carrayptrs(matin);
cout=pyvector_to_Carrayptrs(matout);
对于其他实用函数，请注意，我们使用不同的秩、维度和 NumPy 参数（例如 NPY_LONG）来告诉我们调用的例程数据类型是什么。
Makefile¶
Makefile 非常简单。它是为 Mac OS X 10.4（作为 BSD Unix）编写的。
# ---- Link ---------------------------
_C_arraytest.so: C_arraytest.o C_arraytest.mak gcc -bundle -flat_namespace -undefined suppress -o _C_arraytest.so C_arraytest.o
# ---- gcc C compile ------------------
C_arraytest.o: C_arraytest.c C_arraytest.h C_arraytest.mak gcc -c C_arraytest.c
-I/Library/Frameworks/Python.framework/Versions/2.4/include/python2.4
-I/Library/Frameworks/Python.framework/Versions/2.4/lib/python2.4/site-packages/numpy/core/include/numpy/
编译步骤非常标准。您需要添加 Python 头文件的路径
-I/Library/Frameworks/Python.framework/Versions/2.4/include/python2.4
以及 NumPy 头文件的路径
-I/Library/Frameworks/Python.framework/Versions/2.4/lib/python2.4/site-packages/numpy/core/include/numpy/
这些路径适用于 Mac OS X 上的 Framework Python 2.4 安装。您需要提供计算机上安装的头文件的路径。它们可能不同。我猜 gcc 标志对于编译来说是相同的。
链接步骤生成实际的模块（_C_arraytest.so），该模块可以导入到 Python 代码中。这是特定于 Mac OS X 系统的。在 Linux 或 Windows 上，您将需要其他东西。我一直在寻找通用的示例，但我不知道我找到的东西是否适合大多数人，所以我选择不显示那里的发现。我无法判断代码是否适合这些系统。
请注意，生成的共享库的名称必须与 C 扩展源代码中初始化和方法定义调用中的名称匹配。因此，_C_arraytest.so 名称中的前导下划线。
这是我在 Linux 下编译此代码的修改后的 Makefile（将其保存为与代码相同的目录中的 Makefile，然后运行 'make'）。
-- PaulIvanov
# ---- Link ---------------------------
_C_arraytest.so: C_arraytest.o
gcc -shared C_arraytest.o -o _C_arraytest.so
# ---- gcc C compile ------------------
C_arraytest.o: C_arraytest.c C_arraytest.h
gcc -c C_arraytest.c -I/usr/include/python2.4 -I/usr/lib/python2.4/site-packages/numpy/core/include/numpy
Python 包装代码¶
与 C 代码一样，我将只展示一个包装函数及其用法的详细描述。其他包装函数的代码非常相似，只要理解一个函数，其他函数就很容易理解。我将再次使用 matsq 函数。这是在您导入包装模块 (C_arraytest.py) 后，在您自己的 Python 代码中调用 matsq 函数时首先调用的代码，该模块会自动导入并使用 (对用户隐藏) _C_arraytest.so 中的实际 C 扩展 (注意开头的下划线，它使名称保持独立)。
导入
导入 C 扩展、NumPy 和系统模块 (用于末尾的 exit 语句，该语句是可选的)。
import _C_arraytest
import numpy as NP
import sys
Python matsq 函数的定义
传递一个 NumPy 矩阵 (matin)、一个 Python 整数 (ifac) 和一个 Python 浮点数 (dfac)。检查参数以确保它们具有正确的类型、维度和大小。这在 Python 端更容易更安全，这就是我在这里这样做的原因，即使我已经展示了在 C 扩展中完成其中一些操作的方法。
def matsq(matin, ifac, dfac):
# .... Check arguments, double NumPy matrices?
test=NP.zeros((2,2)) # create a NumPy matrix as a test object to check matin
typetest= type(test) # get its type
datatest=test.dtype # get data type
if type(matin) != typetest:
raise 'In matsq, matrix argument is not *NumPy* array'
if len(NP.shape(matin)) != 2:
raise 'In matsq, matrix argument is not NumPy *matrix*'
if matin.dtype != datatest:
raise 'In matsq, matrix argument is not *Float* NumPy matrix'
if type(ifac) != type(1):
raise 'In matsq, ifac argument is not an int'
if type(dfac) != type(1.0):
raise 'In matsq, dfac argument is not a python float'
最后，调用 C 扩展来对 matin 进行实际计算。
# .... Call C extension function
return _C_arraytest.matsq(matin, ifac, dfac)
您可以看到 Python 部分是最简单的。
使用您的 C 扩展
如果测试函数 mattest2 在另一个模块中 (您正在编写的模块)，那么您可以使用以下方法在脚本中调用包装的 matsq 函数。
from C_arraytest import * # Note, NO underscore since you are importing the wrapper.
import numpy as NP
import sys
def mattest2():
print "\n--- Test matsq ------------------------------"
print " Each 2nd matrix component should = square of 1st matrix component x ifac x dfac \n"
n=4 # Number of columns # Make 2 x n matrices
z=NP.arange(float(n))
x=NP.array([z,3.0*z])
jfac=2
xfac=1.5
y=matsq(x, jfac, xfac)
print "x=",x
print "y=",y
if __name__ == '__main__':
mattest2()
输出如下所示
--- Test matsq ------------------------------
Each 2nd matrix component should = square of 1st matrix component x ifac x dfac
x= [[ 0. 1. 2. 3.]
[ 0. 3. 6. 9.]]
y= [[ 0. 3. 12. 27.]
[ 0. 27. 108. 243.]]
所有测试函数的输出都在 C_arraytest 输出文件中。
总结¶
这是我编写的 C 扩展的第一个草稿 (在帮助下)。如果您对代码、错误等有任何意见，请在 pythonmac 邮件列表中发布。我会看到的。
我编写了 C 扩展来与 NumPy 一起使用，尽管它们最初是为 Numeric 编写的。如果您必须使用 Numeric，您应该测试它们以查看它们是否兼容。我怀疑像 NPY_DOUBLE 这样的名称将不兼容。我强烈建议您升级到 NumPy，因为它将成为 Python 中 Numeric 的未来。这值得付出努力。
评论？！¶
请注意，虽然此行在上面的头文件中，但它在 tar 包中的 .h 文件中缺失。
static PyObject *rowx2_v2(PyObject *self, PyObject *args);
-- PaulIvanov
输出文件名应为 _C_arraytest_d.pyd（调试版本）和 _C_arraytest.pyd（发布版本）。
-- Geoffrey Zhu
ptrvector() 分配了 n*sizeof(double)，但实际上应该分配指向 double 的指针；所以：n*sizeof(double *)
-- Peter Meerwald
在 vecsq() 中，第 86 行，由于您将创建一个一维向量
int i,j,n,m, dims[2];
应该是
int i,j,n,m, dims[1];
在 pyvector_to_Carrayptrs() 中，n 从未被使用。
-- FredSpiessens
章节作者：未知[33]，DavidLinke，TravisOliphant，未知[34]，未知[35]，未知[36]，未知[37]，未知[38]，未知[39]，未知[40]
附件