%% Värmeledningsekvationen u_t = u_xx för 0<=x<=1, t>0

T = 2;            % sluttid

n = 81;             % antal punkter i x
m = 81;             % antal punkter i t

h = 1/(n-1);        % steglängd i x
x = 0:h:1;          % diskretisering av x
t = linspace(0,T,m);% diskretisering av t

%% Skapa finit differensmatris
e = ones(n,1);
D2 = spdiags(1/h^2*[e -2*e e],-1:1,n,n);
D2 = D2(2:n-1,:);   % kapa bort första och sista raden

%% Randvillkor du/dt = g(t) för x = 0 och x = 1
g0 = @(t) -500*exp(-5*t);
g1 = @(t) 0;

%% Begynnelsevillkor (ej kontinuerligt!)
i = find(x<=0.5);
u0 = zeros(n,1);
u0(i) = 100;

%% ODE-system du/dt = F(t,u)
F = @(t,u) [g0(t); D2*u; g1(t)];

%% Lösning med inbyggd lösare
tic
[t,U] = ode45(F,t,u0);
toc

%% Plott av resultatet
mesh(x,t,U)
xlabel('x'), ylabel('t'), zlabel('u')