function Yp = F(Y,Iinj)
% defines the Hodgkin-Huxley model
 
V=Y(1);
n=Y(2);
m=Y(3);
h=Y(4);

Gna=120;
Ena=115;

Gk=36;
Ek=-12;
Cm=1;
Gm=0.3;
Vr=10.613;
 
alpham=(25-V)/(10*(exp((25-V)/10)-1));
betam=4*exp(-V/18);
mp=alpham*(1-m)-betam*m;
 
alphah=0.07*exp(-V/20);
betah=1/(exp((30-V)/10)+1);
hp=alphah*(1-h)-betah*h;

alphan=(10-V)/(100*(exp((10-V)/10)-1));
betan=0.125*exp(-V/80);
ninfinity=alphan/(alphan+betan);
taun=1/(alphan+betan);
np=(ninfinity-n)/taun;

Vp=1/Cm*(Gna*m^3*h*(Ena-V)+Gk*n^4*(Ek-V)+Gm*(Vr-V)+Iinj);

Yp(1)=Vp;
Yp(2)=np;
Yp(3)=mp;
Yp(4)=hp;
