如何直接调用ctypes SciPy被积函数？

Question

我有一些数字代码要在SciPy中运行。它涉及到多项式指数的复杂有理函数，所以它的计算成本相当高；因此，我用C编写了它，并用ctypes调用它。最重要的用例是作为scipy.integrate.quad公司，但我偶尔也需要直接打电话。你知道吗

我函数的“自然”签名是

double f(double x, double y, double z){}

ctypes documentation建议与相应的Python一起使用

import ctypes
so = ctypes.CDLL('/path/to/f.so')
f = so.f
f.restype = ctypes.c_double
f.argtypes = (ctypes.c_double, ctypes.c_double, ctypes.c_double)

无论如何，把它称为被积函数，SciPy requires a certain function signature

double f(int n, double args[n]){}

在Python代码中指定为

import ctypes
so = ctypes.CDLL('/path/to/f.so')
f = so.f
f.restype = ctypes.c_double
f.argtypes = (ctypes.c_int, ctypes.c_double)

要在执行积分时传递参数y和z，它们将作为名为args的元组传递给quad。你知道吗

scipy.integrate.quad(f, lowerlimit, upperlimit, args=(y,z))

这就不清楚如何直接调用f。我天真的尝试是

f(3, (x,y,z))

但这会导致参数2的类型错误。有几个变种也同样失败了。这并不奇怪；ctypes期望函数调用只有一个整数参数，后跟一个双参数。你知道吗

我完全不知道quad如何y和z进入f。我试着查看SciPy源代码，但我必须承认，我在尝试将Python到C到Fortran的调用跟踪时迷失了方向。你知道吗

我可以只编写另一个函数作为直接调用或集成的包装器，但只使用一种形式似乎更为优雅，至少，我想了解SciPy调用是如何工作的。你知道吗

如何直接调用f的SciPy被积函数形式，传递所有三个参数x、y和z？你知道吗

我使用的是python3.4.3、numpy1.11.2和scipy0.18.1。你知道吗

---

编辑：注意f可以通过更改其argtypes来调用：

f.argtypes = (ctypes.c_int, 3*ctypes.c_double)
f(3, (3*ctypes.c_double)(x, y, z))

不过，我还是很好奇西皮在做什么。一直来回设置argtypes充其量也是不雅观和不方便的。你知道吗

---

编辑2:请注意，经过上一次编辑之后，这个问题现在基本上是this one的副本，它在右侧列中弹出了一个有用的窗口。你知道吗

Answer 1

我不知道这是否有助于ctypes的情况，但是当调用Python integrad时，这些scipy函数将自由变量与参数连接起来。换句话说

def bar(x,y,z):
   return np.sin(x*y*z)
In [43]: quad(bar,0,1, args=(1,2))
Out[43]: (0.7080734182735712, 7.861194120923578e-15)

当以0.5计算时，它会(x,)+args：

In [49]: bar(*(.5,)+(1,2))
Out[49]: 0.8414709848078965

In [50]: bar(.5,1,2)
Out[50]: 0.8414709848078965

所以用c签名：

f_for_scipy(int n, double args[n])

n是参数数，args[n]是表示值数组的指针。你知道吗

我正在使用cython及其扩展类型来试验这类调用

https://cython.readthedocs.io/en/latest/src/tutorial/cdef_classes.html

Passing a cython function vs a cython method to scipy.integrate

Answer 2

您需要有正确的函数签名，但是Python argtypes类型应该是POINTER(c_double)。C中的数组衰减到函数参数中的指针：

C示例（Windows）

#include <stdio.h>
__declspec(dllexport) double f(int n, double args[])
{
    double sum = 0;
    int nn;
    for(nn = 0; nn < n; ++nn)
    {
        printf("args[%d] = %ld\n",nn,args[nn]);
        sum += args[nn];
    }
    return sum;
}

c类型示例

>>> from ctypes import *
>>> dll = CDLL('x')
>>> dll.f.restype = c_double
>>> dll.f.argtypes = c_int,POINTER(c_double)
>>> L=[1.1,2.2,3.3,4.4,5.5]
>>> args = (c_double*len(L))(*L)
>>> dll.f(len(args),args)
args[0] = 1.100000
args[1] = 2.200000
args[2] = 3.300000
args[3] = 4.400000
args[4] = 5.500000
16.5

Answer 3

您无法在Python级别修复此问题。文档https://docs.scipy.org/doc/scipy/reference/tutorial/integrate.html#faster-integration-using-ctypes说明

Write an integrand function in C with the function signature double f(int n, double args[n])

在您的情况下，您必须在C级别添加一个函数

double f_for_scipy(int n, double args[n]) {
    return f(args[0], args[1], args[2]);
}

把它喂给四轮车。你知道吗