Why am i getting NaN values in my code?

2 次查看（过去 30 天）
smith
smith 2023-5-6
编辑： Torsten 2023-5-6
Ran in:
I'm writing a code that predicts the surface temperature of an annulus with specific thickness running through it a hot air and outisde it is under a water bath of a specific water temperature, the code should give me the way the temperature of the surface behaves across it's length.
% Input parameters
h_air = 500; % Heat transfer coefficient of hot air [W/m^2*K]
T_air_in = 600; % Inlet temperature of hot air [K]
T_water = 376; % Temperature of surrounding water [K]
k_cylinder = 50; % Thermal conductivity of cylinder material [W/m*K]
r = 0.05; % Cylinder radius [m]
L = 0.1; % Cylinder length [m]
dx = 0.005; % Grid size [m]
dt = 0.1; % Time step [s]
nt = 3600; % Simulation time [s]
d = 0.01; % Cylinder thickness [m]
h_water = 10000; % Heat transfer coefficient of surrounding water [W/m^2*K]
% Calculate grid dimensions
nx = round(L/dx) + 1;
% Set initial and boundary conditions
T = T_water * ones(nx, 1);
T(1) = T_air_in;
T(end) = T_water;
% Calculate constants for Euler's method
alpha = k_cylinder/(dx^2);
beta = h_air*dx/k_cylinder;
% Perform time integration using Euler's method
for i = 1:nt
T_old = T;
for j = 2:nx-1
% Calculate the thermal resistance of the pipe wall
r_pipe = log((r+d)/r)/(2*pi*k_cylinder*L);
% Calculate the thermal resistance of the fluid film
r_film = 1/(h_air*2*pi*r*dx);
% Calculate the heat transfer coefficient
h = 1/(r_pipe+r_film);
% Calculate the temperature gradient
dTdr = (T_old(j+1)-T_old(j-1))/(2*dx);
% Calculate the temperature at the current time step
T(j) = T_old(j) + alpha*dt*(T_old(j+1)-2*T_old(j)+T_old(j-1)) + h*2*pi*r*dt*(T_air_in-T_old(j)) + 2*pi*r*h_water*dt*(T_water-T_old(j));
end
end
% Plot temperature profile
x = linspace(0, L, nx);
plot(x, T);
xlabel('Position [m]');
ylabel('Temperature [K]');
title('Temperature profile of cylinder surface');
Upon running the code, you'll see that it won't show a curve nor values across the cylinder.
  1 个评论
Torsten
Torsten 2023-5-6
编辑：Torsten 2023-5-6
Your update
T(j) = T_old(j) + alpha*dt*(T_old(j+1)-2*T_old(j)+T_old(j-1)) + h*2*pi*r*dt*(T_air_in-T_old(j)) + 2*pi*r*h_water*dt*(T_water-T_old(j));
cannot be correct because the unit of e.g.
alpha*dt*(T_old(j+1)-2*T_old(j)+T_old(j-1))
is
W/(m^3*K) * s * K = J/m^3
but it should be
K
Similar problems with the other terms.

请先登录，再进行评论。

请先登录，再回答此问题。

回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Thermal Analysis 的更多信息

标签

产品


版本

R2019a

Community Treasure Hunt

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

Start Hunting!

Translated by