% Usage: [u x t] = heat(f,T,n,nt) or u = heat(f,T,n,nt)
% Heat equation on [0,1] with finite differencing
%
% Input:
% f - Matlab inline function f(x) for the initial condition
% T - final time moment
% n - number of subintervals (panels)
% nt - number of time intervals
%
% Output:
% u - computed solution
% x - nodes of the spacial grid
% t - time steps
%
% Examples:
% [u x t] = heat(@bump,1,100,100);          here 'bump' is a user-defined function in M-file
% [u x t] = heat(@bump,1,100,80);
% [u x t] = heat(@bump,1,100,99);

function [u x t] = heat(f,T,n,nt);
x = (0:n)/n;
t = (0:nt)*T/nt;
u = zeros(n+1,nt+1);
h = 1/n;
tau = T/nt;
c = tau/h^2;
c1=1-2*c;
for i=1:n+1
    u(i,1) = f((i-1)*h);
end;
b(1)=u(1,1);
b(2)=u(n+1,1);
for k=1:nt-1
    for i=2:n
        u(i,k+1) = c1 * u(i,k) + c * (u(i+1,k)+u(i-1,k));
    end;
    u(1,k+1) = b(1);
    u(n+1,k+1) = b(2);
end;