Sind DAXPY, DCOPY, DSCAL Overkills?

Ich habe CG in FORTRAN implementiert, indem ich es mit Intel MKL verknüpft habe.

Wenn es Aussagen gibt wie: ( Siehe Wikipedia )

 p=r; 
 x=x+alpha*p
 r=r-alpha*Ap;

oder ähnliche in QMR (in viel größerer Menge)

v_tld = r;
y = v_tld;
rho = norm( y );
w_tld = r;
z = w_tld;
xi = norm( z ); (and more)

Ist es sinnvoll, BLAS Level 1-Implementierungen wie DAXPY, DCOPY, DSCAL zu verwenden? Die Motivation für meine Frage ist:

1. Ich habe 2 Implementierungen der Algorithmen. Eine, bei der ich nur Normen und MatVecs mit MKL verknüpft habe; Das Kopieren, Skalieren und Hinzufügen erfolgt über Fortrans intrinsische Funktionen und eine andere, bei der jede mögliche Unterroutine von BLAS ausgeführt wird.

2. Ich war der Meinung, dass nichts schneller werden kann als BLAS. Es stellte sich jedoch heraus, dass mein Code mit Fortrans intrinsischen Funktionen 100% schneller lief als einer mit BLAS Level 1-Subroutinen (FWIW, dies war kein kleines Problem, sondern löste ein dichtes System der Größe 13k x 13k, das meine 4 füllte GB RAM). Ich habe beide auf 2 Threads (auf einem 2-Core-Rechner) ifort QMR.f90 -mklmit ausgeführtMKL_DYNAMIC=TRUE

3. Ich hatte eine SO-Frage zur Erweiterung von BLAS gestellt, aber als ich versuchte, BLAS Level 1 in meinen Code aufzunehmen, wurde mein Code immer langsamer.

Mache ich etwas falsch oder wird das erwartet?

Ist es auch sinnvoll, BLAS zu erweitern, um nicht offensichtliche Operationen wie y = 2.89*xby DCOPY(n,2.89*x,1,y,1) or even DSCAL then DCOPYauszuführen?

Was auch interessant ist, DDOTund die DNRM2Leistung zu verbessern. Ich schrieb es der Tatsache zu, dass es hilfreich sein könnte, sie parallel zu schalten, da sie Multiplikationen mit doppelter Genauigkeit ausführen.

Zusatzfrage: Wann entscheiden Sie, ob eine BLAS Level 1-Operation tatsächlich die Leistung verbessern wird?

Hinzufügen: Derzeit läuft auf einem i3 2,13-GHz-Laptop mit 4 GB RAM und Debian 64-Bit- Proc-Informationen hier . Aber ich bekomme ähnliche Antworten auf einer Intel Xeon 12 Core Workstation mit 24 GB RAM.

Wenn Ihr Ziel wirklich darin besteht, so viel Leistung wie möglich herauszuholen, ist es wichtig, sich daran zu erinnern:

1. Die (BLAS) -Bibliothek wurde möglicherweise nicht genau auf Ihr System / Ihre Konfiguration abgestimmt.
2. Bibliotheksentwickler machen Fehler.

Eine vom Hersteller optimierte BLAS-Bibliothek sollte sicherlich Ihr Standardansatz sein. Wenn Sie sich jedoch die Zeit genommen haben, einzelne Kernel zu verwenden, und festgestellt haben, dass eine andere Implementierung schneller ist, verwenden Sie auf jeden Fall die andere Implementierung. Das Fehlen der Verwendung von Vektor-Intrinsics könnte möglicherweise zu einem großen Leistungsunterschied führen.

Es ist möglich, dass Ihre beste Wahl für einfache Routinen wie daxpy und dscal eine handgeschriebene Schleife ist, die die Vektoreigenschaften ausnutzt.

Angesichts der Zustand der Compiler jetzt Tage zu optimieren, ich glaube nicht , dass es in den linearen BLAS Routinen viel Voodoo ist, zum Beispiel DAXPY, DCOPYund DSCAL, dass Ihr Compiler bereits nicht tun, zB SSE-Vektorisierung und Schleifenentrollen.

Wenn der Code derselbe ist, besteht der einzige Unterschied zwischen Ihrer Routine und einem Aufruf von MKLs BLAS im Overhead des Funktionsaufrufs und in der zusätzlichen Magie, die MKL dort möglicherweise versucht. In diesem Fall sollte der Unterschied zwischen Ihrem Code und dem MKL-Code unabhängig von der Problem- / Vektorgröße konstant sein.

Diese Frage hat interessante Echos dieser Frage , die auch DAXPYals Beispiel dient.

Der BLAS-Standard enthält tatsächlich mehrere Überprüfungen der Richtigkeit der Funktionsargumente, die in vielen Situationen nicht erforderlich sind. Siehe diese Referenzimplementierung von daxpy.f. Darüber hinaus INCXsind Ihnen Konstanten wie diese normalerweise zur Kompilierungszeit bekannt, werden jedoch möglicherweise von der Implementierung nicht angenommen. BLAS ruft Cross Compilation Units auf, und mir sind keine Compiler bekannt, die diese optimieren können, ohne die Optimierung des gesamten Programms zu aktivieren.

• Eine lustige Nebenbemerkung dazu ist, dass der Intel Compiler jetzt BLAS 3-Matrix-Matrix-Multiplikationsschleifen erkennt und diesen Code in den entsprechenden x-gemm Aufruf umwandelt , wobei genügend Optimierung aktiviert ist.

BLAS1-Funktionen stellen eine Reihe von Kerneln dar, deren Bandbreite begrenzt ist, da ihre Rechenintensität gering ist. Insbesondere machen diese Kernel O (1) -Flops pro Speicherzugriff. Dies bedeutet, dass moderne Hardware nur einen kleinen Bruchteil der Spitzenleistung erbringt und Sie im Wesentlichen nichts dagegen tun können. Die beste Implementierung von BLAS1 überprüft die Ausrichtung und moduliert die FPU-Vektorlänge und führt eine Bandbreitenspitze aus, die wahrscheinlich 5-10% der Berechnungsspitze beträgt.

Wenn Sie diese Operationen explizit in die Quelle schreiben, erkennt ein guter Compiler sie sofort und fügt eine optimale Implementierung ein, die der oben genannten BLAS1-Implementierung entspricht. Da der Compiler jedoch mehr über den Kontext weiß, kann er bestimmte Verzweigungen (nicht dass diese so wichtig sind) und den Funktionsaufruf-Overhead vermeiden sowie möglicherweise Transformationen höherer Ordnung im Code durchführen, die durch einen Funktionsaufruf blockiert würden eine undurchsichtige Bibliothek.

Sie können verschiedene Experimente durchführen, um festzustellen, was sich tatsächlich auf die Leistung Ihres Codes auswirkt. Sie sind ziemlich offensichtlich, deshalb werde ich sie hier nicht auflisten.