Benchmarking (python vs c ++ en utilisant BLAS) et (numpy)

Je souhaite écrire un programme qui utilise largement les fonctionnalités d'algèbre linéaire BLAS et LAPACK. Puisque la performance est un problème, j'ai effectué des analyses comparatives et j'aimerais savoir si l'approche que j'ai adoptée est légitime.

J'ai, pour ainsi dire, trois candidats et je souhaite tester leurs performances avec une simple multiplication matrice-matrice. Les concurrents sont:

1. Numpy, en utilisant uniquement la fonctionnalité de dot.
2. Python, appelant les fonctionnalités BLAS via un objet partagé.
3. C ++, appelant les fonctionnalités BLAS via un objet partagé.

Scénario

J'ai implémenté une multiplication matrice-matrice pour différentes dimensions i. iva de 5 à 500 avec un incrément de 5 et les matricies m1et m2sont configurés comme ceci:

m1 = numpy.random.rand(i,i).astype(numpy.float32)
m2 = numpy.random.rand(i,i).astype(numpy.float32)

1. Numpy

Le code utilisé ressemble à ceci:

tNumpy = timeit.Timer("numpy.dot(m1, m2)", "import numpy; from __main__ import m1, m2")
rNumpy.append((i, tNumpy.repeat(20, 1)))

2. Python, appelant BLAS via un objet partagé

Avec la fonction

_blaslib = ctypes.cdll.LoadLibrary("libblas.so")
def Mul(m1, m2, i, r):

    no_trans = c_char("n")
    n = c_int(i)
    one = c_float(1.0)
    zero = c_float(0.0)

    _blaslib.sgemm_(byref(no_trans), byref(no_trans), byref(n), byref(n), byref(n), 
            byref(one), m1.ctypes.data_as(ctypes.c_void_p), byref(n), 
            m2.ctypes.data_as(ctypes.c_void_p), byref(n), byref(zero), 
            r.ctypes.data_as(ctypes.c_void_p), byref(n))

le code de test ressemble à ceci:

r = numpy.zeros((i,i), numpy.float32)
tBlas = timeit.Timer("Mul(m1, m2, i, r)", "import numpy; from __main__ import i, m1, m2, r, Mul")
rBlas.append((i, tBlas.repeat(20, 1)))

3. c ++, appel de BLAS via un objet partagé

Maintenant, le code C ++ est naturellement un peu plus long, donc je réduis les informations au minimum.
Je charge la fonction avec

void* handle = dlopen("libblas.so", RTLD_LAZY);
void* Func = dlsym(handle, "sgemm_");

Je mesure le temps avec gettimeofdaycomme ça:

gettimeofday(&start, NULL);
f(&no_trans, &no_trans, &dim, &dim, &dim, &one, A, &dim, B, &dim, &zero, Return, &dim);
gettimeofday(&end, NULL);
dTimes[j] = CalcTime(start, end);

où jest une boucle exécutée 20 fois. Je calcule le temps passé avec

double CalcTime(timeval start, timeval end)
{
double factor = 1000000;
return (((double)end.tv_sec) * factor + ((double)end.tv_usec) - (((double)start.tv_sec) * factor + ((double)start.tv_usec))) / factor;
}

Résultats

Le résultat est affiché dans le graphique ci-dessous:

Des questions

1. Pensez-vous que mon approche est juste ou y a-t-il des frais généraux inutiles que je peux éviter?
2. Vous attendez-vous à ce que le résultat montre un écart aussi important entre l'approche c ++ et python? Les deux utilisent des objets partagés pour leurs calculs.
3. Puisque je préfère utiliser python pour mon programme, que puis-je faire pour augmenter les performances lors de l'appel des routines BLAS ou LAPACK?

Télécharger

Le benchmark complet peut être téléchargé ici . (JF Sebastian a rendu ce lien possible ^^)

J'ai exécuté votre benchmark . Il n'y a aucune différence entre C ++ et numpy sur ma machine:

