Why do I have a never ending loop ?

2 views (last 30 days)
Jalyn-Rose
Jalyn-Rose on 8 Mar 2024
Commented: Torsten on 9 Mar 2024
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 Comments
Jalyn-Rose
Jalyn-Rose on 8 Mar 2024
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 on 8 Mar 2024
To answer your question:
Your while loop never stops because Tnew(1) is always 1 from the setting Tnew=ones(1,npoints);, thus < Tcondition.

Sign in to comment.

Accepted Answer

Torsten
Torsten on 8 Mar 2024
Edited: Torsten on 8 Mar 2024
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 Comments
Jalyn-Rose
Jalyn-Rose on 8 Mar 2024
Im uncertain of how I would see the time?
Torsten
Torsten on 9 Mar 2024
time = 0;
while ...
time = time + dt;
...
end

Sign in to comment.

More Answers (0)

Categories

Find more on Programming in Help Center and File Exchange

Community Treasure Hunt

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

Start Hunting!