/* * Jafri-Rice-Winslow Ventricular Model 1998 * * Model Status * * This model is known to run in both PCEnv and COR, and has been * curated by Penny Noble of Oxford University. This variant also * contains an embedded CellML description of Niederer, Hunter * and Smith's quantitative model of cardiac myocyte regulation. * The reference for this paper is given below. * * Model Structure * * ABSTRACT: We construct a detailed mathematical model for Ca2+ * regulation in the ventricular myocyte that includes novel descriptions * of subcellular mechanisms based on recent experimental findings: * 1) the Keizer-Levine model for the ryanodine receptor (RyR), * which displays adaptation at elevated Ca2+; 2) a model for the * L-type Ca2+ channel that inactivates by mode switching; and * 3) a restricted subspace into which the RyRs and L-type Ca2+ * channels empty and interact via Ca2+. We add membrane currents * from the Luo-Rudy Phase II ventricular cell model to our description * of Ca2+ handling to formulate a new model for ventricular action * potentials and Ca2+ regulation. The model can simulate Ca2+ * transients during an action potential similar to those seen * experimentally. The subspace [Ca2+] rises more rapidly and reaches * a higher level (10-30 microM) than the bulk myoplasmic Ca2+ * (peak [Ca2+]i approximately 1 microM). Termination of sarcoplasmic * reticulum (SR) Ca2+ release is predominately due to emptying * of the SR, but is influenced by RyR adaptation. Because force * generation is roughly proportional to peak myoplasmic Ca2+, * we use [Ca2+]i in the model to explore the effects of pacing * rate on force generation. The model reproduces transitions seen * in force generation due to changes in pacing that cannot be * simulated by previous models. Simulation of such complex phenomena * requires an interplay of both RyR adaptation and the degree * of SR Ca2+ loading. This model, therefore, shows improved behavior * over existing models that lack detailed descriptions of subcellular * Ca2+ regulatory mechanisms. * * The original paper reference is cited below: * * Cardiac Calcium Dynamics: The Roles of Ryanodine Receptor Adaptation * and Sarcoplasmic Reticulum Load, M. Saleet Jafri, J. Jeremy * Rice and Raimond L. Winslow, 1998, Biophysical Journal, 74, * 1149-1168. PubMed ID: 9512016 * * The reference for the embedded Niederer Hunter Smith model of * cardiac myocyte relaxation is: "A Quantitative Analysis of Cardiac * Myocyte Relaxation: A Simulation Study" Niederer, S.A., Hunter, * P.J., Smith, N.P, Biophysical Journal, Volume 90, March 2006, * pp. 1697-1722. * * cell diagram of the Jafri-Rice-Winslow model showing ionic currents, * pumps and exchangers within the sarcolemma and the sarcoplasmic * reticulum * * [[Image file: jafri_rice_winslow_1998.png]] * * A schematic diagram describing the current flows across the * cell membrane that are captured in the Jafri-Rice-Winslow model. */ import nsrunit; unit conversion on; unit ms=.001 second^1; unit per_ms=1E3 second^(-1); unit per_mm=1E3 meter^(-1); unit mV=.001 kilogram^1*meter^2*second^(-3)*ampere^(-1); unit per_mV=1E3 kilogram^(-1)*meter^(-2)*second^3*ampere^1; unit per_mV_ms=1E6 kilogram^(-1)*meter^(-2)*second^2*ampere^1; unit mS_per_mm2=1E3 kilogram^(-1)*meter^(-4)*second^3*ampere^2; unit mm2=1E-6 meter^2; unit mM=1 meter^(-3)*mole^1; unit mM_per_ms=1E3 meter^(-3)*second^(-1)*mole^1; unit per_mM_per_ms=1E3 meter^3*second^(-1)*mole^(-1); unit per_mM3_per_ms=1E3 meter^9*second^(-1)*mole^(-3); unit per_mM4_per_ms=1E3 meter^12*second^(-1)*mole^(-4); unit uF_per_mm2=1 kilogram^(-1)*meter^(-4)*second^4*ampere^2; unit uA_per_mm2=1 meter^(-2)*ampere^1; unit uA_per_mmcu=1E3 meter^(-3)*ampere^1; unit mm_per_ms=1E3 meter^(-1); unit gas_constant_units=.001 kilogram^1*meter^2*second^(-2)*kelvin^(-1)*mole^(-1); unit faradays_constant_units=1 second^1*ampere^1*mole^(-1); // unit micrometer predefined unit N_per_mm2=1E6 kilogram^1*meter^(-1)*second^(-2); unit balancing=.001 meter^2*second^(-1); math main { realDomain time ms; time.min=0; extern time.max; extern time.delta; real lambda dimensionless; lambda=1; real dlambdadt per_ms; dlambdadt=0; real V(time) mV; when(time=time.min) V=-84.1638; real R gas_constant_units; R=8.3145e3; real T kelvin; T=310; real F faradays_constant_units; F=9.6845e4; real i_Na(time) uA_per_mm2; real i_Ca_L_Ca(time) uA_per_mm2; real i_Ca_L_K(time) uA_per_mm2; real i_K(time) uA_per_mm2; real i_K1(time) uA_per_mm2; real i_NaCa(time) uA_per_mm2; real i_Kp(time) uA_per_mm2; real i_p_Ca(time) uA_per_mm2; real i_Na_b(time) uA_per_mm2; real i_Ca_b(time) uA_per_mm2; real i_NaK(time) uA_per_mm2; real i_ns_Ca(time) uA_per_mm2; real Cm uF_per_mm2; Cm=0.01; real I_stim(time) uA_per_mm2; real stim_start ms; stim_start=100; real stim_end ms; stim_end=9000; real stim_period ms; stim_period=750; real stim_duration ms; stim_duration=1; real stim_amplitude uA_per_mm2; stim_amplitude=-100; real E_Na(time) mV; real g_Na mS_per_mm2; g_Na=0.128; real Nai(time) mM; when(time=time.min) Nai=10.2042; real Nao mM; Nao=140; real m(time) dimensionless; when(time=time.min) m=0.0328302; real h(time) dimensionless; when(time=time.min) h=0.988354; real j(time) dimensionless; when(time=time.min) j=0.99254; real alpha_m(time) per_ms; real beta_m(time) per_ms; real alpha_h(time) per_ms; real beta_h(time) per_ms; real alpha_j(time) per_ms; real beta_j(time) per_ms; real P_Ca mm_per_ms; P_Ca=33.75e-6; real P_K mm_per_ms; P_K=1e-9; real p_k(time) mm_per_ms; real i_Ca_L_Ca_half uA_per_mm2; i_Ca_L_Ca_half=-4.58e-3; real i_Ca_L_Ca_max(time) uA_per_mm2; real O(time) dimensionless; when(time=time.min) O=9.84546e-21; real O_Ca(time) dimensionless; when(time=time.min) O_Ca=0; real alpha(time) per_ms; real beta(time) per_ms; real gamma(time) per_ms; real alpha_a(time) per_ms; real beta_b(time) per_ms; real L_type_Ca_channel.a dimensionless; L_type_Ca_channel.a=2; real b dimensionless; b=2; real g per_ms; g=2; real f per_ms; f=0.3; real g_ per_ms; g_=0; real f_ per_ms; f_=0; real omega per_ms; omega=0.01; real C0(time) dimensionless; when(time=time.min) C0=0.997208; real C1(time) dimensionless; when(time=time.min) C1=6.38897e-5; real C2(time) dimensionless; when(time=time.min) C2=1.535e-9; real C3(time) dimensionless; when(time=time.min) C3=1.63909e-14; real C4(time) dimensionless; when(time=time.min) C4=6.56337e-20; real C_Ca0(time) dimensionless; when(time=time.min) C_Ca0=2.72826e-3; real C_Ca1(time) dimensionless; when(time=time.min) C_Ca1=6.99215e-7; real C_Ca2(time) dimensionless; when(time=time.min) C_Ca2=6.71989e-11; real C_Ca3(time) dimensionless; when(time=time.min) C_Ca3=2.87031e-15; real C_Ca4(time) dimensionless; when(time=time.min) C_Ca4=4.59752e-20; real Ca_SS(time) mM; when(time=time.min) Ca_SS=1.36058e-4; real Cao mM; Cao=1.8; real Ko(time) mM; when(time=time.min) Ko=5.4; real Ki(time) mM; when(time=time.min) Ki=143.727; real y(time) dimensionless; when(time=time.min) y=0.998983; real y_infinity(time) dimensionless; real tau_y(time) ms; real g_K(time) mS_per_mm2; real g_K_max mS_per_mm2; g_K_max=0.001128; real E_K(time) mV; real P_NaK dimensionless; P_NaK=0.01833; real X(time) dimensionless; when(time=time.min) X=0.000928836; real Xi(time) dimensionless; real alpha_X(time) per_ms; real beta_X(time) per_ms; real E_K1(time) mV; real g_K1(time) mS_per_mm2; real g_K1_max mS_per_mm2; g_K1_max=7.5e-3; real K1_infinity(time) dimensionless; real alpha_K1(time) per_ms; real beta_K1(time) per_ms; real E_Kp(time) mV; real g_Kp mS_per_mm2; g_Kp=8.28e-5; real Kp(time) dimensionless; real k_NaCa uA_per_mm2; k_NaCa=50; real K_mNa mM; K_mNa=87.5; real K_mCa mM; K_mCa=1.38; real k_sat dimensionless; k_sat=0.1; real eta dimensionless; eta=0.35; real Cai(time) mM; when(time=time.min) Cai=9.94893e-11; real K_mpCa mM; K_mpCa=0.5e-3; real I_pCa uA_per_mm2; I_pCa=1.15e-2; real g_Nab mS_per_mm2; g_Nab=1.41e-5; real E_NaN(time) mV; real g_Cab mS_per_mm2; g_Cab=6.032e-5; real E_CaN(time) mV; real I_NaK uA_per_mm2; I_NaK=0.013; real f_NaK(time) dimensionless; real K_mNai mM; K_mNai=10; real K_mKo mM; K_mKo=1.5; real sigma dimensionless; real i_ns_Na(time) uA_per_mm2; real i_ns_K(time) uA_per_mm2; real I_ns_Na(time) uA_per_mm2; real I_ns_K(time) uA_per_mm2; real K_m_ns_Ca mM; K_m_ns_Ca=1.2e-3; real P_ns_Ca mm_per_ms; P_ns_Ca=1.75e-9; real EnsCa(time) mV; real VnsCa(time) mV; real Am per_mm; Am=546.69; real V_myo dimensionless; V_myo=0.92; real RyR_open(time) dimensionless; real P_O1(time) dimensionless; when(time=time.min) P_O1=1.19168e-3; real P_O2(time) dimensionless; when(time=time.min) P_O2=6.30613e-9; real P_C1(time) dimensionless; when(time=time.min) P_C1=0.762527; real P_C2(time) dimensionless; when(time=time.min) P_C2=0.236283; real v1 per_ms; v1=1.8; real v2 per_ms; v2=0.58e-4; real v3 mM_per_ms; v3=1.8e-3; // Var below replaced by constant in model eqns to satisfy unit correction // real nCa dimensionless; // nCa=4; // Var below replaced by constant in model eqns to satisfy unit correction // real mCa dimensionless; // mCa=3; real k_a_plus per_mM4_per_ms; k_a_plus=1.215e10; real k_a_minus per_ms; k_a_minus=0.1425; real k_b_plus per_mM3_per_ms; k_b_plus=4.05e7; real k_b_minus per_ms; k_b_minus=1.93; real k_c_plus per_ms; k_c_plus=0.018; real k_c_minus per_ms; k_c_minus=0.0008; real k_htrpn_plus per_mM_per_ms; k_htrpn_plus=20; real k_htrpn_minus per_ms; k_htrpn_minus=0.066e-3; real k_ltrpn_plus per_mM_per_ms; k_ltrpn_plus=40; real k_ltrpn_minus per_ms; k_ltrpn_minus=0.04; real tau_tr ms; tau_tr=34.48; real Ca_JSR(time) mM; when(time=time.min) Ca_JSR=1.17504; real Ca_NSR(time) mM; when(time=time.min) Ca_NSR=1.243891; real V_JSR dimensionless; real V_NSR dimensionless; real V_SS dimensionless; real K_mup mM; K_mup=0.5e-3; real K_mCMDN mM; K_mCMDN=2.38e-3; real K_mCSQN mM; K_mCSQN=0.8; real tau_xfer ms; tau_xfer=3.125; real HTRPN_tot mM; HTRPN_tot=0.14; real LTRPN_tot mM; LTRPN_tot=0.07; real HTRPNCa(time) mM; when(time=time.min) HTRPNCa=0.13598; real LTRPNCa(time) mM; when(time=time.min) LTRPNCa=0.00635; real CSQN_tot mM; CSQN_tot=15; real CMDN_tot mM; CMDN_tot=0.05; real Bi(time) dimensionless; real B_SS(time) dimensionless; real B_JSR(time) dimensionless; real J_rel(time) mM_per_ms; real J_leak(time) mM_per_ms; real J_up(time) mM_per_ms; real J_tr(time) mM_per_ms; real J_xfer(time) mM_per_ms; real J_trpn(time) mM_per_ms; real J_htrpn(time) mM_per_ms; real J_ltrpn(time) mM_per_ms; real Tref N_per_mm2; Tref=56.2; real beta0 dimensionless; beta0=4.9; real Myofilaments.a dimensionless; Myofilaments.a=0.35; real Q1(time) dimensionless; when(time=time.min) Q1=0; real Q2(time) dimensionless; when(time=time.min) Q2=0; real Q3(time) dimensionless; when(time=time.min) Q3=0; real A1 dimensionless; A1=-29; real A2 dimensionless; A2=138; real A3 dimensionless; A3=129; real alpha1 dimensionless; alpha1=0.03; real alpha2 dimensionless; alpha2=0.13; real alpha3 dimensionless; alpha3=0.625; real Ca50ref mM; Ca50ref=1.05e-3; real zp dimensionless; zp=0.85; real beta1 dimensionless; beta1=-4; real alpha0 per_ms; alpha0=8e-3; real alphar1 per_ms; alphar1=2e-3; real alphar2 per_ms; alphar2=1.75e-3; real nRel dimensionless; nRel=3; real Kz dimensionless; Kz=0.15; real nHill dimensionless; nHill=3; real kon per_mM_per_ms; kon=100; real koff per_ms; koff=0.2; real gamma_trpn dimensionless; gamma_trpn=2; real TRPN_tot mM; TRPN_tot=0.07; real T0(time) N_per_mm2; real T0max N_per_mm2; real z(time) dimensionless; when(time=time.min) z=0; real z_max dimensionless; real Q(time) dimensionless; real Cab(time) mM; when(time=time.min) Cab=0; real Ca50 mM; real CaTRPN50 mM; real K_2 dimensionless; real K_1 dimensionless; real Tension(time) N_per_mm2; // // V:time=((i_Na+i_Ca_L_Ca+i_Ca_L_K+i_K+i_NaCa+i_K1+i_Kp+i_p_Ca+i_Na_b+i_Ca_b+i_NaK+i_ns_Ca+I_stim)/Cm); I_stim=(if (((time>=stim_start) and (time<=stim_end)) and ((time-stim_start-floor((time-stim_start)/stim_period)*stim_period)<=stim_duration)) stim_amplitude else (0 uA_per_mm2)); // i_Na=(g_Na*m^3*h*j*(V-E_Na)); E_Na=(R*T/F*ln(Nao/Nai)); // alpha_m=((.32 per_mV_ms)*(V+(47.13 mV))/(1-exp((-1)*(.1 per_mV)*(V+(47.13 mV))))); beta_m=((.08 per_ms)*exp((-1)*V/(11 mV))); m:time=(alpha_m*(1-m)-beta_m*m); // alpha_h=(if (V<((-1)*(40 mV))) (.135 per_ms)*exp(((80 mV)+V)/((-1)*(6.8 mV))) else (0 per_ms)); beta_h=(if (V<((-1)*(40 mV))) (3.56 per_ms)*exp((.079 per_mV)*V)+(3.1E5 per_ms)*exp((.35 per_mV)*V) else 1/((.13 ms)*(1+exp((V+(10.66 mV))/((-1)*(11.1 mV)))))); h:time=(alpha_h*(1-h)-beta_h*h); // alpha_j=(if (V<((-1)*(40 mV))) ((-1)*(127140 per_mV_ms)*exp((.2444 per_mV)*V)-(3.474E-5 per_mV_ms)*exp((-1)*(.04391 per_mV)*V))*(V+(37.78 mV))/(1+exp((.311 per_mV)*(V+(79.23 mV)))) else (0 per_ms)); beta_j=(if (V<((-1)*(40 mV))) (.1212 per_ms)*exp((-1)*(.01052 per_mV)*V)/(1+exp((-1)*(.1378 per_mV)*(V+(40.14 mV)))) else (.3 per_ms)*exp((-1)*(2.535E-7 per_mV)*V)/(1+exp((-1)*(.1 per_mV)*(V+(32 mV))))); j:time=(alpha_j*(1-j)-beta_j*j); // i_Ca_L_Ca=(i_Ca_L_Ca_max*y*(O+O_Ca)); i_Ca_L_K=((1 balancing)*p_k*y*(O+O_Ca)*V*F^2/(R*T)*(Ki*exp(V*F/(R*T))-Ko)/(exp(V*F/(R*T))-1)); p_k=(P_K/(1+i_Ca_L_Ca_max/i_Ca_L_Ca_half)); i_Ca_L_Ca_max=((1 balancing)*(P_Ca*4*V*F^2/(R*T)*((.001 mM)*exp(2*V*F/(R*T))-.341*Cao)/(exp(2*V*F/(R*T))-1))); alpha=((.4 per_ms)*exp((V+(12 mV))/(10 mV))); beta=((.05 per_ms)*exp((V+(12 mV))/((-1)*(13 mV)))); alpha_a=(alpha*L_type_Ca_channel.a); beta_b=(beta/b); gamma=((.1875 per_ms)*Ca_SS/(1 mM)); C0:time=(beta*C1+omega*C_Ca0-(4*alpha+gamma)*C0); C1:time=(4*alpha*C0+2*beta*C2+omega/b*C_Ca1-(beta+3*alpha+gamma*L_type_Ca_channel.a)*C1); C2:time=(3*alpha*C1+3*beta*C3+omega/b^2*C_Ca2-(beta*2+2*alpha+gamma*L_type_Ca_channel.a^2)*C2); C3:time=(2*alpha*C2+4*beta*C4+omega/b^3*C_Ca3-(beta*3+alpha+gamma*L_type_Ca_channel.a^3)*C3); C4:time=(alpha*C3+g*O+omega/b^4*C_Ca4-(beta*4+f+gamma*L_type_Ca_channel.a^4)*C4); O:time=(f*C4-g*O); C_Ca0:time=(beta_b*C_Ca1+gamma*C_Ca0-(4*alpha_a+omega)*C_Ca0); C_Ca1:time=(4*alpha_a*C_Ca0+2*beta_b*C_Ca2+gamma*L_type_Ca_channel.a*C1-(beta_b+3*alpha_a+omega/b)*C_Ca1); C_Ca2:time=(3*alpha_a*C_Ca1+3*beta_b*C_Ca3+gamma*L_type_Ca_channel.a^2*C2-(beta_b*2+2*alpha_a+omega/b^2)*C_Ca2); C_Ca3:time=(2*alpha_a*C_Ca2+4*beta_b*C_Ca4+gamma*L_type_Ca_channel.a^3*C3-(beta_b*3+alpha_a+omega/b^3)*C_Ca3); C_Ca4:time=(alpha_a*C_Ca3+g_*O_Ca+gamma*L_type_Ca_channel.a^4*C4-(beta_b*4+f_+omega/b^4)*C_Ca4); O_Ca:time=(f_*C_Ca4-g_*O_Ca); // y:time=((y_infinity-y)/tau_y); y_infinity=(1/(1+exp((V+(55 mV))/(7.5 mV)))+.1/(1+exp(((-1)*V+(21 mV))/(6 mV)))); tau_y=((20 ms)+(600 ms)/(1+exp((V+(30 mV))/(9.5 mV)))); // g_K=(g_K_max*sqrt(Ko/(5.4 mM))); E_K=(R*T/F*ln((Ko+P_NaK*Nao)/(Ki+P_NaK*Nai))); i_K=(g_K*Xi*X^2*(V-E_K)); // alpha_X=((7.19E-5 per_mV_ms)*(V+(30 mV))/(1-exp((-1)*(.148 per_mV)*(V+(30 mV))))); beta_X=((1.31E-4 per_mV_ms)*(V+(30 mV))/((-1)*1+exp((.0687 per_mV)*(V+(30 mV))))); X:time=(alpha_X*(1-X)-beta_X*X); // Xi=(1/(1+exp((V-(56.26 mV))/(32.1 mV)))); // g_K1=(g_K1_max*sqrt(Ko/(5.4 mM))); E_K1=(R*T/F*ln(Ko/Ki)); i_K1=(g_K1*K1_infinity*(V-E_K1)); // alpha_K1=((1.02 per_ms)/(1+exp((.2385 per_mV)*(V-E_K1-(59.215 mV))))); beta_K1=((.49124 per_ms)*(exp((.08032 per_mV)*(V+(5.476 mV)-E_K1))+exp((.06175 per_mV)*(V-(E_K1+(594.31 mV)))))/(1+exp((-1)*(.5143 per_mV)*(V-E_K1+(4.753 mV))))); K1_infinity=(alpha_K1/(alpha_K1+beta_K1)); // E_Kp=E_K1; Kp=(1/(1+exp(((7.488 mV)-V)/(5.98 mV)))); i_Kp=(g_Kp*Kp*(V-E_Kp)); // i_NaCa=(k_NaCa*1/(K_mNa^3+Nao^3)*1/(K_mCa+Cao)*1/(1+k_sat*exp((eta-1)*V*F/(R*T)))*(exp(eta*V*F/(R*T))*Nai^3*Cao-exp((eta-1)*V*F/(R*T))*Nao^3*Cai)); // i_p_Ca=(I_pCa*Cai/(K_mpCa+Cai)); // E_NaN=E_Na; i_Na_b=(g_Nab*(V-E_NaN)); // E_CaN=(R*T/(2*F)*ln(Cao/Cai)); i_Ca_b=(g_Cab*(V-E_CaN)); // f_NaK=(1/(1+.1245*exp((-1)*.1*V*F/(R*T))+.0365*sigma*exp((-1)*V*F/(R*T)))); sigma=(1/7*(exp(Nao/(67.3 mM))-1)); i_NaK=(I_NaK*f_NaK*1/(1+(K_mNai/Nai)^1.5)*Ko/(Ko+K_mKo)); // EnsCa=(R*T/F*ln((Ko+Nao)/(Ki+Nai))); VnsCa=(V-EnsCa); i_ns_Na=(I_ns_Na*1/(1+(K_m_ns_Ca/Cai)^3)); i_ns_K=(I_ns_K*1/(1+(K_m_ns_Ca/Cai)^3)); i_ns_Ca=(i_ns_Na+i_ns_K); I_ns_Na=((1 balancing)*P_ns_Ca*1^2*VnsCa*F^2/(R*T)*(.75*Nai*exp(VnsCa*F/(R*T))-.75*Nao)/(exp(VnsCa*F/(R*T))-1)); I_ns_K=((1 balancing)*P_ns_Ca*1^2*VnsCa*F^2/(R*T)*(.75*Ki*exp(VnsCa*F/(R*T))-.75*Ko)/(exp(VnsCa*F/(R*T))-1)); // V_SS=(5.828E-5*V_myo); V_NSR=(.081*V_myo); V_JSR=(.00464*V_myo); J_rel=(v1*RyR_open*(Ca_JSR-Ca_SS)); RyR_open=(P_O1+P_O2); P_C1:time=((-1)*k_a_plus*Ca_SS^4*P_C1+k_a_minus*P_O1); P_O1:time=(k_a_plus*Ca_SS^4*P_C1-(k_a_minus*P_O1+k_b_plus*Ca_SS^3*P_O1+k_c_plus*P_O1)+k_b_minus*P_O2+k_c_minus*P_C2); P_O2:time=(k_b_plus*Ca_SS^3*P_O1-k_b_minus*P_O2); P_C2:time=(k_c_plus*P_O1-k_c_minus*P_C2); J_leak=(v2*(Ca_NSR-Cai)); J_up=(v3*Cai^2/(K_mup^2+Cai^2)); J_tr=((Ca_NSR-Ca_JSR)/tau_tr); J_xfer=((Ca_SS-Cai)/tau_xfer); J_htrpn=(k_htrpn_plus*Cai*(HTRPN_tot-HTRPNCa)-k_htrpn_minus*HTRPNCa); J_ltrpn=(k_ltrpn_plus*Cai*(LTRPN_tot-LTRPNCa)-k_ltrpn_minus*LTRPNCa); J_trpn=(J_htrpn+J_ltrpn); HTRPNCa:time=J_htrpn; LTRPNCa:time=J_ltrpn; Bi=(1/(1+CMDN_tot*K_mCMDN/(K_mCMDN+Cai)^2)); B_SS=(1/(1+CMDN_tot*K_mCMDN/(K_mCMDN+Ca_SS)^2)); B_JSR=(1/(1+CSQN_tot*K_mCSQN/(K_mCSQN+Ca_JSR)^2)); Cai:time=(Bi*(J_leak+J_xfer-(J_up+J_trpn+(i_Ca_b-2*i_NaCa+i_p_Ca)*Am/(2*V_myo*F)))); Ca_SS:time=(B_SS*(J_rel*V_JSR/V_SS-J_xfer*V_myo/V_SS-i_Ca_L_Ca*Am/(2*V_SS*F))); Ca_JSR:time=(B_JSR*(J_tr-J_rel)); Ca_NSR:time=((J_up-J_leak)*V_myo/V_NSR-J_tr*V_JSR/V_NSR); // Nai:time=((-1)*(i_Na+i_Na_b+i_ns_Na+i_NaCa*3+i_NaK*3)*Am/(V_myo*F)); Ki:time=((-1)*(i_Ca_L_K+i_K+i_K1+i_Kp+i_ns_K+(-1)*i_NaK*2)*Am/(V_myo*F)); Ko:time=((i_Ca_L_K+i_K+i_K1+i_Kp+i_ns_K+(-1)*i_NaK*2)*Am/(V_myo*F)); // Cab:time=(kon*Cai*(TRPN_tot-Cab)-koff*(1-Tension/(gamma_trpn*Tref))*Cab); K_2=((1 ms)*(alphar2*zp^nRel/(zp^nRel+Kz^nRel)*(1-nRel*Kz^nRel/(zp^nRel+Kz^nRel)))); K_1=((1 ms)*(alphar2*zp^(nRel-1)*nRel*Kz^nRel/(zp^nRel+Kz^nRel)^2)); z_max=(alpha0*(1 ms)/(CaTRPN50/TRPN_tot)^nHill-K_2*(1 per_ms)/(alphar1+K_1*(1 per_ms)+alpha0/(CaTRPN50/TRPN_tot)^nHill)); Ca50=(Ca50ref*(1+beta1*(lambda-1))); CaTRPN50=(Ca50*TRPN_tot/(Ca50+koff/kon*(1-(1+beta0*(lambda-1))*.5/gamma_trpn))); z:time=(alpha0*(Cab/CaTRPN50)^nHill*(1-z)+(-1)*z*alphar1+(-1)*alphar2*z^nRel/(z^nRel+Kz^nRel)); T0max=(Tref*(1+beta0*(lambda-1))); T0=(T0max*z/z_max); Tension=(if (Q<0) T0*(Myofilaments.a*Q+1)/(1-Q) else T0*(1-(Myofilaments.a+2)*Q)/(1+Q)); Q1:time=(A1*dlambdadt-(1 per_ms)*(alpha1*Q1)); Q2:time=(A2*dlambdadt-(1 per_ms)*(alpha2*Q2)); Q3:time=(A3*dlambdadt-(1 per_ms)*(alpha3*Q3)); Q=(Q1+Q2+Q3); }