MATLAB question regarding loop

1 view (last 30 days)
Chin Bei yi
Chin Bei yi on 4 May 2018
Edited: Chin Bei yi on 6 May 2018
if true
clc;
clear;
FL=0.04; %Fin length
BL=0.005; %Base length
TL=BL+FL; %Total length
FW=0.010; %Fin width
BW=0.030; %Base width
dx=0.001;
dy=dx;
%Total Nodes
M=round(((BW/2)/dx+1)); %i Direction
N=round(((TL)/dx+1)); %j Direction
%Initial Temperature
T=zeros(M,N);
loop=10000; %Number of Iterations
err_max = 1e-6; %Maximum Error
n=0; %Number of Iterations
while n<loop
n=n+1;
Tp=T;
for i=1:1:M
for j=1:1:N
%Corner Node 8
if i==1 & j==1 %(i=0, j=0)
T(i,j)=(1/2)*(T(i,j+1)+T(i+1,j))+(7500/89)*dx;
%Corner Node 13
elseif i==M & j==1 %(i=15, j=0)
T(i,j)=(1/2)*(T(i,j+1)+T(i-1,j))+(7500/89)*dx;
%Corner Node 9
elseif i==1 & j==N %(i=0, j=45)
T(i,j)= (89*(T(i+1,j)+T(i,j-1))+4500*dx)/(150*dx+178);
%Constant heat flux boundary node 5
elseif i>1 & i<M & j==1 %(i=1 to 14, j=0)
T(i,j)=(1/4)*(T(i-1,j)+T(i+1,j))+(1/2)*(T(i,j+1))+(7500/89)*dx;
%Symmetrical boundary node 2
elseif i==1 & j>1 & j<N %(i=1 , j=1 to 45)
T(i,j)=(1/4)*(T(i,j+1)+T(i,j-1))+(1/2)*(T(i+1,j));
%Convective boundary node 7
elseif i>1 & i<6 & j==N %(i=1 to 4 , j= 45)
T(i,j)=((89)*(T(i-1,j)+T(i+1,j))+(178)*(T(i,j-1))+9000*dx)/(300*dx+356);
%Convective Corner Node 10
elseif i==6 & j==N %(i=5 , j= 45)
T(i,j)=((89)*(T(i-1,j)+T(i,j-1))+9000*dx)/(300*dx+178);
%Convective Boundary Node 3
elseif i==6 & j>5 & j<N %(i=5, j= 5 to 44)
T(i,j)=((89)*(T(i,j+1)+T(i,j-1))+178*T(i-1,j)+9000*dx)/(9000*dy+356);
%Corner Node 11
elseif i==6 & j==6 %(i=5, j=5)
T(i,j)=((178)*(T(i-1,j)+T(i,j-1)+T(i+1,j)+T(i,j+1))+4500*dy)/(150*dy+534);
%Corner Node 12
elseif i==M & j==6 %(i=15, j=5)
T(i,j)=((89)*(T(i-1,j)+T(i,j-1))+4500*dx)/(150*dx+178);
%Convective Boundary Node 6
elseif i>5 & i<M & j==6 %(i=5 to 15, j=5)
T(i,j)=((178)*(T(i,j-1))+(89)*(T(i+1,j)+T(i-1,j))+9000*dx)/(9000*dx+356);
%Symmetrical Boundary Node 4
elseif i==M & j>1 &j<6 %(i=15, j=1 to 4)
T(i,j)=(1/4)*(T(i,j-1)+T(i,j+1))+(1/2)*(T(i-1,j));
%Internal Node 1 (SOLID)
elseif i>1 & i<M & j>1 & j<6
T(i,j)=(1/4)*(T(i,j-1)+T(i,j+1)+T(i+1,j)+T(i-1,j));
%Internal Node 1 (SOLID)
elseif i>1 & i<6 & j>5 & j<N
T(i,j)=(1/4)*(T(i,j-1)+T(i,j+1)+T(i+1,j)+T(i-1,j));
end;
err = max(abs(T - Tp))
end
end
if err < err_max
break
end
end
end

Accepted Answer

