/* * Model Status * * This CellML model represents the ventricular M cell. The model * has been checked in both COR and PCEnv and it runs to recreate * the published results. The units have been checked and they * balance. As with the CellML version of the Hund and Rudy 2004 * model, on which this particular model is based, a stimulus protocol * has been added to allow the model to simulate action potentials * for 5 seconds. The parameter 'tissue' has been added to switch * between the original (default value 0) and 'tissue' (any other * value, for example 1) models. * * Model Structure * * ABSTRACT: We have constructed computational models of canine * ventricular cells and tissues, ultimately combining detailed * tissue architecture and heterogeneous transmural electrophysiology. * The heterogeneity is introduced by modifying the Hund-Rudy canine * cell model in order to reproduce experimentally reported electrophysiological * properties of endocardial, midmyocardial (M) and epicardial * cells. These models are validated against experimental data * for individual ionic current and action potential characteristics, * and their rate dependencies. 1D and 3D heterogeneous virtual * tissues are constructed, with detailed tissue architecture (anisotropy * and orthotropy, due to fibre orientation and sheet structure) * of the left ventricular wall wedge extracted from a diffusion * tensor imaging data set. The models are used to study the effects * of tissue heterogeneity and class III drugs on transmural propagation * and tissue vulnerability to re-entry. We have determined relationships * between the transmural dispersion of action potential duration * (APD) and the vulnerable window in the 1D virtual ventricular * wall, and demonstrated how changes in the transmural heterogeneity, * and hence tissue vulnerability, can lead to generation of re-entry * in the 3D ventricular wedge. Two class III drugs with opposite * qualitative effects on transmural APD heterogeneity are considered: * d-sotalol that increases transmural APD dispersion, and amiodarone * that decreases it. Simulations with the 1D virtual ventricular * wall show that under d-sotalol conditions the vulnerable window * is substantially wider compared to amiodarone conditions, primarily * in the epicardial region where unidirectional conduction block * persists until the adjacent M cells are fully repolarised. Further * simulations with the 3D ventricular wedge have shown that ectopic * stimulation of the epicardial region results in generation of * sustained re-entry under d-sotalol conditions, but not under * amiodarone conditions or in control. Again, APD increase in * M cells was identified as the major contributor to tissue vulnerability--re-entry * was initiated primarily due to ectopic excitation propagating * around the unidirectional conduction block in the M cell region. * This suggests an electrophysiological mechanism for the anti- * and proarrhythmic effects of the class III drugs: the relative * safety of amiodarone in comparison to d-sotalol can be explained * by relatively low transmural APD dispersion, and hence, a narrow * vulnerable window and low probability of re-entry in the tissue. * * model diagram * * [[Image file: hund_2004.png]] * * Schematic diagram of the Hund and Rudy 2004 canine ventricular * cell model (some of the model parameters and I_to equations * have been updated in the current Benson et al. 2008 model). * * The original paper reference is cited below: * * The canine virtual ventricular wall: a platform for dissecting * pharmacological effects on propagation and arrhythmogenesis, * Alan P. Benson, Oleg V. Aslanidi, Henggui Zhang and Arun V. * Holden, 2008 Progress in Biophysics and Molecular Biology, 96, * 187-208. PubMed ID: 17915298 */ import nsrunit; unit conversion on; unit m2u=.001 dimensionless; unit ms=.001 second^1; unit per_ms=1E3 second^(-1); unit uA_per_uF=1 kilogram^1*meter^2*second^(-4)*ampere^(-1); unit mS_per_uF=1E3 second^(-1); unit uL=1E-9 meter^3; unit mV=.001 kilogram^1*meter^2*second^(-3)*ampere^(-1); unit per_mV=1E3 kilogram^(-1)*meter^(-2)*second^3*ampere^1; unit uF=1E-6 kilogram^(-1)*meter^(-2)*second^4*ampere^2; unit cm=.01 meter^1; unit cm2=1E-4 meter^2; unit mJ_per_mole_K=.001 kilogram^1*meter^2*second^(-2)*kelvin^(-1)*mole^(-1); unit C_per_mole=1 second^1*ampere^1*mole^(-1); unit mM=1 meter^(-3)*mole^1; unit mM_per_ms=1E3 meter^(-3)*second^(-1)*mole^1; unit mM4=1 meter^(-12)*mole^4; unit uF_per_cm2=.01 kilogram^(-1)*meter^(-4)*second^4*ampere^2; unit uF_mole_per_C=1E-6 kilogram^(-1)*meter^(-2)*second^3*ampere^1*mole^1; unit L_per_F_ms=1 kilogram^1*meter^5*second^(-5)*ampere^(-2); math main { realDomain time ms; time.min=0; extern time.max; extern time.delta; real F C_per_mole; F=96485; real K_o mM; K_o=5.4; real Ca_o mM; Ca_o=1.8; real Na_o mM; Na_o=140; real Cl_o mM; Cl_o=100; real FonRT per_mV; FonRT=0.0374358835078; real tissue dimensionless; tissue=0; real V(time) mV; when(time=time.min) V=-85.781844107117; real INa(time) uA_per_uF; real ICaL(time) uA_per_uF; real IK1(time) uA_per_uF; real IKp(time) uA_per_uF; real IKs(time) uA_per_uF; real IKr(time) uA_per_uF; real IpCa(time) uA_per_uF; real ICab(time) uA_per_uF; real INaCa(time) uA_per_uF; real INaK(time) uA_per_uF; real Ito(time) uA_per_uF; real Ito2(time) uA_per_uF; real IClb(time) uA_per_uF; real INal(time) uA_per_uF; real caiont(time) uA_per_uF; real naiont(time) uA_per_uF; real kiont(time) uA_per_uF; real clont(time) uA_per_uF; real l cm; l=0.01; real a cm; a=0.0011; real vcell uL; real ageo cm2; real Acap uF; real vmyo uL; real vmito uL; real vsr uL; real vnsr uL; real vjsr uL; real vss uL; real AF uF_mole_per_C; real stim_offset ms; stim_offset=0; real stim_period ms; stim_period=1e3; real stim_duration ms; stim_duration=3; real stim_amplitude uA_per_uF; stim_amplitude=-20; real i_Stim(time) uA_per_uF; real past(time) ms; real ENa(time) mV; real GNa mS_per_uF; real gNa(time) mS_per_uF; real H(time) dimensionless; when(time=time.min) H=0.987317750543; real m(time) dimensionless; when(time=time.min) m=0.001356538159; real J(time) dimensionless; when(time=time.min) J=0.991924983076; real am(time) per_ms; real bm(time) per_ms; real ah(time) per_ms; real bh(time) per_ms; real aj(time) per_ms; real bj(time) per_ms; real Ca_ss(time) mM; when(time=time.min) Ca_ss=0.00012271265; real d(time) dimensionless; when(time=time.min) d=0.00000164013; real dp(time) dimensionless; when(time=time.min) dp=8.98230672628; real f(time) dimensionless; when(time=time.min) f=0.999961508634; real fca(time) dimensionless; when(time=time.min) fca=0.97836624923; real fca2(time) dimensionless; when(time=time.min) fca2=0.893052931249; real f2(time) dimensionless; when(time=time.min) f2=0.992234519148; real pca L_per_F_ms; pca=2.43e-4; real gacai dimensionless; gacai=1; real gacao dimensionless; gacao=0.341; real CaMKactive(time) dimensionless; real ibarca(time) uA_per_uF; real dss(time) dimensionless; real taud(time) ms; real fss(time) dimensionless; real f2ss(time) dimensionless; real tauf(time) ms; real tauf2(time) ms; real dpss(time) dimensionless; real fcass(time) dimensionless; real fca2ss(time) dimensionless; real taufca(time) ms; real taufca2(time) ms; real EK(time) mV; real ak1(time) per_ms; real bk1(time) per_ms; real gkr mS_per_uF; real r(time) dimensionless; real xr(time) dimensionless; when(time=time.min) xr=0.00000724074; real xrss(time) dimensionless; real tauxr(time) ms; real Ca_i(time) mM; when(time=time.min) Ca_i=0.00012131666; real gks(time) mS_per_uF; real EKs(time) mV; real xss(time) dimensionless; real tauxs(time) ms; real xs1(time) dimensionless; when(time=time.min) xs1=0.019883138161; real xs2(time) dimensionless; when(time=time.min) xs2=0.019890650554; real gitodv mS_per_uF; gitodv=0.1805; real rv(time) dimensionless; real ay(time) per_ms; real by(time) per_ms; real ay2(time) per_ms; real by2(time) per_ms; real ay3(time) per_ms; real by3(time) per_ms; real ydv(time) dimensionless; when(time=time.min) ydv=0.013970786703; real ydv2(time) dimensionless; when(time=time.min) ydv2=0.99996472752; real zdv(time) dimensionless; when(time=time.min) zdv=0.829206149767; real Na_i(time) mM; when(time=time.min) Na_i=12.972433387269; real kmnai mM; kmnai=10; real kmko mM; kmko=1.5; real ibarnak uA_per_uF; ibarnak=0.61875; real sigma dimensionless; real fnak(time) dimensionless; real ca_i_NaCa(time) mM; real KmCa mM; KmCa=1.25e-4; real allo(time) dimensionless; real NCXmax uA_per_uF; NCXmax=5.85; real ksat dimensionless; ksat=0.27; real eta dimensionless; eta=0.35; real KmNai mM; KmNai=12.3; real KmNao mM; KmNao=87.5; real KmCai mM; KmCai=0.0036; real KmCao mM; KmCao=1.3; real num(time) mM4; real denom1(time) dimensionless; real denom2(time) mM4; real denom3(time) mM4; real ibarpca uA_per_uF; ibarpca=0.0575; real kmpca mM; kmpca=0.5e-3; real Cl_i(time) mM; when(time=time.min) Cl_i=15.59207157178; real PCl L_per_F_ms; PCl=4e-7; real AA(time) dimensionless; when(time=time.min) AA=0.000816605172; real Ito2_max(time) uA_per_uF; real AAss(time) dimensionless; real Kmto2 mM; Kmto2=0.1502; real ECl(time) mV; real GClb mS_per_uF; GClb=2.25e-4; real GNaL mS_per_uF; GNaL=65e-4; real mL(time) dimensionless; when(time=time.min) mL=0.001356538159; real hL(time) dimensionless; when(time=time.min) hL=0.26130711759; real amL(time) per_ms; real bmL(time) per_ms; real hLss(time) dimensionless; real K_i(time) mM; when(time=time.min) K_i=135.469546216758; real prnak dimensionless; prnak=0.01833; real Ca_jsr(time) mM; when(time=time.min) Ca_jsr=1.737580994071; real Grel(time) per_ms; real dro_inf(time) dimensionless; real dtau_rel_max ms; dtau_rel_max=10; real dtau_rel(time) ms; real ross(time) dimensionless; real riss(time) dimensionless; real tauri(time) ms; real irelcicr(time) mM_per_ms; real CaMK0 dimensionless; CaMK0=0.05; real Km mM; Km=0.0015; real KmCaMK dimensionless; KmCaMK=0.15; real CaMKbound(time) dimensionless; real CaMKtrap(time) dimensionless; when(time=time.min) CaMKtrap=0.021123704774; real ro(time) dimensionless; when(time=time.min) ro=0; real ri(time) dimensionless; when(time=time.min) ri=0.862666650318; real vg(time) dimensionless; real cafac(time) dimensionless; real dKmPLBmax mM; dKmPLBmax=0.00017; real dJupmax dimensionless; dJupmax=0.75; real dKmPLB(time) mM; real dJup(time) dimensionless; real iupmax mM_per_ms; iupmax=0.004375; real Kmup mM; Kmup=0.00092; real nsrmax mM; nsrmax=15; real iup(time) mM_per_ms; real ileak(time) mM_per_ms; real Ca_nsr(time) mM; when(time=time.min) Ca_nsr=1.832822335168; real idiff(time) mM_per_ms; real itr(time) mM_per_ms; real CTNaCl(time) mM_per_ms; real CTNaClmax mM_per_ms; CTNaClmax=9.8443e-6; real CTKCl(time) mM_per_ms; real CTKClmax mM_per_ms; CTKClmax=7.0756e-6; real kmt mM; kmt=0.5e-3; real kmc mM; kmc=2.38e-3; real tbar mM; tbar=70e-3; real cbar mM; cbar=50e-3; real kmcsqn mM; kmcsqn=0.8; real csqnbar mM; csqnbar=10; real bcsqn(time) dimensionless; real bmyo(time) dimensionless; real BSRmax mM; BSRmax=0.047; real KmBSR mM; KmBSR=0.00087; real BSLmax mM; BSLmax=1.124; real KmBSL mM; KmBSL=0.0087; real bss(time) dimensionless; // // vcell=((1E3 m2u)*3.141592653589793*a*a*l); ageo=(2*3.141592653589793*a*a+2*3.141592653589793*a*l); Acap=(ageo*(2 uF_per_cm2)); vmyo=(vcell*.68); vmito=(vcell*.26); vsr=(vcell*.06); vnsr=(vcell*.0552); vjsr=(vcell*.0048); vss=(vcell*.02); AF=(Acap/F); past=(floor(time/stim_period)*stim_period); i_Stim=(if (((time-past)>=stim_offset) and ((time-past)<=(stim_offset+stim_duration))) stim_amplitude else (0 uA_per_uF)); caiont=(ICaL+ICab+IpCa-2*INaCa); naiont=(INa+3*INaCa+3*INaK+INal); kiont=(IKr+IKs+IK1+IKp-2*INaK+Ito+.5*i_Stim); clont=(IClb+Ito2+.5*i_Stim); V:time=((-1)*(naiont+kiont+caiont+clont)); // GNa=(if (tissue=0) (11.5 mS_per_uF) else 4*(11.5 mS_per_uF)); gNa=(GNa*m*m*m*H*J); INa=(gNa*(V-ENa)); ah=(if (V>=((-1)*(40 mV))) (0 per_ms) else (.135 per_ms)*exp(((80 mV)+V)/((-1)*(6.8 mV)))); aj=(if (V>=((-1)*(40 mV))) (0 per_ms) else ((-1)*(127140.00000000001 per_ms)*exp((.2444 per_mV)*V)-(3.4739999999999996E-5 per_ms)*exp((-1)*(.04391 per_mV)*V))*(1 per_mV)*(V+(37.78 mV))/(1+exp((.311 per_mV)*(V+(79.23 mV))))); bh=(if (V>=((-1)*(40 mV))) (1 per_ms)/(.13*(1+exp((V+(10.66 mV))/((-1)*(11.1 mV))))) else (3.56 per_ms)*exp((.079 per_mV)*V)+(3.1E5 per_ms)*exp((.35 per_mV)*V)); bj=(if (V>=((-1)*(40 mV))) (.3 per_ms)*exp((-1)*(2.535E-7 per_mV)*V)/(1+exp((-1)*(.1 per_mV)*(V+(32 mV)))) else (.1212 per_ms)*exp((-1)*(.01052 per_mV)*V)/(1+exp((-1)*(.1378 per_mV)*(V+(40.14 mV))))); am=((.32 per_ms)*(1 per_mV)*(V+(47.13 mV))/(1-exp((-1)*(.1 per_mV)*(V+(47.13 mV))))); bm=((.08 per_ms)*exp((-1)*V/(11 mV))); H:time=(ah*(1-H)-bh*H); m:time=(am*(1-m)-bm*m); J:time=(aj*(1-J)-bj*J); // ibarca=(pca*4*(V-(15 mV))*F*FonRT*(gacai*Ca_ss*exp(2*(V-(15 mV))*FonRT)-gacao*Ca_o)/(exp(2*(V-(15 mV))*FonRT)-1)); ICaL=(if (tissue=0) d^dp*f*f2*fca*fca2*ibarca else d*f*f2*fca*fca2*ibarca); dss=(1/(1+exp((-1)*(V-(4 mV))/(6.74 mV)))); taud=((.59 ms)+(.8 ms)*exp((.052 per_mV)*(V+(13 mV)))/(1+exp((.132 per_mV)*(V+(13 mV))))); fss=(.7/(1+exp((V+(17.12 mV))/(7 mV)))+.3); f2ss=(.77/(1+exp((V+(17.12 mV))/(7 mV)))+.23); tauf=((1 ms)/(.2411*exp((-1)*((.045 per_mV)*(V-(9.6914 mV)))^2)+.0529)); tauf2=((1 ms)/(.0423*exp((-1)*((.059 per_mV)*(V-(18.5726 mV)))^2)+.0054)); dpss=(9-8/(1+exp((-1)*(V+(65 mV))/(3.4 mV)))); fcass=(.3/(1-ICaL/(.05 uA_per_uF))+.55/(1+Ca_ss/(.003 mM))+.15); fca2ss=(1/(1-ICaL/(.01 uA_per_uF))); taufca=((10 ms)*CaMKactive/(.15+CaMKactive)+(1 ms)/(1+Ca_ss/(.003 mM))+(.5 ms)); taufca2=((300 ms)/(1+exp(((-1)*ICaL-(.175 uA_per_uF))/(.04 uA_per_uF)))+(125 ms)); d:time=((dss-d)/taud); dp:time=((dpss-dp)/(10 ms)); f:time=((fss-f)/tauf); f2:time=((f2ss-f2)/tauf2); fca:time=((fcass-fca)/taufca); fca2:time=((fca2ss-fca2)/taufca2); // ak1=((1.02 per_ms)/(1+exp((.2385 per_mV)*(V-EK-(59.215 mV))))); bk1=(((.49124 per_ms)*exp((.08032 per_mV)*(V-EK+(5.476 mV)))+(1 per_ms)*exp((.06175 per_mV)*(V-EK-(594.31 mV))))/(1+exp((-1)*(.5143 per_mV)*(V-EK+(4.753 mV))))); IK1=((.5 mS_per_uF)*sqrt(K_o/(5.4 mM))*ak1/(ak1+bk1)*(V-EK)); // gkr=((.0138542 mS_per_uF)*sqrt(K_o/(5.4 mM))); r=(1/(1+exp((V+(10 mV))/(15.4 mV)))); IKr=(gkr*xr*r*(V-EK)); xrss=(1/(1+exp((-1)*(V+(10.085 mV))/(4.25 mV)))); tauxr=((1 ms)/((6.000000000000001E-4 per_mV)*(V-(1.7384 mV))/(1-exp((-1)*(.136 per_mV)*(V-(1.7384 mV))))+(3.0000000000000003E-4 per_mV)*(V+(38.3608 mV))/(exp((.1522 per_mV)*(V+(38.3608 mV)))-1))); xr:time=((xrss-xr)/tauxr); // gks=((.0074692 mS_per_uF)*(1+.6/(1+((3.7999999999999995E-5 mM)/Ca_i)^1.4))); xss=(1/(1+exp((-1)*(V-(10.5 mV))/(24.7 mV)))); tauxs=((1 ms)/((7.609999999999999E-5 per_mV)*(V+(44.6 mV))/(1-exp((-1)*(9.97 per_mV)*(V+(44.6 mV))))+(3.6E-4 per_mV)*(V-(.55 mV))/(exp((.128 per_mV)*(V-(.55 mV)))-1))); IKs=(gks*xs1*xs2*(V-EKs)); xs1:time=((xss-xs1)/tauxs); xs2:time=((xss-xs2)/tauxs/2); // rv=exp(V/(300 mV)); ay=((25 per_ms)*exp((V-(40 mV))/(25 mV))/(1+exp((V-(40 mV))/(25 mV)))); by=((25 per_ms)*exp((-1)*(V+(90 mV))/(25 mV))/(1+exp((-1)*(V+(90 mV))/(25 mV)))); ay2=((.03 per_ms)/(1+exp((V+(60 mV))/(5 mV)))); by2=((.2 per_ms)*exp((V+(25 mV))/(5 mV))/(1+exp((V+(25 mV))/(5 mV)))); ay3=((.00225 per_ms)/(1+exp((V+(60 mV))/(5 mV)))); by3=((.1 per_ms)*exp((V+(25 mV))/(5 mV))/(1+exp((V+(25 mV))/(5 mV)))); Ito=(gitodv*ydv^3*ydv2*zdv*rv*(V-EK)); ydv:time=(ay*(1-ydv)-by*ydv); ydv2:time=(ay2*(1-ydv2)-by2*ydv2); zdv:time=(ay3*(1-zdv)-by3*zdv); // sigma=((exp(Na_o/(67.3 mM))-1)/7); fnak=(1/(1+.1245*exp((-1)*.1*V*FonRT)+.0365*sigma*exp((-1)*V*FonRT))); INaK=(ibarnak*fnak*1/(1+(kmnai/Na_i)^2)*K_o/(K_o+kmko)); // ca_i_NaCa=(1.5*Ca_i); allo=(1/(1+(KmCa/ca_i_NaCa)^2)); num=(Na_i^3*Ca_o*exp(eta*V*FonRT)-Na_o^3*ca_i_NaCa*exp((eta-1)*V*FonRT)); denom1=(1+ksat*exp((eta-1)*V*FonRT)); denom2=(KmCao*Na_i^3+KmNao^3*ca_i_NaCa+KmNai^3*Ca_o*(1+ca_i_NaCa/KmCai)); denom3=(KmCai*Na_o^3*(1+(Na_i/KmNai)^3)+Na_i^3*Ca_o+Na_o^3*ca_i_NaCa); INaCa=(NCXmax*allo*num/(denom1*(denom2+denom3))); // IKp=((.00276 mS_per_uF)*(V-EK)/(1+exp(((7.488 mV)-V)/(5.98 mV)))); // IpCa=(ibarpca*Ca_i/(kmpca+Ca_i)); // ICab=((1.995084E-7 L_per_F_ms)*4*V*F*FonRT*(Ca_i*exp(2*V*FonRT)-.341*Ca_o)/(exp(2*V*FonRT)-1)); // Ito2_max=(PCl*V*F*FonRT*(Cl_i-Cl_o*exp(V*FonRT))/(1-exp(V*FonRT))); AAss=(1/(1+Kmto2/Ca_ss)); Ito2=(Ito2_max*AA); AA:time=((AAss-AA)/(1 ms)); // IClb=(GClb*(V-ECl)); // INal=(GNaL*mL^3*hL*(V-ENa)); amL=((.32 per_ms)*(1 per_mV)*(V+(47.13 mV))/(1-exp((-1)*(.1 per_mV)*(V+(47.13 mV))))); bmL=((.08 per_ms)*exp((-1)*V/(11 mV))); hLss=(1/(1+exp((V+(91 mV))/(6.1 mV)))); mL:time=(amL*(1-mL)-bmL*mL); hL:time=((hLss-hL)/(600 ms)); // ENa=(ln(Na_o/Na_i)/FonRT); EK=(ln(K_o/K_i)/FonRT); EKs=(ln((K_o+prnak*Na_o)/(K_i+prnak*Na_i))/FonRT); ECl=((-1)*ln(Cl_o/Cl_i)/FonRT); // CaMKbound=(CaMK0*(1-CaMKtrap)/(1+Km/Ca_ss)); CaMKactive=(CaMKbound+CaMKtrap); Grel=((3E3 per_ms)*vg); vg=(if (tissue=0) 1/(1+exp((ibarca+(13 uA_per_uF))/(5 uA_per_uF))) else 1); cafac=(1/(1+exp((ICaL+(.05 uA_per_uF))/(.015 uA_per_uF)))); dro_inf=(Ca_jsr^1.9/(Ca_jsr^1.9+((49.28 mM)*Ca_ss/(Ca_ss+(.0028 mM)))^1.9)); dtau_rel=(dtau_rel_max*CaMKactive/(KmCaMK+CaMKactive)); ross=(dro_inf/(((1 uA_per_uF)/ICaL)^2+1)); riss=(1/(1+exp((Ca_ss-(4E-4 mM)+(.002 mM)*cafac)/(2.4999999999999998E-5 mM)))); tauri=((3 ms)+dtau_rel+((350 ms)-dtau_rel)/(1+exp((Ca_ss-(.003 mM)+(.003 mM)*cafac)/(2E-4 mM)))); irelcicr=(Grel*ro*ri*(Ca_jsr-Ca_ss)); ro:time=((ross-ro)/(3 ms)); ri:time=((riss-ri)/tauri); CaMKtrap:time=((.05 per_ms)*CaMKactive*(CaMKactive-CaMKtrap)-(6.8E-4 per_ms)*CaMKtrap); // dKmPLB=(dKmPLBmax*CaMKactive/(KmCaMK+CaMKactive)); dJup=(dJupmax*CaMKactive/(KmCaMK+CaMKactive)); iup=((dJup+1)*iupmax*Ca_i/(Ca_i+Kmup-dKmPLB)); ileak=(iupmax*Ca_nsr/nsrmax); // idiff=((Ca_ss-Ca_i)/(.2 ms)); itr=((Ca_nsr-Ca_jsr)/(120 ms)); // CTNaCl=(CTNaClmax*(ENa-ECl)^4/((ENa-ECl)^4+(87.8251 mV)^4)); Na_i:time=((-1)*naiont*AF/vmyo+CTNaCl); // CTKCl=(CTKClmax*(EK-ECl)/(EK-ECl+(87.8251 mV))); K_i:time=((-1)*kiont*AF/vmyo+CTKCl); // Cl_i:time=(clont*AF/vmyo+CTNaCl+CTKCl); // Ca_i:time=(bmyo*((-1)*(ICab+IpCa-2*INaCa)*AF/(vmyo*2)+(ileak-iup)*vnsr/vmyo+idiff*vss/vmyo)); Ca_ss:time=(bss*((-1)*ICaL*AF/(vss*2)+irelcicr*vjsr/vss-idiff)); Ca_nsr:time=(iup-itr*vjsr/vnsr-ileak); Ca_jsr:time=(bcsqn*(itr-irelcicr)); bcsqn=(1/(1+kmcsqn*csqnbar/(Ca_jsr+kmcsqn)^2)); bmyo=(1/(1+cbar*kmc/(Ca_i+kmc)^2+kmt*tbar/(Ca_i+kmt)^2)); bss=(1/(1+BSRmax*KmBSR/(KmBSR+Ca_ss)^2+BSLmax*KmBSL/(KmBSL+Ca_ss)^2)); // // // }