NumPy Matrixmultiplikation Effizienz für Matrix mit bekannter Struktur

Question

Ich habe zwei NxN-Matrizen, die ich zusammen multiplizieren will: A und B. In NumPy, habe ich:

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

Allerdings habe ich vor, dass für die Matrix wissen B nur Zeile n und Spalte n sind ungleich Null (dies kommt direkt von der analytischen Formel, die die Matrix erzeugt hat und ist ohne Zweifel immer der Fall).

der Hoffnung, sich diese Tatsache zunutze zu machen und die Anzahl der Multiplikationen zu produzieren C erforderlich reduzieren, ersetzte ich das oben mit:

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]) 

Analytisch dies die Gesamtkomplexität verringern, sollten Sie wie folgt vor: Im allgemeinen Fall (nicht irgendwelche fancy Tricks, nur grundlegende Matrix-Multiplikation) C = AB, wo A und B beide NxN sind, sollte O (N^3) sein. Das heißt, alle N Zeilen müssen alle N Spalten multiplizieren, und jedes dieser Skalarprodukte enthält N Multiplikationen => 0 (N N N) = O (N^3) #

Ausnutzen der Struktur von B als Ich habe oben aber sollte gehen als O (N^2 + N^2) = O (2N^2) = O (N^2). Das heißt, alle N Zeilen müssen alle N Spalten multiplizieren, für alle (außer denen mit B [:, n]) ist jedoch nur eine Skalarmultiplikation erforderlich: nur ein Element von 'B [:, m]' ist nicht Null für m! = n. Wenn n == m, was N mal vorkommen wird (einmal für jede Zeile von A, die die Spalte n von B multiplizieren muss), müssen N skalare Multiplikationen auftreten. #

Der erste Code-Block (mit np.dot (A, B)) ist wesentlich schneller. Ich bin mir bewusst (über Informationen wie: Why is matrix multiplication faster with numpy than with ctypes in Python?), dass die Low-Level-Implementierungsdetails von np.dot wahrscheinlich dafür verantwortlich sind. Meine Frage lautet also: Wie kann ich die Struktur von Matrix B ausnutzen, um die Multiplikationseffizienz zu verbessern, ohne die Implementierungseffizienz von NumPy, , zu opfern, ohne meine eigene Low-Level-Matrixmultiplikation in c zu erstellen?

Diese Methode ist Teil einer numerischen Optimierung über viele Variablen, daher ist O (N^3) nicht praktikabel, während O (N^2) wahrscheinlich die Aufgabe erfüllt.

Vielen Dank für jede Hilfe. Außerdem bin ich neu in SO, also bitte entschuldigen Sie alle Anfängerfehler.

Answer 1

Wenn ich A und B richtig verstanden habe, dann verstehe ich nicht die for Loops und warum Sie vermehren sich nicht nur durch die beiden Nicht-Null-Vektoren:

# 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) 

nehmen Sie die Nicht-Null-Zeile und Spalte von B:

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

multiply A von diesen beiden vector:

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

überprüfen Kontrolle:

X - np.dot(A, B) 

mit 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

I timed es, und sparse mit schneller:

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) 

Ergebnis:

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

Da die Anzahl der Zeilen von B ansteigt, sollte auch der Zeitvorteil der Multiplikation mit Sparse-Matrix zunehmen.