function nagumo_pde
%

   
m = 0;
x = linspace(0,100,200);
t = linspace(0,50,100);

sol = pdepe(m,@pr,@ic,@bc,x,t);
u = sol(:,:,1);
v = sol(:,:,2);
% 
% figure(1);
% surf(x,t,u);
% title('u(x,t)');
% xlabel('Distance x');
% ylabel('Time t');
% 
% figure(2);
% surf(x,t,v);
% title('v(x,t)');
% xlabel('Distance x');
% ylabel('Time t');

figure(3)
for i=1:length(t)
    plot(x,u(i,:),'g',x,v(i,:),'r');
    title(t(i))
    legend('u','v')
  %  ylim([-1.2 1.2])
    drawnow;
    
end
% --------------------------------------------------------------------------

function [c,f,s] = pr(x,t,u,DuDx)
%   The PDEs are
%   
%      D(u)/Dt = D^2(u)/Dx^2 +u(a-u)(u-1)-v
%      D(v)/Dt = bu - gv 
% 
%
%   In the form expected by PDEPE, the equations are
%
%   |1|         |u|      | D(u)/Dx |    | u(a-u)(u-1)-v |
%   | | .*  D_  | | = D_ |         | +  |               |
%   |1|     Dt  |v|   Dx |   0     |    | bu - gv       |
%
%   ---         ---   ----------------  ---------------
%    c           u      f(x,t,u,Du/Dx)   s(x,t,u,Du/Dx)
%
a=6;
b=12;
g=2;
c = [1; 1];                                  
f = [1; 0] .* DuDx;  
s1= u(1)*(a-u(1))*(u(1)-1)  -   u(2);
s2= b*u(1) - g*u(2);
s = [s1; s2]; %[0 ;0];%                                

% --------------------------------------------------------------------------

function u0 = ic(x)
%   The initial condition is u(x,0) = 0 and v(x,0) = 0 for 0 <= x <= 1.
u0 = [0; 0]; 
if (x>49) & (x<51)
u0=[5;0];
end
% --------------------------------------------------------------------------

function [pl,ql,pr,qr] = bc(xl,ul,xr,ur,t)
%   The left boundary condition is 
%      |0   |      |1|     | D(u)/Dx |   |0|
%      |    |  +   | | .*  |         | = | |
%      |0   |      |1|     | D(v)/Dx |   |0|
%
%      ---        ---     ------------------   ---
%    p(0,t,u)    q(0,t)     f(0,t,u,Du/Dx)      0
%

pl = [ur(1);ur(2)];  

ql = [-10; -10];  


%   The right boundary condition is
%
%      |u |       |0|     | D(u)/Dx |   |0|
%      |  |   +   | | .*  |         | = | |
%      |v |       |0|     | D(v)/Dx |   |0|
%
%      -------       -----   ----------------   ---
%      p(1,t,u)      q(1,t)  f(1,t,u,Du/Dx)      0
%
                                
pr = [ur(1); ur(2)];                            
qr = [0; 0];                                  

