arrayfun kann erheblich langsamer sein als eine explizite Schleife in matlab. Warum?

Betrachten Sie den folgenden einfachen Geschwindigkeitstest für arrayfun:

T = 4000;
N = 500;
x = randn(T, N);
Func1 = @(a) (3*a^2 + 2*a - 1);

tic
Soln1 = ones(T, N);
for t = 1:T
    for n = 1:N
        Soln1(t, n) = Func1(x(t, n));
    end
end
toc

tic
Soln2 = arrayfun(Func1, x);
toc

Auf meinem Computer (Matlab 2011b unter Linux Mint 12) lautet die Ausgabe dieses Tests:

Elapsed time is 1.020689 seconds.
Elapsed time is 9.248388 seconds.

Was zum?!? arrayfunObwohl es zugegebenermaßen eine sauberere Lösung ist, ist es eine Größenordnung langsamer. Was geht hier vor sich?

Außerdem habe ich einen ähnlichen Teststil für durchgeführt cellfunund festgestellt, dass er ungefähr dreimal langsamer ist als eine explizite Schleife. Auch dieses Ergebnis ist das Gegenteil von dem, was ich erwartet hatte.

Meine Frage ist: Warum arrayfunund cellfunso viel langsamer? Und gibt es vor diesem Hintergrund gute Gründe, sie zu verwenden (außer den Code gut aussehen zu lassen)?

Hinweis: Ich spreche hier von der Standardversion arrayfun, NICHT von der GPU-Version aus der Parallelverarbeitungs-Toolbox.

EDIT: Um ganz klar zu sein, ich bin mir bewusst, dass Func1oben, wie von Oli hervorgehoben, vektorisiert werden kann. Ich habe es nur gewählt, weil es einen einfachen Geschwindigkeitstest für die Zwecke der eigentlichen Frage liefert.

EDIT: Auf Vorschlag von Grungetta habe ich den Test mit erneut durchgeführt feature accel off. Die Ergebnisse sind:

Elapsed time is 28.183422 seconds.
Elapsed time is 23.525251 seconds.

Mit anderen Worten, es scheint, dass ein großer Teil des Unterschieds darin besteht, dass der JIT-Beschleuniger die explizite forSchleife viel besser beschleunigt als er arrayfun. Dies erscheint mir seltsam, da es arrayfuntatsächlich mehr Informationen liefert, dh seine Verwendung zeigt, dass die Reihenfolge der Anrufe Func1keine Rolle spielt. Außerdem habe ich festgestellt, dass mein System unabhängig davon, ob der JIT-Beschleuniger ein- oder ausgeschaltet ist, immer nur eine CPU verwendet ...

Sie können sich die Idee einfallen lassen, indem Sie andere Versionen Ihres Codes ausführen. Erwägen Sie, die Berechnungen explizit aufzuschreiben, anstatt eine Funktion in Ihrer Schleife zu verwenden

tic
Soln3 = ones(T, N);
for t = 1:T
    for n = 1:N
        Soln3(t, n) = 3*x(t, n)^2 + 2*x(t, n) - 1;
    end
end
toc

Zeit zum Rechnen auf meinem Computer:

Soln1  1.158446 seconds.
Soln2  10.392475 seconds.
Soln3  0.239023 seconds.
Oli    0.010672 seconds.

Während die vollständig "vektorisierte" Lösung eindeutig die schnellste ist, können Sie sehen, dass das Definieren einer Funktion, die für jeden x-Eintrag aufgerufen werden soll, einen enormen Aufwand darstellt. Nur das explizite Ausschreiben der Berechnung hat uns um den Faktor 5 beschleunigt. Ich denke, dies zeigt, dass der MATLABs JIT-Compiler keine Inline-Funktionen unterstützt . Nach der Antwort von gnovice ist es eigentlich besser, eine normale Funktion zu schreiben, als eine anonyme. Versuch es.

Nächster Schritt - Entfernen (Vektorisieren) der inneren Schleife:

tic
Soln4 = ones(T, N);
for t = 1:T
    Soln4(t, :) = 3*x(t, :).^2 + 2*x(t, :) - 1;
end
toc

Soln4  0.053926 seconds.

Ein weiterer Faktor 5 Beschleunigung: Diese Aussagen enthalten etwas, das besagt, dass Sie Schleifen in MATLAB vermeiden sollten ... Oder gibt es das wirklich? Schauen Sie sich das dann an

tic
Soln5 = ones(T, N);
for n = 1:N
    Soln5(:, n) = 3*x(:, n).^2 + 2*x(:, n) - 1;
end
toc

