Why do I have a never ending loop ?

3 ビュー (過去 30 日間)
Jalyn-Rose
Jalyn-Rose 2024 年 3 月 8 日
コメント済み: Torsten 2024 年 3 月 9 日
Not sure why my loop continues? Solving a diffusion equation for context. Continues to print out time but should stop when conditions are met.
% Define constraints
alpha = 0.001; % cm^2/s
height = 2; % cm
dx = 0.1; % cm
dt = 0.1; % s
% Boundary conditions
Tbottom = 130; % C degrees
Tcondition=60; %degrees C
% Grid spacing
npoints = int32(height/dx)+1;
% Mesh
T=ones(1,npoints)*4.0; %initial temp is 4
disp(T)
time=0;
condition=true;
while condition
time=time+dt;
fprintf('Current Time is: %f\n', time);
Tnew=ones(1,npoints);
T(1)=Tbottom;
for i=2:npoints-1
Tnew(i)=(T(i))+alpha*dt/dx/dx * (T(i+1)+T(i-1)-2*T(i));
end
i=npoints-1;
Tnew(i)=T(i)-alpha*dt/dx/dx * (T(i)-T(i-1));
Tmin=inf;
for i=1:npoints
T(i)=Tnew(i);
if T(i)<Tmin
Tmin=T(i);
end
end
if Tmin>Tcondition
condition=false;
end
end
disp(T);
  3 件のコメント
Jalyn-Rose
Jalyn-Rose 2024 年 3 月 8 日
I used the x=height for my grid spacing not for the boundary. My only boundaries are the T=130 and 60. T=4 is the initial temp across right side of my diffusion.
Torsten
Torsten 2024 年 3 月 8 日
To answer your question:
Your while loop never stops because Tnew(1) is always 1 from the setting Tnew=ones(1,npoints);, thus < Tcondition.

サインインしてコメントする。

採用された回答

Torsten
Torsten 2024 年 3 月 8 日
編集済み: Torsten 2024 年 3 月 8 日
The initial condition at t = 0 seems to be T = 4 for all x in your code.
But you must set boundary conditions for T at x = 0 and x = 2.
At x = 0, this seems to be T = 130 for all t. But what about x = 2 ?
% Define constraints
alpha = 0.001; % cm^2/s
height = 2; % cm
dx = 0.1; % cm
nx = height/dx+1;
% Boundary conditions
Tbottom = 130; % C degrees
Tcondition=60; %degrees C
dt = 0.1; % s
nt = 5000;
T = 4*ones(nx,1);
T(1) = Tbottom;
Tnew = zeros(size(T));
condition = true;
counter = 0;
hold on
while condition
counter = counter + 1;
Tnew(1) = T(1); % T = 130 at x = 0
Tnew(2:nx-1) = T(2:nx-1)+dt*alpha/dx^2*(T(1:nx-2)-2*T(2:nx-1)+T(3:nx));
Tnew(nx) = Tnew(nx-1); %dT/dx = 0 at x = 2
T = Tnew;
condition = min(T) < Tcondition;
if counter == 500
counter = 0;
plot((0:nx-1)*dx,T)
end
end
hold off
grid on
  7 件のコメント
Jalyn-Rose
Jalyn-Rose 2024 年 3 月 8 日
Im uncertain of how I would see the time?
Torsten
Torsten 2024 年 3 月 9 日
time = 0;
while ...
time = time + dt;
...
end

サインインしてコメントする。

その他の回答 (0 件)

カテゴリ

Help Center および File ExchangeLoops and Conditional Statements についてさらに検索

Community Treasure Hunt

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

Start Hunting!

Translated by