SciComp@UNC Courseware

Interactive demo:

Initialization code

(* Ion names *)
ions={"(SuperPlus[Na])","(SuperPlus[K])","(SuperMinus[Cl])"};

(* Physical constants.[J/K/mol,K,C] *)
{R,T,F}={8.314,280.,96485.};

(* Intracellular steady-state concentrations *)
ci={15.,150.,7.5};

(* Extracellular steady-state concentrations *)
co={145.,5.,110.};

(* Transmembrane potential differences.[mV] *)
V0=R T/F Log[ci/co] 10^3;
{ENa,EK,El}=V0;

(* Maximum conductances [m-mhos/cm^2] *)
{gNa,gK,gl}={70.7,24.34,0.3};

(* Sodium channel (m) *)
fa[V_,Vr_]=0.1(25-(V-Vr));am[V_,Vr_]=fa[V,Vr]/(Exp[fa[V,Vr]]-1);
bm[V_,Vr_]=4. Exp[-(V-Vr)/18.];
iNa[V_,m_,h_]=gNa m^3 h (V-ENa);
msol = DSolve[{m'[t]==am[V,Vr] (1-m[t]) - bm[V,Vr] m[t],m[0]==0},m[t],t][[1,1]];

(* Potassium channel (n) *)
ga[V_,Vr_]=0.01(10-(V-Vr));an[V_,Vr_]=ga[V,Vr]/(Exp[ga[V,Vr]]-1);
bn[V_,Vr_]=0.125 Exp[-(V-Vr)/80.];
iK[V_,n_]=gK n^4 (V-EK);
nsol = DSolve[{n'[t]==an[V,Vr] (1-n[t]) - bn[V,Vr] n[t],n[0]==0},n[t],t][[1,1]];

(* Chloride channel (h) *)
ah[V_,Vr_]=0.07 Exp[-0.05(V-Vr)];
bh[V_,Vr_]=1./(Exp[0.1(30-(V-Vr))]+1.);
il[V_]=gl (V-El);hsol = DSolve[{h'[t]==ah[V,Vr] (1-h[t]) - bh[V,Vr] h[t],h[0]==0},h[t],t][[1,1]];

(* Currents for imposed voltage, resting voltage *)
HHcurrents[V_,Vr_]:=Module[{mt,nt,ht},
{mt,nt,ht}={m[t] /. msol, n[t] /. nsol, h[t] /. hsol};
Return[{iNa[V,mt,ht],iK[V,nt],il[V]}]
]

Manipulate

(* Manipulate is set up for Plot. Modify as needed *)
tf=1.;
Manipulate[
Plot[ 
Evaluate[ HHcurrents[V,Vr] /. {V->U+P HeavisidePi[(t-t0)/δt],Vr->U}],
{t,0,tf},
Frame->True,FrameLabel->{"t [ms]","i [mA]"},
PlotLabel->"Current response to V(t)=P Π((t-t0)/δt)",
GridLines->Automatic,
PlotLegends->ions],
{{P,-20},-40,0},
{{t0,0.02},0,tf},
{{δt,0.01},(10^(-6)),.1},
{{U,-70},-75,-65},
SaveDefinitions->True]