function [out1, out2, out3] = FHN2(t,x,flag)
global eps R E I_in
V_C=x(1);
I_L=x(2);
if nargin < 3 | isempty(flag)
   out1=[-(F(E+V_C)+I_L+I_in)/eps; 
        (V_C-R*I_L)];% inward forcing current
else
   switch(flag)
    case 'init'
      out1 = [0; 50];
      out2 = [.5,.5];
      out3 = []; 
     case 'events' 
     sec1 = -(F(E+V_C)+I_L+I_in);
     out1 = [sec1];
     out2 = [0];
     out3 = [0];
   otherwise
      error(['Unknown flag ''' flag ''.']);
   end
end; 
