k = log(4/3); Te = 22;
yp = -k*(y - Te);
% This is the RHS of the DE
end
Note that MatLab uses % beginning any comments that are not executed.
The next step is to use MatLab’s numerical ODE solver to create the solution of the differential
equation from the time of death,
t
d
, found earlier for about 12 hours after the murder. Below we
show the solving of the model and plotting of the graph.
[t1,y1] = ode23(@murder,[0,td],30);
[t2,y2] = ode23(@murder,[0,12],30);
plot(t1,y1,’b-’,t2,y2,’b-’);grid;
This produces a blue solution from the time of death until 12 hours after
t
= 0. This is quick way
to produce a graph, but we want to illustrate the power of MatLab to produce a very good looking,
professional graph.
Below we present an M-file,
Script
, which executes a series of commands to produce a good
graph. There are comments interior to the program to help understand what is happening, and
clearly one could explore other options to improve the graph.
clear
% Clear previous definitions
figure(1)
% Assign figure number
clf
% Clear previous figures
hold off
% Start with fresh graph
mytitle = ’Body Temperature’;
% Title
xlab = ’$t$’;
% X-label
ylab = ’$T(t)$’;
% Y-label
f = @(k) 22+8*exp(-k)-28;
% Function for finding heat coefficient
k = fzero(f,0.3);
% Find the heat coefficient
ft = @(t) 22+8*exp(-k*t)-37; % Function for finding the time of death
td = fzero(ft,-5);
% Find the time of death
[t1,y1] = ode23(@murder,[0,td],30); % Simulate the heat equation
[t2,y2] = ode23(@murder,[0,12],30); % Simulate the heat equation
plot(t1,y1,’k-’,’LineWidth’,1.5);
% Plot from death to t = 0
hold on
% Plots Multiple graphs
plot(t2,y2,’k-’,’LineWidth’,1.5);
% Plot from t = 0 to 12
plot([-9,td],[37,37],’r-’,’LineWidth’,1.5);
% Plot from t = -9 to td
plot([-9,12],[22,22],’b:’,’LineWidth’,1.5);
% Plot room temperature
plot([td,td],[20,40],’k:’,’LineWidth’,1.5);
% Plot time of murder
grid
% Adds Gridlines
text(-2.7,28,’Time of death, $t_d$’,’rot’,90,’FontSize’,14,...