Hello!
I have some difficulties with Kernik model generated from the .cellml file (source: «A computational model of induced pluripotent stem-cell derived cardiomyocytes incorporating experimental variability from multiple data sources»). File .model is generated on the openCARP v12.0.
The resulting action potential does not similar to the shape of potential in the research paper. Also, the excitation does not propagate beyond the stimulation area on a 2d or 3d object (regardless of the strength of the stimulation or conductivities).
Could you tell me what the mistake might be?
The kernik_2019.model is attached below
#generated from: kernik_2019.cellml
V; .external(Vm); .nodal();
Iion; .external(); .nodal();
# Constants
Cm = 60.0;
#(There are constants here, they are missing due to message size restrictions.)
# Initial values
#Analogically
V_init = -7.55966016388546791e+01;
#Cell
V_tot_tenT = (Vc_tenT+VSR_tenT);
Vc = (V_tot*(Vc_tenT/V_tot_tenT));
V_SR = (V_tot*(VSR_tenT/V_tot_tenT));
#erev
E_Ca = 0.5*RTF*log(Cao/Cai);
E_K = RTF*log(Ko/Ki);
E_Na = RTF*log(Nao/Nai);
#nai
diff_Nai = (-Cm/(F*Vc))*(i_Na+i_b_Na+i_fNa+3.0*i_NaK+3.0*i_NaCa+i_CaL_Na);
#casr
diff_Ca_SR = (((Ca_SR_bufSR*Vc)/V_SR)*(i_up - i_rel - i_leak));
Ca_SR_bufSR = (1.0/(1.0+((Buf_SR*Kbuf_SR)/pow((Ca_SR+Kbuf_SR),2.0))));
#ibca
g_b_Ca = (0.000592*0.62);
i_b_Ca = (g_b_Ca*(V - E_Ca));
#ibna
g_b_Na = (0.00029*1.5);
i_b_Na = (g_b_Na*(V - E_Na));
#cai
Cai_bufc = 1.0/(1.0+(Buf_C*Kbuf_C)/pow((Cai+Kbuf_C),2.0));
diff_Cai = Cai_bufc*(-i_up+i_leak+i_rel - (Cm/(2.0*Vc*F))*(i_CaL_Ca+i_CaT+i_b_Ca+i_PCa - 2.0*i_NaCa));
#ical
d3 = (d5*d1);
d4 = (1.0/((1.0/d2)+(1.0/d6)));
alpha_ical_d = (d1*exp((V/d2)));
beta_ical_d= (d3*exp((V/d4)));
ical_d_inf = (alpha_ical_d/(alpha_ical_d+beta_ical_d));
tau_ical_d = ((1.0/(alpha_ical_d+beta_ical_d))+taud_const);
diff_d = ((ical_d_inf - d)/tau_ical_d);
f3 = (f5*f1);
f4 = (1.0/((1.0/f2)+(1.0/f6)));
alpha_ical_f = (f1*exp((V/f2)));
beta_ical_f = (f3*exp((V/f4)));
ical_f_inf = (alpha_ical_f/(alpha_ical_f+beta_ical_f));
tau_ical_f = ((1.0/(alpha_ical_f+beta_ical_f))+tauf_const);
diff_f = ((ical_f_inf - f)/tau_ical_f);
i_cai_fCa_gate_alpha_fCa = (1.0/(1.0+pow(((scale*Cai)/0.000325),8.0)));
i_cai_fCa_gate_beta_fCa = (0.1/(1.0+exp((((scale*Cai) - 0.0005)/0.0001))));
i_cai_fCa_gate_gamma_fCa = (0.2/(1.0+exp((((scale*Cai) - 0.00075)/0.0008))));
fCa_inf = ((i_cai_fCa_gate_alpha_fCa+i_cai_fCa_gate_beta_fCa+i_cai_fCa_gate_gamma_fCa+0.23)/1.46);
tau_fCa= 2;
k_fca = ((fCa_inf>fCa and V>-60.0) ? 0.0 : 1.0);
diff_fCa = (k_fca * ( fCa_inf - fCa ) / tau_fCa);
i_CaL = (i_CaL_Ca+i_CaL_Na+i_CaL_K);
i_CaL_Ca = (ibarca*d*f*fCa);
i_CaL_K = (ibark*d*f*fCa);
i_CaL_Na = (ibarna*d*f*fCa);
ibarca = (((p_CaL_Ca*4.0*V*FFRT)*(((0.341*Cai)*exp(((2.0*V)*FRT))) - (0.341*Cao)))/(exp(((2.0*V)*FRT)) - 1.0));
ibark = (((p_CaL_K*V*FFRT)*(((0.75*Ki)*exp((V*FRT))) - (0.75*Ko)))/(exp((V*FRT)) - 1.0));
ibarna = (((p_CaL_Na*V*FFRT)*(((0.75*Nai)*exp((V*FRT))) - (0.75*Nao)))/(exp((V*FRT)) - 1.0));
p_CaL_Ca = (p_CaL_shannonCap*p_CaL);
p_CaL_K = (p_CaL_shannonKp*p_CaL);
p_CaL_Na = (p_CaL_shannonNap*p_CaL);
p_CaL_shannonCap = (p_CaL_shannonCa/p_CaL_shannonTot);
p_CaL_shannonKp = (p_CaL_shannonK/p_CaL_shannonTot);
p_CaL_shannonNap = (p_CaL_shannonNa/p_CaL_shannonTot);
p_CaL_shannonTot = ((p_CaL_shannonCa+p_CaL_shannonNa)+p_CaL_shannonK);
#icat
d_cat_inf = (1.0/(1.0+exp(((V+26.3)/-6.0))));
tau_d_cat = (1.0/((1.068*exp(((V+26.3)/30.0)))+(1.068*exp(((V+26.3)/-30.0)))));
diff_d_cat = ((d_cat_inf - d_cat)/tau_d_cat);
f_cat_inf = (1.0/(1.0+exp(((V+61.7)/5.6))));
tau_f_cat = (1.0/((0.0153*exp((-(V+61.7)/83.3)))+(0.015*exp(((V+61.7)/15.38)))));
diff_f_cat = ((f_cat_inf - f_cat)/tau_f_cat);
i_CaT = ((g_CaT*d_cat*f_cat)*(V - E_Ca));
#ifunny
xF3 = (xF5*xF1);
xF4 = (1.0/((1.0/xF2)+(1.0/xF6)));
alpha_ifunny_Xf = (xF1*exp((V/xF2)));
beta_ifunny_Xf = (xF3*exp((V/xF4)));
ifunny_Xf_inf = (alpha_ifunny_Xf/(alpha_ifunny_Xf+beta_ifunny_Xf));
tau_ifunny_Xf = ((1.0/(alpha_ifunny_Xf+beta_ifunny_Xf))+xF_const);
diff_Xf = ((ifunny_Xf_inf - Xf)/tau_ifunny_Xf);
Na_frac = (NatoK_ratio/(NatoK_ratio+1.0));
i_fNa = (((Na_frac*g_f)*Xf)*(V - E_Na));
i_fK = ((((1.0 - Na_frac)*g_f)*Xf)*(V - E_K));
i_f = (i_fNa+i_fK);
#ik1
alpha_xK1 = (xK11*exp(((V+xK13)/xK12)));
beta_xK1 = (1.0*exp(((V+xK15)/xK14)));
xK1_inf = (alpha_xK1/(alpha_xK1+beta_xK1));
i_K1 = ((g_K1*sqrt(Ko/5.4)*xK1_inf)*(V - E_K));
#ikr
Xr1_3 = (Xr1_5*Xr1_1);
Xr2_3 = (Xr2_5*Xr2_1);
Xr1_4 = (1.0/((1.0/Xr1_2)+(1.0/Xr1_6)));
Xr2_4 = (1.0/((1.0/Xr2_2)+(1.0/Xr2_6)));
alpha_ikr_Xr1 = (Xr1_1*exp((V/Xr1_2)));
beta_ikr_Xr1 = (Xr1_3*exp((V/Xr1_4)));
ikr_Xr1_inf = (alpha_ikr_Xr1/(alpha_ikr_Xr1+beta_ikr_Xr1));
tau_ikr_Xr1 = ((1.0/(alpha_ikr_Xr1+beta_ikr_Xr1))+tau_1_offset);
diff_Xr1 = ((ikr_Xr1_inf - Xr1)/tau_ikr_Xr1);
alpha_ikr_Xr2 = (Xr2_1*exp((V/Xr2_2)));
beta_ikr_Xr2 = (Xr2_3*exp((V/Xr2_4)));
ikr_Xr2_inf = (alpha_ikr_Xr2/(alpha_ikr_Xr2+beta_ikr_Xr2));
tau_ikr_Xr2 = ((1.0/(alpha_ikr_Xr2+beta_ikr_Xr2))+tau_2_offset);
diff_Xr2 = ((ikr_Xr2_inf - Xr2)/tau_ikr_Xr2);
i_Kr = ((g_Kr*sqrt(Ko/5.4)*Xr1*Xr2)*(V - E_K));
#iks
ks3 = (ks5*ks1);
ks4 = (1.0/((1.0/ks2)+(1.0/ks6)));
alpha_iks_Xs = (ks1*exp((V/ks2)));
beta_iks_Xs = (ks3*exp((V/ks4)));
iks_Xs_inf = (alpha_iks_Xs/(alpha_iks_Xs+beta_iks_Xs));
tau_iks_Xs = ((1.0/(alpha_iks_Xs+beta_iks_Xs))+tauks_const);
diff_Xs = ((iks_Xs_inf - Xs)/tau_iks_Xs);
i_Ks = ((g_Ks*pow(Xs,2.0))*(V - E_K));
#ileak
i_leak = ((Ca_SR - Cai)*V_leak);
#ina
m3 = (m5*m1);
m4 = (1.0/((1.0/m2)+(1.0/m6)));
h3 = (h5*h1);
h4 = (1.0/((1.0/h2)+(1.0/h6)));
j5 = h5;
j6 = h6;
j3 = (j5*j1);
j4 = (1.0/((1.0/j2)+(1.0/j6)));
alpha_ina_h = (h1*exp((V/h2)));
beta_ina_h = (h3*exp((V/h4)));
ina_h_inf = (alpha_ina_h/(alpha_ina_h+beta_ina_h));
tau_ina_h = ((1.0/(alpha_ina_h+beta_ina_h))+tau_h_const);
alpha_ina_j = (j1*exp((V/j2)));
beta_ina_j = (j3*exp((V/j4)));
ina_j_inf = (alpha_ina_j/(alpha_ina_j+beta_ina_j));
tau_ina_j = ((1.0/(alpha_ina_j+beta_ina_j))+tau_j_const);
alpha_ina_m = (m1*exp((V/m2)));
beta_ina_m = (m3*exp((V/m4)));
ina_m_inf = (alpha_ina_m/(alpha_ina_m+beta_ina_m));
tau_ina_m = ((1.0/(alpha_ina_m+beta_ina_m))+tau_m_const);
diff_h = ((ina_h_inf - h)/tau_ina_h);
diff_j = ((ina_j_inf - j)/tau_ina_j);
diff_m = ((ina_m_inf - m)/tau_ina_m);
i_Na = (g_Na*m*m*m*h*j)*(V - E_Na);
#inaca
i_NaCa = ((kNaCa*(((exp(gamma*V*FRT)*pow(Nai,3.0))*Cao) - (((exp((((gamma - 1.0)*V)*FRT))*pow(Nao,3.0))*Cai)*alpha)))/(((pow(KmNai,3.0)+pow(Nao,3.0))*(KmCa+Cao))*(1.0+(Ksat*exp((((gamma - 1.0)*V)*FRT))))));
#inak
i_NaK = ((PNaK*Ko*Nai)/(((Ko+Km_K)*(Nai+Km_Na))*((1.0+(0.1245*exp(((-0.1*V)*FRT))))+(0.0353*exp((-V*FRT))))));
#ipca
g_PCa = (0.025*10.5);
KPCa = 0.0005;
i_PCa = ((g_PCa*Cai)/(Cai+KPCa));
#irel
diff_R = ((kim*RI - kiSRCa*Cai*R - ((koSRCa*pow(Cai,2.0))*R))+(kom*O));
diff_O = (((((koSRCa*pow(Cai,2.0))*R) - (kom*O)) - ((kiSRCa*Cai)*O))+(kim*I));
diff_I = (((kiSRCa*Cai*O - kim*I) - (kom*I))+((koSRCa*pow(Cai,2.0))*RI));
i_rel = (((ks*O)*(Ca_SR - Cai))*(V_SR/Vc));
kCaSR = (MaxSR - ((MaxSR - MinSR)/(1.0+pow((ec50SR/Ca_SR),2.5))));
kiSRCa = (kiCa*kCaSR);
koSRCa = (koCa/kCaSR);
RI = 1.0 - R - O - I;
#ito
r3 = (r5*r1);
r4 = (1.0/((1.0/r2)+(1.0/r6)));
s3 = (s5*s1);
s4 = (1.0/((1.0/s2)+(1.0/s6)));
alpha_ito_s = (s1*exp((V/s2)));
beta_ito_s = (s3*exp((V/s4)));
ito_s_inf = (alpha_ito_s/(alpha_ito_s+beta_ito_s));
tau_ito_s = ((1.0/(alpha_ito_s+beta_ito_s))+tau_s_const);
diff_s = ((ito_s_inf - s)/tau_ito_s);
alpha_ito_r = (r1*exp((V/r2)));
beta_ito_r = (r3*exp((V/r4)));
ito_r_inf = (alpha_ito_r/(alpha_ito_r+beta_ito_r));
tau_ito_r = ((1.0/(alpha_ito_r+beta_ito_r))+tau_r_const);
diff_r = ((ito_r_inf - r)/tau_ito_r);
i_to = ((g_to*r*s)*(V - E_K));
#iup
i_up = (VmaxUp/(1.0+(pow(Kup,2.0)/pow(Cai,2.0))));
#ki
diff_Ki = (-Cm/(F*Vc))*(i_K1+i_to+i_Kr+i_Ks+i_fK +i_CaL_K- 2.0*i_NaK);
#membrane
Iion = i_K1+i_Na+i_to+i_Kr+i_Ks+i_CaL+i_CaT+i_NaK+i_NaCa+i_PCa+i_f+i_b_Na+i_b_Ca;
#phys
FFRT = (F*FRT);
FRT = (F/(R_phys*T));
RTF = ((R_phys*T)/F);
group {
V;
Iion;
}.trace();
group {
Ki;
Nai;
Cai;
Ca_SR;
d;
f;
fCa;
d_cat;
f_cat;
Xf;
Xr1;
Xr2;
Xs;
h;
j;
m;
I;
O;
R;
r;
s;
}.method(rush_larsen);
