编写您自己的ufunc

I have the Power!

— He-Man

创建一个新的ufunc

在阅读本文之前,通过阅读/略读扩展和嵌入Python解释器的第1部分中的教程以及如何扩展NumPy,可以帮助您熟悉Python的C扩展基础知识。

umath模块是一个计算机生成的C模块,可以创建许多ufunc。它提供了许多如何创建通用函数的示例。使用ufunc机制创建自己的ufunc也不困难。假设您有一个函数,您想要在其输入上逐个元素地操作。通过创建一个新的ufunc,您将获得一个处理的函数

  • 广播
  • N维循环
  • 自动类型转换,内存使用量最少
  • 可选的输出数组

创建自己的ufunc并不困难。所需要的只是您想要支持的每种数据类型的1-d循环。每个1-d循环必须具有特定签名,并且只能使用固定大小数据类型的ufunc。下面给出了用于创建新的ufunc以处理内置数据类型的函数调用。使用不同的机制为用户定义的数据类型注册ufunc。

在接下来的几节中,我们提供了可以轻松修改的示例代码,以创建自己的ufunc。这些示例是logit函数的连续更完整或复杂版本,这是统计建模中的常见功能。Logit也很有趣,因为由于IEEE标准(特别是IEEE 754)的神奇之处,下面创建的所有logit函数都自动具有以下行为。

  1. >>> logit(0)
  2. -inf
  3. >>> logit(1)
  4. inf
  5. >>> logit(2)
  6. nan
  7. >>> logit(-2)
  8. nan

这很好,因为函数编写器不必手动传播infs或nans。

示例非ufunc扩展名

为了比较和阅读器的一般启发,我们提供了一个简单的logit C扩展实现,它没有使用numpy。

为此,我们需要两个文件。第一个是包含实际代码的C文件,第二个是用于创建模块的setup.py文件。

  1. #include <Python.h>
  2. #include <math.h>
  3. /*
  4. * spammodule.c
  5. * This is the C code for a non-numpy Python extension to
  6. * define the logit function, where logit(p) = log(p/(1-p)).
  7. * This function will not work on numpy arrays automatically.
  8. * numpy.vectorize must be called in python to generate
  9. * a numpy-friendly function.
  10. *
  11. * Details explaining the Python-C API can be found under
  12. * 'Extending and Embedding' and 'Python/C API' at
  13. * docs.python.org .
  14. */
  15. /* This declares the logit function */
  16. static PyObject* spam_logit(PyObject *self, PyObject *args);
  17. /*
  18. * This tells Python what methods this module has.
  19. * See the Python-C API for more information.
  20. */
  21. static PyMethodDef SpamMethods[] = {
  22. {"logit",
  23. spam_logit,
  24. METH_VARARGS, "compute logit"},
  25. {NULL, NULL, 0, NULL}
  26. };
  27. /*
  28. * This actually defines the logit function for
  29. * input args from Python.
  30. */
  31. static PyObject* spam_logit(PyObject *self, PyObject *args)
  32. {
  33. double p;
  34. /* This parses the Python argument into a double */
  35. if(!PyArg_ParseTuple(args, "d", &p)) {
  36. return NULL;
  37. }
  38. /* THE ACTUAL LOGIT FUNCTION */
  39. p = p/(1-p);
  40. p = log(p);
  41. /*This builds the answer back into a python object */
  42. return Py_BuildValue("d", p);
  43. }
  44. /* This initiates the module using the above definitions. */
  45. #if PY_VERSION_HEX >= 0x03000000
  46. static struct PyModuleDef moduledef = {
  47. PyModuleDef_HEAD_INIT,
  48. "spam",
  49. NULL,
  50. -1,
  51. SpamMethods,
  52. NULL,
  53. NULL,
  54. NULL,
  55. NULL
  56. };
  57. PyMODINIT_FUNC PyInit_spam(void)
  58. {
  59. PyObject *m;
  60. m = PyModule_Create(&moduledef);
  61. if (!m) {
  62. return NULL;
  63. }
  64. return m;
  65. }
  66. #else
  67. PyMODINIT_FUNC initspam(void)
  68. {
  69. PyObject *m;
  70. m = Py_InitModule("spam", SpamMethods);
  71. if (m == NULL) {
  72. return;
  73. }
  74. }
  75. #endif