Pensez-vous que mon approche est juste ou y a-t-il des frais généraux inutiles que je peux éviter?

Cela semble juste car il n'y a pas de différence dans les résultats.

Vous attendez-vous à ce que le résultat montre un écart aussi important entre l'approche c ++ et python? Les deux utilisent des objets partagés pour leurs calculs.

Non.

Puisque je préfère utiliser python pour mon programme, que puis-je faire pour augmenter les performances lors de l'appel des routines BLAS ou LAPACK?

Assurez-vous que numpy utilise une version optimisée des bibliothèques BLAS / LAPACK sur votre système.

MISE À JOUR (30.07.2014):

Je relance le benchmark sur notre nouveau HPC. Le matériel ainsi que la pile logicielle ont changé par rapport à la configuration dans la réponse d'origine.

Je mets les résultats dans une feuille de calcul google (contient également les résultats de la réponse originale).

Matériel

Notre HPC possède deux nœuds différents, l'un avec les processeurs Intel Sandy Bridge et l'autre avec les processeurs Ivy Bridge les plus récents:

Sandy (MKL, OpenBLAS, ATLAS):

• Processeur : 2 x 16 Intel (R) Xeon (R) E2560 Sandy Bridge à 2,00 GHz (16 cœurs)
• Mémoire RAM : 64 Go

Ivy (MKL, OpenBLAS, ATLAS):

• Processeur : 2 x 20 Intel (R) Xeon (R) E2680 V2 Ivy Bridge à 2,80 GHz (20 cœurs, avec HT = 40 cœurs)
• Mémoire RAM : 256 Go

Logiciel

La pile logicielle est pour les deux nœuds le sam. Au lieu de GotoBLAS2 , OpenBLAS est utilisé et il existe également un ATLAS BLAS multithread défini sur 8 threads (codé en dur).

• OS : Suse
• Compilateur Intel : ictce-5.3.0
• Numpy: 1.8.0
• OpenBLAS: 0.2.6
• ATLAS :: 3.8.4

Benchmark des produits scalaires

Le code de référence est le même que ci-dessous. Cependant, pour les nouvelles machines, j'ai également utilisé le benchmark pour les tailles de matrice 5000 et 8000 .
Le tableau ci-dessous comprend les résultats de référence de la réponse originale (renommée: MKL -> Nehalem MKL, Netlib Blas -> Nehalem Netlib BLAS, etc.)

Performances à filetage unique:

Performances multi-threads (8 threads):

Threads vs taille de la matrice (Ivy Bridge MKL) :

Suite Benchmark

Performances à filetage unique:

Performances multi-threads (8 threads):

Conclusion

Les nouveaux résultats de référence sont similaires à ceux de la réponse originale. OpenBLAS et MKL fonctionnent au même niveau, à l'exception du test Eigenvalue . Le test Eigenvalue ne fonctionne que raisonnablement bien sur OpenBLAS en mode mono-thread . En mode multi-thread, les performances sont pires.

Le "tableau de la taille de la matrice par rapport aux threads" montre également que bien que MKL ainsi qu'OpenBLAS s'adaptent généralement bien avec le nombre de cœurs / threads, cela dépend de la taille de la matrice. Pour les petites matrices, l'ajout de cœurs supplémentaires n'améliorera pas beaucoup les performances.

Il y a également une augmentation des performances d'environ 30% de Sandy Bridge à Ivy Bridge, ce qui pourrait être dû à une fréquence d'horloge plus élevée (+ 0,8 GHz) et / ou à une meilleure architecture.

Réponse originale (04.10.2011):

Il y a quelque temps, j'ai dû optimiser certains calculs / algorithmes d'algèbre linéaire écrits en python en utilisant numpy et BLAS, j'ai donc évalué / testé différentes configurations numpy / BLAS.

Plus précisément, j'ai testé:

• Numpy avec ATLAS
• Numpy avec GotoBlas2 (1.13)
• Numpy avec MKL (11.1 / 073)
• Numpy avec Accelerate Framework (Mac OS X)

