Archiv verlassen und diese Seite im Standarddesign anzeigen : MATLAB parfor
Naitsabes
2016-11-27, 22:50:59
Hey,
meine Matlabkünste sind doch eher bescheiden, sodass ich eine kleine Frage zu parfor habe.
Ich habe eine 9x9 Matrix, bei der ich einige Werte je Z mal ändern muss. Da Z in diesem Fall 100*10^6 ist, dauert das elendig lange und ich würde es gerne so weit wie möglich parallelisieren.
Aktuell sieht der Code so aus:
for i=1:length(f)
w=2*pi*f(i);
M(2,2)=1/R1+1j*w*C1;
M(3,3)=1/R2+1j*w*C3;
M(4,4)=1/R3+1j*w*C3;
M(5,5)=1/R4+1j*w*C2;
M(3,4)=-1j*w*C3;
M(4,3)=-1j*w*C3;
M(6,6)=1/R5+1j*w*C2;
M(6,5)=-1j*w*C2;
M(5,6)=-1j*w*C2;
M(9,9)=-1j*w*L1;
b=[0;0;0;0;0;0;v1;0;0];
x=inv(M)*b;
e3(i)=x(3);
e6(i)=x(6);
i_vcvs(i)=x(8);
end
Wenn ich anstelle des for ein parfor eingebe,dann meckert MATLAB
"The PARFOR loop cannot run due to the way variable 'M' is used"
Klar, M gibt es nur einmal und darum dachte ich mir, dass ich aus M eine 9x9xZ Matrix mache und dann die Teifenebenen in der parfor-Schleife fülle. Aber auch da gibt es den gleichen Fehler. :confused:
Wäre super, wenn ihr mich auf meinen Denkfehler aufmerksam machen könntet :freak:
edit.
Ich habe die Erzeugung der Matrix in eine Funktion gepackt, die jetzt in der Schleife aufgerufen wird. Damit läuft es dann :)
Ist das eine Übungsaufgabe oder ein echtes Problem?
Falls letzteres:
Wenn ich raten müsste dann versuchst du Impedanzen eines elektrischen Netzes zu berechnen. Schau Dir folgende Befehle in der Doku an:
profile (on|viewer) -> Profiler zeigt wo genau deine Rechenzeit verbraten wird
inv -> Es wird in der Doku davon abgeraten. Besser M\b verwenden. Genauer und schneller.
Willst (musst!) Du wirklich von 1 Hz bis 100 MHz in 1 Hz Schritten rechnen? Üblich wäre eine logarithmische Abstufung (siehe logspace).
Alles außer das Lösen des Gleichungssystems lässt sich vektorisiert schreiben und kann daher außerhalb der for-Schleife passieren.
Alle Ergebnisse sollten vor der Schleife alloziert werden, damit das nicht auch noch in der Schleife passieren muss.
Naitsabes
2016-11-28, 15:27:29
Ist eine Übungsaufgabe zur Knotenpotentialanalyse, genau.
Es funktioniert aber bereits und ich möchte das jetzt nur eleganter/schneller lösen. Auch 1Hz Schritte sind natürlich absolut nicht nötig, aber ich wollte eben mal rumspielen ;)
Wenn ich statt inv(M)*b M/b nutze, dann spuckt MATLAB eine Fehler raus, dass die Dimensionen nicht übereinstimmen würden. M ist 9x9 und b ist ein zeilenvektor mit 9 Elementen. Das sollte also eigentlich passen.
Ist eine Übungsaufgabe zur Knotenpotentialanalyse, genau.
Es funktioniert aber bereits und ich möchte das jetzt nur eleganter/schneller lösen. Auch 1Hz Schritte sind natürlich absolut nicht nötig, aber ich wollte eben mal rumspielen ;)
Wenn ich statt inv(M)*b M/b nutze, dann spuckt MATLAB eine Fehler raus, dass die Dimensionen nicht übereinstimmen würden. M ist 9x9 und b ist ein zeilenvektor mit 9 Elementen. Das sollte also eigentlich passen.
Nicht M/b sondern M\b.
Bevor du dich mit parfor abgibst solltest du leben deinen Code in Vektoren zu schreiben. Das bringt eigentlich immer das meiste.
Naitsabes
2016-11-28, 20:00:19
Nicht M/b sondern M\b.
Daaaaanke, ich bin echt ein Blindfisch :freak:
Wenn ich demnächst wieder etwas Zeit habe mache ich mir dann mal Gedanken über eine Vektorisierung des Codes.
Parfor erschien mir nur auf die schnelle als sinnvoll.
Danke für deine Hilfe :smile:
Daaaaanke, ich bin echt ein Blindfisch :freak:
Wenn ich demnächst wieder etwas Zeit habe mache ich mir dann mal Gedanken über eine Vektorisierung des Codes.
Parfor erschien mir nur auf die schnelle als sinnvoll.
Danke für deine Hilfe :smile:
Musste auch zweimal schauen ;-)
parfor ist sicher sinnvoll aber auch wieder Teil einer extra zu bezahlenden Toolbox. Im universitären Umfeld sicher nicht das Problem. Wenn Du aber die nächsten 5k€ für die 1123zigste Toolbox abdrücken darfst, überlegst Du ganz schnell ob es nicht auch ohne geht...
Flyinglosi
2016-12-21, 10:06:51
Wie dicht ist denn M überhaupt besetzt? Eventuell fährst du hier mit sparse-Matrizen weit besser.
Nicht M/b sondern M\b.
Bevor du dich mit parfor abgibst solltest du leben deinen Code in Vektoren zu schreiben. Das bringt eigentlich immer das meiste.
Wie meinst du das genau?
mfg
Stephan
raumfahrer
2016-12-22, 00:32:02
Ich bin etwas spät auf der Party...
Deine hundert Millionen sind etwas exzessiv um es mal vorsichtig
auszudrücken. Wie dem auch sei, ein paar Hinweise.
1) Vektorisieren
2) Zum Lösen von linearen Gleichungssystemen (A*X = B) '\' oder linsolve
verwenden
3) Matlab parallelisiert viele Vektoroperationen automatisch, sofern die Anzahl
der Elemente es sinnvoll macht [1]
Zu deiner Frage kann ich leider nichts sagen, da ich nur Octave habe und da ist
parfor noch nicht implementiert.
Was ich dir aber zeigen kann ist, was vektorisierter Code mit deinem macht.
f = [1:1e6];
e3 = zeros(size(f));
e6 = zeros(size(f));
i_vcvs = zeros(size(f));
R1 = 1;
R2 = 1;
R3 = 1;
R4 = 1;
R5 = 1;
C1 = 1;
C2 = 1;
C3 = 1;
L1 = 1;
v1 = 1;
tic
for i=1:length(f)
w=2*pi*f(i);
M = zeros(9);
M(2,2)=1/R1+1j*w*C1;
M(3,3)=1/R2+1j*w*C3;
M(4,4)=1/R3+1j*w*C3;
M(5,5)=1/R4+1j*w*C2;
M(3,4)=-1j*w*C3;
M(4,3)=-1j*w*C3;
M(6,6)=1/R5+1j*w*C2;
M(6,5)=-1j*w*C2;
M(5,6)=-1j*w*C2;
M(9,9)=-1j*w*L1;
M = eye(9);
b=[0;0;0;0;0;0;v1;0;0];
x = M \ b;
e3(i) = x(3);
e6(i) = x(6);
i_vcvs(i) = x(8);
end
toc
Ergebnis:
Elapsed time is 112.289 seconds.
Ohne 'M = eye(9)' und 'x = M \ b' => 'x = b' braucht die Sache noch etwa 104
Sekunden.
% Werte aus den Haaren, deswegen unten die Verrenkung mit eye
R1 = 1;
R2 = 1;
R3 = 1;
R4 = 1;
R5 = 1;
C1 = 1;
C2 = 1;
C3 = 1;
L1 = 1;
v1 = 1;
% Eine Millionen Matrizen aufstellen und Lösungsvektor allozieren
tic
w = 2 * pi * [1:1e6];
M = zeros(9, 9, length(w));
M(2,2,:) = 1j*w*C1 + 1/R1;
M(3,3,:) = 1j*w*C3 + 1/R2;
M(4,4,:) = 1j*w*C3 + 1/R3;
M(5,5,:) = 1j*w*C2 + 1/R4;
M(6,6,:) = 1j*w*C2 + 1/R5;
M(3,4,:) = -1j*w*C3;
M(4,3,:) = M(3,4,:);
M(5,6,:) = -1j*w*C2;
M(6,5,:) = M(5,6,:);
M(9,9,:) = -1j*w*L1;
b = zeros(9, 1); b(7) = v1;
x = zeros(9, length(w));
toc
% Eine Millionen Einheitsmatrizen aufstellen
tic
M = zeros(9, 9, length(w));
for i=1:9
M(i,i,:) = 1;
end
toc
% Eine Millionen Gleichungssysteme lösen
tic
for i = 1:length(w)
x(:,i) = M(:,:,i) \ b;
end
toc
% Relevanten Werte rausziehen
tic
e3 = x(3,:);
e6 = x(6,:);
i_vcvs = x(8,:);
toc
Ergebnis:
Elapsed time is 0.609573 seconds.
Elapsed time is 0.14601 seconds.
Elapsed time is 16.0487 seconds.
Elapsed time is 0.0159659 seconds.
Am meisten Zeit kostet der Zugriff auf die "Submatrix" 'M(:,:,i)' in der
Schleife. Wäre das eine fixe 9x9 Matrix (i.e. 'x(:,i) = M \ b') würde die
Schleife nur etwa 6.4 Sekunden brauchen.
Der Code lässt sich jetzt nicht auf deine 100e6 hochziehen (Speicherbelegung
und Indizierung der Einträge), aber da solltest du dir auch nochmal überlegen,
ob das wirklich notwendig ist.
Wie dicht ist denn M überhaupt besetzt? Eventuell fährst du hier mit sparse-Matrizen weit besser.
Wie meinst du das genau?
mfg
Stephan
12 % nicht-null Elemente, da kann man schon über sparse nachdenken. Allerdings
ist die Prämisse von 100 Millionen Auswertepunkten nicht wirklich zielführend.
Zur Frage, welchen Teil seines Posts willst du näher erläutert haben?
[1] https://www.mathworks.com/matlabcentral/answers/95958-which-matlab-functions-benefit-from-multithreaded-computation
Flyinglosi
2016-12-23, 12:58:48
@raumfahrer: Ich steht aktuell total am Schlauch, was genau man in diesem Fall unter "Vektorisieren" versteht. Und da ich auch nicht selten Matlab-Code fabriziere, würde ich gern wissen worums geht ;)
raumfahrer
2016-12-23, 22:07:15
Matlab benutzt für Vektor- und Matrixoperationen die BLAS-Routinen (IIRC
genauer MKL von Intel). Die sind stark optimiert auf Eigenschaften des
Prozessors (Cachegröße, Vektoreinheiten usw.) und auch sonst sehr hardwarenah
(in Fortran/C/Assembly) geschrieben. Kann man Matlab also dazu bringen, diese
gewrappten Routinen zu nutzen hat man a) Fortran/C-Geschwindigkeit und viel
wichtiger b) ziemlich gut auf die Hardware abgestimmte Algorithmen.
Ein anderer Punkt ist der Speicherzugriff. Matlab legt die Einträge von
Vektoren/Matrizen spaltenweise am Stück in den Speicher. Folgendes Beispiel:
a = ones(1e6,1);
b = zeros(1e6,1);
for i=1:length(a)
b(i) = 2*a(i);
end
I.d.R. ist Matlabs JIT relativ gut bei sowas, muss aber nicht. Was passiert?
Der Matlab Interpreter zieht sich das i-te Element aus a, multipliziert das mit
2 und schiebt das Ergebnis an i-te Stelle von b. Soweit so gut, i ist monoton
steigend, also einfach zu erkennen. Matlab kann also prinzipiell wissen,
welches Element als nächstes gebraucht wird. Problematisch kann aber immer noch
das b(i) sein. Jedes Ergebnis wird eins nach dem anderen 'irgendwo'
reingeschrieben. Es muss auch insbesondere nicht immer sicher sein, dass der
ganze Vektor genutzt wird (komplizierte verschachtelte Schleifen mit
Abbruchkriterien etc.).
Machen wir jetzt folgendes:
a = ones(1e6,1);
b = 2*a;
Hier kann Matlab jetzt aus der Syntax sicher sein, dass der gesamte Vektor
genutzt wird und das Ergebnis als ganzes gespeichert wird. Matlab alloziert
automatisch den benötigten Zielspeicher und stopft a, b und 2 in die
entsprechende BLAS-Routine. Die berechnet dann optimal und in "einem Rutsch"
das Ergebnis und speichert es auch in "einem Rutsch" weg. Anführungszeichen,
weil bei großen Vektoren nicht alles in eine Cacheline passt und es so gesehen
dann nicht ein "Rutsch" ist.
Genau hier hat der OP nicht aufgepasst bzw. war sich dessen nicht bewusst. Sein
Code ist Variante 1, die Schleife. Er zieht sich jedes mal ein neues Element
aus f und berechnet sein omega (omega = 2*pi*f(i)). Besser ist direkt alle
omegas zu berechnen (omega = 2*pi*f).
Wie geschrieben kann der Matlab JIT mittlerweile vieles erkennen und
optimieren. Es ist aber eigentlich besser, zu überlegen welche Operationen man
direkt auf ganze Vektoren/Matrizen anwenden kann anstatt sie Elementweise
durchzuführen.
Wenn du schreibst, dass du Matlab nutzt, war das hoffentlich nicht zu basic.
Timbaloo
2016-12-23, 23:20:28
Wie geschrieben kann der Matlab JIT mittlerweile vieles erkennen
Naaaja. Er wurde vor einer Weile (R2010a/b?) deutlich überarbeitet und ist dann bei Schleifen und ähnlichem nicht mehr unendlich langsam sondern nur noch gähnend langsam, aber der JIT (wenn man es so überhaupt nennen kann) ist ein Witz, denn bereits dein Schleifenbeispiel könnte ein guter JIT direkt sauber in die zweite Variante überführen wenn er doch nur die komplette (und hier triviale) Schleife betrachten würde. Machen wir uns nichts vor, der Interpreter ist scheisse, deswegen sind Vektor/Matrixausdrücke nötig. Mal ganz davon abgesehen, dass es idR auch deutlich ausdruckstärker ist.
Wie es richtig geht zeigt Julia.
Naitsabes
2016-12-24, 12:12:37
Sehr interessant mal genauere Einblicke in die Programmierung mit matlab zu bekommen.
Bei mir ist das aktuell immer nur Trial & Error und ich versuche mir das irgendwie selber beizubringen. Das ist bei uns an der Uni echt chaotisch - im Master wird in fast jedem Fach selbstverständlich davon ausgegangen, dass man Matlab drauf hat. Aber es gibt keinen einzigen Kurs, in dem es einem näher gebracht wird (klar, Selbststudium gehört zum Studium dazu, aber es ist nervig wenn alle Beispiele nur noch in Form von Matlab präsentiert werden (und der Code nicht verfügbar gemacht wird). Dann sitzt man in der Übung da und versucht halbwegs nachzuvollziehen was da wohl gemacht wurde und hofft dann die Lösung der Aufgabe zu verstehen :uponder:)
Flyinglosi
2016-12-24, 12:24:01
@raumfahrer: Danke für die ausführliche Erklärung. Mich amüsiert diesbezüglich, dass ich selbst zu Gunsten von kompaktem Code ohnehin wenn möglich stets die "vektorisierte" Variante genutzt habe.
@Naitsabes: Du wirst im beruflichen Alltag häufig neue Tools erlernen müssen, ohne dass dir jemand dabei hilft. Ich würde dir aber diesbezüglich denn Tipp geben, dich jetzt schon mit opensource-Lösungen zu beschäftigen. Im Studium ist der Reiz gross kommerzielle Lösugen wie Matlab einzuseten, aber später darf dein Arbeitgeber dafür tief in die Tasche greifen.
Mit freundlichen Grüßen
Stephan
Mit freundlichen Grüßen
Stephan
Naitsabes
2016-12-24, 12:40:20
Ja, das ist mir auch klar ;)
Dieses Semester ist es nur wirklich sehr extrem - ich versenke mehr Zeit damit mich in die jeweiligen Programme einzuarbeiten, als mit den eigentlichen Themen der Vorlesungen und musste mich gerade etwas auskotzen :D
edit.
also gefühlt mehr Zeit in die Programme investiert. Ist natürlich Quatsch, aber ihr kennt das ja, dass man manchmal einfach nur genervt ist :ufinger:
vBulletin®, Copyright ©2000-2024, Jelsoft Enterprises Ltd.