% ROD2.M solves PDE of heat transmission along insulated rod.
% This M-file uses sparse-matrix facilities.
% See book, Section 17.3.
k = 0.835; L = 10; Nx = 50;
x = linspace(0, L, Nx);
dt = 0.01;
dx = L/(Nx - 1);
la = k*dt/(dx*dx); % coefficient lambda
% build matrix T
e = ones(Nx, 1);
A = la*e; D = (1 - 2*la)*e; B = la*e;
T = spdiags([A, D, B], -1:1, Nx, Nx);
T(1, 1) = 1; T(1, 2) = 0;
T(Nx, Nx - 1) = 0;
% define initial conditions
T(Nx, Nx) = 1;
U = zeros(Nx, 1);
U(1) = 50; % temperature at hot end, deg C
% simulate evolution of rod temperature
for k = 1:501
U(:, k+1) = T*U(:, k);
end
% plot results
plot(x, U(:, 101), x, U(:, 201), x, U(:, 301),...
x, U(:, 401), x, U(:, 501))
xlabel('x')
ylabel('Temperature, deg C')
for m = 1:5
text(x(10+2*m), U(10+2*m, (1+100*m)), num2str(m))
end
