2017-03-10 20 views
2

Ich habe Zugriff auf eine Reihe von Matrix-Bibliotheken, aber für dieses Projekt verwende ich Eigen, aufgrund seiner Definition der Kompilierungszeit und seiner Einbeziehung von SVD.Effiziente Matrix Transponierte Matrix Multiplikation in Eigen

Nun, ich mache die folgende Operation:

Eigen::Matrix<double,M,N> A;  // populated in the code 

Eigen::Matrix<double,N,N> B = A.transpose() * A; 

Wie ich verstehe, dies macht eine Kopie von A und bildet die transponieren, die wiederum von A multipliziert wird. Diese Operation wird auf relativ kleinen Matrizen durchgeführt (M = 20-30, N = 3), aber viele Millionen Male pro Sekunde, was bedeutet, dass sie so schnell wie möglich sein muss.

ich gelesen, dass die folgende Verwendung schneller ist:

B.noalias() = A.transpose() * A; 

Ich konnte mein eigenes Unterprogramm schreiben, die A als Eingabe akzeptiert und füllt B, aber ich frage mich, ob es eine effiziente, bestehende Implementierung, die verwendet die geringste Anzahl von Zyklen.

+0

Betrachten Sie dies: http://scicomp.stackexchange.com/questions/25283/beating-typical-blas-libraries-matrix-multiplication-performance –

+0

Hilft das? http://stackoverflow.com/questions/39606224/does-eigen-have-self-transpoin-multiply-optimization-like-h-transpoheh – kennytm

Antwort

1

Erstens, da Eigen auf Vorlagenausdrücke beruht, wird A.transpose() nicht in eine temporäre ausgewertet.

Zweitens, in:

Matrix<double,N,N> B = A.transpose() * A; 

Eigen weiß, dass B nicht auf der rechten Seite des Ausdrucks erscheinen kann (weil hier der Compiler den Konstruktor von B nennt), und daher wird keine temporäre überhaupt erstellt . Dies entspricht:

Matrix<double,N,N> B;    // declare first 
B.noalias() = A.transpose() * A; // eval later 

schließlich für solche kleine Matrizen, ich erwarte nicht, dass die Verwendung von B.selfadjointView().rankUpdate(A) helfen wird (wie in kennytm Kommentar vorgeschlagen).

Auf der otherhand, mit N = 3, könnte es versuchen, die faul Implementierung wert sein:

B = A.transpose().lazyProduct(A)

nur sicher zu sein, zu. Eigens verfügt über integrierte Heuristiken, um die beste Produktimplementierung auszuwählen, aber da die Heuristik einfach und schnell zu bewerten sein muss, ist sie möglicherweise nicht zu 100% richtig.

+0

Vielen Dank. Der faule Projekttipp ist sehr wie. Jetzt habe ich etwas komplett anderes gemacht, seit ich herausgefunden habe, dass Eigen nicht in der cuda auf der GPU arbeitet. Ich mag die Bibliothek. Außerdem ist es am effizientesten, überhaupt kein A zu bauen, was ich getan habe. –