要使用setup.py文件,请将setup.py和spammodule.c放在同一文件夹中。然后python setup.py build将构建要导入的模块,或者setup.py install将模块安装到您的site-packages目录。

  1. '''
  2. setup.py file for spammodule.c
  3. Calling
  4. $python setup.py build_ext --inplace
  5. will build the extension library in the current file.
  6. Calling
  7. $python setup.py build
  8. will build a file that looks like ./build/lib*, where
  9. lib* is a file that begins with lib. The library will
  10. be in this file and end with a C library extension,
  11. such as .so
  12. Calling
  13. $python setup.py install
  14. will install the module in your site-packages file.
  15. See the distutils section of
  16. 'Extending and Embedding the Python Interpreter'
  17. at docs.python.org for more information.
  18. '''
  19. from distutils.core import setup, Extension
  20. module1 = Extension('spam', sources=['spammodule.c'],
  21. include_dirs=['/usr/local/lib'])
  22. setup(name = 'spam',
  23. version='1.0',
  24. description='This is my spam package',
  25. ext_modules = [module1])

将垃圾邮件模块导入python后,您可以通过spam.logit调用logit。请注意,上面使用的函数不能按原样应用于numpy数组。为此,我们必须在其上调用numpy.vectorize。例如,如果在包含垃圾邮件库或垃圾邮件的文件中打开了python解释器,则可以执行以下命令:

  1. >>> import numpy as np
  2. >>> import spam
  3. >>> spam.logit(0)
  4. -inf
  5. >>> spam.logit(1)
  6. inf
  7. >>> spam.logit(0.5)
  8. 0.0
  9. >>> x = np.linspace(0,1,10)
  10. >>> spam.logit(x)
  11. TypeError: only length-1 arrays can be converted to Python scalars
  12. >>> f = np.vectorize(spam.logit)
  13. >>> f(x)
  14. array([ -inf, -2.07944154, -1.25276297, -0.69314718, -0.22314355,
  15. 0.22314355, 0.69314718, 1.25276297, 2.07944154, inf])

结果编辑功能并不快!numpy.vectorize只是循环遍历spam.logit。循环在C级完成,但numpy数组不断被解析并重新构建。这很贵。当作者将numpy.vectorize(spam.logit)与下面构造的logit ufuncs进行比较时,logit ufuncs几乎快4倍。当然,取决于功能的性质,可以实现更大或更小的加速。

一种dtype的NumPy ufunc示例

为简单起见,我们为单个dtype提供了一个ufunc,即’f8’双精度型。与前一节一样,我们首先给出.c文件,然后是用于创建包含ufunc的模块的setup.py文件。

