¿Por qué Cython es mucho más lento que Numba cuando se repite en matrices NumPy?

Al iterar sobre matrices NumPy, Numba parece dramáticamente más rápido que Cython.
¿Qué optimizaciones de Cython me faltan?

Aquí hay un ejemplo simple:

Código de Python puro:

import numpy as np def f(arr): res=np.zeros(len(arr)) for i in range(len(arr)): res[i]=(arr[i])**2 return res arr=np.random.rand(10000) %timeit f(arr) 

fuera: 4.81 ms ± 72.2 µs por bucle (media ± desviación estándar de 7 carreras, 100 bucles cada una)


Código de Cython (dentro de Jupyter):

 %load_ext cython %%cython import numpy as np cimport numpy as np cimport cython from libc.math cimport pow #@cython.boundscheck(False) #@cython.wraparound(False) cpdef f(double[:] arr): cdef np.ndarray[dtype=np.double_t, ndim=1] res res=np.zeros(len(arr),dtype=np.double) cdef double[:] res_view=res cdef int i for i in range(len(arr)): res_view[i]=pow(arr[i],2) return res arr=np.random.rand(10000) %timeit f(arr) 

Fuera: 445 µs ± 5.49 µs por bucle (media ± desviación estándar de 7 carreras, 1000 bucles cada una)


Numba código:

 import numpy as np import numba as nb @nb.jit(nb.float64[:](nb.float64[:])) def f(arr): res=np.zeros(len(arr)) for i in range(len(arr)): res[i]=(arr[i])**2 return res arr=np.random.rand(10000) %timeit f(arr) 

Fuera: 9.59 µs ± 98.8 ns por bucle (media ± desviación estándar de 7 ejecuciones, 100000 bucles cada una)


En este ejemplo, Numba es casi 50 veces más rápido que Cython.
Siendo un principiante de Cython, supongo que me estoy perdiendo algo.

Por supuesto, en este caso simple, el uso de la función vectorizada NumPy square habría sido mucho más adecuado:

 %timeit np.square(arr) 

Fuera: 5.75 µs ± 78.9 ns por bucle (media ± desviación estándar de 7 carreras, 100000 bucles cada una)

Como @Antonio ha señalado, usar pow para una simple multiplicación no es muy sabio y conduce a una sobrecarga:

Por lo tanto, la sustitución de pow(arr[i], 2) través de arr[i]*arr[i] conduce a una aceleración bastante grande:

 cython-pow-version 356 µs numba-version 11 µs cython-mult-version 14 µs 

La diferencia restante probablemente se deba a la diferencia entre los comstackdores y los niveles de optimizaciones (llvm vs MSVC en mi caso). Es posible que desee utilizar Clang para que coincida con el rendimiento numba (consulte, por ejemplo, esta respuesta SO )

Para facilitar la optimización para el comstackdor, debe declarar la entrada como matriz continua, es decir, double[::1] arr (vea esta pregunta por qué es importante para la vectorización), use @cython.boundscheck(False) (use opción -a para ver que hay menos amarillo) y también agregue indicadores de comstackción (es decir, -march=native , -march=native o similar, dependiendo de su comstackdor para habilitar la vectorización, tenga cuidado con los indicadores de comstackción que se utilizan de forma predeterminada, lo que puede impedir cierta optimización) , por ejemplo -fwrapv ). Al final, es posible que desee escribir el bucle de trabajo en C, comstackr con la combinación correcta de banderas / comstackdor y usar Cython para envolverlo.

Por cierto, al escribir los parámetros de la función como nb.float64[:](nb.float64[:]) , disminuye el rendimiento de numba; ya no se puede suponer que la matriz de entrada es continua, lo que descarta la vectorización. Permita que numba detecte los tipos (o nb.float64[::1](nb.float64[::1] como continuo, es decir, nb.float64[::1](nb.float64[::1] ), y obtendrá un mejor rendimiento:

 @nb.jit(nopython=True) def nb_vec_f(arr): res=np.zeros(len(arr)) for i in range(len(arr)): res[i]=(arr[i])**2 return res 

Conduce a la siguiente mejora:

 %timeit f(arr) # numba version # 11.4 µs ± 137 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each) %timeit nb_vec_f(arr) # 7.03 µs ± 48.9 ns per loop (mean ± std. dev. of 7 runs, 100000 loops each) 

Y como lo señaló @ max9111, no tenemos que inicializar la matriz resultante con ceros, pero podemos usar np.empty(...) lugar de np.zeros(...) : esta versión incluso supera el np.square() del numpy np.square()

Las actuaciones de diferentes enfoques en mi máquina son:

 numba+vectorization+empty 3µs np.square 4µs numba+vectorization 7µs numba missed vectorization 11µs cython+mult 14µs cython+pow 356µs