MATLAB Answers

How to replace two for loops with matrix expression

11 views (last 30 days)
Fei Lai
Fei Lai on 18 May 2016
Answered: Nobel Mondal on 18 May 2016
the original code is like this:
%
for i=2:nhx-1
for j=2:nhy-1
Unew(i,j)= U(i,j)-dt*(P(i+1,j)-P(i-1,j))/(2*hx)...
+nu*dt*(1/(hx*hx)*(U(i+1,j)-2.*U(i,j)+U(i-1,j))...
+1/(hy*hy)*(U(i,j+1)-2.*U(i,j)+U(i,j-1)))...
-dt*U(i,j)/(hx)*(U(i,j)-U(i-1,j))...
-dt*V(i,j)/(2*hy)*(U(i,j+1)-U(i,j-1));
end
end
how can I replace two for loops to increase the calculation speed?
  1 Comment
Todd Leonhardt
Todd Leonhardt on 18 May 2016
Edit your post to format the code by indenting all code by at least 2 spaces. That will make it more readable.

Sign in to comment.

Answers (1)

Nobel Mondal
Nobel Mondal on 18 May 2016
You may directly use "vectorization" and "element-wise operation to solve this :
Here is an example code -
%%Assumptions - as I don't have the actual data
U = magic(100);
P = magic(length(U));
V = magic(length(U));
dt = 0.1;
hx = rand;
hy = rand;
nu = rand;
nhx = size(U,1);
nhy = size(U,2);
Unew = zeros(nhx, nhy);
%%Old algo
for i=2:nhx-1
for j=2:nhy-1
Unew(i,j)= U(i,j)-dt*(P(i+1,j)-P(i-1,j))/(2*hx)...
+nu*dt*(1/(hx*hx)*(U(i+1,j)-2.*U(i,j)+U(i-1,j))...
+1/(hy*hy)*(U(i,j+1)-2.*U(i,j)+U(i,j-1)))...
-dt*U(i,j)/(hx)*(U(i,j)-U(i-1,j))...
-dt*V(i,j)/(2*hy)*(U(i,j+1)-U(i,j-1));
end
end
%%Vectorized algo
myUnew = zeros(nhx, nhy);
myUnew(2:nhx-1, 2:nhy-1)= U(2:nhx-1, 2:nhy-1)-dt.*(P(3:nhx,2:nhy-1)-P(1:nhx-2,2:nhy-1))./(2*hx)...
+nu.*dt.*(1/(hx*hx).*(U(3:nhx,2:nhy-1)-2.*U(2:nhx-1, 2:nhy-1)+U(1:nhx-2,2:nhy-1))...
+1/(hy*hy).*(U(2:nhx-1,3:nhy)-2.*U(2:nhx-1, 2:nhy-1)+U(2:nhx-1,1:nhy-2)))...
-dt.*U(2:nhx-1, 2:nhy-1)/(hx).*(U(2:nhx-1, 2:nhy-1)-U(1:nhx-2,2:nhy-1))...
-dt.*V(2:nhx-1, 2:nhy-1)/(2*hy).*(U(2:nhx-1,3:nhy)-U(2:nhx-1,1:nhy-2));
%%Now verify functionality
success = isequal(Unew, myUnew);
if success
disp('Test passed')
else
disp('Test failed. Please check the algorithm.')
end

Tags

No tags entered yet.

Products

Community Treasure Hunt

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

Start Hunting!