How to make a multiple for-loop run faster, the function of which is to find the steady state?

1 view (last 30 days)
Hi. I am programing a double for-loop to find the steady state of a iterative physical process (include fft/ifft of matrix). The program works well with default single thread calculation but takes to much time. Whereas the 'parfor' I've tried seems to be useless here, and even spent more time. I attach a part of the codes below. Much appreciate if anyone could give me advice to help improve the speed of my codes, or a better way to find the steady state.
%% input parameters and create variables
nt = 2^12; % sampling points
N = 2000; % max number of iteration
A = zeros(N,nt); % create initial A
A(1,:) = randn(1,nt)+1i*randn(1,nt);
% ...and other parameters
%% Main loop
rmse = ones(N,1); % Create Root-mean-square error variable
B = fftshift(fft(ifftshift(A)));
C = abs(B).^2;
flag = 0; % initial steady-state mark == 0
for m = 1:10 % to find the minimum m value which leads to steady state (linear, when m>=threshold, become steady)
for n = 1:N-1 % iterative process
A(n+1,:) = Function(A(n,:),m);% a custom function, input A(n) and m, output A(n+1)
B(n+1,:) = fftshift(fft(ifftshift(A(n+1,:))));
C(n+1,:) = abs(B(n+1,:)).^2;
rmse(n+1) = sqrt(mean((C(n+1,:)-C(n,:)).^2)); % calculate Root-mean-square error between C(n) and C(n+1)
fin = n+1; % store n value
if (rmse(fin) < 0.02) && (fin > 2) % check out steady state
flag = 1; % set mark == 1 when at steady state
break
end
end
if flag == 1 % junmp out the loop
break
end
end
%% Visualization
figure;
subplot(2,1,1);
mesh(A);
subplot(2,1,2);
mesh(C);

Accepted Answer

Raymond Norris
Raymond Norris on 9 Nov 2021
Yuu ought to be able to change the last several lines as such
end
temp(m,1) = rmsep;
end
rmse = temp;
With that said, I'm not going to go line by line, other than to say, I don't think what you originally posted maps to your parfor suggestion. Specifically, the original example might exit early, whereas your parfor example will not (well, only partially).

More Answers (1)

Raymond Norris
Raymond Norris on 8 Nov 2021
Edited: Raymond Norris on 8 Nov 2021
Can you elaborate a bit on the following:
  • how are you measuring the time (walltime, tic/toc, profiler)
  • how long it took to run serially and what speedup you're looking for (e.g. takes 2 days, lookng for 3 hours, etc.)
  • what compute do you have access to that you're running your code on (e.g. 4 cores, 16 GB RAM, etc.)
  • I'm assuming you're trying to parallelize the outter for-loop since rmse has a dependency on C(n+1) and C(n). Regardless of which one you chose, you'd need to re-write the for-loop as they're currently written. Can you post what your parfor loops look like.
steady-state problems (i.e. "needle in the haystack") won't work with parfor, where you can't preemptively end the loop. Instead, you're better off with parfeval, where you can see the results and cancel "futures" when you've met a threashold.
For instance
% Set the size of the pool how you choose (can be less then the number of
% iterations
p = parpool(10);
for m = 10:-1:1
f(m) = p.parfeval(@calc_steady_state,3,m);
end
for m = 1:10
[indx, flag, A, C] = f(m).fetchNext;
if flag==true
% Found our steady-state early
f.cancel
end
end
%% Visualization
figure;
subplot(2,1,1);
mesh(A);
subplot(2,1,2);
mesh(C);
function [flag, A, C] = calc_steady_state
%% input parameters and create variables
nt = 2^12; % sampling points
N = 2000; % max number of iteration
A = zeros(N,nt); % create initial A
A(1,:) = randn(1,nt)+1i*randn(1,nt);
% ...and other parameters
%% Main loop
rmse = ones(N,1); % Create Root-mean-square error variable
B = fftshift(fft(ifftshift(A)));
C = abs(B).^2;
flag = 0;
for n = 1:N-1 % iterative process
A(n+1,:) = Function(A(n,:),m);% a custom function, input A(n) and m, output A(n+1)
B(n+1,:) = fftshift(fft(ifftshift(A(n+1,:))));
C(n+1,:) = abs(B(n+1,:)).^2;
rmse(n+1) = sqrt(mean((C(n+1,:)-C(n,:)).^2)); % calculate Root-mean-square error between C(n) and C(n+1)
fin = n+1; % store n value
if (rmse(fin) < 0.02) && (fin > 2) % check out steady state
flag = 1;
break
end
end
end
  1 Comment
ChiQAQ
ChiQAQ on 8 Nov 2021
Hi Norris. Thanks for your quick and professional answer. I will try parfeval as the example you demonstrated. Frankly speaking I'm not very familiar with parallel pool function. Let me learn more about it via Help Center first.
FYI,
  • Use profiler to measure time.
  • It takes 4 min to complete one internal for-loop consists of 2000 iterations. So the whole program would cost 40 min, counting the external loop. However, this is a simplified test code. The one I am handling includes at least m=1:380, and even need triple for-loop to address more parameters. Just considering the case of m=1:10, I prefer no more than 10 min, if possible.
  • i7-6700@3.4GHz, 4 cores, 8 threads, 16GB RAM. (Also have GPU but didn't use it)
  • You're right. I've tried to parallelize the outter loop. Below is the parfor loop I imitate from another case.
% create A, B, C, rmse, fin
A = zeros(N,nt,10);
A(1,:,:) = repmat(randn(1,nt)+1i*randn(1,nt),1,1,10);
B = fftshift(fft(ifftshift(A)));
C = abs(B).^2;
rmse = ones(10,1);
fin = zeros(10,1);
% parfor loop
parfor m = 1:10
Ap = A; % Cannot directly use A, cuz A is not recognized as sliced variable.
% Or error happens "Valid indices for 'A' are restricted in PARFOR loops"
Bp = B; % B, C, rmse, fin have same reason
Cp = C;
rmsep = rmse;
finp = fin;
for n = 1:N-1
Ap(n+1,:,m) = Function(Ap(n,:,m),m);
Bp(n+1,:,m) = fftshift(fft(ifftshift(Ap(n+1,:,m))));
Cp(n+1,:,m) = abs(Bp(n+1,:,m)).^2;
rmsep(m) = sqrt(mean((Cp(n+1,:,m)-Cp(n,:,m)).^2));
finp(m) = n+1;
if (rmsep(m) < 0.02) && (finp > 2)
break
end
end
temp{m} = rmsep;
end
for m=1:10
rmse(m) = temp{m};
end

Sign in to comment.

Community Treasure Hunt

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

Start Hunting!