代码中与ufunc的实际计算相对应的位置标有/ BEGIN main ufunc computation /和/ END main ufunc computation /。这些行之间的代码是必须更改以创建自己的ufunc的主要事物。

  1. #include "Python.h"
  2. #include "math.h"
  3. #include "numpy/ndarraytypes.h"
  4. #include "numpy/ufuncobject.h"
  5. #include "numpy/npy_3kcompat.h"
  6. /*
  7. * single_type_logit.c
  8. * This is the C code for creating your own
  9. * NumPy ufunc for a logit function.
  10. *
  11. * In this code we only define the ufunc for
  12. * a single dtype. The computations that must
  13. * be replaced to create a ufunc for
  14. * a different function are marked with BEGIN
  15. * and END.
  16. *
  17. * Details explaining the Python-C API can be found under
  18. * 'Extending and Embedding' and 'Python/C API' at
  19. * docs.python.org .
  20. */
  21. static PyMethodDef LogitMethods[] = {
  22. {NULL, NULL, 0, NULL}
  23. };
  24. /* The loop definition must precede the PyMODINIT_FUNC. */
  25. static void double_logit(char **args, npy_intp *dimensions,
  26. npy_intp* steps, void* data)
  27. {
  28. npy_intp i;
  29. npy_intp n = dimensions[0];
  30. char *in = args[0], *out = args[1];
  31. npy_intp in_step = steps[0], out_step = steps[1];
  32. double tmp;
  33. for (i = 0; i < n; i++) {
  34. /*BEGIN main ufunc computation*/
  35. tmp = *(double *)in;
  36. tmp /= 1-tmp;
  37. *((double *)out) = log(tmp);
  38. /*END main ufunc computation*/
  39. in += in_step;
  40. out += out_step;
  41. }
  42. }
  43. /*This a pointer to the above function*/
  44. PyUFuncGenericFunction funcs[1] = {&double_logit};
  45. /* These are the input and return dtypes of logit.*/
  46. static char types[2] = {NPY_DOUBLE, NPY_DOUBLE};
  47. static void *data[1] = {NULL};
  48. #if PY_VERSION_HEX >= 0x03000000
  49. static struct PyModuleDef moduledef = {
  50. PyModuleDef_HEAD_INIT,
  51. "npufunc",
  52. NULL,
  53. -1,
  54. LogitMethods,
  55. NULL,
  56. NULL,
  57. NULL,
  58. NULL
  59. };
  60. PyMODINIT_FUNC PyInit_npufunc(void)
  61. {
  62. PyObject *m, *logit, *d;
  63. m = PyModule_Create(&moduledef);
  64. if (!m) {
  65. return NULL;
  66. }
  67. import_array();
  68. import_umath();
  69. logit = PyUFunc_FromFuncAndData(funcs, data, types, 1, 1, 1,
  70. PyUFunc_None, "logit",
  71. "logit_docstring", 0);
  72. d = PyModule_GetDict(m);
  73. PyDict_SetItemString(d, "logit", logit);
  74. Py_DECREF(logit);
  75. return m;
  76. }
  77. #else
  78. PyMODINIT_FUNC initnpufunc(void)
  79. {
  80. PyObject *m, *logit, *d;
  81. m = Py_InitModule("npufunc", LogitMethods);
  82. if (m == NULL) {
  83. return;
  84. }
  85. import_array();
  86. import_umath();
  87. logit = PyUFunc_FromFuncAndData(funcs, data, types, 1, 1, 1,
  88. PyUFunc_None, "logit",
  89. "logit_docstring", 0);
  90. d = PyModule_GetDict(m);
  91. PyDict_SetItemString(d, "logit", logit);
  92. Py_DECREF(logit);
  93. }
  94. #endif

这是上面代码的setup.py文件。和以前一样,可以通过在命令提示符下调用python setup.py build来构建模块,也可以通过python setup.py install将其安装到site-packages。

  1. '''
  2. setup.py file for logit.c
  3. Note that since this is a numpy extension
  4. we use numpy.distutils instead of
  5. distutils from the python standard library.
  6. Calling
  7. $python setup.py build_ext --inplace
  8. will build the extension library in the current file.
  9. Calling
  10. $python setup.py build
  11. will build a file that looks like ./build/lib*, where
  12. lib* is a file that begins with lib. The library will
  13. be in this file and end with a C library extension,
  14. such as .so
  15. Calling
  16. $python setup.py install
  17. will install the module in your site-packages file.
  18. See the distutils section of
  19. 'Extending and Embedding the Python Interpreter'
  20. at docs.python.org and the documentation
  21. on numpy.distutils for more information.
  22. '''
  23. def configuration(parent_package='', top_path=None):
  24. import numpy
  25. from numpy.distutils.misc_util import Configuration
  26. config = Configuration('npufunc_directory',
  27. parent_package,
  28. top_path)
  29. config.add_extension('npufunc', ['single_type_logit.c'])
  30. return config
  31. if __name__ == "__main__":
  32. from numpy.distutils.core import setup
  33. setup(configuration=configuration)

