
clear;                  % always a good idea in Matlab 

%
%   Traveling wave in the Fitzhugh-Nagumo system
%

a=.7;  b=.8;  eps=0.08; % parameters of the FN-system
 
dt=.05;  dx=1;  D=3;    % time step, step in space, diffusion constant

v0=-1.2;  w0=-0.624;    % fixed point of the FN-system for I=0

id=200;                 % number of points in space

v(1:id)=v0; w(1:id)=w0; % all points in space at fixed point 

x=linspace(0,id,id);    % needed for plotting 

I=[0.3*ones(15,1)];     % current at first 10 points is 0.3
I(1:2)=0;  I(16:id)=0;  % current at points 1,2 and all others is 0

figure('Position',[50 450 900 250],'Color',[1 1 1]);  % make a figure
axes('Position',[.03 .1 .95 .88]);  % define the plot region

for t=0:dt:120                      % we calculate until t=120 
  if t>5 I(1:id)=0; end             % at t=5 current is switched off
  for ix=1:id                       % for every point in space
    if ix>1 & ix<id 
      vp(ix)=v(ix)-v(ix)^3/3-w(ix)+I(ix) ...    % the dynamics for v
           +D*(v(ix+1)+v(ix-1)-2*v(ix))/dx^2;   % with diffusion
    end
    wp(ix)=eps*(v(ix)+a-b*w(ix));   % the dynamics for w, no diffusion
  end;  

  vp(1)=0;  vp(id)=0;               % keep resting state at boundaries
                                    % comment out for periodic boundaries
                                    
%   vp(1)=v(1)-v(1)^3/3-w(1)+I(1) ...       % periodic boundaries
%       +D*(v(2)+v(id)-2*v(1))/dx^2;        % uncomment to see how 
%   vp(id)=v(id)-v(id)^3/3-w(id)+I(id) ...  % action potentials annihilate
%       +D*(v(1)+v(id-1)-2*v(id))/dx^2;     % comment out the line above

  v=v+vp*dt;  w=w+wp*dt;                    % update v and w
  plot(x,v,'r',x,w,'b','Linewidth',2);      % and plot both curves
  
  set(gca,'Linewidth',2,'Fontsize',14, ...  % make the graphics a
        'Fontweight','demi');               % little more facy
    
  axis([0 200 -2.2 2.2]);                   % same axis for all plots 
  drawnow;  pause(.01);                     % draw it and wait a moment
end                                         % that's it :)
