Stuck at a For Loop
4 views (last 30 days)
Show older comments
Hi, good day! I am having a problem with my coding.
I want to develop a code for Neumann boundary condition, where one edge of a heat plate is insulated, using Liebmann (or Jacobi?) Method.
However, my final answer returns the original matrix T (representing the temperatures), which means my iteration failed, but I don't know why...
Appreciate much for any guidance! Thanks!
% Overrelaxation, lambda
lam = input('Please enter the overrelaxation value, lambda: ');
e_s = input('Please enter the stopping criterion in %: ');
err = e_s/100;
disp('The unit of k value is suggested to be in cal/s.cm.Celcius')
disp('As long as the unit for dimensions is consistent')
k = input('Please enter the k value of the plate material: ');
disp('The dimensions are suggested to be in cm')
l = input('Please enter the length of the heated plate: ');
w = input('Please enter the width of the heated plate: ');
pt = input('Please enter the desired number of points per row within the plate: ');
% Calculating the temperatures
delx = w/(pt + 1);
dely = l/(pt + 1);
T = zeros(pt + 2);
T(:,1) = 75;
T(:, (pt + 2)) = 50;
T(pt + 2, (2:pt + 1)) = 100;
disp('T_up + T_down + T_left + T_right = 4*T_current')
disp('At insulated edge, T_up is not available, use fictitious point')
disp('Hence the equation becomes:')
disp('2*T_down + 2*dely*dTdy + T_left + T_right = 4*T_current')
disp('Since dTdy = 0 at insulated edge, the equation becomes:')
disp('T_current = (2*T_down + T_left + T_right)/4')
A = T;
e1 = ones((pt + 1), pt);
while max(max(e1)) > err
% At insulated edge
for r = 1
for c = 2:(pt + 1)
A(r,c) = (T(r,(c-1)) + T(r, (c+1)) + (2*T((r+1),c)))/4;
end
end
% At other points
for r = 2:(pt + 1)
for c = 2:(pt + 1)
A(r,c) = (T(r,(c-1)) + T(r,(c+1)) + T((r-1),c) + T((r+1),c))/4;
end
end
T(r,c) = (lam*A(r,c)) + ((1-lam)*T(r,c));
e1(r,(c-1)) = abs((T(r,c) - A(r,c))/T(r,c));
end
disp('The temperature distribution on the plate is:')
T_plate = flipud(T);
disp(T_plate)
2 Comments
Answers (1)
per isakson
on 27 May 2019
Edited: per isakson
on 27 May 2019
Notes from my debugging session:
- added your assignments of input values to the top of the script
- commented out the input() statements
- ran the script
- the while-loop looped forever
- the following two line seems to be in the wrong place
T(r,c) = (lam*A(r,c)) + ((1-lam)*T(r,c));
e1(r,(c-1)) = abs((T(r,c) - A(r,c))/T(r,c));
- using a loop counter, in this case r and c, outside the loop is allowed in Matlab, but it smells
- the values of r and c, don't change. ( pt doesn't change value inside the while-loop.) At these two lines r and c, holds the end-values of the previous for-loops, respectively.
- thus only the lower right values of T and e are changed, which explains why the while condition never takes the value false
Your turn
0 Comments
See Also
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!