安装完上述内容后,可以按如下方式导入和使用。

  1. >>> import numpy as np
  2. >>> import npufunc
  3. >>> npufunc.logit(0.5)
  4. 0.0
  5. >>> a = np.linspace(0,1,5)
  6. >>> npufunc.logit(a)
  7. array([ -inf, -1.09861229, 0. , 1.09861229, inf])

示例具有多个dtypes的NumPy ufunc

我们最后给出了一个完整的ufunc示例,内部循环用于半浮点数,浮点数,双精度数和长双精度数。与前面的部分一样,我们首先给出.c文件,然后是相应的setup.py文件。

代码中与ufunc的实际计算相对应的位置标有/ BEGIN main ufunc computation /和/ END main ufunc computation /。这些行之间的代码是必须更改以创建自己的ufunc的主要事物。

  1. #include "Python.h"
  2. #include "math.h"
  3. #include "numpy/ndarraytypes.h"
  4. #include "numpy/ufuncobject.h"
  5. #include "numpy/halffloat.h"
  6. /*
  7. * multi_type_logit.c
  8. * This is the C code for creating your own
  9. * NumPy ufunc for a logit function.
  10. *
  11. * Each function of the form type_logit defines the
  12. * logit function for a different numpy dtype. Each
  13. * of these functions must be modified when you
  14. * create your own ufunc. The computations that must
  15. * be replaced to create a ufunc for
  16. * a different function are marked with BEGIN
  17. * and END.
  18. *
  19. * Details explaining the Python-C API can be found under
  20. * 'Extending and Embedding' and 'Python/C API' at
  21. * docs.python.org .
  22. *
  23. */
  24. static PyMethodDef LogitMethods[] = {
  25. {NULL, NULL, 0, NULL}
  26. };
  27. /* The loop definitions must precede the PyMODINIT_FUNC. */
  28. static void long_double_logit(char **args, npy_intp *dimensions,
  29. npy_intp* steps, void* data)
  30. {
  31. npy_intp i;
  32. npy_intp n = dimensions[0];
  33. char *in = args[0], *out=args[1];
  34. npy_intp in_step = steps[0], out_step = steps[1];
  35. long double tmp;
  36. for (i = 0; i < n; i++) {
  37. /*BEGIN main ufunc computation*/
  38. tmp = *(long double *)in;
  39. tmp /= 1-tmp;
  40. *((long double *)out) = logl(tmp);
  41. /*END main ufunc computation*/
  42. in += in_step;
  43. out += out_step;
  44. }
  45. }
  46. static void double_logit(char **args, npy_intp *dimensions,
  47. npy_intp* steps, void* data)
  48. {
  49. npy_intp i;
  50. npy_intp n = dimensions[0];
  51. char *in = args[0], *out = args[1];
  52. npy_intp in_step = steps[0], out_step = steps[1];
  53. double tmp;
  54. for (i = 0; i < n; i++) {
  55. /*BEGIN main ufunc computation*/
  56. tmp = *(double *)in;
  57. tmp /= 1-tmp;
  58. *((double *)out) = log(tmp);
  59. /*END main ufunc computation*/
  60. in += in_step;
  61. out += out_step;
  62. }
  63. }
  64. static void float_logit(char **args, npy_intp *dimensions,
  65. npy_intp* steps, void* data)
  66. {
  67. npy_intp i;
  68. npy_intp n = dimensions[0];
  69. char *in=args[0], *out = args[1];
  70. npy_intp in_step = steps[0], out_step = steps[1];
  71. float tmp;
  72. for (i = 0; i < n; i++) {
  73. /*BEGIN main ufunc computation*/
  74. tmp = *(float *)in;
  75. tmp /= 1-tmp;
  76. *((float *)out) = logf(tmp);
  77. /*END main ufunc computation*/
  78. in += in_step;
  79. out += out_step;
  80. }
  81. }
  82. static void half_float_logit(char **args, npy_intp *dimensions,
  83. npy_intp* steps, void* data)
  84. {
  85. npy_intp i;
  86. npy_intp n = dimensions[0];
  87. char *in = args[0], *out = args[1];
  88. npy_intp in_step = steps[0], out_step = steps[1];
  89. float tmp;
  90. for (i = 0; i < n; i++) {
  91. /*BEGIN main ufunc computation*/
  92. tmp = *(npy_half *)in;
  93. tmp = npy_half_to_float(tmp);
  94. tmp /= 1-tmp;
  95. tmp = logf(tmp);
  96. *((npy_half *)out) = npy_float_to_half(tmp);
  97. /*END main ufunc computation*/
  98. in += in_step;
  99. out += out_step;
  100. }
  101. }
  102. /*This gives pointers to the above functions*/
  103. PyUFuncGenericFunction funcs[4] = {&half_float_logit,
  104. &float_logit,
  105. &double_logit,
  106. &long_double_logit};
  107. static char types[8] = {NPY_HALF, NPY_HALF,
  108. NPY_FLOAT, NPY_FLOAT,
  109. NPY_DOUBLE,NPY_DOUBLE,
  110. NPY_LONGDOUBLE, NPY_LONGDOUBLE};
  111. static void *data[4] = {NULL, NULL, NULL, NULL};
  112. #if PY_VERSION_HEX >= 0x03000000
  113. static struct PyModuleDef moduledef = {
  114. PyModuleDef_HEAD_INIT,
  115. "npufunc",
  116. NULL,
  117. -1,
  118. LogitMethods,
  119. NULL,
  120. NULL,
  121. NULL,
  122. NULL
  123. };
  124. PyMODINIT_FUNC PyInit_npufunc(void)
  125. {
  126. PyObject *m, *logit, *d;
  127. m = PyModule_Create(&moduledef);
  128. if (!m) {
  129. return NULL;
  130. }
  131. import_array();
  132. import_umath();
  133. logit = PyUFunc_FromFuncAndData(funcs, data, types, 4, 1, 1,
  134. PyUFunc_None, "logit",
  135. "logit_docstring", 0);
  136. d = PyModule_GetDict(m);
  137. PyDict_SetItemString(d, "logit", logit);
  138. Py_DECREF(logit);
  139. return m;
  140. }
  141. #else
  142. PyMODINIT_FUNC initnpufunc(void)
  143. {
  144. PyObject *m, *logit, *d;
  145. m = Py_InitModule("npufunc", LogitMethods);
  146. if (m == NULL) {
  147. return;
  148. }
  149. import_array();
  150. import_umath();
  151. logit = PyUFunc_FromFuncAndData(funcs, data, types, 4, 1, 1,
  152. PyUFunc_None, "logit",
  153. "logit_docstring", 0);
  154. d = PyModule_GetDict(m);
  155. PyDict_SetItemString(d, "logit", logit);
  156. Py_DECREF(logit);
  157. }
  158. #endif

