Having trouble modelling 2-dimensional transient heat transfer
5 views (last 30 days)
Show older comments
I'm having trouble modelling 2-dimensional transient heat transfer in MATLAB. I've reduced the problem down to a very simple model, simply an insulated box split into a 2x2 matrix with a heat generation term inside, but the values I'm getting out are way off what should be correct. The code makes use of the ODE45 command in MATLAB. This is the code I'm using, any help is much appreciated! The version of MATLAB that I'm using is 7.10 (R2010a)
This is the code I'm using:
Main file
%%Basic 2 dimensional heat transfer code
clear all
close
clc
global l h k mass spec_heat T_amb m n z
% Box parameters
l = 1; % length
h = 1; % height
k = 0.2; % W/m.K
mass = 100; % kg
spec_heat = 4000; % J/kg
% Initial conditions
T_amb = 288; % Kelvin
% Number of elements
m = 2; % x direction
n = 2; % y direction
z = m*n; % total number of elements
T0 = zeros(m*n, 1);
for i = 1:(m*n)
T0(i,1) = T_amb; % T0(i,1) is initial temp of the material in K
end
tdur = 3456000; % Duration of simulation in seconds (s).
tstep = 3600; % Time step (of 1 hour) in seconds.
tspan = [0:tstep:tdur];
[t, T] = ode45(@twod_function, tspan, T0);
plot (t,T(:,1),'-r',t,T(:,2),'b',t,T(:,3),'-y',t,T(:,4),'--')
And for the function file
function dTdt = twod_function(t,T)
global l h k mass spec_heat T_amb m n z
R = (l/m)/(k*l*h); % resistance of material in box
Q_gen = 50;
Q_cond = zeros(1,(m*n));
for i = 1:(m*n)
if i == 1
Q_cond(i) = (T(i) - 288)/R;
else
Q_cond(i) = (T(i)-T(i-1))/R; % heat transfer by conduction
end
end
dTdt = zeros(4,1);
for i = 1:(m*n)
if i == 1
dTdt(i,1) = (Q_gen - Q_cond(i+1) - Q_cond(i+2))/((1/m)*mass*spec_heat);
elseif i == 2
dTdt(i,1) = (Q_cond(i-1) - Q_cond(i+2))/((1/m)*mass*spec_heat);
elseif i == 3
dTdt(i,1) = (Q_cond(i-2) - Q_cond(i+1))/((1/m)*mass*spec_heat);
else
dTdt(i,1) = (Q_cond(i-2) + Q_cond(i-1))/((1/m)*mass*spec_heat);
end
end
end
2 Comments
Youssef Khmou
on 12 Feb 2013
Hi , i tried your code, but it gives the same oscillatory Temp in the end, i think you should work on a 2d matrix , can you ?
Youssef Khmou
on 12 Feb 2013
try this as preliminary test :
surf(10*log10(abs(T)))
>> shading interp
Answers (1)
Youssef Khmou
on 12 Feb 2013
Edited: Youssef Khmou
on 12 Feb 2013
hi, please, change the duration as :
tspan=[0:4:7200];
and after running the file, type in the command prmpt :
surf(T);
Does that mesh make sens as Decreasing Gradient of Temperature ??
3 Comments
Youssef Khmou
on 13 Feb 2013
hi, Mark,
yes you try per example initializing a 2D T matrix , and put the equation you want to solve as Laplace(T) = constant ( just as example, you put those convection and conduction equations in general form ) , next it will be easy to write a function that differentiate the 2D T .
See Also
Products
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!