MATLAB question regarding loop

1 visualización (últimos 30 días)
Chin Bei yi
Chin Bei yi el 4 de Mayo de 2018
Editada: Chin Bei yi el 6 de Mayo de 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

Respuesta aceptada

Image Analyst
Image Analyst el 4 de Mayo de 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);

Más respuestas (1)

Wick
Wick el 4 de Mayo de 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.

Categorías

Más información sobre Creating and Concatenating Matrices en Help Center y File Exchange.

Etiquetas

Aún no se han introducido etiquetas.

Community Treasure Hunt

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

Start Hunting!

Translated by