How to replace a "for loop" with vectorization?
219 views (last 30 days)
Show older comments
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)
%======= Defining the variables=====================
%===End of defining the variables===============
%=== Construct the equations ================
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?
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));
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.
Stephan on 6 May 2018
Edited: Stephan on 6 May 2018
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.
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
More Answers (1)
Stephan on 1 May 2018
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:
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.
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 ;-)
Find more on Loops and Conditional Statements in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!Start Hunting!