Purtroppo, a causa della mia svista, ho avuto una versione precedente di MKL (11.1) collegato contro NumPy. La versione più recente di MKL (11.3.1) offre le stesse prestazioni in C e quando viene chiamata da python.
Ciò che oscurava le cose, era anche se collegava le librerie condivise compilate esplicitamente con il più recente MKL, e indicando le variabili LD_ * a loro, e poi in python facendo import numpy, stava facendo in modo che Python chiamasse le vecchie librerie MKL. Solo rimpiazzando nella cartella lib python tutti libmkl _ *. Quindi con MKL più recente ero in grado di eguagliare le prestazioni nelle chiamate python e C.
Informazioni su sfondo/libreria.
La moltiplicazione della matrice è stata effettuata tramite chiamate della libreria MKL di Intel (single-precision) e dgemm (doppia precisione), tramite la funzione numpy.dot. L'effettiva chiamata delle funzioni della libreria può essere verificata con ad es. oprof.
Utilizzando qui CPU core 2x18 E5-2699 v3, quindi un totale di 36 core fisici. KMP_AFFINITY = scatter. Funzionando su Linux.
TL; DR
1) Perché numpy.dot, anche se si sta chiamando le stesse funzioni di libreria MKL, due volte più lento nella migliore rispetto a C codice compilato?
2) Perché tramite numpy.dot si ottiene una riduzione delle prestazioni con l'aumento del numero di core, mentre lo stesso effetto non si osserva nel codice C (chiamando le stesse funzioni della libreria).
Il problema
ho osservato che facendo la moltiplicazione di matrici di singola doppia precisione/galleggia in numpy.dot, così come chiamare cblas_sgemm/dgemm direttamente da un C libreria condivisa compilato dare notevolmente peggiore prestazioni rispetto a chiamare le stesse funzioni MKL cblas_sgemm/dgemm dall'interno del puro codice C.
import numpy as np
import mkl
n = 10000
A = np.random.randn(n,n).astype('float32')
B = np.random.randn(n,n).astype('float32')
C = np.zeros((n,n)).astype('float32')
mkl.set_num_threads(3); %time np.dot(A, B, out=C)
11.5 seconds
mkl.set_num_threads(6); %time np.dot(A, B, out=C)
6 seconds
mkl.set_num_threads(12); %time np.dot(A, B, out=C)
3 seconds
mkl.set_num_threads(18); %time np.dot(A, B, out=C)
2.4 seconds
mkl.set_num_threads(24); %time np.dot(A, B, out=C)
3.6 seconds
mkl.set_num_threads(30); %time np.dot(A, B, out=C)
5 seconds
mkl.set_num_threads(36); %time np.dot(A, B, out=C)
5.5 seconds
facendo esattamente la stessa di cui sopra, ma con doppia precisione A, B e C, si ottiene: 3 nuclei: 20s, 6 core: 10s, 12 core: 5s, 18 core: 4.3s, 24 core: 3, 30 core: 2,8 s, 36 core: 2,8 s.
Il rabbocco di velocità per punti flottanti a precisione singola sembra essere associato a errori di cache. Per 28 core run, ecco l'output di perf. Per precisione singola:
perf stat -e task-clock,cycles,instructions,cache-references,cache-misses ./ptestf.py
631,301,854 cache-misses # 31.478 % of all cache refs
E doppia precisione:
93,087,703 cache-misses # 5.164 % of all cache refs
C libreria condivisa, compilato con
/opt/intel/bin/icc -o comp_sgemm_mkl.so -openmp -mkl sgem_lib.c -lm -lirc -O3 -fPIC -shared -std=c99 -vec-report1 -xhost -I/opt/intel/composer/mkl/include
#include <stdio.h>
#include <stdlib.h>
#include "mkl.h"
void comp_sgemm_mkl(int m, int n, int k, float *A, float *B, float *C);
void comp_sgemm_mkl(int m, int n, int k, float *A, float *B, float *C)
{
int i, j;
float alpha, beta;
alpha = 1.0; beta = 0.0;
cblas_sgemm(CblasRowMajor, CblasNoTrans, CblasNoTrans,
m, n, k, alpha, A, k, B, n, beta, C, n);
}
funzione wrapper Python, chiamando il sopra libreria compilata:
def comp_sgemm_mkl(A, B, out=None):
lib = CDLL(omplib)
lib.cblas_sgemm_mkl.argtypes = [c_int, c_int, c_int,
np.ctypeslib.ndpointer(dtype=np.float32, ndim=2),
np.ctypeslib.ndpointer(dtype=np.float32, ndim=2),
np.ctypeslib.ndpointer(dtype=np.float32, ndim=2)]
lib.comp_sgemm_mkl.restype = c_void_p
m = A.shape[0]
n = B.shape[0]
k = B.shape[1]
if np.isfortran(A):
raise ValueError('Fortran array')
if m != n:
raise ValueError('Wrong matrix dimensions')
if out is None:
out = np.empty((m,k), np.float32)
lib.comp_sgemm_mkl(m, n, k, A, B, out)
Tuttavia, le chiamate esplicite da un binario C-compilato che chiama MKL's cblas_sgemm/cblas_dgemm, con le matrici allocate tramite malloc in C, offrono prestazioni quasi 2x rispetto al codice python, ovvero la chiamata numpy.dot. Inoltre, l'effetto del degrado delle prestazioni con un numero crescente di core NON viene osservato. La migliore prestazione è stata di 900 ms per la moltiplicazione di matrice a precisione singola ed è stata ottenuta utilizzando tutti i 36 core fisici tramite mkl_set_num_cores ed eseguendo il codice C con numactl --interleave = all.
Forse qualche strumento o consiglio di fantasia per profilare/ispezionare/comprendere ulteriormente questa situazione? Anche qualsiasi materiale di lettura è molto apprezzato.
UPDATE Seguendo il consiglio @Hristo Iliev, in esecuzione numactl --interleave = tutti ./ipython non ha modificato i tempi (entro rumore), ma migliora i puri C runtime binari.
Probabilmente non stai raggiungendo il limite di scalabilità con il doppio poiché è più di 2 volte più lavoro che in precisione singola. Se ridurrai le dimensioni della matrice potresti osservare lo stesso comportamento anche con la doppia precisione. – Elalfer
Ho dovuto ridurre la dimensione della matrice a n = 1000 per la doppia precisione, in modo che la degradazione delle prestazioni diventasse osservabile con l'aggiunta di più core. Con le taglie più alte, è solo il massimo. Inoltre, non è solo 2x più lavoro (a causa della vettorizzazione), ma 2 volte più memoria da trasferire. –
Provare a eseguire l'interprete Python come 'numactl --interleave = nodes python' e ripetere di nuovo i benchmark. –