clear, close all;

%% Problem 3, Homework 3
% Solve
% u_t = C_0 u_{xx} +x  \qquad 0<x<L, t>0 \\
% u_x(0,t) = 0, \quad u(L,t) = 0, \quad t>0 \\
% u(x,0) = f(x) \quad 0<x<L
%
% Let u(x,t) = v(x,t) + U(x), where
% \frac{d^2 U}{dx^2} = -x/C_0, \\
% U'(0) = U(L) = 0
% which, integrating twice and applying the boundary conditions leads to U(x) = \frac{1}{6 C_0} (x^3-L^3). The system
% v_t = C_0 v_{xx}, \\
% v_x(0,t) = v(L,t) = 0, \\
% v(x,0) = f(x),
% is solved by separation of variables: v(x,t) = H(t)G(x),
% then
% H_t = -C_0 k^2 H \\
% G^{''} + k^2 G = 0, \qquad
% G'(0) = G(L) = 0.
% the G equation is a Sturm Liouville problem with eigenfuntions and eigenvalues
% \phi_n = \cos k_n x \\
% k_n = \frac{1}{L} (n+\frac{1}{2}) \pi, \quad n= 0, 1, 2, ...
% The H equation is a simple separable equation with a decaying solution. The general solution for v is
% v(x,t) = \sum_{n=0}^\infty a_n e^{-C_0 k_n^2 t} \phi_n(x)
% when we apply the initial conditions on v(x,t) and use orthogonality we obtain
% a_n = \frac{2}{L} \int_0^L f(x) \phi_n(x) dx, \quad n=0,1,2,...
%
% Suppose f(x) = 1-(x/L)^2
% , the solution
% u=v+U
% is


syms N An F Phi X L KN Cint
KN = (N+1/2)*pi/L
F = 1-(X/L)^2
Phi = cos(KN*X)
Cint = 2*F*Phi/L
An = int(Cint,X)

% a_n =\frac{2  \cos(n \pi)}{\pi^3 (8 n^3+12 n^2+6 n+1)}
% hence
% u = \sum_{n=0}^\infty \frac{2 \cos(n \pi)}{\pi^3 ( 8 n^3 + 12 n^2 + 6 n + 1)} \cos(k_n x) e^{-C_0 k_n^2 t} + \frac{1}{6 C_0} (x^3-L^3)


LL=1; dx = 0.05; x=0:dx:LL; 
t_0 = 0; t_f = 0.5; dt = t_f/100;
t=(t_0:dt:t_f)';
N = 10;

C_0 = 25;

U = (x.^3 - LL^3)/(6*C_0);

v = zeros(length(x),length(t));
for ii = 1:N+1
    nn=ii-1;
    kn = (nn+0.5)*pi/LL;
    an = 2*cos(pi*nn)/(pi^3*(8*nn^3+12*nn^2+6*nn+1));
    for jj = 1:length(t)
        for kk = 1:length(x)
    v(kk,jj) = v(kk,jj) + an*cos(kn*x(1,kk)) ...
        .*exp(-C_0*kn*t(jj,1))+U(1,kk);
        end
    end
end
[TT,XX]=meshgrid(t,x);
subplot(3,1,1), mesh(XX,TT,v)
subplot(3,1,2), plot(x,U)
        for kk = 1:length(x)
    v(kk,:) = v(kk,:) +U(1,kk);
        end
subplot(3,1,3), surf(XX,TT,v)