使用回调将 C 库 (GSL) 包装在 cython 代码中

Posted 2023-02-23

【中文标题】使用回调将 C 库 (GSL) 包装在 cython 代码中

【英文标题】：wrapping a C library (GSL) in a cython code by using callback

【发布时间】：2014-07-01 14:49:28

【问题描述】：

我是cython 和c 的新手。我想使用 cython 来加快我的代码的性能。我想在我的代码中使用gsl_integration 库进行集成。 更新：test_gsl.pyx

cdef extern from "math.h":
    double log(double x) nogil
cdef extern from "gsl/gsl_math.h":
    ctypedef struct gsl_function:
        double (* function) (double x, void * params)
        void * params

cdef extern from "gsl/gsl_integration.h":
    ctypedef struct gsl_integration_workspace
    gsl_integration_workspace *  gsl_integration_workspace_alloc(size_t n)
    void  gsl_integration_workspace_free(gsl_integration_workspace * w)

    int  gsl_integration_qags(const gsl_function * f, double a, double b, double epsabs, double epsrel, size_t limit, gsl_integration_workspace * workspace, double *result, double *abserr)


cdef double do_callback(double x, void* params):
    return (<MyCallback>params).eval(x)

cdef class MyCallback:
    cdef double a
    def __init__(self, a):
        self.a = a
    cpdef double eval(self, double x):
        return self.a * log(x+1) * x
    def call_gsl(self):
        cdef gsl_integration_workspace* w =gsl_integration_workspace_alloc (1000)

        cdef gsl_function F
        F.function = &do_callback
        F.params = <void*>self

        cdef double result = 3, error = 5
        gsl_integration_qags(&F, 0, 1, 0, 1e-7, 1000, w, &result, &error)
        print result, error
        gsl_integration_workspace_free(w)

此.pyx 代码使用以下setup.py 文件编译，不会引发任何错误消息：

from distutils.core import setup
from distutils.extension import Extension
from Cython.Distutils import build_ext
import numpy
import sys 
ext = Extension("test_gsl", ["test_gsl.pyx"],
    include_dirs=[numpy.get_include(),
                  "/usr/include/"],
    library_dirs=["/usr/lib/"],
    libraries=["gsl"])

setup(ext_modules=[ext],
    cmdclass = 'build_ext': build_ext)

甚至使用命令行：

cython test_gsl.pyx
gcc -pthread -fno-strict-aliasing -g -O2 -DNDEBUG -g -fwrapv -O3 -Wall -Wstrict-prototypes -fPIC -I/vol/anaconda/lib/python2.7/site-packages -I/usr/include -I/vol/anaconda/include/python2.7 -c test_gsl.c `pkg-config --cflags gsl` 
gcc -pthread -shared test_gsl.o -L/usr/lib -L/vol/anaconda/lib -lgsl -lgslcblas -lpython2.7 `pkg-config --libs gsl` -o test_gsl.so

当它在 python 中如下导入时，它确实会引发错误：

>>> import pyximport; pyximport.install()
(None, <pyximport.pyximport.PyxImporter object at 0x7f0c7e888150>)
>>> import gsl
Traceback (most recent call last):
  File "<stdin>", line 1, in <module>
  File "/anaconda/lib/python2.7/site-packages/Cython-0.20.1-py2.7-linux-x86_64.egg/pyximport/pyximport.py", line 431, in load_module
    language_level=self.language_level)
  File "/anaconda/lib/python2.7/site-packages/Cython-0.20.1-py2.7-linux-x86_64.egg/pyximport/pyximport.py", line 210, in load_module
    mod = imp.load_dynamic(name, so_path)
ImportError: Building module gsl failed: ['ImportError: /users/dalek/.pyxbld/lib.linux-x86_64-2.7/gsl.so: undefined symbol: gsl_integration_qags\n']

gsl_integration_qags已经被正确定义了，我不明白为什么我又得到了这个错误？

【问题讨论】：