J'ai exécuté deux benchmarks différents:

1. produit scalaire simple de matrices de différentes tailles
2. Benchmark suite qui peut être trouvée ici .

Voici mes résultats:

Machines

Linux (MKL, ATLAS, No-MKL, GotoBlas2):

• Système d'exploitation : Ubuntu Lucid 10.4 64 bits.
• Processeur : 2 x 4 Intel (R) Xeon (R) E5504 à 2,00 GHz (8 cœurs)
• Mémoire RAM : 24 Go
• Compilateur Intel : 11.1 / 073
• Scipy : 0,8
• Numpy : 1,5

Mac Book Pro (Accelerate Framework):

• Système d'exploitation : Mac OS X Snow Leopard (10.6)
• CPU : 1 Intel Core 2 Duo 2,93 Ghz (2 cœurs)
• Mémoire RAM : 4 Go
• Scipy : 0,7
• Numpy : 1,3

Mac Server (Accelerate Framework):

• Système d'exploitation : Mac OS X Snow Leopard Server (10.6)
• Processeur : 4 X Intel (R) Xeon (R) E5520 @ 2,26 Ghz (8 cœurs)
• Mémoire RAM : 4 Go
• Scipy : 0,8
• Numpy : 1.5.1

Benchmark produit dot

Code :

import numpy as np
a = np.random.random_sample((size,size))
b = np.random.random_sample((size,size))
%timeit np.dot(a,b)

Résultats :

    Système | taille = 1000 | taille = 2000 | taille = 3000 |
netlib BLAS | 1350 ms | 10900 ms | 39200 ms |    
ATLAS (1 processeur) | 314 ms | 2560 ms | 8700 ms |     
MKL (1 processeurs) | 268 ms | 2110 ms | 7120 ms |
MKL (2 processeurs) | - | - | 3660 ms |
MKL (8 processeurs) | 39 ms | 319 ms | 1000 ms |
GotoBlas2 (1 processeur) | 266 ms | 2100 ms | 7280 ms |
GotoBlas2 (2 processeurs) | 139 ms | 1009 ms | 3690 ms |
GotoBlas2 (8 processeurs) | 54 ms | 389 ms | 1250 ms |
Mac OS X (1 processeur) | 143 ms | 1060 ms | 3605 ms |
Serveur Mac (1 processeur) | 92 ms | 714 ms | 2130 ms |

Suite Benchmark

Code :
Pour plus d'informations sur la suite de référence, cliquez ici .

Résultats :

    Système | valeurs propres | svd | det | inv | point |
netlib BLAS | 1688 ms | 13102 ms | 438 ms | 2155 ms | 3522 ms |
ATLAS (1 processeur) | 1210 ms | 5897 ms | 170 ms | 560 ms | 893 ms |
MKL (1 processeurs) | 691 ms | 4475 ms | 141 ms | 450 ms | 736 ms |
MKL (2 processeurs) | 552 ms | 2718 ms | 96 ms | 267 ms | 423 ms |
MKL (8 processeurs) | 525 ms | 1679 ms | 60 ms | 137 ms | 197 ms |  
GotoBlas2 (1 processeur) | 2124 ms | 4636 ms | 147 ms | 456 ms | 743 ms |
GotoBlas2 (2 processeurs) | 1560 ms | 3278 ms | 116 ms | 295 ms | 460 ms |
GotoBlas2 (8 processeurs) | 741 ms | 2914 ms | 82 ms | 262 ms | 192 ms |
Mac OS X (1 processeur) | 948 ms | 4339 ms | 151 ms | 318 ms | 566 ms |
Serveur Mac (1 processeur) | 1033 ms | 3645 ms | 99 ms | 232 ms | 342 ms |

Installation

L'installation de MKL comprenait l'installation de la suite complète de compilateurs Intel, ce qui est assez simple. Cependant, à cause de certains bogues / problèmes, la configuration et la compilation de numpy avec le support MKL étaient un peu compliquées.

