//The following is C code that was used in a MatLab mex file to implement the neuron model.
//The differential equations can be solved by various means,
// where the input is STIM and the output is Y
//3D simplified HodgkinHuxley neuron model with slow Na inactivation gate
//Parameters to implement model in Figure 11 of Lundstrom et al., J Comp Neurosci
// Parameters:
cm = 1; // uF / cm2
Gl = 15; // mS/cm2
Gk = 36;
Gna = 50;
El = 54; // mV
Ek = 77;
Ena = 50;
//Differential equations
mVh = 40;
mk = 7;
nVh = 30;
nk = 5;
Cbase = 50;
Camp = 2000;
Vmax = 40;
sig = 5;
ntau= 3;
ninf = 1.0 / (1+exp( (nVh  y[0])/nk) );
minf = 1.0 / (1+exp( (mVh  y[0])/mk) );
h = .89  1.1*y[1];
s = .89  1.1*y[2];
//time constant for slow Na inactivation
stau= Cbase + Camp*exp(pow(Vmax  y[0],2)/pow(sig,2));
//stim is the variable that passes the stimulus value for each time t
dy[0] = 1*(Gl*(y[0]  El) + Gk*pow(y[1],4)*(y[0]  Ek) + Gna*pow(minf,3)*h*s*(y[0]  Ena))/cm + stim/cm;
dy[1] = (ninf  y[1]) / ntau;
dy[2] = (ninf  y[2]) / stau;