这是上面代码的setup.py文件。和以前一样,可以通过在命令提示符下调用python setup.py build来构建模块,也可以通过python setup.py install将其安装到site-packages。

  1. '''
  2. setup.py file for logit.c
  3. Note that since this is a numpy extension
  4. we use numpy.distutils instead of
  5. distutils from the python standard library.
  6. Calling
  7. $python setup.py build_ext --inplace
  8. will build the extension library in the current file.
  9. Calling
  10. $python setup.py build
  11. will build a file that looks like ./build/lib*, where
  12. lib* is a file that begins with lib. The library will
  13. be in this file and end with a C library extension,
  14. such as .so
  15. Calling
  16. $python setup.py install
  17. will install the module in your site-packages file.
  18. See the distutils section of
  19. 'Extending and Embedding the Python Interpreter'
  20. at docs.python.org and the documentation
  21. on numpy.distutils for more information.
  22. '''
  23. def configuration(parent_package='', top_path=None):
  24. import numpy
  25. from numpy.distutils.misc_util import Configuration
  26. from numpy.distutils.misc_util import get_info
  27. #Necessary for the half-float d-type.
  28. info = get_info('npymath')
  29. config = Configuration('npufunc_directory',
  30. parent_package,
  31. top_path)
  32. config.add_extension('npufunc',
  33. ['multi_type_logit.c'],
  34. extra_info=info)
  35. return config
  36. if __name__ == "__main__":
  37. from numpy.distutils.core import setup
  38. setup(configuration=configuration)