Image Analyst
Image Analyst on 4 May 2018
Lots wrong with that code. Gladd to see that you're using a failsafe (limit to number of iterations though). Here, just try this. I fixed several things.
FL=0.04; %Fin length
BL=0.005; %Base length
TL=BL+FL; %Total length
FW=0.010; %Fin width
BW=0.030; %Base width
dx=0.001;
dy=dx;
%Total Nodes
M=round(((BW/2)/dx+1)); % i Direction
N=round(((TL)/dx+1)); % j Direction
%Initial Temperature
T=zeros(M,N);
maxLoopIterations = 10000; % Number of Iterations
err_max = 1e-6; % Maximum Error
n=0; %Number of Iterations
err = zeros(1, maxLoopIterations);
while n < maxLoopIterations
n=n+1;
Tp=T;
for i=1:1:M
for j=1:1:N
%Corner Node 8
if i==1 && j==1 %(i=0, j=0)
T(i,j)=(1/2)*(T(i,j+1)+T(i+1,j))+(7500/89)*dx;
%Corner Node 13
elseif i==M && j==1 %(i=15, j=0)
T(i,j)=(1/2)*(T(i,j+1)+T(i-1,j))+(7500/89)*dx;
%Corner Node 9
elseif i==1 && j==N %(i=0, j=45)
T(i,j)= (89*(T(i+1,j)+T(i,j-1))+4500*dx)/(150*dx+178);
%Constant heat flux boundary node 5
elseif i>1 && i<M && j==1 %(i=1 to 14, j=0)
T(i,j)=(1/4)*(T(i-1,j)+T(i+1,j))+(1/2)*(T(i,j+1))+(7500/89)*dx;
%Symmetrical boundary node 2
elseif i==1 && j>1 && j<N %(i=1 , j=1 to 45)
T(i,j)=(1/4)*(T(i,j+1)+T(i,j-1))+(1/2)*(T(i+1,j));
%Convective boundary node 7
elseif i>1 && i<6 && j==N %(i=1 to 4 , j= 45)
T(i,j)=((89)*(T(i-1,j)+T(i+1,j))+(178)*(T(i,j-1))+9000*dx)/(300*dx+356);
%Convective Corner Node 10
elseif i==6 && j==N %(i=5 , j= 45)
T(i,j)=((89)*(T(i-1,j)+T(i,j-1))+9000*dx)/(300*dx+178);
%Convective Boundary Node 3
elseif i==6 && j>5 && j<N %(i=5, j= 5 to 44)
T(i,j)=((89)*(T(i,j+1)+T(i,j-1))+178*T(i-1,j)+9000*dx)/(9000*dy+356);
%Corner Node 11
elseif i==6 && j==6 %(i=5, j=5)
T(i,j)=((178)*(T(i-1,j)+T(i,j-1)+T(i+1,j)+T(i,j+1))+4500*dy)/(150*dy+534);
%Corner Node 12
elseif i==M && j==6 %(i=15, j=5)
T(i,j)=((89)*(T(i-1,j)+T(i,j-1))+4500*dx)/(150*dx+178);
%Convective Boundary Node 6
elseif i>5 && i<M && j==6 %(i=5 to 15, j=5)
T(i,j)=((178)*(T(i,j-1))+(89)*(T(i+1,j)+T(i-1,j))+9000*dx)/(9000*dx+356);
%Symmetrical Boundary Node 4
elseif i==M && j>1 &&j<6 %(i=15, j=1 to 4)
T(i,j)=(1/4)*(T(i,j-1)+T(i,j+1))+(1/2)*(T(i-1,j));
%Internal Node 1 (SOLID)
elseif i>1 && i<M && j>1 && j<6
T(i,j)=(1/4)*(T(i,j-1)+T(i,j+1)+T(i+1,j)+T(i-1,j));
%Internal Node 1 (SOLID)
elseif i>1 && i<6 && j>5 && j<N
T(i,j)=(1/4)*(T(i,j-1)+T(i,j+1)+T(i+1,j)+T(i-1,j));
end
end
end
diffMatrix = abs(T - Tp);
err(n) = max(diffMatrix(:));
if err(n) < err_max
break
end
end
% Crop off unused part of array
err = err(1:n);
plot(err, 'b-', 'LineWidth', 2);
grid on;
drawnow;
fontSize = 20;
title('Error vs. Iteration Number', 'FontSize', fontSize);
xlabel('Iteration Number', 'FontSize', fontSize);
ylabel('Error', 'FontSize', fontSize);
xticks(0:100:n);
yticks(0:.01:.15);

More Answers (1)

Wick
Wick on 4 May 2018
In general, when inside a loop you can keep track of "something" simply by assigning that value to an element in a vector using the index of the loop to place the value. Example:
loop_length = 10000;
error_vector = zeros(1,loop_length)
for jj = 1:loop_length
error_vector(jj) = abs(10000 - jj);
end
Obviously, this is an inefficient way to go about defining the error_vector as I've written but it gives you an idea of how you might use the loop index to assign a value to a position in a vector.

Categories

Find more on Creating and Concatenating Matrices in Help Center and File Exchange

Tags

No tags entered yet.

Community Treasure Hunt

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

Start Hunting!