55 integer,
parameter ::
jna = 1
56 integer,
parameter ::
js = 2
57 integer,
parameter ::
jx1 = 3
58 integer,
parameter ::
jk1 = 4
61 integer,
parameter ::
ym = 1
62 integer,
parameter ::
yh = 2
63 integer,
parameter ::
yj = 3
64 integer,
parameter ::
yd = 4
65 integer,
parameter ::
yf = 5
66 integer,
parameter ::
yx1 = 6
67 integer,
parameter ::
ycai = 7
68 integer,
parameter ::
yv = 8
75 real(RP),
parameter ::
rtsf = 26.713760659695648_rp
77 real(RP),
parameter ::
cm = 1.0_rp
78 real(RP),
parameter ::
g_na = 4.0_rp
79 real(RP),
parameter ::
g_nac = 0.003_rp
80 real(RP),
parameter ::
g_s = 0.09_rp
81 real(RP),
parameter ::
ena = 50.0_rp
82 real(RP),
parameter ::
c_ix1 = exp(0.04_rp*(77._rp-35._rp))
89 real(RP),
dimension(BR_NY),
intent(out) :: Y
91 y(
ym) = 1.09819687232595410246167255529228788
e-0002_rp
92 y(
yh) = 0.987721175487584100912586864929795661_rp
93 y(
yj) = 0.974838138981649082196488362602052421_rp
94 y(
yd) = 2.97072466320784053758100507598202674
e-0003_rp
95 y(
yf) = 0.999981333893687946493220122499403334_rp
96 y(
yx1) = 5.62865057052210353702214502232095934
e-0003_rp
97 y(
ycai) = 1.78200721561980600333356651988516442
e-0007_rp
98 y(
yv) = -84.5737561222608724445622259686443887_rp
119 real(RP),
dimension(BR_NW),
intent(out) :: a, b
120 real(RP) ,
intent(in) :: V
122 b(
ym) = -(v + 47._rp)/(exp(-0.1_rp * (v + 47._rp)) - 1._rp)
123 a(
ym) = 40._rp * exp(-0.056_rp*(v + 72._rp))
125 b(
yh) = 0.126_rp * exp(-0.25_rp*(v + 77._rp))
126 a(
yh) = 1.7_rp / (exp(-0.082_rp*(v+22.5_rp))+1._rp)
128 b(
yj) = 0.055_rp*exp(-0.25_rp*(v + 78._rp))/&
129 & (exp(-0.2_rp*(v+78._rp))+1._rp)
130 a(
yj) = 0.3_rp / (exp( -0.1_rp*(v+32._rp))+1._rp)
132 b(
yd) = 0.095_rp*exp(-0.01_rp*(v-5._rp))/&
133 & (1._rp+exp(-0.072_rp*(v-5._rp)))
134 a(
yd) = 0.07_rp*exp(-0.017_rp*(v+44._rp))/&
135 & (1._rp+exp(0.05_rp*(v+44._rp)))
137 b(
yf) = 0.012_rp*exp(-0.008_rp*(v+28._rp))/&
138 & (1._rp+exp(0.15_rp*(v+28._rp)))
139 a(
yf) = 0.0065_rp*exp(-0.02_rp*(v+30._rp))/&
140 & (1._rp+exp(-0.2_rp*(v+30._rp)))
142 b(
yx1) = 0.0005_rp*exp(0.083_rp*(v+50._rp))/&
143 & (1._rp+exp(0.057_rp*(v+50._rp)))
144 a(
yx1) = 0.0013_rp*exp(-0.06_rp*(v+20._rp))/&
145 & (1._rp+exp(-0.04_rp*(v+20._rp)))
158 subroutine br_ab_0(a, b, I, Y, N, Na)
159 real(RP),
dimension(Na),
intent(out) :: a
160 real(RP),
dimension(N) ,
intent(out) :: b
161 real(RP) ,
intent(in) :: I
162 real(RP),
dimension(N) ,
intent(in) :: Y
163 integer ,
intent(in) :: N, Na
165 real(RP),
dimension(BR_NW) :: alpha
171 b(
ycai) = -1.
e-7_rp*b(
js) + 0.07_rp*(1.
e-7_rp - y(
ycai))
188 subroutine br_ab_w(a, b, I, Y, N, Na)
189 real(RP),
dimension(Na),
intent(out) :: a
190 real(RP),
dimension(N) ,
intent(out) :: b
191 real(RP) ,
intent(in) :: I
192 real(RP),
dimension(N) ,
intent(in) :: Y
193 integer ,
intent(in) :: N, Na
199 b(
ycai) = -1.
e-7_rp*b(
js) + 0.07_rp*(1.
e-7_rp - y(
ycai))
216 subroutine br_ab_wv(a, b, I, Y, N, Na)
217 real(RP),
dimension(Na),
intent(out) :: a
218 real(RP),
dimension(N) ,
intent(out) :: b
219 real(RP) ,
intent(in) :: I
220 real(RP),
dimension(N) ,
intent(in) :: Y
221 integer ,
intent(in) :: N, Na
226 b(
ycai) = -1.
e-7_rp*b(1) + 0.07_rp*(1.
e-7_rp - y(
ycai))
232 b(
yv) = (- b(
yv) + i ) /
cm 247 real(RP),
dimension(BR_NI),
intent(out) :: IList
248 real(RP),
dimension(BR_NY),
intent(in) :: Y
261 real(RP),
intent(in),
dimension(BR_NY) :: Y
273 real(RP),
dimension(BR_NY),
intent(in) :: Y
282 f = exp(0.04_rp * (y(
yv) + 85._rp) )
283 i_k1 = 4._rp*( f - 1._rp )
285 f = exp(0.04_rp * (y(
yv) + 53._rp) )
286 i_k1 = i_k1 / ( f + f**2 )
288 f = exp(-0.04_rp * (y(
yv) + 23._rp) )
289 i_k1 = i_k1 + 0.2_rp * ( y(
yv) + 23._rp) / ( 1._rp - f )
291 i_k1 = i_k1 * 0.35_rp
298 real(RP),
dimension(BR_NY),
intent(in) :: Y
301 i_s = -82.3_rp - 13.0287_rp*log( y(
ycai) )
302 i_s =
g_s * y(
yd) * y(
yf) * (y(
yv) - i_s)
310 real(RP),
dimension(BR_NY),
intent(in) :: Y
315 i_x1 =
c_ix1 - exp( -0.04_rp * (y(
yv) + 35._rp) )
316 i_x1 = i_x1 * y(
yx1) * 0.8_rp
324 real(RP),
dimension(BR_NY),
intent(in) :: Y
328 i_na = i_na * ( y(
yv) -
ena )
336 real(RP),
intent(in),
dimension(BR_NY) :: Y
345 real(RP),
intent(in) :: x, a, b
361 function i_s_sp(Y)
result(I_s)
362 real(RP),
dimension(BR_NY),
intent(in) :: Y
365 real(RP),
parameter :: g_S_modified = 0.01_rp
367 i_s = -82.3_rp - 13.0287_rp*log( y(
ycai) )
368 i_s = g_s_modified * y(
yd) * y(
yf) * (y(
yv) - i_s)
375 real(RP),
dimension(BR_NI),
intent(out) :: IList
376 real(RP),
dimension(BR_NY),
intent(in) :: Y
389 real(RP),
dimension(BR_NY) :: Y
391 y(
ym) = 1.06278410016131944945954637794827772
e-0002_rp
392 y(
yh) = 0.988691698520646794656617541431759221_rp
393 y(
yj) = 0.975989197046760783619339482293370610_rp
394 y(
yd) = 2.90835646511938874988339682604803243
e-0003_rp
395 y(
yf) = 0.999982195505595856678047191499808037_rp
396 y(
yx1) = 5.49126107877021382608497656435722875
e-0003_rp
397 y(
ycai) = 1.17485016892666208390185442837442507
e-0007_rp
398 y(
yv) =-84.8250996515381622468412766225817842_rp
408 real(RP),
dimension(Na),
intent(out) :: a
409 real(RP),
dimension(N) ,
intent(out) :: b
410 real(RP) ,
intent(in) :: I
411 real(RP),
dimension(N) ,
intent(in) :: Y
412 integer ,
intent(in) :: N, Na
418 b(
ycai) = -1.
e-7_rp*b(
js) + 0.07_rp*(1.
e-7_rp - y(
ycai))
subroutine, public br_sp_y0(Y)
Model IONIC_BR_SP: rest state
subroutine, public br_y0(Y)
re-computed rest state (solving F(Y0)=0)
real(rp), parameter c_ix1
integer, parameter, public rp
real(kind=RP) = real precision in the code REAL_TOL = epsilon to test real equality ...
real(rp) function i_na(Y)
integer, parameter, public br_nw
REAL NUMBERS PRECISION IN CHORAL: selects simple/double/quad
integer, parameter, public br_ni
Beeler Reuter ionic model
subroutine, public br_ab_w(a, b, I, Y, N, Na)
Model IONIC_BR: reaction term
real(rp) function i_x1(Y)
real(rp) function e(x, v1, v2)
subroutine, public br_sp_ab_w(a, b, I, Y, N, Na)
Model IONIC_BR_SP: reaction term
subroutine, public br_sp_ilist(IList, Y)
Model IONIC_BR_SP: ionic current list
real(rp) function i_s_sp(Y)
subroutine, public br_ilist(IList, Y)
integer, parameter, public br_ny
subroutine i_na_dv(a, b, Y)
subroutine, public br_ab_0(a, b, I, Y, N, Na)
Model IONIC_BR: reaction term
real(rp), parameter g_nac
subroutine itot_dv(a, b, Y)
real(rp) function ff(x, a, b)
real(rp) function i_k1(Y)
subroutine br_gating_variables(a, b, V)
Model IONIC_BR: Gating variables dw_i/dt = (w_i-w_{i,~})/{i} =1,,{ NW} 1/ -w_{i,~})/
subroutine, public br_ab_wv(a, b, I, Y, N, Na)
Model IONIC_BR: reaction term