安装完上述内容后,可以按如下方式导入和使用。

  1. >>> import numpy as np
  2. >>> import npufunc
  3. >>> npufunc.logit(0.5)
  4. 0.0
  5. >>> a = np.linspace(0,1,5)
  6. >>> npufunc.logit(a)
  7. array([ -inf, -1.09861229, 0. , 1.09861229, inf])

示例具有多个参数/返回值的NumPy ufunc

我们的最后一个例子是一个带有多个参数的ufunc。它是对具有单个dtype的数据的logit ufunc的代码的修改。我们计算 (A*B, logit(A*B))。

我们只给出 C 代码,因为setup.py文件与一种dtype的NumPy ufunc示例中的setup.py文件完全相同,只是一行

  1. config.add_extension('npufunc', ['single_type_logit.c'])

被替换为

  1. config.add_extension('npufunc', ['multi_arg_logit.c'])

C文件如下。生成的ufunc接受两个参数A和B.它返回一个元组,其第一个元素是A B,第二个元素是logit(A B)。请注意,它会自动支持广播以及ufunc的所有其他属性。

  1. #include "Python.h"
  2. #include "math.h"
  3. #include "numpy/ndarraytypes.h"
  4. #include "numpy/ufuncobject.h"
  5. #include "numpy/halffloat.h"
  6. /*
  7. * multi_arg_logit.c
  8. * This is the C code for creating your own
  9. * NumPy ufunc for a multiple argument, multiple
  10. * return value ufunc. The places where the
  11. * ufunc computation is carried out are marked
  12. * with comments.
  13. *
  14. * Details explaining the Python-C API can be found under
  15. * 'Extending and Embedding' and 'Python/C API' at
  16. * docs.python.org .
  17. *
  18. */
  19. static PyMethodDef LogitMethods[] = {
  20. {NULL, NULL, 0, NULL}
  21. };
  22. /* The loop definition must precede the PyMODINIT_FUNC. */
  23. static void double_logitprod(char **args, npy_intp *dimensions,
  24. npy_intp* steps, void* data)
  25. {
  26. npy_intp i;
  27. npy_intp n = dimensions[0];
  28. char *in1 = args[0], *in2 = args[1];
  29. char *out1 = args[2], *out2 = args[3];
  30. npy_intp in1_step = steps[0], in2_step = steps[1];
  31. npy_intp out1_step = steps[2], out2_step = steps[3];
  32. double tmp;
  33. for (i = 0; i < n; i++) {
  34. /*BEGIN main ufunc computation*/
  35. tmp = *(double *)in1;
  36. tmp *= *(double *)in2;
  37. *((double *)out1) = tmp;
  38. *((double *)out2) = log(tmp/(1-tmp));
  39. /*END main ufunc computation*/
  40. in1 += in1_step;
  41. in2 += in2_step;
  42. out1 += out1_step;
  43. out2 += out2_step;
  44. }
  45. }
  46. /*This a pointer to the above function*/
  47. PyUFuncGenericFunction funcs[1] = {&double_logitprod};
  48. /* These are the input and return dtypes of logit.*/
  49. static char types[4] = {NPY_DOUBLE, NPY_DOUBLE,
  50. NPY_DOUBLE, NPY_DOUBLE};
  51. static void *data[1] = {NULL};
  52. #if PY_VERSION_HEX >= 0x03000000
  53. static struct PyModuleDef moduledef = {
  54. PyModuleDef_HEAD_INIT,
  55. "npufunc",
  56. NULL,
  57. -1,
  58. LogitMethods,
  59. NULL,
  60. NULL,
  61. NULL,
  62. NULL
  63. };
  64. PyMODINIT_FUNC PyInit_npufunc(void)
  65. {
  66. PyObject *m, *logit, *d;
  67. m = PyModule_Create(&moduledef);
  68. if (!m) {
  69. return NULL;
  70. }
  71. import_array();
  72. import_umath();
  73. logit = PyUFunc_FromFuncAndData(funcs, data, types, 1, 2, 2,
  74. PyUFunc_None, "logit",
  75. "logit_docstring", 0);
  76. d = PyModule_GetDict(m);
  77. PyDict_SetItemString(d, "logit", logit);
  78. Py_DECREF(logit);
  79. return m;
  80. }
  81. #else
  82. PyMODINIT_FUNC initnpufunc(void)
  83. {
  84. PyObject *m, *logit, *d;
  85. m = Py_InitModule("npufunc", LogitMethods);
  86. if (m == NULL) {
  87. return;
  88. }
  89. import_array();
  90. import_umath();
  91. logit = PyUFunc_FromFuncAndData(funcs, data, types, 1, 2, 2,
  92. PyUFunc_None, "logit",
  93. "logit_docstring", 0);
  94. d = PyModule_GetDict(m);
  95. PyDict_SetItemString(d, "logit", logit);
  96. Py_DECREF(logit);
  97. }
  98. #endif

