// Ref 1: Mechanisms of altered excitation-contraction coupling // in Canine Tachycardia-induced heart Failure, II // Raimond Winslow, et al, Circ Res, 84:571-586, 1999. // Ref 2: Role of the Calcium-independent transient outward current // Ito1 in shaping action potential morphology and duration // Joseph Greenstein et al., Circ Res, 1026-1033, 2000. // JSim v1.1 import nsrunit; unit conversion off; math example1 { // real domain t[0,500,1001] ms; realDomain t ms; t.min=0.0; t.max=500.0; t.delta=0.5; extern real I(t) uA/uF; extern real Vclamp(t) mV; real Csc = 1 uF/cm^2, // (Table 1) Acap = 1.534*10^(-4) cm^2, Vmyo = 25.84*10^(-6) uL, VJSR = 0.16 *10^(-6) uL, VNSR = 2.10 *10^(-6) uL, Vss = 1.2 *10^(-9) uL, K_o = 4.0 mM, // (Table 2) Na_o = 138.0 mM, Ca_o = 2.0 mM, F = 96.5 coulomb/mmol, // (Table 3) Temp = 310 K, R = 8.314 J/(mol*K), RT_over_F = R*Temp/F, GKr = 0.00339 mS/uF, GKs = 0.00271 mS/uF, GK1 = 2.8 mS/uF, GKp = 0.002216 mS/uF, GNa = 12.8 mS/uF, kNaCa = 0.30 uA/uF, KmNa = 87.5 mM, KmCa = 1.38 mM, KmK1 = 13.0 mM, ksat = 0.2, eta = 0.35, INaKmax = 0.693 uA/uF, KmNa_i = 10.0 mM, KmK_o = 1.5 mM, IpCamax = 0.05 uA/uF, KmpCa = 0.0005 mM, GCa_b = 0.00007684 mS/uF, GNa_b = 0.0031 mS/uF, v1 = 1.8 1/ms, // (Table 4) Kfb = 0.168e-3 mM, Krb = 3.29 mM, KSR = 1.0, Nfb = 1.2, Nrb = 1.0, vmaxf = 0.748*10^(-4) mM/ms, vmaxr = 0.318*10^(-3) mM/ms, tautr = 0.5747 ms, tauxfer= 26.7 ms, kap = 12.15e9 mM^(-4)*ms^(-1), kam = 0.576 1/ms, kbp = 4.05e6 mM^(-3)*ms^(-1), kbm = 1.930 1/ms, kcp = 0.1 1/ms, kcm = 0.0008 1/ms, ncoop = 4, mcoop = 3, fL = 0.3 1/ms, // (Table 5) gL = 2.0 1/ms, // fp = 0.005 1/ms, // gp = 7 1/ms, fp = 0 1/ms, gp = 0 1/ms, b = 1.9356, a = 2.0, omega = 0.01 1/ms, PCa = 2.469*10^(-4) cm/ms, PK = 4.574*10^(-7) cm/ms, ICahalf = -0.265 uA/uF, GKv43 = 0.1194 mS/uF, // (Table 6) PKv14 = 1.709*10^(-7) cm/s, alpha_NaK = 0.02, alpha_a0Kv43 = 0.5437 1/ms, a_aKv43 = 0.02898 1/mV, beta_a0Kv43 = 0.08019 1/ms, b_aKv43 = 0.04684 1/mV, alpha_i0Kv43 = 0.04984 1/ms, a_iKv43 = 3.7302*10^(-4) 1/mV, beta_i0Kv43 = 8.1948*10^(-4) 1/ms, b_iKv43 = 5.374*10^(-8) 1/mV, f1Kv43 = 1.8936, f2Kv43 = 14.225, f3Kv43 = 158.574, f4Kv43 = 142.937, b1Kv43 = 6.7735, b2Kv43 = 15.621, b3Kv43 = 28.753, b4Kv43 = 524.576, alpha_a0Kv14 = 1.8931 1/ms, a_aKv14 = 0.006950 1/mV, beta_a0Kv14 = 0.01179 1/ms, b_aKv14 = 0.08527 1/mV, alpha_i0Kv14 = 0.002963 1/ms, a_iKv14 = 0.0 1/mV, beta_i0Kv14 = 1.0571*10^(-4) 1/ms, b_iKv14 = 0.0 1/mV, f1Kv14 = 0.2001, f2Kv14 = 0.3203, f3Kv14 = 13.509, f4Kv14 = 1151.745, b1Kv14 = 2.2300, b2Kv14 = 12.000, b3Kv14 = 5.370, b4Kv14 = 5.240, LTRPNtot = 70*10^(-3) mM, // (Table 7) HTRPNtot = 140*10^(-3) mM, Khtrpnp = 20 1/mM/ms, Khtrpnm = 66*10^(-6) 1/ms, Kltrpnp = 40 1/mM/ms, Kltrpnm = 0.04 1/ms, CMDNtot = 50*10^(-3) mM, CSQNtot = 15 mM, KmCMDN = 2.38*10^(-3) mM, KmCSQN = 0.8 mM, Vini = -95.431 mV; real clamp_flag = 0, Iion(t) uA/uF; real chf_flag = 0, // CHF Parameters chfsc_IKv43 = 0.143, chfsc_IK1 = 0.68, chfsc_Jup = 0.38, chfsc_INaCa = 1.75; //Variables for all the algebraic equations real INa(t) uA/uF, ENa(t) mV, alpham(t), betam(t), alphah(t), betah(t), alphaj(t), betaj(t); real IKr(t) uA/uF, EK(t) mV, RV(t), fK_o(t), K12(t), K21(t), XKr_inf(t), tauXKr(t); real IKs(t) uA/uF, EKs(t) mV, XKs_inf(t), tauXKs(t); real Ito1(t) uA/uF; real IK1(t) uA/uF, K1_inf(t); real IKp(t) uA/uF, Kp(t); real INaCa(t) uA/uF; real INaK(t) uA/uF, fNaK(t), sigma(t); real IpCa(t) uA/uF, ICa_b(t) uA/uF, ECa(t) mV; real INa_b(t) uA/uF; real alpha(t), beta(t), alphap(t), betap(t), gamma(t), ICamax(t) uA/uF, ICa(t) uA/uF, ICa_K(t) uA/uF, PKp(t), y_inf(t), tauy(t); real Jrel(t) mM/ms, fb(t), rb(t), Jup(t) mM/ms, Jtr(t) mM/ms, Jxfer(t) mM/ms, Jtrpn(t) mM/ms, betai(t), betass(t), betaJSR(t); real alpha_aKv43(t), beta_aKv43(t), alpha_iKv43(t), beta_iKv43(t), alpha_aKv14(t), beta_aKv14(t), alpha_iKv14(t), beta_iKv14(t); //State variables for all the ODEs real m(t), h(t), j(t), VV(t) mV, V(t) mV, Na_i(t) mM, K_i(t) mM, Ca_i(t) mM, Ca_NSR(t) mM, Ca_ss(t) mM, Ca_JSR(t) mM, PC1(t), PO1(t), PO2(t), PC2(t), XKr(t), XKs(t), C0(t), C1(t), C2(t), C3(t), C4(t), O(t), CCa0(t), CCa1(t), CCa2(t), CCa3(t), CCa4(t), OCa(t), y(t), LTRPNCa(t) mM, HTRPNCa(t) mM, IKv43(t) uA/uF, IKv14K(t) uA/uF, IKv14Na(t) uA/uF, C0Kv43(t), C1Kv43(t), C2Kv43(t), C3Kv43(t), OKv43(t), CI0Kv43(t), CI1Kv43(t), CI2Kv43(t), CI3Kv43(t), OIKv43(t), C0Kv14(t), C1Kv14(t), C2Kv14(t), C3Kv14(t), OKv14(t), CI0Kv14(t), CI1Kv14(t), CI2Kv14(t), CI3Kv14(t), OIKv14(t); // Temporary variables private real VFsqr_over_RT(t), VF_over_RT(t), tmp1(t), dHTRPNCadt(t), dLTRPNCadt(t), dPC1dt(t), dPC2dt(t), dPO2dt(t); when (t=t.min) { VV = if (clamp_flag = 0) Vini else Vclamp; // (Table 8) m = 2.65431*10^(-4); h = 0.998552; j = 0.99877; XKr = 0.26747; XKs = 1.5017*10^(-4); Na_i = 10; K_i = 155.464; Ca_i = 7.9547*10^(-5); Ca_NSR = 0.2291; Ca_ss = 1.2031*10^(-4); Ca_JSR = 0.22875; PC1 = 0.4952836; PO1 = 5.99064*10^(-4); PO2 = 9.0*10^(-7); PC2 = 0.504120; C0 = 0.99838; C1 = 7.0871*10^(-6); C2 = 0; C3 = 0; C4 = 0; O = 0; CCa0 = 1.606*10^(-3); CCa1 = 2.263*10^(-6); CCa2 = 0; CCa3 = 0; CCa4 = 0; OCa = 0; y = 1.0; LTRPNCa = 5.2396*10^(-3); HTRPNCa = 136.484*10^(-3); C0Kv43 = 0.590957; C1Kv43 = 0.155217; C2Kv43 = 0.015288; C3Kv43 = 6.6924*10^(-4); OKv43 = 1.0785*10^(-5); CI0Kv43 = 0.221; CI1Kv43 = 0.014351; CI2Kv43 = 5.9681*10^(-4); CI3Kv43 = 6.4201*10^(-5); OIKv43 = 1.84528*10^(-3); C0Kv14 = 0.764631; C1Kv14 = 0.073937; C2Kv14 = 2.68108*10^(-3); C3Kv14 = 4.32122*10^(-5); OKv14 = 2.6204*10^(-7); CI0Kv14 = 0.149331; CI1Kv14 = 6.4425*10^(-3); CI2Kv14 = 2.0126*10^(-3); CI3Kv14 = 6.1286*10^(-4); OIKv14 = 3.0843*10^(-4); } VFsqr_over_RT = V*F^2*1000/(R*Temp); VF_over_RT = V*F/(R*Temp); INa = GNa*m^3*h*j*(V-ENa); //(A.1) ENa = RT_over_F*ln(Na_o/Na_i); //(A.2) m:t = alpham*(1-m)-betam*m; //(A.3) h:t = alphah*(1-h)-betah*h; //(A.4) j:t = alphaj*(1-j)-betaj*j; //(A.5) alpham = if (abs(V+47.13) > 0.0001) 0.32*(V+47.13)/(1-exp(-0.1*(V+47.13))) else 0.32/(0.1-0.005*(V+47.13)); //(A.6) betam = 0.08*exp(-V/11); //(A.7) alphah = if (V >= -40) 0 else 0.135*exp(-(80+V)/6.8); //(A.8 & 12) betah = if (V >= -40) 1/(0.13*(1+exp(-(V+10.66)/11.1))) else (3.56*exp(0.079*V)+310000*exp(0.35*V)); //(A.9 & 13) alphaj = if (V >= -40) 0 else (-127140*exp(0.2444*V)-3.474*10^(-5)* exp(-0.04391*V))*(V+37.78) /(1+exp(0.311*(V+79.23))); //(A.10 & 14) betaj = if (V >= -40) 0.3*exp(-2.535*10^(-7)*V)/(1+exp(-0.1*(V+32))) else 0.1212*exp(-0.01052*V) /(1+exp(-0.1378*(V+40.14))); //(A.11 & 15) IKr = GKr*fK_o*RV*XKr*(V-EK); //(A.16) EK = RT_over_F*ln(K_o/K_i); //(A.17) RV = 1/(1+1.4945*exp(0.0446*V)); //(A.18) fK_o = sqrt(K_o/4); //(A.19) XKr:t = (XKr_inf-XKr)/tauXKr; //not in the paper K12 = exp(-5.495+0.1691*V); //(A.21) K21 = exp(-7.677-0.0128*V); //(A.22) XKr_inf= K12/(K12+K21); tauXKr = 1/(K12+K21) + 27; IKs = GKs*XKs^2*(V-EKs); //(A.23) EKs = RT_over_F*ln((K_o+0.01833*Na_o) /(K_i+0.01833*Na_i)); //(A.24) XKs:t = (XKs_inf-XKs)/tauXKs; //(A.25) XKs_inf = 1/(1+exp(-(V-24.7)/13.6)); //(A.26) tauXKs = if (V=10) 1/(0.0000719/0.148 +0.000131/0.0687) else 1/(0.0000719*(V-10)/ (1-exp(-0.148*(V-10))) +0.000131*(V-10)/ (exp(0.0687*(V-10))-1)); //(A.27) Ito1 = if (chf_flag = 0) IKv43 + IKv14K + IKv14Na else chfsc_IKv43*IKv43 + IKv14K + IKv14Na; //(A.28) IKv43 = GKv43 * OKv43 * (V-EK); //(A.29) IKv14K = if (V=0) PKv14/Csc*OKv14 * F * 1000 * (K_i*exp(VF_over_RT)-K_o) else PKv14/Csc*OKv14 * VFsqr_over_RT * (K_i*exp(VF_over_RT)-K_o) / (exp(VF_over_RT)-1); //(A.30) IKv14Na = if (V=0) alpha_NaK*PKv14/Csc*OKv14 * F * 1000 * (Na_i*exp(VF_over_RT)-Na_o) else alpha_NaK*PKv14/Csc*OKv14 * VFsqr_over_RT * (Na_i*exp(VF_over_RT)-Na_o) / (exp(VF_over_RT)-1); //(A.31) alpha_aKv43 = alpha_a0Kv43 * exp(a_aKv43*V); //(A.32) beta_aKv43 = beta_a0Kv43 * exp(-b_aKv43*V); //(A.33) alpha_iKv43 = alpha_i0Kv43 * exp(-a_iKv43*V); //(A.34) beta_iKv43 = beta_i0Kv43 * exp(b_iKv43*V); //(A.35) C0Kv43:t = alpha_iKv43 * CI0Kv43 + beta_aKv43 * C1Kv43 - (beta_iKv43 + 4*alpha_aKv43)*C0Kv43; //(A.36) C1Kv43:t = 4*alpha_aKv43 * C0Kv43 + alpha_iKv43/b1Kv43 * CI1Kv43 + 2*beta_aKv43 * C2Kv43 - (beta_aKv43 + f1Kv43*beta_iKv43 + 3*alpha_aKv43)*C1Kv43; //(A.37) C2Kv43:t = 3*alpha_aKv43 * C1Kv43 + alpha_iKv43/b2Kv43 * CI2Kv43 + 3*beta_aKv43 * C3Kv43 - (2*beta_aKv43 + f2Kv43*beta_iKv43 + 2*alpha_aKv43) * C2Kv43; //(A.38) C3Kv43:t = 2*alpha_aKv43 * C2Kv43 + alpha_iKv43/b3Kv43 * CI3Kv43 + 4*beta_aKv43 * OKv43 - (3*beta_aKv43 + f3Kv43*beta_iKv43 + alpha_aKv43) * C3Kv43; //(A.39) OKv43:t = alpha_aKv43 * C3Kv43 + alpha_iKv43/b4Kv43 * OIKv43 - (4*beta_aKv43 + f4Kv43*beta_iKv43)* OKv43; //(A.40) CI0Kv43:t = beta_iKv43 * C0Kv43 + beta_aKv43/f1Kv43 * CI1Kv43 - (alpha_iKv43 + b1Kv43*4*alpha_aKv43)* CI0Kv43; //(A.41) CI1Kv43:t = b1Kv43*4*alpha_aKv43 * CI0Kv43 + f1Kv43*beta_iKv43 * C1Kv43 + f1Kv43*2*beta_aKv43/f2Kv43 * CI2Kv43 - (beta_aKv43/f1Kv43 + alpha_iKv43/b1Kv43 + b2Kv43*3*alpha_aKv43/b1Kv43)* CI1Kv43; //(A.42) CI2Kv43:t = b2Kv43*3*alpha_aKv43/b1Kv43 * CI1Kv43 + f2Kv43*beta_iKv43 * C2Kv43 + f2Kv43*3*beta_aKv43/f3Kv43 * CI3Kv43 - (f1Kv43*2*beta_aKv43/f2Kv43 + alpha_iKv43/b2Kv43 + b3Kv43*2*alpha_aKv43/b2Kv43) * CI2Kv43; //(A.43) CI3Kv43:t = b3Kv43*2*alpha_aKv43/b2Kv43 * CI2Kv43 + f3Kv43*beta_iKv43 * C3Kv43 + f3Kv43*4*beta_aKv43/f4Kv43 * OIKv43 - (f2Kv43*3*beta_aKv43/f3Kv43 + alpha_iKv43/b3Kv43 + b4Kv43*alpha_aKv43/b3Kv43) * CI3Kv43; //(A.44) OIKv43:t = b4Kv43*alpha_aKv43/b3Kv43 * CI3Kv43 + f4Kv43*beta_iKv43 * OKv43 - (f3Kv43*4*beta_aKv43/f4Kv43 + alpha_iKv43/b4Kv43) * OIKv43; //(A.45) alpha_aKv14 = alpha_a0Kv14 * exp(a_aKv14*V); //(A.32') beta_aKv14 = beta_a0Kv14 * exp(-b_aKv14*V); //(A.33') alpha_iKv14 = alpha_i0Kv14 * exp(-a_iKv14*V); //(A.34') beta_iKv14 = beta_i0Kv14 * exp(b_iKv14*V); //(A.35') C0Kv14:t = alpha_iKv14 * CI0Kv14 + beta_aKv14 * C1Kv14 - (beta_iKv14 + 4*alpha_aKv14)*C0Kv14; //(A.36') C1Kv14:t = 4*alpha_aKv14 * C0Kv14 + alpha_iKv14/b1Kv14 * CI1Kv14 + 2*beta_aKv14 * C2Kv14 - (beta_aKv14 + f1Kv14*beta_iKv14 + 3*alpha_aKv14)*C1Kv14; //(A.37') C2Kv14:t = 3*alpha_aKv14 * C1Kv14 + alpha_iKv14/b2Kv14 * CI2Kv14 + 3*beta_aKv14 * C3Kv14 - (2*beta_aKv14 + f2Kv14*beta_iKv14 + 2*alpha_aKv14) * C2Kv14; //(A.38') C3Kv14:t = 2*alpha_aKv14 * C2Kv14 + alpha_iKv14/b3Kv14 * CI3Kv14 + 4*beta_aKv14 * OKv14 - (3*beta_aKv14 + f3Kv14*beta_iKv14 + alpha_aKv14) * C3Kv14; //(A.39') OKv14:t = alpha_aKv14 * C3Kv14 + alpha_iKv14/b4Kv14 * OIKv14 - (4*beta_aKv14 + f4Kv14*beta_iKv14)* OKv14; //(A.40') CI0Kv14:t = beta_iKv14 * C0Kv14 + beta_aKv14/f1Kv14 * CI1Kv14 - (alpha_iKv14 + b1Kv14*4*alpha_aKv14)* CI0Kv14; //(A.41') CI1Kv14:t = b1Kv14*4*alpha_aKv14 * CI0Kv14 + f1Kv14*beta_iKv14 * C1Kv14 + f1Kv14*2*beta_aKv14/f2Kv14 * CI2Kv14 - (beta_aKv14/f1Kv14 + alpha_iKv14/b1Kv14 + b2Kv14*3*alpha_aKv14/b1Kv14)* CI1Kv14; //(A.42') CI2Kv14:t = b2Kv14*3*alpha_aKv14/b1Kv14 * CI1Kv14 + f2Kv14*beta_iKv14 * C2Kv14 + f2Kv14*3*beta_aKv14/f3Kv14 * CI3Kv14 - (f1Kv14*2*beta_aKv14/f2Kv14 + alpha_iKv14/b2Kv14 + b3Kv14*2*alpha_aKv14/b2Kv14) * CI2Kv14; //(A.43') CI3Kv14:t = b3Kv14*2*alpha_aKv14/b2Kv14 * CI2Kv14 + f3Kv14*beta_iKv14 * C3Kv14 + f3Kv14*4*beta_aKv14/f4Kv14 * OIKv14 - (f2Kv14*3*beta_aKv14/f3Kv14 + alpha_iKv14/b3Kv14 + b4Kv14*alpha_aKv14/b3Kv14) * CI3Kv14; //(A.44') OIKv14:t = b4Kv14*alpha_aKv14/b3Kv14 * CI3Kv14 + f4Kv14*beta_iKv14 * OKv14 - (f3Kv14*4*beta_aKv14/f4Kv14 + alpha_iKv14/b4Kv14) * OIKv14; //(A.45') IK1 = if (chf_flag = 0) GK1*K1_inf*(K_o/(K_o+KmK1))*(V-EK) else chfsc_IK1 * GK1*K1_inf*(K_o/(K_o+KmK1))*(V-EK); //(A.46) K1_inf = 1/(2+exp(1.5/RT_over_F*(V-EK))); //(A.47) IKp = GKp*Kp*(V-EK); //(A.48) Kp = 1/(1+exp((7.488-V)/5.98)); //(A.49) INaCa = if (chf_flag = 0) kNaCa*5000/(KmNa^3+Na_o^3)/(KmCa+Ca_o) / (1+ksat*exp((eta-1)*VF_over_RT)) * (exp(eta*VF_over_RT)*Na_i^3*Ca_o -exp((eta-1)*VF_over_RT)*Na_o^3*Ca_i) else chfsc_INaCa * kNaCa*5000/(KmNa^3+Na_o^3)/(KmCa+Ca_o) / (1+ksat*exp((eta-1)*VF_over_RT)) * (exp(eta*VF_over_RT)*Na_i^3*Ca_o -exp((eta-1)*VF_over_RT)*Na_o^3*Ca_i); //(A.50) INaK = INaKmax*fNaK/(1+(KmNa_i/Na_i)^1.5) * K_o/(K_o+KmK_o); //(A.51) fNaK = 1/(1+0.1245*exp(-0.1*VF_over_RT) +0.0365*sigma*exp(-VF_over_RT)); //(A.52) sigma = (exp(Na_o/67.3)-1)/7; //(A.53) IpCa = IpCamax*Ca_i/(KmpCa+Ca_i); //(A.54) ICa_b = GCa_b*(V-ECa); //(A.55) ECa = 0.5*RT_over_F*ln(Ca_o/Ca_i); //(A.56) INa_b = GNa_b*(V-ENa); //(A.57) VV:t = if (clamp_flag = 0) I-(INa+ICa+ICa_K+IKr+IKs+Ito1+IK1+IKp +INaCa+INaK+IpCa+ICa_b+INa_b) else 0; //(A.58) V = if (clamp_flag = 0) VV else Vclamp; alpha = 0.416*exp((V+2)/10); //(A.59) beta = 0.049*exp(-(V+2)/13); //(A.60) alphap = alpha*a; //(A.61) betap = beta/b; //(A.62) gamma = 0.09233*Ca_ss; //(A.63) C0:t = beta*C1+omega*CCa0-(4*alpha+gamma)*C0; //(A.64) C1:t = 4*alpha*C0+2*beta*C2+omega/b*CCa1 - (beta+3*alpha+gamma*a)*C1; //(A.65) C2:t = 3*alpha*C1+3*beta*C3+omega/b^2*CCa2 - (2*beta+2*alpha+gamma*a^2)*C2; //(A.66) C3:t = 2*alpha*C2+4*beta*C4+omega/b^3*CCa3 - (3*beta+alpha+gamma*a^3)*C3; //(A.67) C4:t = alpha*C3+gL*O+omega/b^4*CCa4 - (4*beta+fL+gamma*a^4)*C4; //(A.68) O:t = fL*C4-gL*O; //(A.69) CCa0:t = betap*CCa1+gamma*C0-(4*alphap+omega)*CCa0; //(A.70) CCa1:t = 4*alphap*CCa0+2*betap*CCa2+gamma*a*C1 - (betap+3*alphap+omega/b)*CCa1; //(A.71) CCa2:t = 3*alphap*CCa1+3*betap*CCa3+gamma*a^2*C2 - (2*betap+2*alphap+omega/b^2)*CCa2; //(A.72) CCa3:t = 2*alphap*CCa2+4*betap*CCa4+gamma*a^3*C3 - (3*betap+alphap+omega/b^3)*CCa3; //(A.73) CCa4:t = alphap*CCa3+gp*OCa+gamma*a^4*C4 - (4*betap+fp+omega/b^4)*CCa4; //(A.74) OCa:t = fp*CCa4-gp*OCa; //(A.75) ICamax = if (V=0) PCa/Csc*4*F*1000 * (0.001*exp(2*VF_over_RT)-0.341*Ca_o) / 2 else PCa/Csc*4*VFsqr_over_RT * (0.001*exp(2*VF_over_RT)-0.341*Ca_o) / (exp(2*VF_over_RT)-1); //(A.76) ICa = ICamax*y*(O+OCa); //(A.77) ICa_K = if (V=0) PKp/Csc*y*(O+OCa)*F*1000 * (K_i*exp(VF_over_RT)-K_o) else PKp/Csc*y*(O+OCa)*VFsqr_over_RT * (K_i*exp(VF_over_RT)-K_o) /(exp(VF_over_RT)-1); //(A.78) // PKp = PK/(1+ICamax/ICahalf); //(A.79) PKp = if (ICamax > 0) PK else PK/(1+ICamax/ICahalf); // Version 2 y:t = (y_inf-y)/tauy; //(A.80) y_inf = 0.745/(1+exp((V+12.5)/5))+0.255; //(A.81) tauy = 1/(0.012/(0.5+exp(-V/4.5)) +0.0035*exp(-V/36.3)); //(A.82) dPC1dt = -kap*Ca_ss^ncoop*PC1+kam*PO1; PC1:t = dPC1dt; //(A.83) PO1:t = -(dPC1dt + dPO2dt + dPC2dt); //(A.84) dPO2dt = kbp*Ca_ss^mcoop*PO1-kbm*PO2; PO2:t = dPO2dt; //(A.85) dPC2dt = kcp*PO1-kcm*PC2; PC2:t = dPC2dt; //(A.86) Jrel = v1*(PO1+PO2)*(Ca_JSR-Ca_ss); //(A.87) fb = (Ca_i/Kfb)^Nfb; //(A.88) rb = (Ca_NSR/Krb)^Nrb; //(A.89) Jup = if (chf_flag = 0) KSR*(vmaxf*fb-vmaxr*rb)/(1+fb+rb) else chfsc_Jup * KSR*(vmaxf*fb-vmaxr*rb)/(1+fb+rb); //(A.90) Jtr = (Ca_NSR-Ca_JSR)/tautr; //(A.91) Jxfer = (Ca_ss-Ca_i)/tauxfer; //(A.92) Jtrpn = dHTRPNCadt + dLTRPNCadt; //(A.93) HTRPNCa:t = dHTRPNCadt; dHTRPNCadt = Khtrpnp*Ca_i*(HTRPNtot-HTRPNCa) - Khtrpnm*HTRPNCa; //(A.94) LTRPNCa:t = dLTRPNCadt; dLTRPNCadt = Kltrpnp*Ca_i*(LTRPNtot-LTRPNCa) - Kltrpnm*LTRPNCa; //(A.95) // Na_i:t = -(INa+INa_b+3*(INaK+INaCa)) // * Acap*Csc/(Vmyo*F*1000); //(A.96) Na_i:t = 0; // Version 2 tmp1 = Acap*Csc/(F*1000); K_i:t = -(IKr+IKs+Ito1+IK1+IKp+ICa_K-2*INaK) * tmp1/Vmyo; //(A.97) Ca_i:t = betai*(Jxfer-Jup-Jtrpn- (ICa_b-2*INaCa+IpCa) *tmp1/(2*Vmyo)); //(A.98) betai = 1/(1+CMDNtot*KmCMDN/(KmCMDN+Ca_i)^2); //(A.99) betass = 1/(1+CMDNtot*KmCMDN/(KmCMDN+Ca_ss)^2); //(A.100) betaJSR= 1/(1+CSQNtot*KmCSQN/(KmCSQN+Ca_JSR)^2); //(A.101) Ca_ss:t= betass*(Jrel*VJSR-Jxfer*Vmyo -ICa*tmp1/2)/Vss; //(A.102) Ca_JSR:t = betaJSR*(Jtr-Jrel); //(A.103) Ca_NSR:t = (Jup*Vmyo-Jtr*VJSR)/VNSR; //(A.104) Iion = INa+ICa+ICa_K+IKr+IKs+Ito1+IK1+IKp + INaCa+INaK+IpCa+ICa_b+INa_b; }