41
Cercavo di capire il modo più veloce per fare moltiplicazione matriciale e provato 3 modi diversi:Perché la moltiplicazione della matrice è più veloce con numpy che con ctypes in Python?
- attuazione Pure pitone: sorprese qui.
- Implementazione di Numpy con
numpy.dot(a, b) - Interfaccia con C utilizzando il modulo
ctypesin Python.
Questo è il codice C che si trasforma in una libreria condivisa:
#include <stdio.h>
#include <stdlib.h>
void matmult(float* a, float* b, float* c, int n) {
int i = 0;
int j = 0;
int k = 0;
/*float* c = malloc(nay * sizeof(float));*/
for (i = 0; i < n; i++) {
for (j = 0; j < n; j++) {
int sub = 0;
for (k = 0; k < n; k++) {
sub = sub + a[i * n + k] * b[k * n + j];
}
c[i * n + j] = sub;
}
}
return ;
}
E il codice Python che lo chiama:
def C_mat_mult(a, b):
libmatmult = ctypes.CDLL("./matmult.so")
dima = len(a) * len(a)
dimb = len(b) * len(b)
array_a = ctypes.c_float * dima
array_b = ctypes.c_float * dimb
array_c = ctypes.c_float * dima
suma = array_a()
sumb = array_b()
sumc = array_c()
inda = 0
for i in range(0, len(a)):
for j in range(0, len(a[i])):
suma[inda] = a[i][j]
inda = inda + 1
indb = 0
for i in range(0, len(b)):
for j in range(0, len(b[i])):
sumb[indb] = b[i][j]
indb = indb + 1
libmatmult.matmult(ctypes.byref(suma), ctypes.byref(sumb), ctypes.byref(sumc), 2);
res = numpy.zeros([len(a), len(a)])
indc = 0
for i in range(0, len(sumc)):
res[indc][i % len(a)] = sumc[i]
if i % len(a) == len(a) - 1:
indc = indc + 1
return res
avrei scommesso che la versione con C sarebbe stato più veloce ... e avrei perso! Qui di seguito è il mio punto di riferimento che sembra dimostrare che io sia fatto in modo errato, o che numpy è stupidamente veloce:
mi piacerebbe capire il motivo per cui la versione numpy è più veloce rispetto alla versione ctypes, io' Non parlo nemmeno della pura implementazione di Python poiché è abbastanza ovvia.
Bella domanda: risulta np.dot() è anche più veloce di un'implementazione ingenua della GPU in C. – user2398029
Una delle cose più importanti che rende il tuo lento C matmul lento è il modello di accesso alla memoria. 'b [k * n + j];' all'interno del ciclo interno (su 'k') ha un passo di' n', quindi tocca una linea cache diversa su ogni accesso. E il tuo ciclo non può auto-vectorize con SSE/AVX. ** Risolvilo trasposando 'b' in avanti, il che costa O (n^2) tempo e si ripaga da solo in carenze di cache ridotte mentre fai carichi O (n^3) da' b'. ** Ciò farebbe comunque essere un'implementazione ingenua senza il blocco della cache (aka loop tiling), però. –
Poiché si utilizza un 'int sum' (per qualche motivo ...), il ciclo potrebbe effettivamente vectorize senza' -ffast-math' se il ciclo interno stava accedendo a due array sequenziali. FP math non è associativo, quindi i compilatori non possono riordinare le operazioni senza '-ffast-math', ma l'intero matematico è associativo (e ha una latenza inferiore rispetto a FP, il che aiuta se non hai intenzione di ottimizzare il tuo ciclo con accumulatori multipli o altre cose che nascondono la latenza). I costi di conversione 'float' ->' int' sono gli stessi di un 'add' FP (in realtà si usa FP aggiungendo ALU su CPU Intel), quindi non ne vale la pena nel codice ottimizzato. –