示例带有结构化数组dtype参数的NumPy ufunc

此示例显示如何为结构化数组dtype创建ufunc。在这个例子中,我们展示了一个简单的ufunc,用于添加两个带有dtype’u8,u8,u8’的数组。该过程与其他示例略有不同,因为调用PyUFunc_FromFuncAndData不会为自定义dtypes和结构化数组dtypes完全注册ufunc。我们还需要调用 PyUFunc_RegisterLoopForDescr完成设置ufunc。

我们只提供C代码,因为setup.py文件与一种dtype的NumPy ufunc示例的setup.py文件完全相同,只有一行。

  1. config.add_extension('npufunc', ['single_type_logit.c'])

被替换为

  1. config.add_extension('npufunc', ['add_triplet.c'])

C文件如下。

  1. #include "Python.h"
  2. #include "math.h"
  3. #include "numpy/ndarraytypes.h"
  4. #include "numpy/ufuncobject.h"
  5. #include "numpy/npy_3kcompat.h"
  6. /*
  7. * add_triplet.c
  8. * This is the C code for creating your own
  9. * NumPy ufunc for a structured array dtype.
  10. *
  11. * Details explaining the Python-C API can be found under
  12. * 'Extending and Embedding' and 'Python/C API' at
  13. * docs.python.org .
  14. */
  15. static PyMethodDef StructUfuncTestMethods[] = {
  16. {NULL, NULL, 0, NULL}
  17. };
  18. /* The loop definition must precede the PyMODINIT_FUNC. */
  19. static void add_uint64_triplet(char **args, npy_intp *dimensions,
  20. npy_intp* steps, void* data)
  21. {
  22. npy_intp i;
  23. npy_intp is1=steps[0];
  24. npy_intp is2=steps[1];
  25. npy_intp os=steps[2];
  26. npy_intp n=dimensions[0];
  27. uint64_t *x, *y, *z;
  28. char *i1=args[0];
  29. char *i2=args[1];
  30. char *op=args[2];
  31. for (i = 0; i < n; i++) {
  32. x = (uint64_t*)i1;
  33. y = (uint64_t*)i2;
  34. z = (uint64_t*)op;
  35. z[0] = x[0] + y[0];
  36. z[1] = x[1] + y[1];
  37. z[2] = x[2] + y[2];
  38. i1 += is1;
  39. i2 += is2;
  40. op += os;
  41. }
  42. }
  43. /* This a pointer to the above function */
  44. PyUFuncGenericFunction funcs[1] = {&add_uint64_triplet};
  45. /* These are the input and return dtypes of add_uint64_triplet. */
  46. static char types[3] = {NPY_UINT64, NPY_UINT64, NPY_UINT64};
  47. static void *data[1] = {NULL};
  48. #if defined(NPY_PY3K)
  49. static struct PyModuleDef moduledef = {
  50. PyModuleDef_HEAD_INIT,
  51. "struct_ufunc_test",
  52. NULL,
  53. -1,
  54. StructUfuncTestMethods,
  55. NULL,
  56. NULL,
  57. NULL,
  58. NULL
  59. };
  60. #endif
  61. #if defined(NPY_PY3K)
  62. PyMODINIT_FUNC PyInit_struct_ufunc_test(void)
  63. #else
  64. PyMODINIT_FUNC initstruct_ufunc_test(void)
  65. #endif
  66. {
  67. PyObject *m, *add_triplet, *d;
  68. PyObject *dtype_dict;
  69. PyArray_Descr *dtype;
  70. PyArray_Descr *dtypes[3];
  71. #if defined(NPY_PY3K)
  72. m = PyModule_Create(&moduledef);
  73. #else
  74. m = Py_InitModule("struct_ufunc_test", StructUfuncTestMethods);
  75. #endif
  76. if (m == NULL) {
  77. #if defined(NPY_PY3K)
  78. return NULL;
  79. #else
  80. return;
  81. #endif
  82. }
  83. import_array();
  84. import_umath();
  85. /* Create a new ufunc object */
  86. add_triplet = PyUFunc_FromFuncAndData(NULL, NULL, NULL, 0, 2, 1,
  87. PyUFunc_None, "add_triplet",
  88. "add_triplet_docstring", 0);
  89. dtype_dict = Py_BuildValue("[(s, s), (s, s), (s, s)]",
  90. "f0", "u8", "f1", "u8", "f2", "u8");
  91. PyArray_DescrConverter(dtype_dict, &dtype);
  92. Py_DECREF(dtype_dict);
  93. dtypes[0] = dtype;
  94. dtypes[1] = dtype;
  95. dtypes[2] = dtype;
  96. /* Register ufunc for structured dtype */
  97. PyUFunc_RegisterLoopForDescr(add_triplet,
  98. dtype,
  99. &add_uint64_triplet,
  100. dtypes,
  101. NULL);
  102. d = PyModule_GetDict(m);
  103. PyDict_SetItemString(d, "add_triplet", add_triplet);
  104. Py_DECREF(add_triplet);
  105. #if defined(NPY_PY3K)
  106. return m;
  107. #endif
  108. }

