Azzera filtri
Azzera filtri

How to replace a "for loop" with vectorization?

65 visualizzazioni (ultimi 30 giorni)
Hi everyone, I have built in Simulink a user-defined function block for solving a set of algebraic nonlinear equations (with sine and cosine). It works alright. However, I used for loop to construct the equations, which is very slow when the system size increases. It takes up to one hour to solve a system consisting of 136 equations with a reasonable accuracy and step size. I have read that vectorizing can enormously speed up the simulation. I was not able to do so, however. Any help is greatly appreciated. This is my code (part of it):
function outt = DAE_Full01(x)
n=68;m=16;
%======= Defining the variables=====================
S=1;
for i1=1:m
Id(i1)=x(S);S=S+1;
end
for i2=1:m
Iq(i2)=x(S);S=S+1;
end
for i3=1:n
V(i3)=x(S);S=S+1;
end
for i4=1:n
TH(i4)=x(S);S=S+1;
end
%===End of defining the variables===============
%=== Construct the equations ================
for i6=1:m
SE1(i6,:)=Edp(i6)-V(i6)*sin(D(i6)-TH(i6))-Rs(i6)*Id(i6)+Xqp(i6)*Iq(i6);
SE2(i6,:)=Eqp(i6)-V(i6)*cos(D(i6)-TH(i6))-Rs(i6)*Iq(i6)-Xdp(i6)*Id(i6);
end
.
.
.
[xx ]=solve...
Where Rs, Xdp, Xqp are constant values. Id, Iq, Edp, Eqp, and D are inputs to the block (variables) with size m. x0 is a vector of initial values already defined. How one can replace the above for loop with vectorization form?

Risposta accettata

Stephan
Stephan il 3 Mag 2018
Modificato: Stephan il 3 Mag 2018
Hi,
I finished ;-)
You could change the code as follows:
1.) To replace the nested loop this is needed:
[b, a] = ndgrid(1:m, 1:n);
2.) I made a function for this calculation because the Performance of a function against a script is usually better. I tested both of it and found the function faster. Call the function this way:
[PV1, PV2] = vectorized_run(a, b, m,n, V, TH, Iq, Id, IC, Yabs, Yang, D);
3.) Due to matlab wants defintions of a function behind all other code, define the function as follows:
function [PV1, PV2] = vectorized_run(a, b, m, n, V, TH, Iq, Id, IC, Yabs, Yang, D)
% The first step saves calculation time, because the first part of both nested
% loops is the same for calculating sum1 and sum2 - same for PV_common
sum_common = V(b).*V(a) .* Yabs(1:m,1:n);
sum1 = sum((sum_common .* cos(TH(b)-TH(a)-Yang(1:m,1:n))),2);
sum2 = sum((sum_common .* sin(TH(b)-TH(a)-Yang(1:m,1:n))),2);
PV_common = (Id(1:m) .* V(1:m));
PV1 = diag(PV_common .* (sin(D(1:m)-TH(1:m))) + Iq(1:m).*V(1:m) .* (cos(D(1:m)-TH(1:m))) - IC(1:m,5) - sum1(1:m));
PV2 = diag(PV_common .* (cos(D(1:m)-TH(1:m))) - Iq(1:m).*V(1:m) .* (sin(D(1:m)-TH(1:m))) - IC(1:m,6) - sum2(1:m));
end
Thats it ;-)
When i tested with different sizes of the problem the vectorized code took 33...40% less time to calculate PV1 and PV2. I guess that's a pretty good result.
Please test for your purposes and give me a Feedback to this.
Best regards
Stephan
  11 Commenti
Stephan
Stephan il 4 Mag 2018
Maybe you hear from me next time again. Got some ideas but i want to try first and think about it
Best regards
Stephan
Stephan
Stephan il 6 Mag 2018
Modificato: Stephan il 6 Mag 2018
Hi Ismael,
i have some pretty good news for you. i could solve the problems regarding the global variables and improve the performance of sum3 and sum4 so that we got this result:
To reach this all funtions are now in only .m-file. The overall execution time reaced <1s and the part for solving DAE is now clearly below 600ms - in the attached test it is 576ms.
You should still keep both variants, because my computer today obviously has a good day and the variant of previously also achieved very similar results. I can not say which of the two variants works better in the end for larger problems - please test. Therefore, I will send you 3 files by email:
1.) Main_vec.m & DAE_vec.m
These now additonaly include the improved code for sum3 and sum4 for optimal performance, but still work with global variables.
2.) solve_DAE_fast.m
This file contains all the improvements that have been incorporated into the other variant, but eliminates the global variables. If you want to know more about how this works, have a look at this link:
So my suggestions for the improvement of your code are implemented and now I really have nothing more to improve ... (or maybe even later? ;-)
Best regards and keep me up to date please
Stephan

Accedi per commentare.

Più risposte (1)

Stephan
Stephan il 1 Mag 2018
Modificato: Stephan il 3 Mag 2018
Hi,
You can find a good start to this topic in the link below.
You assign the decision variables x(1) ... x(168) in the 4 loops by increasing S successively. Instead you could omit S and take advantage of the known quantities of m and n like this:
Id(1:m)=x(1:m); % x(1) ... x(16)
Iq(1:m)=x(m+1:2*m); % x(17) ... x(32)
V(1:n)=x(2*m+1:2*m+n); % x(33) ... x(100)
TH(1:n)=x(2*m+1+n:2*m+2*n); % x(101) ... x(168))
I hope that's what you wanted to achieve in your code.
The equations below in the code should also work without the for loop:
SE1(1:m,:)=Edp(1:m)-V(1:m)*sin(D(1:m)-TH(1:m))-Rs(1:m)*Id(1:m)+Xqp(1:m)*Iq(1:m);
SE2(1:m,:)=Eqp(1:m)-V(1:m)*cos(D(1:m)-TH(1:m))-Rs(1:m)*Iq(1:m)-Xdp(1:m)*Id(1:m);
Since I could not test the code, I would be very grateful for any feedback as to whether it works. I have found that with every help that I offer here in the forum, I am learning and would be happy if you contributed with your feedback.
I would also be interested to know if you have seen an improvement in the performance of your code and how large it is.
Best regards
Stephan
  10 Commenti
Ismaeel
Ismaeel il 3 Mag 2018
Thanks, Stephan. You did your best, and sorry for taking your time. Please just delete the files, appreciate it. I don't think generating matrices is precise or practical. If you want to ask about the problem, please don't mention this post or a piece of the code. Thanks once again.
Stephan
Stephan il 3 Mag 2018
Modificato: Stephan il 3 Mag 2018
Hi,
i could find the issue... the sums sum1 and sum2 work fine now. i saved about 50% calculation time for the loops for the m=68 and n=16 dimensions. Excited to hear about the runtime behavior for larger problems.
This works for both sums, but still i have to look at the PV values.
The matrices and vectors i have made are only for comparision of the results with your formulas to eliminate possible errors that i could have made. If my calculation brings the same values as your calculation i made no errors. I guess that should work for this purpose.
I keep you updated ;-)
Best regards
Stephan

Accedi per commentare.

Categorie

Scopri di più su Loops and Conditional Statements in Help Center e File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by