function xprime = FHN(t,x)
global eps R E I_in
xprime=zeros(2,1);
V_C=x(1);I_L=x(2);
xprime(1)=-(F(E+V_C)+I_L+I_in)/eps; % inward forcing current
xprime(2)=(V_C-R*I_L);