返回的ufunc对象是一个可调用的Python对象。它应该放在一个(模块)字典中,其名称与ufunc-creation例程的name参数中使用的字典相同。以下示例是从umath模块改编而来的

  1. static PyUFuncGenericFunction atan2_functions[] = {
  2. PyUFunc_ff_f, PyUFunc_dd_d,
  3. PyUFunc_gg_g, PyUFunc_OO_O_method};
  4. static void* atan2_data[] = {
  5. (void *)atan2f,(void *) atan2,
  6. (void *)atan2l,(void *)"arctan2"};
  7. static char atan2_signatures[] = {
  8. NPY_FLOAT, NPY_FLOAT, NPY_FLOAT,
  9. NPY_DOUBLE, NPY_DOUBLE, NPY_DOUBLE,
  10. NPY_LONGDOUBLE, NPY_LONGDOUBLE, NPY_LONGDOUBLE
  11. NPY_OBJECT, NPY_OBJECT, NPY_OBJECT};
  12. ...
  13. /* in the module initialization code */
  14. PyObject *f, *dict, *module;
  15. ...
  16. dict = PyModule_GetDict(module);
  17. ...
  18. f = PyUFunc_FromFuncAndData(atan2_functions,
  19. atan2_data, atan2_signatures, 4, 2, 1,
  20. PyUFunc_None, "arctan2",
  21. "a safe and correct arctan(x1/x2)", 0);
  22. PyDict_SetItemString(dict, "arctan2", f);
  23. Py_DECREF(f);
  24. ...