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õesCython versão 0.22 (versão mais recente do Pypi quando a pergunta foi feita) e GCC 4.9.2.
RepositórioCriei 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.