1D Heat Equation Explicit Scheme (fixed)

9 次查看(过去 30 天)
clear all
p=2050; % density(kg/m^3)
k=0.52; % conductivity(W/m-K)
c=1840; % specific heat capacity(J/kg-K)
T_ini=20; % initial temperature(K)
T_inf=-15; % external temperature(K)
alpha=k/(p*c); % Let thermal diffusivity be alpha
t=60*60*60*24; %time in seconds
DELTA_t=100; % change in time per timestep
M=t/DELTA_t; % number of timesteps
N=300; % number of steps of depth
L=3; % maximum depth being analyzed
DELTA_x=L/(N); % change in depth per step
t_diff=t/M;
time=0:t_diff:t;
%explicit scheme
% Boundary Conditions
T(1:N,1:M+1)=T_ini;
T(1,1:M+1)=T_inf;
CFL=(alpha*DELTA_t)/(DELTA_x^2); % Courant number, which should be less than 0.5
% Finding temperature for next timestep for each depth
for j=1:M
for i=2:N-1
T(i,j+1)=T(i,j)+((alpha*DELTA_t)/((DELTA_x)^2))*(T(i+1,j)-2*T(i,j)+T(i-1,j));
end
end
Z=N*(0.1/L); % t=0.1 seconds in terms of the time step
xm=(69/300)*L;
plot(time,T(Z,:), 'b',time,T(69,:), 'r')
xlabel('Time')
ylabel('Temperature')
grid on
set(gca,'fontsize',16)
legend('0.1','xm')

采纳的回答

Torsten
Torsten 2022-4-19
You use
for i=0:0.1:50
...
T(i,j) = ...
But only positive integer values are possible as array indices. Thus something like T(0,j), T(0.1,j) etc is not allowed.
Additionally, you need a boundary condition at x=xstart and x=xend. What are these conditions ?
Delta_x and Delta_t are scalar values. Thus addressing them in the if-statement
if (alpha*DELTA_t(j))/(DELTA_x(i)^2)<=0.5
as vector is wrong.

更多回答(1 个)

Alan Stevens
Alan Stevens 2022-4-19
With
for i=0:0.1:50
you are trying to index variables with zero. However, Matlab's indexing starts at 1.

类别

Help CenterFile Exchange 中查找有关 Heat and Mass Transfer 的更多信息

产品


版本

R2021b

Community Treasure Hunt

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

Start Hunting!

Translated by