%% problem description
% Determine the temperature vs time for hot iron with convection and
% radiation.
HW3_prob4_driver()
%% driver function (do not change!)
function HW3_prob4_driver
% to run: click the green "Run" button under the "EDITOR" tab
% to publish: click the "Publish" button under the "PUBLISH" tab
clc
close all
h = 25; %step size on time (seconds)
t = 0:h:400; % time over constant (h)
T = HW3_prob4_SOL(t,h);
% plot the numerical solution
figure
plot(t,T,'-o','MarkerSize',2)
xlabel('time [sec]')
ylabel('temperature [K]')
end
%% solution function(s)
function T = HW3_prob4_SOL(t,h)
% Your entire solution should be contained in this m-file. If you want to
% call additional functions, add them to the bottom of this script. hint:
% "Euler_ODE" and "the_function" will be VERY similar to what we wrote in
% class to solve ODEs using Euler's method.
T = Euler_ODE(t,h);
end
function T = Euler_ODE(t,h)
T = zeros(1,length(t)); %initialize the solution vector
T(1) = 500; % initial temperature at t=0 (K)
for i = 1:length(t)-1
T(i+1) = T(i)+the_function(t(i),T(i))*h; %Euler's method formula
end
end
%% y' or f(x,y) for diff eqns of form dy/dx = f(x,y)
function dTdt = the_function(t,T)
Cp = 4.*(10.^5); % Volumetric Heat Capacity (J/m^3K)
d = 0.005; % iron plate's thickness (m)
H = 20; % Convection coefficient (W/m^2K)
T_inf = 300; % Surrounding Air temperature (K)
eps = 0.5; % iron's emissivity
sig = 5.67.*(10.^-8); % constant
T_surr = 300; % Surrounding Air Temperature (K)
dTdt = (-H.*(T-T_inf)-eps.*sig.*(T.^4-(T_surr.^4)))/(Cp.*d);
end
