Cython 将双复数返回到浮点复数导致表达式不在纯 C 中

Posted 2023-04-14

【中文标题】Cython 将双复数返回到浮点复数导致表达式不在纯 C 中

【英文标题】：Cython returning a double complex to a float complex causes the expression not to be in pure C

【发布时间】：2018-11-14 11:47:26

【问题描述】：

我在 Cython 中尝试使用 complex64_t 时遇到问题。这是我的简单 cython 示例。

cimport numpy as cnp

cdef extern from "complex.h":
    double complex cexp(double complex)

cpdef example():
    cdef float b = 2.0
    cdef cnp.complex64_t temp1
    cdef cnp.complex128_t temp2

    temp1 = cexp(1j * b)
    temp2 = cexp(1j * b)

当我使用以下 setup.py 对文件进行 cythonize 处理时

from distutils.core import setup
from Cython.Build import cythonize
from distutils.extension import Extension
import numpy as np


ext_modules = [
    Extension(
        "bug_example",
        ["bug_example.pyx"],
        include_dirs=[np.get_include()],
    )
]


setup(
    name='bug_example',
    ext_modules=cythonize(ext_modules, annotate=True,
                          compiler_directives='boundscheck': False)
)

一切编译都没有问题，但我在包含的行上得到黄色（不是纯 C）

temp1 = cexp(1j * b)

但不在

temp2 = cexp(1j * b)

这似乎是一个将双复数返回到浮点复数的问题。我尝试将其显式转换为浮动复杂，例如：

temp1 = <float complex>(cexp(1j * b))

但这并没有什么不同。

谁能帮我修复我的代码，使带有 temp1 的行不再有黄色并且是纯 C。这将允许我在 cython 中使用 openmp。

【问题讨论】：

它实际上是否使用任何 Python API 调用（或以任何其他方式阻止 openmp）？正如它所说，注释着色是一个提示；如果单击+ 符号并查看C 代码，它在做什么？我认为它应该调用一些宏或内联函数进一步定义在 C 文件中，如果它不使用 Python（加上宏调用以获取 Cython complex128 的实部和图像部分），它只会执行类似 (x) + (y)*(_CYTHON_COMPLEX_64_I_CONST_HERE) 的操作，您应该能够很容易地验证。

……虽然它可能会很混乱，因为 IIRC，Cython 为 C++98、C11、C99 和 C89 编写了每个宏的替代版本，并通过检查预处理器标志在它们之间切换。

另外，您没有指定语言。如果你为 C89 构建，它将使用 complex64 和 complex128 的结构，这并不像使用 Python 或其他任何东西那么慢，但仍然不如内置的 C 复杂类型好（并且对于案例可能不正确涉及无穷大，IIRC）。

您使用的是哪个编译器？如果是gcc，那么我怀疑这条线会阻止openmp的使用，因为宏__Pyx_CREAL扩展为__real__扩展。

谢谢大家的帮助。

【参考方案1】：

黄色是由__Pyx_CREAL 和__Pyx_CIMAG 引起的，这应该不是问题，但谁知道...

为了避免这种情况，您必须避免从 double 转换到 float 并返回。

例如：

cimport numpy as cnp

#take the float version (cexpf) instead of double-version (cexp)
cdef extern from "complex.h":
     float complex cexpf(float complex)

#1j maps to double complex, so create a float version
cdef float complex float_1j = 1j

cpdef example():
    cdef float b_float = 2.0                              #use float not double
    cdef cnp.complex64_t temp1 = cexpf(float_1j*b_float)  #everything is float