Vectorizing multiple for loops in Matlab

1 view (last 30 days)
I am trying to code a mathematical model, and it involves computing a particular quantity over a grid of values thousands of times, with some changing model parameters. Currently, this is far too slow and I am looking for advice on vectorizing the most intensive part of my model.
The equation I am evaluating is:
I've currently got a basic implementation of it for ease of reading, but now want to vectorize the entire code segment below if possible. A minimal example of the code segment is:
% Setup grid to evaluate and results vector
T_max = 10000;
eval_points = linspace(0, T_max, 1000);
results = zeros(size(eval_points));
% Function that is used in computation
Z_func = @(x, omega) (1./(omega.*sqrt(2*pi))).*exp( -(x.^2)./(2.*omega.*omega) );
% Random data for now, known in full problem
historic_weights = rand(1,100);
historic_times = rand(1,100);
% Fixed single parameter omega
omega = 0.5;
% Time evaluation
tic()
for eval_counter = 1:size(eval_points,2)
temp_result = 0;
for historic_counter = 1:size(historic_weights,2)
for k = 0:1:T_max
temp_result = temp_result + Z_func( eval_points(eval_counter) - historic_times(historic_counter) + 1440*floor(historic_times(historic_counter)/1440) - 1440*k, omega );
end % End of looping over k
results(eval_counter) = results(eval_counter) + historic_weights(historic_counter)*temp_result;
end % End of looping over weights
end % End of looping over evaluation points
toc()
On my computer, this took just over 60 seconds to evaluate. I do not wish to use the parallel toolbox, as I am already using that elsewhere, and the shown segment of code is called on every process.

Accepted Answer

darova
darova on 20 Sep 2019
I changed your code a bit
See attached script
  2 Comments
Kieran Kalair
Kieran Kalair on 20 Sep 2019
Thanks, there was a typo that you spotted so your solution gives the correct results. RAM is not too much of an issuie so it even scales to quite large arrays for my avaliable hardware.
Rik
Rik on 20 Sep 2019
For those not wanting to download the script to see the solution:
clc,clear
% Setup grid to evaluate and results vector
T_max = 5;
t = linspace(0, T_max, 100);
res = zeros(size(t));
% Function that is used in computation
Z_func = @(x, omega) (1./(omega.*sqrt(2*pi))).*exp( -(x.^2)./(2.*omega.*omega) );
% Random data for now, known in full problem
hw = rand(1,10);
ht = rand(1,10);
% Fixed single parameter omega
omega = 0.5;
% Time evaluation
tic()
for i = 1:size(t,2)
for j = 1:size(hw,2)
tres = 0;
for k = 0:1:T_max
x = t(i) - ht(j) + 1440*floor(ht(j)/1440) - 1440*k;
tres = tres + Z_func( x, omega );
end % End of looping over k
res(i) = res(i) + hw(j)*tres;
end % End of looping over weights
end % End of looping over evaluation points
toc()
% VECTORIZE VERSION
tic
[T, HT, K] = ndgrid(t,ht, 0:T_max);
[~, HW] = ndgrid(t,hw);
X = T - HT + 1440*floor(HT/1440) - 1440*K;
tres = Z_func(X, omega);
tres = HW.*sum(tres,3);
res1 = sum(tres,2);
toc
res' - res1

Sign in to comment.

More Answers (0)

Products


Release

R2019a

Community Treasure Hunt

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

Start Hunting!