Soln5   0.013875 seconds.

Viel näher an der "vollständig" vektorisierten Version. Matlab speichert Matrizen spaltenweise. Sie sollten Ihre Berechnungen immer (wenn möglich) so strukturieren, dass sie "spaltenweise" vektorisiert werden.

Wir können jetzt zu Soln3 zurückkehren. Die Schleifenreihenfolge dort ist "zeilenweise". Lass es uns ändern

tic
Soln6 = ones(T, N);
for n = 1:N
    for t = 1:T
        Soln6(t, n) = 3*x(t, n)^2 + 2*x(t, n) - 1;
    end
end
toc

Soln6  0.201661 seconds.

Besser, aber immer noch sehr schlecht. Einzelschleife - gut. Doppelschleife - schlecht. Ich denke, MATLAB hat einige anständige Arbeit geleistet, um die Leistung von Loops zu verbessern, aber der Loop-Overhead ist immer noch da. Wenn Sie etwas schwerere Arbeit im Inneren hätten, würden Sie es nicht bemerken. Da diese Berechnung jedoch an die Speicherbandbreite gebunden ist, sehen Sie den Schleifen-Overhead. Und Sie werden noch deutlicher sehen, wie viel Aufwand es kostet, dort Func1 aufzurufen.

Also, was ist los mit Arrayfun? Auch dort gibt es keine Funktion, also viel Overhead. Aber warum so viel schlimmer als eine doppelt verschachtelte Schleife? Tatsächlich wurde das Thema der Verwendung von cellfun / arrayfun viele Male ausführlich diskutiert (z. B. hier , hier , hier und hier ). Diese Funktionen sind einfach langsam, Sie können sie nicht für solche feinkörnigen Berechnungen verwenden. Sie können sie für Code-Kürze und ausgefallene Konvertierungen zwischen Zellen und Arrays verwenden. Die Funktion muss jedoch schwerer sein als das, was Sie geschrieben haben:

tic
Soln7 = arrayfun(@(a)(3*x(:,a).^2 + 2*x(:,a) - 1), 1:N, 'UniformOutput', false);
toc

Soln7  0.016786 seconds.

Beachten Sie, dass Soln7 jetzt eine Zelle ist. Manchmal ist das nützlich. Die Codeleistung ist jetzt recht gut. Wenn Sie eine Zelle als Ausgabe benötigen, müssen Sie Ihre Matrix nicht konvertieren, nachdem Sie die vollständig vektorisierte Lösung verwendet haben.

Warum ist Arrayfun langsamer als eine einfache Schleifenstruktur? Leider können wir das nicht mit Sicherheit sagen, da kein Quellcode verfügbar ist. Sie können nur vermuten, dass Arrayfun eine Allzweckfunktion ist, die alle Arten von unterschiedlichen Datenstrukturen und Argumenten verarbeitet. In einfachen Fällen, die Sie direkt als Schleifennester ausdrücken können, ist es nicht unbedingt sehr schnell. Woher der Overhead kommt, können wir nicht wissen. Könnte der Overhead durch eine bessere Implementierung vermieden werden? Vielleicht nicht. Leider können wir nur die Leistung untersuchen, um die Fälle zu identifizieren, in denen es gut funktioniert, und diejenigen, in denen dies nicht der Fall ist.

Update Da die Ausführungszeit dieses Tests kurz ist, habe ich jetzt eine Schleife um die Tests hinzugefügt, um zuverlässige Ergebnisse zu erhalten:

for i=1:1000
   % compute
end

Einige Zeiten unten angegeben:

Soln5   8.192912 seconds.
Soln7  13.419675 seconds.
Oli     8.089113 seconds.

Sie sehen, dass der Arrayfun immer noch schlecht ist, aber mindestens drei Größenordnungen schlechter als die vektorisierte Lösung. Auf der anderen Seite ist eine einzelne Schleife mit spaltenweisen Berechnungen so schnell wie die vollständig vektorisierte Version ... Das alles wurde auf einer einzelnen CPU durchgeführt. Die Ergebnisse für Soln5 und Soln7 ändern sich nicht, wenn ich zu 2 Kernen wechsle. In Soln5 müsste ich ein Parfor verwenden, um es parallel zu machen. Vergessen Sie die Beschleunigung ... Soln7 läuft nicht parallel, weil arrayfun nicht parallel läuft. Olis vektorisierte Version auf der anderen Seite:

Oli  5.508085 seconds.

Das, weil!!!!

x = randn(T, N); 

ist kein gpuarrayTyp;

Alles was Sie tun müssen, ist

x = randn(T, N,'gpuArray');