NumPy moltiplicazione di matrici di efficienza per matrice struttura nota

Question

Ho due matrici NxN che voglio moltiplicare insieme: A e B. In NumPy, ho usato:

import numpy as np 
C = np.dot(A, B) 

Tuttavia, mi capita di sapere che per Matrix B solo la riga n e la colonna n sono diversi da zero (questo deriva direttamente dalla formula analitica che ha prodotto la matrice ed è senza dubbio sempre il caso).

speranza di trarre vantaggio da questo fatto e ridurre il numero di moltiplicazioni necessarie per produrre C, ho sostituito il precedente con:

import numpy as np 
for row in range(0, N): 
    for col in range(0, N): 
     if col != n: 
      C[row, col] = A[row, n]*B[n, col] #Just one scalar multiplication 
     else: 
      C[row, col] = np.dot(A[row], B[:, n]) 

Analiticamente, questo dovrebbe ridurre la complessità totale come segue: Nel caso generale (non usando trucchi fantasiosi, solo moltiplicazione matrice di base) C = AB, dove A e B sono entrambi NxN, dovrebbe essere O (N^3). Cioè, ogni N righe devono moltiplicare tutte le colonne N, e ciascuno di questi prodotti dot contiene N moltiplicazioni => O (N N N) = O (N^3). #

sfruttando la struttura di B come Ho fatto sopra tuttavia dovrebbe andare come O (N^2 + N^2) = O (2N^2) = O (N^2). Ovvero, tutte le righe N devono moltiplicare tutte le colonne N, tuttavia, per tutte queste (eccetto quelle che riguardano 'B [:, n]') è richiesta una sola moltiplicazione scalare: solo un elemento di 'B [:, m]' è diverso da zero per m! = n. Quando n == m, che si verifica N volte (una volta per ogni riga di A che deve moltiplicare la colonna n di B), devono verificarsi N moltiplicazioni scalari. #

Tuttavia, il primo blocco di codice (utilizzando np.dot (A, B)) è sostanzialmente più veloce. Sono a conoscenza (tramite informazioni come: Why is matrix multiplication faster with numpy than with ctypes in Python?) che i dettagli di implementazione di basso livello di np.dot sono probabilmente da incolpare per questo. Quindi la mia domanda è questa: come posso sfruttare la struttura della matrice B per migliorare l'efficienza di moltiplicazione senza sacrificare l'efficienza di implementazione di NumPy, senza creare la mia propria moltiplicazione di matrice di basso livello in c?

Questo metodo fa parte di un'ottimizzazione numerica su molte variabili, quindi O (N^3) è intrattabile mentre O (N^2) probabilmente completerà il lavoro.

Grazie per qualsiasi aiuto. Inoltre, sono nuovo di SO, quindi per favore perdonare eventuali errori newbie.

Answer 1

Se ho capito A e B correttamente, quindi non capisco le for loop e perché non sono solo moltiplicando per i due vettori non nulli:

# say A & B are like this: 
n, N = 3, 5 
A = np.array(np.random.randn(N, N)) 

B = np.zeros_like(A) 
B[ n ] = np.random.randn(N) 
B[:, n] = np.random.randn(N) 

prendere la non-zero riga e colonna di B:

rowb, colb = B[n,:], np.copy(B[:,n]) 
colb[ n ] = 0 

moltiplicare A da questi due vettore:

X = np.outer(A[:,n], rowb) 
X[:,n] += np.dot(A, colb) 

per verificare di controllo:

X - np.dot(A, B) 

con N=100:

%timeit np.dot(A, B) 
1000 loops, best of 3: 1.39 ms per loop 

%timeit colb = np.copy(B[:,n]); colb[ n ] = 0; X = np.outer(A[:,n], B[n,:]); X[:,n] += np.dot(A, colb) 
10000 loops, best of 3: 98.5 µs per loop 

Answer 2

ho cronometrato, e l'utilizzo di sparse è più veloce:

import numpy as np 
from scipy import sparse 

from timeit import timeit 

A = np.random.rand(100,100) 
B = np.zeros(A.shape, dtype=np.float) 

B[3] = np.random.rand(100) 
B[:,3] = np.random.rand(100) 

sparse_B = sparse.csr_matrix(B) 

n = 1000 

t1 = timeit('np.dot(A, B)', 'from __main__ import np, A, B', number=n) 
print 'dense way : {}'.format(t1) 
t2 = timeit('A * sparse_B', 'from __main__ import A, sparse_B',number=n) 
print 'sparse way : {}'.format(t2) 

Risultato:

dense way : 1.15117192268 
sparse way : 0.113152980804 
>>> np.allclose(np.dot(A, B), A * sparse_B) 
True 

Come numero di righe di B aumenta, così dovrebbe il vantaggio tempo di moltiplicazione utilizzando matrice sparsa.