Simulation of the Hodgkin Huxley Model

The system of four differential equations in the Hodgkin-Huxley model is: $$\begin{aligned}C_m\frac{\dd V_m}{\dd t} &= - {x_{\rm{K}}}\overline {{g_{\rm{K}}}} {n^4}({V_m} - {E_{\rm{K}}}) - {x_{{\rm{Na}}}}\overline {{g_{{\rm{Na}}}}} {m^3}h({V_m} - {E_{{\rm{Na}}}}) - \overline {{g_{\rm{L}}}} ({V_m} - {E_{\rm{L}}}) + I\\[0.25cm] \frac{\dd n}{\dd t}&=\alpha_n(V_m)(1-n)-\beta_n(V_m)n\\[0.25cm] \frac{\dd m}{\dd t}&=\alpha_m(V_m)(1-m)-\beta_m(V_m)m \\[0.25cm] \frac{\dd h}{\dd t}&=\alpha_h(V_m)(1-h)-\beta_h(V_m)h \end{aligned}$$ where $V_m$ is the membrane potential, $E_X$ is the Nernst potential for a given ion $X$ or the leakage potential (for $X=L$), and the gate functions $n$, $m$, $h$ satisfy the following initial value problems of exponentially restricted growth equations with voltage-dependent coefficients:
$$\begin{aligned} \alpha_n(v) &= \frac{0.01\bigl(10-(v-V_r)\bigr)}{\exp\left(\frac{10-(v-V_r)}{10}\right)-1}\\[0.25cm] \beta_n(v)&=0.125\exp\left(\frac{-(v-V_r)}{80}\right)\\[0.25cm] \alpha_m(v)&= \frac{0.01\bigl(25-(v-V_r)\bigr)}{\exp\left(\frac{25-(v-V_r)}{10}\right)-1}\\[0.25cm] \beta_m(v)&=4\exp\left(\frac{-(v-V_r)}{18}\right)\\[0.25cm] \alpha_h(v)&=0.07\exp\left(\frac{-(v-V_r)}{20}\right)\\[0.25cm] \beta_h(v) &= \frac{1}{\exp\left(\frac{30-(v-V_r)}{10}\right)+1}\\[0.25cm] n(0) &= 0.32\\[0.25cm] m(0)&=0.06\\[0.25cm] h(0) &= 0.6 \end{aligned}$$ We have also introduced blocking of potassium and sodium channels by adding to the original Hodgkin-Huxley model parameters $x_\mathrm{K}$ and $x_\mathrm{Na}$, defined as the fractions of active ion channels (in terms of percentages $b_\mathrm{K}$ and $b_\mathrm{Na}$ of blocked ion channels): $$x_\mathrm{K}=1-\frac{b_\mathrm{K}}{100}\qquad\text{and}\qquad x_\mathrm{Na}=1-\frac{b_\mathrm{Na}}{100}$$
Blocking of ion channels often is a result of neurotoxines: Tetrodoxin (TTX), isolated from the Japanese pufferfish fugu, blocks sodium channels, and tetraethylammonium (TEA) blocks potassium channels.