: \$Id: boltz_cvode.inc,v 1.1.1.1 2005/12/15 15:16:39 hines Exp \$ TITLE Boltzmann eqn definition channel COMMENT Each channel has activation and inactivation particles as in the original Hodgkin Huxley formulation. The activation particle mm and inactivation particle hh go from on to off states according to kinetic variables alpha and beta which are voltage dependent. The form of the alpha and beta functions were dissimilar in the HH study. The formulae are: Inf = 1./(1.+exp((v+vhalf)/kconst)) Tau must be set separately by the function tauset() chanexp : number for exponentiating the state variable; e.g. original HH Na channel use m^3, note that chanexp = 0 will turn off this state variable erev : reversal potential for the channel celsius : temperature at which experiment was done (Tau will will be adjusted using a q10 of 3.0) vhalf : (a voltage) determines the voltage at which the value of the sigmoid function for Inf is 1/2 vrest : voltage shift for vhalf kconst: the Boltzmann K : determines steepness of the sigmoid ENDCOMMENT INDEPENDENT {t FROM 0 TO 1 WITH 1 (ms)} NEURON { RANGE gmax, g, i GLOBAL erev, Inf, Tau, vrest } : end NEURON CONSTANT { FARADAY = 96489.0 : Faraday's constant R= 8.31441 : Gas constant } : end CONSTANT UNITS { (mA) = (milliamp) (mV) = (millivolt) (umho) = (micromho) } : end UNITS COMMENT ** parameter values should come from files specific to particular channels ASSIGNED { ina ena } : eg for Na PARAMETER { erev = 0 (mV) gmax = 0 (S/cm2) vrest = 0 (mV) mvhalf = 0 mkconst = 0 exptemp = 0 mq10 = 3 mexp = 0 hvhalf = 0 hkconst = 0 hq10 = 3 hexp = 0 } : end PARAMETER ENDCOMMENT PARAMETER { cao (mM) cai (mM) celsius (degC) dt (ms) v (mV) } ASSIGNED { i (mA/cm^2) g (mho/cm^2) Inf[2] : 0 = m and 1 = h Tau[2] : 0 = m and 1 = h mexp_val hexp_val } : end ASSIGNED STATE { m h } INITIAL { mh(v) m = Inf[0] h = Inf[1] } BREAKPOINT { LOCAL hexp_val, index, mexp_val SOLVE states METHOD cnexp hexp_val = 1 mexp_val = 1 : Determining h's exponent value if (hexp > 0) { FROM index=1 TO hexp { hexp_val = h * hexp_val } } : Determining m's exponent value if (mexp > 0) { FROM index = 1 TO mexp { mexp_val = m * mexp_val } } : mexp hexp : Note that mexp_val is now = m and hexp_val is now = h g = gmax * mexp_val * hexp_val iassign() } : end BREAKPOINT : ASSIGNMENT PROCEDURES : Must be overwritten by user routines in parameters.multi : PROCEDURE iassign () { i = g*(v-erev) ina=i } : PROCEDURE iassign () { i = g*ghkca(v) ica=i } :------------------------------------------------------------------- : I suppose we have 2 choices, to use the DERIVATIVE function or : to explicitly state m+ and h+. If you were to use the DERIVATIVE : function, then you will do as follows: : DERIVATIVE deriv { : m' = (-m + minf) / mtau : h' = (-h + hinf) / htau : } DERIVATIVE states { mh(v) m' = (-m + Inf[0]) / Tau[0] h' = (-h + Inf[1]) / Tau[1] } :------------------------------------------------------------------- : NOTE : 0 = m and 1 = h PROCEDURE mh (v) { LOCAL a, b, j, mqq10, hqq10, mv0, hv0 mv0 = mvhalf + vrest hv0 = hvhalf + vrest mqq10 = mq10^((celsius-exptemp)/10.) : Calculater Inf and Tau values for h and m Inf[0]=1./(1.+exp((v+mv0)/mkconst)) if (hexp == 0) { Inf[1] = 1 Tau[1]=1 } else { Inf[1]=1./(1.+exp((v+hv0)/hkconst)) hqq10 = hq10^((celsius-exptemp)/10.) Tau[1]=settau(1,v)/hqq10 } Tau[0]=settau(0,v)/mqq10 } : end PROCEDURE mh (v) :------------------------------------------------------------------- FUNCTION FRT(temperature) { FRT = FARADAY * 0.001 / R / (temperature + 273.15) } : end FUNCTION FRT (temperature) :------------------------------------------------------------------- FUNCTION ghkca (v) { : Goldman-Hodgkin-Katz eqn LOCAL nu, efun nu = v*2*FRT(celsius) if(fabs(nu) < 1.e-6) { efun = 1.- nu/2. } else { efun = nu/(exp(nu)-1.) } ghkca = -FARADAY*2.e-3*efun*(cao - cai*exp(nu)) } : end FUNCTION ghkca() :________________________________________________________________ FUNCTION interp (v,ystart,yend,vmin,vmax) { interp = ystart + ((v-vmin)/(vmax-vmin))*(yend-ystart) } : end FUNCTION interp() COMMENT : j==0 m; j==1 h FUNCTION settau(j,v) { if (j==0) { if (v <= -110) { settau = 15 } else if (v < -40) { settau = interp(v,15., 80., -110.,-40.) } else if (v < 40) { settau = interp(v,80., 20, -40.,40.) } else { settau = 20 } } else { if (v <= -110) { settau = 200 } else if (v < -75) { settau = interp(v,200., 1000., -110.,-75.) } else { settau = 1000 } } } : end FUNCTION settau (j) ENDCOMMENT