function xprime = ring40(t,x)
%
% Compute diffusion reaction described in Turing's paper "The Chemical Basis
% if Morphogenesis" (p. 61)
%
%  Call:
%        load /home/vision/v1/matlab/neuron-modeling/ring_Values.mat
%        [t,x] = ode45('ring40',0,200,[ring_Values]);
%

fprintf(' %25.21f\n',t);
my_Eps = 0.01;                            % for robustness
cx = 5*10^(-4);                           % X diffusion constant
cy = 2.5*10^(-4);                         % Y diffusion constant
gamma = 0.1;                              % an arbitrary value in [-0.5, 0.5]
num_cells = length(x)/2;
X = x(1:2:2*num_cells-1);
Y = x(2:2:2*num_cells);

for index = 1 : num_cells      % 40 values in x, 20 cells; for goes from 1:20
           %Xi and Yi in cell i
  rate = (1/32)*(50*X(index)*Y(index) + 7*X(index)^2 - 55*(1 + gamma));
        %take diffusion from both of your neighboring cells
  indexf = index+1;
  if indexf>num_cells    % set forward index
    indexf = 1;
  end;
  indexb = index-1;      % set backward index
  if indexb==0
    indexb = 20;
  end;

  Xdiffuse = cx * ( X(indexf) + X(indexb) - 2*X(index) );
  Xprime(index) = -rate + (1/16) + Xdiffuse;
  if abs(Y(index)) > my_Eps
    Ydiffuse = cy * ( Y(indexf) + Y(indexb) - 2*Y(index) );
    Yprime(index) = rate - Y(index)/16 + Ydiffuse;
  end                %if x(yind) is 0, keep it at 0

end            %end for loop
xprime = zeros(1,2*num_cells);
xprime(1,1:2:2*num_cells-1) = X;
xprime(1,2:2:2*num_cells) = Y;
