Numpy vs Cython speed

j'ai un code d'analyse qui fait de lourdes opérations numériques en utilisant numpy. Juste par curiosité, j'ai essayé de le compiler avec cython avec de petits changements, puis je l'ai réécrit en utilisant des boucles pour la partie pépère.

a ma grande surprise, le code basé sur les boucles était beaucoup plus rapide (8x). Je ne peux pas poster le code complet, mais j'ai assemblé un calcul indépendant très simple qui montre un comportement similaire (bien que la différence de timing n'est pas si grande):

Version 1 (sans cython)

import numpy as np

def _process(array):

    rows = array.shape[0]
    cols = array.shape[1]

    out = np.zeros((rows, cols))

    for row in range(0, rows):
        out[row, :] = np.sum(array - array[row, :], axis=0)

    return out

def main():
    data = np.load('data.npy')
    out = _process(data)
    np.save('vianumpy.npy', out)

Version 2 (création d'un module avec cython)

import cython
cimport cython

import numpy as np
cimport numpy as np

DTYPE = np.float64
ctypedef np.float64_t DTYPE_t

@cython.boundscheck(False)
@cython.wraparound(False)
@cython.nonecheck(False)
cdef _process(np.ndarray[DTYPE_t, ndim=2] array):

    cdef unsigned int rows = array.shape[0]
    cdef unsigned int cols = array.shape[1]
    cdef unsigned int row
    cdef np.ndarray[DTYPE_t, ndim=2] out = np.zeros((rows, cols))

    for row in range(0, rows):
        out[row, :] = np.sum(array - array[row, :], axis=0)

    return out

def main():
    cdef np.ndarray[DTYPE_t, ndim=2] data
    cdef np.ndarray[DTYPE_t, ndim=2] out
    data = np.load('data.npy')
    out = _process(data)
    np.save('viacynpy.npy', out)

Version 3 (Construire un module avec cython)

import cython
cimport cython

import numpy as np
cimport numpy as np

DTYPE = np.float64
ctypedef np.float64_t DTYPE_t

@cython.boundscheck(False)
@cython.wraparound(False)
@cython.nonecheck(False)
cdef _process(np.ndarray[DTYPE_t, ndim=2] array):

    cdef unsigned int rows = array.shape[0]
    cdef unsigned int cols = array.shape[1]
    cdef unsigned int row
    cdef np.ndarray[DTYPE_t, ndim=2] out = np.zeros((rows, cols))

    for row in range(0, rows):
        for col in range(0, cols):
            for row2 in range(0, rows):
                out[row, col] += array[row2, col] - array[row, col]

    return out

def main():
    cdef np.ndarray[DTYPE_t, ndim=2] data
    cdef np.ndarray[DTYPE_t, ndim=2] out
    data = np.load('data.npy')
    out = _process(data)
    np.save('vialoop.npy', out)

avec une matrice 10000x10 sauvegardée dans les données.npy, les temps sont:

$ python -m timeit -c "from version1 import main;main()"
10 loops, best of 3: 4.56 sec per loop

$ python -m timeit -c "from version2 import main;main()"
10 loops, best of 3: 4.57 sec per loop

$ python -m timeit -c "from version3 import main;main()"
10 loops, best of 3: 2.96 sec per loop

Est-ce attendu ou y a-t-il une optimisation qui me manque? Le fait que les versions 1 et 2 donnent le même résultat est en quelque sorte attendu, mais pourquoi la version 3 est plus rapide?

Ps.- Ce n'est PAS le calcul que je dois faire, juste un exemple simple qui montre la même chose.

33
demandé sur Hernan 2011-10-18 01:46:35

5 réponses

comme mentionné dans les autres réponses, la version 2 est essentiellement la même que la version 1 puisque cython n'est pas en mesure de creuser dans l'opérateur array access pour l'optimiser. Il y a 2 raisons à cela

  • tout d'abord, il y a un certain montant de frais généraux dans chaque appel à une fonction numpy, par rapport au code C optimisé. Toutefois, ces frais généraux deviendront moins importants si chaque opération traite de grands tableaux

  • deuxièmement, il y a le création de tableaux intermédiaires. Ceci est plus clair si vous considérez une opération plus complexe comme out[row, :] = A[row, :] + B[row, :]*C[row, :]. Dans ce cas, de toute une panoplie B*C doit être créé en mémoire, puis ajouté à A. Cela signifie que le cache CPU est en train d'être écrasé, car les données sont lues et écrites en mémoire plutôt que d'être conservées dans le CPU et utilisées immédiatement. Fait important, ce problème devient pire si vous avez affaire à de grands tableaux.

D'autant plus que vous dites que votre code réel est plus complexe que l'exemple, et il montre beaucoup plus rapide, je pense que la deuxième raison est sans doute le principal facteur dans votre cas.

a part cela, si vos calculs sont suffisamment simples, vous pouvez surmonter cet effet en utilisant numexpr, bien que cython soit bien sûr utile dans bien d'autres situations, donc c'est peut-être la meilleure approche pour vous.

34
répondu DaveP 2018-04-03 06:46:47

avec une légère modification, la version 3 devient deux fois plus rapide:

@cython.boundscheck(False)
@cython.wraparound(False)
@cython.nonecheck(False)
def process2(np.ndarray[DTYPE_t, ndim=2] array):

    cdef unsigned int rows = array.shape[0]
    cdef unsigned int cols = array.shape[1]
    cdef unsigned int row, col, row2
    cdef np.ndarray[DTYPE_t, ndim=2] out = np.empty((rows, cols))

    for row in range(rows):
        for row2 in range(rows):
            for col in range(cols):
                out[row, col] += array[row2, col] - array[row, col]

    return out

Le goulot d'étranglement dans votre calcul d'accès à la mémoire. Votre tableau d'entrées est commandé par C, ce qui signifie que se déplacer le long du dernier axe fait le plus petit saut en mémoire. Par conséquent, votre boucle interne doit être le long de l'axe 1, et non de l'axe 0. Ce changement réduit de moitié le temps de course.

si vous devez utiliser cette fonction sur de petits tableaux d'entrées, alors vous pouvez réduire les frais généraux en utilisant np.empty au lieu np.ones. Pour réduire les frais généraux utilisation ultérieure PyArray_EMPTY à partir de la numpy C API.

si vous utilisez cette fonction sur de très grands tableaux d'entrées (2**31) alors les entiers utilisés pour l'indexation (et dans le range fonction) va déborder. Pour être sûr d'utilisation:

cdef Py_ssize_t rows = array.shape[0]
cdef Py_ssize_t cols = array.shape[1]
cdef Py_ssize_t row, col, row2

au lieu de

cdef unsigned int rows = array.shape[0]
cdef unsigned int cols = array.shape[1]
cdef unsigned int row, col, row2

Durée:

In [2]: a = np.random.rand(10000, 10)
In [3]: timeit process(a)
1 loops, best of 3: 3.53 s per loop
In [4]: timeit process2(a)
1 loops, best of 3: 1.84 s per loop

process est votre version 3.

39
répondu kwgoodman 2012-05-07 22:41:39

je recommande d'utiliser l'option-a pour que cython génère le fichier html qui montre ce qui est traduit en C pur vs appelant L'API python:

http://docs.cython.org/src/quickstart/cythonize.html

la Version 2 donne à peu près le même résultat que la Version 1, parce que tout le levage lourd est fait par L'API Python (via numpy) et cython ne fait rien pour vous. En fait, sur ma machine, numpy est construit contre MKL, donc quand je compiler le code c généré par cython en utilisant gcc, la Version 3 est en fait un peu plus lente que les deux autres.

Cython brille quand vous faites une manipulation de tableau que numpy ne peut pas faire d'une manière' vectorisée', ou quand vous faites quelque chose d'intensif en mémoire qui vous permet d'éviter de créer un grand tableau temporaire. J'ai eu 115x speed-ups en utilisant cython vs numpy pour mon propre code.:

https://github.com/synapticarbors/pylangevin-integrator

une partie de cela était d'appeler le répertoire randomkit au niveau du code c au lieu de l'appeler par numpy.random, mais la plupart de cela était cython traduisant l'intensif computationnellement pour les boucles en c PUR sans appels à python.

7
répondu JoshAdel 2011-10-18 00:21:55

la différence peut être due aux versions 1 et 2 faisant un appel au niveau de Python à np.sum() pour chaque ligne, alors que la version 3 compile probablement vers une boucle c serrée et pure.

L'étude de la différence entre la version 2 et la version 3 de Cython-generated C source devrait être instructive.

3
répondu Pi Delport 2011-10-17 22:17:04

je suppose que le principal overhead que vous sauvegardez est les tableaux temporaires créés. Vous créez un grand tableau array - array[row, :], puis le réduire dans un plus petit tableau en utilisant sum. Mais construire ce grand tableau temporaire ne sera pas libre, surtout si vous avez besoin d'allouer la mémoire.

1
répondu wisty 2011-10-17 23:14:58