Números complexos em Cython

Qual é a maneira correta de trabalhar com números complexos no Cython?

Gostaria de escrever um loop C puro usando um numpy.ndarray do dtype np.complex128. No Cython, o tipo C associado é definido emCython/Includes/numpy/__init__.pxd Como

ctypedef double complex complex128_t

parece que este é apenas um simples complexo duplo C.

No entanto, é fácil obter comportamentos estranhos. Em particular, com essas definições

cimport numpy as np
import numpy as np
np.import_array()

cdef extern from "complex.h":
    pass

cdef:
    np.complex128_t varc128 = 1j
    np.float64_t varf64 = 1.
    double complex vardc = 1j
    double vard = 1.

a linha

varc128 = varc128 * varf64

pode ser compilado pelo Cython, mas o gcc não pode compilar o código C produzido (o erro é "testcplx.c: 663: 25: erro: dois ou mais tipos de dados nos especificadores de declaração" e parece estar relacionado à linhatypedef npy_float64 _Complex __pyx_t_npy_float64_complex;) Este erro já foi relatado (por exemplo,aqui), mas não encontrei nenhuma boa explicação e / ou solução limpa.

Sem inclusão decomplex.h, não há erro (acho que porque otypedef então não está incluído).

No entanto, ainda existe um problema, pois no arquivo html produzido pelocython -a testcplx.pyx, a linhavarc128 = varc128 * varf64 é amarelo, o que significa que não foi traduzido para C. puro. O código C correspondente é:

__pyx_t_2 = __Pyx_c_prod_npy_float64(__pyx_t_npy_float64_complex_from_parts(__Pyx_CREAL(__pyx_v_8testcplx_varc128), __Pyx_CIMAG(__pyx_v_8testcplx_varc128)), __pyx_t_npy_float64_complex_from_parts(__pyx_v_8testcplx_varf64, 0));
__pyx_v_8testcplx_varc128 = __pyx_t_double_complex_from_parts(__Pyx_CREAL(__pyx_t_2), __Pyx_CIMAG(__pyx_t_2));

e a__Pyx_CREAL e__Pyx_CIMAG são laranja (chamadas em Python).

Curiosamente, a linha

vardc = vardc * vard

não produz nenhum erro e é traduzido para C puro (apenas__pyx_v_8testcplx_vardc = __Pyx_c_prod(__pyx_v_8testcplx_vardc, __pyx_t_double_complex_from_parts(__pyx_v_8testcplx_vard, 0));), embora seja muito parecido com o primeiro.

Eu posso evitar o erro usando variáveis intermediárias (e isso se traduz em C puro):

vardc = varc128
vard = varf64
varc128 = vardc * vard

ou simplesmente lançando (mas não se traduz em C puro):

vardc = <double complex>varc128 * <double>varf64

Então o que acontece? Qual é o significado do erro de compilação? Existe uma maneira limpa de evitá-lo? Por que a multiplicação de um np.complex128_t e np.float64_t parece envolver chamadas de Python?

Versões

Cython versão 0.22 (versão mais recente do Pypi quando a pergunta foi feita) e GCC 4.9.2.

Repositório

Criei um pequeno repositório com o exemplo (hg clone https://bitbucket.org/paugier/test_cython_complex) e um pequeno Makefile com 3 destinos (make clean, make build, make html), por isso é fácil testar qualquer coisa.

questionAnswers(0)

yourAnswerToTheQuestion