GotoBlas2 est un petit paquet qui peut être facilement compilé en tant que bibliothèque partagée. Cependant, à cause d'un bogue, vous devez recréer la bibliothèque partagée après l'avoir construite afin de l'utiliser avec numpy.
En plus de cette construction, il pour plusieurs plates-formes cibles n'a pas fonctionné pour une raison quelconque. J'ai donc dû créer un fichier .so pour chaque plateforme pour laquelle je souhaite avoir un fichier libgoto2.so optimisé .

Si vous installez numpy à partir du référentiel d'Ubuntu, il installera et configurera automatiquement numpy pour utiliser ATLAS . L'installation d' ATLAS à partir des sources peut prendre un certain temps et nécessite des étapes supplémentaires (fortran, etc.).

Si vous installez numpy sur une machine Mac OS X avec des ports Fink ou Mac, il configurera numpy pour utiliser ATLAS ou Accelerate Framework d'Apple . Vous pouvez vérifier en exécutant ldd sur le fichier numpy.core._dotblas ou en appelant numpy.show_config () .

Conclusions

MKL fonctionne mieux suivi de près par GotoBlas2 .
Dans le test de valeur propre , GotoBlas2 fonctionne étonnamment moins bien que prévu. Je ne sais pas pourquoi c'est le cas.
L'Accelerate Framework d'Apple fonctionne vraiment bien, en particulier en mode mono-thread (par rapport aux autres implémentations BLAS).

Les deux GotoBlas2 et MKL s'adaptent très bien avec le nombre de threads. Donc, si vous devez gérer de grandes matrices, l'exécuter sur plusieurs threads vous aidera beaucoup.

Dans tous les cas, n'utilisez pas l' implémentation par défaut de netlib blas car elle est bien trop lente pour tout travail de calcul sérieux.

Sur notre cluster, j'ai également installé l'ACML d'AMD et les performances étaient similaires à MKL et GotoBlas2 . Je n'ai pas de chiffres difficiles.

Personnellement, je recommanderais d'utiliser GotoBlas2 car il est plus facile à installer et gratuit.

Si vous voulez coder en C ++ / C, consultez également Eigen3 qui est censé surpasser MKL / GotoBlas2 dans certains cas et est également assez facile à utiliser.

Voici un autre benchmark (sous Linux, tapez simplement make): http://dl.dropbox.com/u/5453551/blas_call_benchmark.zip

http://dl.dropbox.com/u/5453551/blas_call_benchmark.png

Je ne vois essentiellement aucune différence entre les différentes méthodes pour les grandes matrices, entre Numpy, Ctypes et Fortran. (Fortran au lieu de C ++ --- et si cela compte, votre benchmark est probablement cassé.)

Votre CalcTimefonction en C ++ semble avoir une erreur de signe. ... + ((double)start.tv_usec))devrait être à la place ... - ((double)start.tv_usec)). Peut-être que votre benchmark a également d'autres bogues, par exemple, la comparaison entre différentes bibliothèques BLAS, ou différents paramètres BLAS tels que le nombre de threads, ou entre le temps réel et le temps CPU?

EDIT : impossible de compter les accolades dans la CalcTimefonction - c'est OK.

À titre indicatif: si vous faites un benchmark, veuillez toujours poster tout le code quelque part. Commenter les benchmarks, surtout quand il est surprenant, sans avoir le code complet n'est généralement pas productif.

Pour savoir à quel BLAS Numpy est lié, procédez comme suit:

$ python
Python 2.7.2+ (par défaut, 16 août 2011, 07:24:41) 
[GCC 4.6.1] sur linux2
Tapez "aide", "copyright", "crédits" ou "licence" pour plus d'informations.
>>> importer numpy.core._dotblas
>>> numpy.core._dotblas .__ fichier__
"/usr/lib/pymodules/python2.7/numpy/core/_dotblas.so"
>>> 
$ ldd /usr/lib/pymodules/python2.7/numpy/core/_dotblas.so
    linux-vdso.so.1 => (0x00007fff5ebff000)
    libblas.so.3gf => /usr/lib/libblas.so.3gf (0x00007fbe618b3000)
    libc.so.6 => /lib/x86_64-linux-gnu/libc.so.6 (0x00007fbe61514000)