我发现值得我花时间观看this video lecture。我就是这样开始的。我在您的代码中看不到任何无法转换的内容（乍一看）。首先尝试将 .h 类和结构转换为 cython 类。一些指针可以找到on the cython website。但是，如果您先看演示文稿，这一切都会更有意义！

【参考方案1】：

第一法则：过早的优化是万恶之源。

第二条规则：不惜一切代价遵循第一条规则。

第三条规则：不要使用 C++ 复杂特性（与 C 相比复杂——这包括类），如果不需要的话（即使你像我一样是 C++ ***者）。如果您将 C++ 与 C 库混合使用，则尤其如此。

我看不出在您的示例中为什么需要 C++ 类的任何原因，尤其是因为您这样做创建了不必要的间接（包装器）！如果您使用编译语言进行编码以提高性能，那么您确实希望避免不必要的步骤和间接操作！你无缘无故地让你的生活变得困难，特别是因为 C 语言中的 GSL 例程将在你的程序中完成 99.9% 的计算。为什么不直接使用 cython-gsl 之类的东西并将您的代码恢复为类似的东西（取自 cython-gsl 示例文件夹）。这要短得多，更干净，鉴于python没有做任何繁重的工作，我看不出它表现不佳的原因（假设函数 foo() 将被转换为 C 似乎是这种情况）！

from cython_gsl cimport *

ctypedef double * double_ptr ctypedef void * void_ptr

cdef double foo(double x, void * params) nogil:
    cdef double alpha, f
    alpha = (<double_ptr> params)[0]
    f = log(alpha*x) / sqrt(x)
    return f


def main():
    cdef gsl_integration_workspace * w
    cdef double result, error, expected, alpha
    w = gsl_integration_workspace_alloc (1000)

    expected = -4.0
    alpha = 1

    cdef gsl_function F
    F.function = &foo
    F.params = &alpha

    gsl_integration_qags (&F, 0, 1, 0, 1e-7, 1000, w, &result, &error)
    print "result          = % .18f\n" % result
    print "estimated error          = % .18f\n" % error

【讨论】：

我也不喜欢复杂性。但我发布的是我更大问题的一个例子。我有一个 MCMC 代码，它每次都调用一个具有集成的类，其中一个瓶颈是集成。所以我决定用cython重写这个类。恐怕我无法避免我的代码中的类结构。

顺便说一句，我不知道我应该如何定义 gsl_integration_workspace 以避免该消息？！！

python 中的 MCMC？ :O。考虑到速度要求，我会完全用编译语言对其进行编码。在尝试 cython/c 之前，您应该测试更高级别的优化。自适应集成方案很慢，因为它们是通用的。他们不知道函数的形状，这通常意味着对给定精度的评估比必要的评估多。您可能对要集成的函数的形状有一个很好的了解。所以自己编写非自适应快速方案（NR 有很多例子）。

如果你真的想要 Python 中的 MCMC，我会使用 Cosmosis，暗能量调查软件 (bitbucket.org/joezuntz/cosmosis/wiki/Home)，因为他们开发了一种非常好的方法来创建许多低级模块并让他们说话彼此（使用他们创建的 C++ 抽象数据类型）

谢谢。这似乎很有趣，特别是因为我做了弱镜头。我已经调查过了！但是您的一般态度是首先删除问题！

【参考方案2】：

使用以下命令编译后：

cython test_gsl.pyx
gcc -m64 -pthread -fno-strict-aliasing -Wstrict-prototypes -DNDEBUG -g -fwrapv -O3 -Wall -fPIC -I/usr/include -I/vol/dalek/anaconda/include/python2.7 -c test_gsl.c -o build/temp.linux-x86_64-2.7/test_gsl.o
gcc -pthread -shared -L/usr/lib/ -L/vol/dalek/anaconda/lib -o test_gsl.so  build/temp.linux-x86_64-2.7/test_gsl.o -lpython2.7  -lgsl -lgslcblas -lm

如果我在没有pyximport 的情况下在python 中导入test_gsl，它会完美运行。由于pyximport 不支持链接外部库。