% This file solves the Fisher's equation
% on a bounded one-dimensional domain with homogeneous 
%Dirichlet boundary conditions

set(0,                           ...
   'defaultaxesfontsize', 20,   ...
   'defaultaxeslinewidth', 1.0, ...
   'defaultlinelinewidth', 1.2, ...
   'defaultpatchlinewidth', 0.7);
clear
 

Llist = [2;5];
d = 1; % diffusion coefficient = 1 wlog
 
N=199; % number of spatial grid points
M = 50; % number of time steps

for j = 1:2
L = Llist(j);
h = L/(N+1);
dt = h;
 scal = d*dt/h^2;
 X = h*(0:N+1)';
X0 = h*(1:N);
    cse = j;
    switch cse % specify initial data
        case 1

            
        V =  ( tanh((X0' )/(L/25))).*( tanh(-(X0'- L)/(L/25)));
        case 2
        V =  (1+tanh((X0'-L/2)/(L/25))).*(1+tanh(-(X0'-L/2)/(L/25)))/16;
 
    end
    


 figure(j)
  plot(X,[0;V;0],'--','linewidth',3)
  axis([0 L  0 1])
  hold on
kstep = 10;
t=0;

% uses Crank Nicolson for the diffusion part
% and explicit Euler for reaction part
Atm=(1+scal)*ones(N,1);
%Atm(1,1) = 1+scal/2;  these give Neumann conditions
%Atm(N,1) = 1+scal/2;

Am=diag(Atm)-diag(scal/2*ones(N-1,1),1)-diag(scal/2*ones(N-1,1),-1);

% The algorithm
  
for n=2:M
for k = 1:kstep
 %F=V.*(1-V).*(V-0.1);  % bistable rhs
 %F =  V.*(1-V)-V.^2./(R*(ep^2+V.^2)); % spruce budworm rhs
 F = V.*(1-V);  % Fisher equation
  
%V_t =  scal/2 * ([V(2:N );V(N )] - 2*V + [V(1 );V(1:N-1 )])  ;  Neumann
%conditions
Vt =  scal/2 * ([V(2:N );0] - 2*V + [0;V(1:N-1 )]) +dt*F  ;
V = Am\(V+Vt);
t=t + dt;
end
figure(j)
  plot(X,[0;V;0])
  axis([0 L  0 1])
 % axis([0 L -.5 1])
xlabel('\xi','fontsize',20)
ylabel('v(\xi,\tau)','fontsize',20)
   hold on

end

end

 