MISE À JOUR : Si vous ne pouvez pas importer numpy.core._dotblas, votre Numpy utilise sa copie de secours interne de BLAS, qui est plus lente et n'est pas destinée à être utilisée dans le calcul des performances! La réponse de @Woltan ci-dessous indique que c'est l'explication de la différence qu'il voit dans Numpy vs Ctypes + BLAS.

Pour corriger la situation, vous avez besoin d'ATLAS ou de MKL --- vérifiez ces instructions: http://scipy.org/Installing_SciPy/Linux La plupart des distributions Linux sont livrées avec ATLAS, donc la meilleure option est d'installer leur libatlas-devpackage (le nom peut varier) .

Compte tenu de la rigueur dont vous avez fait preuve dans votre analyse, je suis surpris des résultats obtenus jusqu'à présent. Je mets cela comme une «réponse», mais seulement parce que c'est trop long pour un commentaire et offre une possibilité (même si je suppose que vous l'avez envisagée).

J'aurais pensé que l'approche numpy / python n'ajouterait pas beaucoup de frais généraux pour une matrice de complexité raisonnable, car à mesure que la complexité augmente, la proportion à laquelle python participe devrait être petite. Je suis plus intéressé par les résultats sur le côté droit du graphique, mais les écarts d'ordres de grandeur montrés là-bas seraient inquiétants.

Je me demande si vous utilisez les meilleurs algorithmes que Numpy peut exploiter. À partir du guide de compilation pour Linux:

"Build FFTW (3.1.2): Versions de SciPy> = 0.7 et Numpy> = 1.2: En raison de problèmes de licence, de configuration et de maintenance, la prise en charge de FFTW a été supprimée dans les versions de SciPy> = 0.7 et NumPy> = 1.2. une version intégrée de fftpack. Il existe plusieurs façons de tirer parti de la vitesse de FFTW si nécessaire pour votre analyse. Revenir à une version Numpy / Scipy qui inclut la prise en charge. Installez ou créez votre propre wrapper de FFTW. Voir http: //developer.berlios.de/projects/pyfftw/ comme exemple non approuvé. "

Avez-vous compilé numpy avec mkl? ( http://software.intel.com/en-us/articles/intel-mkl/ ). Si vous utilisez Linux, les instructions pour compiler numpy avec mkl sont ici: http://www.scipy.org/Installing_SciPy/Linux#head-7ce43956a69ec51c6f2cedd894a4715d5bfff974 (malgré l'url). La partie clé est:

[mkl]
library_dirs = /opt/intel/composer_xe_2011_sp1.6.233/mkl/lib/intel64
include_dirs = /opt/intel/composer_xe_2011_sp1.6.233/mkl/include
mkl_libs = mkl_intel_lp64,mkl_intel_thread,mkl_core 

Si vous êtes sous Windows, vous pouvez obtenir un binaire compilé avec mkl, (et également obtenir pyfftw, et de nombreux autres algorithmes connexes) à: http://www.lfd.uci.edu/~gohlke/pythonlibs/ , avec un Remerciements à Christoph Gohlke du Laboratoire de Dynamique de la Fluorescence, UC Irvine.

Attention, dans les deux cas, il y a de nombreux problèmes de licence et ainsi de suite dont il faut être conscient, mais la page Intel les explique. Encore une fois, j'imagine que vous avez envisagé cela, mais si vous remplissez les conditions de licence (ce qui est très facile à faire sous Linux), cela accélérerait considérablement la partie numpy par rapport à l'utilisation d'une simple construction automatique, sans même FFTW. Je serai intéressé de suivre ce fil et de voir ce que les autres en pensent. Quoi qu'il en soit, excellente rigueur et excellente question. Merci de l'avoir publié.