Converted Kernik 2019 ionic model does not work properly

Question

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);

Answer 1

The problem has been solved recently. The thing is that units of measurement were initially annotated after each line, for example:

Buf_C = 0.06; .units(mol/m^3);

Deleting them fixed the problem for some reason.