OpenCMISS-Iron Internal API Documentation
electrophysiology_cell_routines.f90
Go to the documentation of this file.
1 
43 
45  USE base_routines
46  use field_routines
47  use kinds
48  use strings
49  use types
50 
51 #include "macros.h"
52 
53  implicit none
54  private
55 
56  !interfaces
57 
60 
61  interface pow
62  module procedure powr
63  module procedure powi
64  end interface pow
65 
66 contains
67  ! for auto generated code
68 
69  real(dp) function powr(a,b)
70  real(dp), intent(in) :: a,b
71  powr = a**b
72  end function powr
73  real(dp) function powi(a,b)
74  real(dp), intent(in) :: a
75  integer(intg), intent(in) :: b
76  powi = a**b
77  end function powi
78 
79  ! Bueno-Orovio (2008) minimal cell model. Membrane potential rescaled, partial evaluation applied.
80 
81  ! initialize field
82  subroutine bueno_orovio_initialize(field,err,error,*)
83  type(field_type), intent(inout), pointer :: field
84  INTEGER(INTG), INTENT(OUT) :: ERR
85  TYPE(varying_string), INTENT(OUT) :: ERROR
86 
87  integer(intg) :: i
88  REAL(DP), PARAMETER, DIMENSION(1:4) :: y0 = (/ -8.09242e+01, 9.99993e-01, 2.16366e-02, 9.84320e-01 /) ! paced initial condition
89 
90  enters('bueno_orovio_init',err,error,*999)
91  DO i=1,4
92  CALL field_component_values_initialise(field,field_u_variable_type,field_values_set_type,i,y0(i),err,error,*999)
93  END DO
94 
95  exits('bueno_orovio_init')
96  return
97 999 errorsexits('bueno_orovio_init',err,error)
98  return 1
99 
100  end subroutine bueno_orovio_initialize
101 
102  subroutine bueno_orovio_dydt(t,y,dydt,activ)
103  real(dp), dimension(:), intent(out) :: dydt(:)
104  real(dp), dimension(:), intent(in) :: y(:)
105  real(dp), intent(in) :: t
106  real(dp), intent(in) :: activ
107 
108  real(dp),parameter :: u_m=3.0e-01, u_p=1.3e-01, u_q=6.0e-03, u_r=6.0e-03, vscale=81.0, period = 1000.0
109  real(dp) :: J_fi, J_si, J_so, Vm, bv, m, p, q, r, tau_o, tau_s, tau_so, tau_v_minus, tau_w_minus, w_inf, J_stim
110 
111  if (mod(t,period)>=0 .and. mod(t,period)<=0.999999) then ! apply stimulus in first ms
112  j_stim = -0.65 * activ ! stimulus: ~52 mV for 1 ms
113  else
114  j_stim = 0
115  end if
116 
117  bv = (1+(y(1)/vscale));
118  if (bv < u_m) then
119  m = 0
120  else
121  m = 1
122  end if
123  if (bv < u_q) then
124  q = 0
125  else
126  q = 1
127  end if
128  if (bv < u_p) then
129  p = 0
130  else
131  p = 1
132  end if
133 
134  if (bv < u_r) then
135  r = 0
136  else
137  r = 1
138  end if
139 
140  j_fi = (9.09090909090909e+00*(-m)*y(2)*(-3.0e-01+bv)*(1.55e+00-bv))
141  tau_v_minus = ((1150*q)+(60*(1-q)))
142  dydt(2) = (((1-m)*(0-y(2))/tau_v_minus)-(6.89655172413793e-01*m*y(2)))
143  tau_o = ((400*(1-r))+(6*r))
144  tau_so = (3.002e+01+(-1.4512e+01*(1+tanh((2.046e+00*(-6.5e-01+bv))))))
145  j_so = ((bv*(1-p)/tau_o)+(p/tau_so))
146  j_si = (5.29801324503311e-01*(-p)*y(4)*y(3))
147  dydt(1) = ((-j_fi-j_so-j_si-j_stim)*vscale)
148  vm = (-83+(8.57e+01*bv))
149  tau_s = ((2.7342e+00*(1-p))+(16*p))
150  dydt(3) = (((5.0e-01*(1+tanh((2.0994e+00*(-9.087e-01+bv)))))-y(3))/tau_s)
151  w_inf = (((1-r)*(1-(1.42857142857143e+01*bv)))+(9.4e-01*r))
152  tau_w_minus = (60+(-2.25e+01*(1+tanh((65*(-3.0e-02+bv))))))
153  dydt(4) = (((1-r)*(w_inf-y(4))/tau_w_minus)-(5.0e-03*r*y(4)))
154  end subroutine bueno_orovio_dydt
155 
156  ! integrates all cell models
157  subroutine bueno_orovio_integrate(cells,materials,t0,t1,err,error,*)
158  type(field_type), intent(inout), pointer :: cells, materials
159  real(dp), intent(in) :: t0, t1
160  INTEGER(INTG), INTENT(OUT) :: ERR
161  TYPE(varying_string), INTENT(OUT) :: ERROR
162 
163  integer, parameter :: celldim = 4
164  real(dp), dimension(1:celldim) :: y, dydt
165  real(dp) :: t, dt, activ
166  real(dp), dimension(:), pointer :: celldata, activdata
167 
168  integer(intg) :: ncells, i, d
169  type(domain_ptr_type), pointer :: domain
170  TYPE(field_variable_type), POINTER :: CELLS_VARIABLE, ACTIV_VARIABLE
171 
172  enters('bueno_orovio_integrate',err,error,*999)
173 
174  domain=>cells%decomposition%domain(1)
175  ncells = domain%ptr%topology%nodes%number_of_nodes ! local
176 
177  cells_variable=>cells%VARIABLE_TYPE_MAP(field_u_variable_type)%PTR
178  activ_variable=>materials%VARIABLE_TYPE_MAP(field_u_variable_type)%PTR
179 
180  CALL field_parameter_set_data_get(cells,field_u_variable_type,field_values_set_type,celldata,err,error,*999)
181  CALL field_parameter_set_data_get(materials,field_u_variable_type,field_values_set_type,activdata,err,error,*999)
182 
183  do i=1,ncells
184  ! field -> y
185  do d=1,celldim
186  !Default to version 1 of each node derivative
187  y(d) = celldata(cells_variable%COMPONENTS(d)%PARAM_TO_DOF_MAP%NODE_PARAM2DOF_MAP%NODES(i)%DERIVATIVES(1)%VERSIONS(1))
188  end do
189  !Default to version 1 of each node derivative
190  activ = activdata(activ_variable%COMPONENTS(1)%PARAM_TO_DOF_MAP%NODE_PARAM2DOF_MAP%NODES(i)%DERIVATIVES(1)%VERSIONS(1))
191 
192  t = t0
193  do while (t < t1 - 1e-6)
194  ! integrate one cell, one time step.
195  ! just use adaptive forward euler, tested in matlab
196  call bueno_orovio_dydt(t,y,dydt,activ)
197  dt = min(t1-t,0.5) ! default max
198  dt = min(dt,1 / abs(dydt(1))) ! maximum increase in Vm
199 ! IF(DT < 1e-6) WRITE(*,*) 'CELL ',i, 'DVDT = ',dydt(1),'dt = ',dt,' activ=',activ,'y=',y
200  y = y + dt * dydt ! FWE
201  t = t + dt
202  end do
203  ! y -> field
204  do d=1,celldim
205  ! call field_parameter_set_update_local_node(cells,field_u_variable_type,field_values_set_type,1,i,d,y(d), err,error,*999)
206  !Default to version 1 of each node derivative
207  celldata(cells_variable%COMPONENTS(d)%PARAM_TO_DOF_MAP%NODE_PARAM2DOF_MAP%NODES(i)%DERIVATIVES(1)%VERSIONS(1)) = y(d)
208  end do
209  end do
210  CALL field_parameter_set_data_restore(cells,field_u_variable_type,field_values_set_type,celldata,err,error,*999)
211  CALL field_parameter_set_data_restore(materials,field_u_variable_type,field_values_set_type,activdata,err,error,*999)
212 
213 
214  exits('bueno_orovio_integrate')
215  return
216 999 errorsexits('bueno_orovio_integrate',err,error)
217  return 1
218  end subroutine bueno_orovio_integrate
219 
220 
221 
222  ! Ten Tusscher and Panfilov (2006) human ventricular cell model
223 
224  ! initialize field
225  subroutine tentusscher06_initialize(field,err,error,*)
226  type(field_type), intent(inout), pointer :: field
227  INTEGER(INTG), INTENT(OUT) :: ERR
228  TYPE(varying_string), INTENT(OUT) :: ERROR
229 
230  integer(intg) :: i
231  REAL(DP), PARAMETER, DIMENSION(1:19) :: y0 = (/ -85.23, 0.9755, 0.9953, 0.7888, 3.64, 0.000126, 0.00036, 0.9073, 0.7444,&
232  & 0.7045, 0.00172, 3.373e-5, 136.89, 0.00621, 0.4712, 0.0095, 8.604, 2.42e-8, 0.999998 /) ! paced initial condition
233 
234  enters('tentusscher06_init',err,error,*999)
235  DO i=1,19
236  CALL field_component_values_initialise(field,field_u_variable_type,field_values_set_type,i,y0(i),err,error,*999)
237  END DO
238 
239  exits('tentusscher06_init')
240  return
241 999 errorsexits('tentusscher06_init',err,error)
242  return 1
243 
244  end subroutine tentusscher06_initialize
245 
246  subroutine tentusscher06_dydt(t,y,dydt,activ)
247  real(dp), dimension(:), intent(out) :: dydt(:)
248  real(dp), dimension(:), intent(in) :: y(:)
249  real(dp), intent(in) :: t
250  real(dp), intent(in) :: activ
251 
252  real(dp),parameter :: period=1000
253  real(dp) :: Ca_i_bufc, Ca_sr_bufsr, Ca_ss_bufss, E_Ca, E_K, E_Ks, E_Na, O, alpha_K1, alpha_d, alpha_m, alpha_xr1, alpha_xr2,&
254  & alpha_xs, beta_K1, beta_d, beta_m, beta_xr1, beta_xr2, beta_xs, d_inf, f2_inf, fCass_inf, f_inf, gamma_d, h_inf, i_CaL,&
255  & i_K1, i_Kr, i_Ks, i_Na, i_NaCa, i_NaK, i_Stim, i_b_Ca, i_b_Na, i_leak, i_p_Ca, i_p_K, i_rel, i_to, i_up, i_xfer, j_inf, k1,&
256  & k2, kcasr, m_inf, r_inf, s_inf, tau_d, tau_f, tau_f2, tau_fCass, tau_h, tau_j, tau_m, tau_r, tau_s, tau_xr1, tau_xr2,tau_xs,&
257  & xK1_inf, xr1_inf, xr2_inf, xs_inf
258 
259  if (mod(t,period)>=0 .and. mod(t,period)<=1.999999) then ! apply stimulus in 2 ms
260  i_stim = -3.57142857142857e+01 * activ
261  else
262  i_stim = 0
263  end if
264 
265  i_cal = (5.75002020961245e-01*y(12)*y(4)*y(2)*y(3)*(-15+y(1))*(-2+(2.5e-01*y(7)*exp((7.48677816454909e-02*(-15+y(1))))))/&
266  & (-1+exp((7.48677816454909e-02*(-15+y(1))))))
267 
268  d_inf = (1/(1+exp((1.33333333333333e-01*(-8-y(1))))))
269  alpha_d = (2.5e-01+(1.4e+00/(1+exp((7.69230769230769e-02*(-35-y(1)))))))
270  beta_d = (1.4e+00/(1+exp((2.0e-01*(5+y(1))))))
271  gamma_d = (1/(1+exp((5.0e-02*(50-y(1))))))
272  tau_d = ((1*alpha_d*beta_d)+gamma_d)
273  f2_inf = (3.3e-01+(6.7e-01/(1+exp((1.42857142857143e-01*(35+y(1)))))))
274  tau_f2 = ((562*exp((4.16666666666667e-03*(-pow((27+y(1)),2)))))+(31/(1+exp((1.0e-01*(25-y(1))))))+&
275  & (80/(1+exp((1.0e-01*(30+y(1)))))))
276  dydt(2) = ((f2_inf-y(2))/tau_f2)
277  fcass_inf = (4.0e-01+(6.0e-01/(1+pow((20*y(7)),2))))
278  tau_fcass = (2+(80/(1+pow((20*y(7)),2))))
279  dydt(3) = ((fcass_inf-y(3))/tau_fcass)
280  f_inf = (1/(1+exp((1.42857142857143e-01*(20+y(1))))))
281  tau_f = (20+(1.1025e+03*exp((4.44444444444444e-03*(-pow((27+y(1)),2)))))+(200/(1+exp((1.0e-01*(13-y(1))))))+&
282  &(180/(1+exp((1.0e-01*(30+y(1)))))))
283  dydt(4) = ((f_inf-y(4))/tau_f)
284  e_ca = (1.33568803298478e+01*log((2/y(6))))
285  i_b_ca = (5.92e-04*(y(1)-e_ca))
286  kcasr = (2.5e+00-(1.5e+00/(1+pow((1.5e+00/y(5)),2))))
287  k1 = (1.5e-01/kcasr)
288  o = (k1*pow(y(7),2)*y(8)/(6.0e-02+(k1*pow(y(7),2))))
289  i_rel = (1.02e-01*o*(y(5)-y(7)))
290  i_up = (6.375e-03/(1+(6.25e-08/pow(y(6),2))))
291  i_leak = (3.6e-04*(y(5)-y(6)))
292  i_xfer = (3.8e-03*(y(7)-y(6)))
293  k2 = (4.5e-02*kcasr)
294  dydt(8) = (((-k2)*y(7)*y(8))+(5.0e-03*(1-y(8))))
295  ca_i_bufc = (1/(1+(2.0e-04/pow((1.0e-03+y(6)),2))))
296  ca_sr_bufsr = (1/(1+(3/pow((3.0e-01+y(5)),2))))
297  ca_ss_bufss = (1/(1+(1.0e-04/pow((2.5e-04+y(7)),2))))
298  i_p_ca = (1.238e-01*y(6)/(5.0e-04+y(6)))
299  i_naca = (8.66622022994244e-05*((2*exp((1.31018617879609e-02*y(1)))*pow(y(17),3))-(6860000*&
300  &exp((-2.43320290347846e-02*y(1)))*y(6)))/(1+(1.0e-01*exp((-2.43320290347846e-02*y(1))))))
301  dydt(6) = (ca_i_bufc*((6.66910509631797e-02*(i_leak-i_up))+i_xfer-(5.84427487220097e-05*(i_b_ca+i_p_ca-(2*i_naca)))))
302  dydt(5) = (ca_sr_bufsr*(i_up-i_rel-i_leak))
303  dydt(7) = (ca_ss_bufss*((-1.75328246166029e-02*i_cal)+(2.00073152889539e+01*i_rel)-(300*i_xfer)))
304  e_na = (2.67137606596956e+01*log((140/y(17))))
305  i_na = (1.4838e+01*pow(y(11),3)*y(9)*y(10)*(y(1)-e_na))
306  h_inf = (1/pow((1+exp((1.34589502018843e-01*(7.155e+01+y(1))))),2))
307  if (y(1) < -40.0) then
308  tau_h = (1/((5.7e-02*exp((1.47058823529412e-01*(-80-y(1)))))+(2.7e+00*exp((7.9e-02*y(1))))+(310000*exp((3.485e-01*y(1))))))
309  else
310  tau_h = (1.68831168831169e-01*(1+exp((-9.00900900900901e-02*(1.066e+01+y(1))))))
311  end if
312  dydt(9) = ((h_inf-y(9))/tau_h)
313  j_inf = (1/pow((1+exp((1.34589502018843e-01*(7.155e+01+y(1))))),2))
314  if (y(1) < -40.0) then
315  tau_j = (1/((1*((-25428*exp((2.444e-01*y(1))))-(6.948e-06*exp((-4.391e-02*y(1)))))*(3.778e+01+y(1))/&
316  &(1+exp((3.11e-01*(7.923e+01+y(1))))))+(2.424e-02*exp((-1.052e-02*y(1)))/(1+exp((-1.378e-01*(4.014e+01+y(1))))))))
317  else
318  tau_j = (1.66666666666667e+00/exp((5.7e-02*y(1)))*(1+exp((-1.0e-01*(32+y(1))))))
319  end if
320  dydt(10) = ((j_inf-y(10))/tau_j)
321  m_inf = (1/pow((1+exp((1.10741971207087e-01*(-5.686e+01-y(1))))),2))
322  alpha_m = (1/(1+exp((2.0e-01*(-60-y(1))))))
323  beta_m = ((1.0e-01/(1+exp((2.0e-01*(35+y(1))))))+(1.0e-01/(1+exp((5.0e-03*(-50+y(1)))))))
324  tau_m = (1*alpha_m*beta_m)
325  dydt(11) = ((m_inf-y(11))/tau_m)
326  e_k = (2.67137606596956e+01*log((5.4e+00/y(13))))
327  alpha_k1 = (1.0e-01/(1+exp((6.0e-02*(-200+y(1)-e_k)))))
328  beta_k1 = (((3*exp((2.0e-04*(100+y(1)-e_k))))+exp((1.0e-01*(-10+y(1)-e_k))))/(1+exp((-5.0e-01*(y(1)-e_k)))))
329  xk1_inf = (alpha_k1/(alpha_k1+beta_k1))
330  i_k1 = (5.405e+00*xk1_inf*(y(1)-e_k))
331  i_to = (2.94e-01*y(18)*y(19)*(y(1)-e_k))
332  i_kr = (1.53e-01*y(14)*y(15)*(y(1)-e_k))
333  e_ks = (2.67137606596956e+01*log((9.6e+00/(y(13)+(3.0e-02*y(17))))))
334  i_ks = (3.92e-01*pow(y(16),2)*(y(1)-e_ks))
335  i_nak = (2.298375e+00*y(17)/(40+y(17))/(1+(1.245e-01*exp((-3.74338908227455e-03*y(1))))+&
336  &(3.53e-02*exp((3.74338908227455e-02*(-y(1)))))))
337  i_b_na = (2.9e-04*(y(1)-e_na))
338  i_p_k = (1.46e-02*(y(1)-e_k)/(1+exp((1.67224080267559e-01*(25-y(1))))))
339  dydt(12) = ((d_inf-y(12))/tau_d)
340  dydt(1) = (-1*(i_k1+i_to+i_kr+i_ks+i_cal+i_nak+i_na+i_b_na+i_naca+i_b_ca+i_p_k+i_p_ca+i_stim))
341  dydt(13) = (-1.16885497444019e-04*(i_k1+i_to+i_kr+i_ks+i_p_k+i_stim-(2*i_nak)))
342  xr1_inf = (1/(1+exp((1.42857142857143e-01*(-26-y(1))))))
343  alpha_xr1 = (450/(1+exp((1.0e-01*(-45-y(1))))))
344  beta_xr1 = (6/(1+exp((8.69565217391304e-02*(30+y(1))))))
345  tau_xr1 = (1*alpha_xr1*beta_xr1)
346  dydt(14) = ((xr1_inf-y(14))/tau_xr1)
347  xr2_inf = (1/(1+exp((4.16666666666667e-02*(88+y(1))))))
348  alpha_xr2 = (3/(1+exp((5.0e-02*(-60-y(1))))))
349  beta_xr2 = (1.12e+00/(1+exp((5.0e-02*(-60+y(1))))))
350  tau_xr2 = (1*alpha_xr2*beta_xr2)
351  dydt(15) = ((xr2_inf-y(15))/tau_xr2)
352  xs_inf = (1/(1+exp((7.14285714285714e-02*(-5-y(1))))))
353  alpha_xs = (1400/sqrt((1+exp((1.66666666666667e-01*(5-y(1)))))))
354  beta_xs = (1/(1+exp((6.66666666666667e-02*(-35+y(1))))))
355  tau_xs = (80+(1*alpha_xs*beta_xs))
356  dydt(16) = ((xs_inf-y(16))/tau_xs)
357  dydt(17) = (-1.16885497444019e-04*(i_na+i_b_na+(3*i_nak)+(3*i_naca)))
358  r_inf = (1/(1+exp((1.66666666666667e-01*(20-y(1))))))
359  tau_r = (8.0e-01+(9.5e+00*exp((5.55555555555556e-04*(-pow((40+y(1)),2))))))
360  dydt(18) = ((r_inf-y(18))/tau_r)
361  s_inf = (1/(1+exp((2.0e-01*(20+y(1))))))
362  tau_s = (3+(85*exp((3.125e-03*(-pow((45+y(1)),2)))))+(5/(1+exp((2.0e-01*(-20+y(1)))))))
363  dydt(19) = ((s_inf-y(19))/tau_s)
364  end subroutine tentusscher06_dydt
365 
366  ! integrates all cell models
367  subroutine tentusscher06_integrate(cells,materials,t0,t1,err,error,*)
368  type(field_type), intent(inout), pointer :: cells, materials
369  real(dp), intent(in) :: t0, t1
370  INTEGER(INTG), INTENT(OUT) :: ERR
371  TYPE(varying_string), INTENT(OUT) :: ERROR
372 
373  integer, parameter :: celldim = 19
374  real(dp), dimension(1:celldim) :: dydt
375  real(dp), dimension(:), pointer :: y
376  real(dp) :: t, dt, activ, m_inf, d_inf, m_inf0, d_inf0, prev_v
377  real(dp), dimension(:), pointer :: celldata, activdata
378 
379  integer(intg) :: ncells, i, d
380  type(domain_ptr_type), pointer :: domain
381  TYPE(field_variable_type), POINTER :: CELLS_VARIABLE, ACTIV_VARIABLE
382 
383  enters('tentusscher06_integrate',err,error,*999)
384 
385  domain=>cells%decomposition%domain(1)
386  ncells = domain%ptr%topology%nodes%number_of_nodes ! local
387 
388  cells_variable=>cells%VARIABLE_TYPE_MAP(field_u_variable_type)%PTR
389  activ_variable=>materials%VARIABLE_TYPE_MAP(field_u_variable_type)%PTR
390 
391  CALL field_parameter_set_data_get(cells,field_u_variable_type,field_values_set_type,celldata,err,error,*999)
392  CALL field_parameter_set_data_get(materials,field_u_variable_type,field_values_set_type,activdata,err,error,*999)
393 
394  do i=1,ncells
395  !Default to version 1 of each node derivative
396  d = cells_variable%COMPONENTS(1)%PARAM_TO_DOF_MAP%NODE_PARAM2DOF_MAP%NODES(i)%DERIVATIVES(1)%VERSIONS(1)
397  y => celldata(d:d+celldim-1)
398  prev_v = y(1)
399  !Default to version 1 of each node derivative
400  activ = activdata(activ_variable%COMPONENTS(1)%PARAM_TO_DOF_MAP%NODE_PARAM2DOF_MAP%NODES(i)%DERIVATIVES(1)%VERSIONS(1))
401 
402 
403  t = t0
404  do while (t < t1 - 1e-6)
405  ! integrate one cell, one time step.
406  ! just use adaptive forward euler, tested in matlab
407  call tentusscher06_dydt(t,y,dydt,activ)
408  dt = min(t1-t,1.0/3) ! default max
409  dt = min(dt,1 / abs(dydt(1))) ! maximum increase in Vm
410 
411  m_inf0 = (1/pow((1+exp((1.10741971207087e-01*(-5.686e+01-y(1))))),2))
412  d_inf0 = (1/(1+exp((1.33333333333333e-01*(-8-y(1))))))
413 
414  y = y + dt * dydt ! FWE
415 
416  m_inf = (1/pow((1+exp((1.10741971207087e-01*(-5.686e+01-y(1))))),2))
417  d_inf = (1/(1+exp((1.33333333333333e-01*(-8-y(1))))))
418 
419  ! FWEoo on m,d gates : prevent gates from crossing midpoint steady state
420 
421  m_inf = (m_inf0 + m_inf)/2
422  d_inf = (d_inf0 + d_inf)/2
423  if(dydt(11) * (y(11) - m_inf) > 0) then
424  y(11) = m_inf
425  end if
426  if(dydt(12) * (y(12) - d_inf) > 0) then
427  y(12) = d_inf
428  end if
429 
430  t = t + dt
431  end do
432  if(prev_v < 0 .and. y(1) > 0) then ! store activation times, where else?
433  !Default to version 1 of each node derivative
434  activdata(activ_variable%COMPONENTS(7)%PARAM_TO_DOF_MAP%NODE_PARAM2DOF_MAP%NODES(i)%DERIVATIVES(1)%VERSIONS(1))=t1
435  end if
436  end do ! for cells
437  CALL field_parameter_set_data_restore(cells,field_u_variable_type,field_values_set_type,celldata,err,error,*999)
438  CALL field_parameter_set_data_restore(materials,field_u_variable_type,field_values_set_type,activdata,err,error,*999)
439 
440  exits('tentusscher06_integrate')
441  return
442 999 errorsexits('tentusscher06_integrate',err,error)
443  return 1
444  end subroutine tentusscher06_integrate
446 
447 
subroutine, public enters(NAME, ERR, ERROR,)
Records the entry into the named procedure and initialises the error code.
subroutine, public bueno_orovio_integrate(cells, materials, t0, t1, err, error,)
A buffer type to allow for an array of pointers to a DOMAIN_TYPE.
Definition: types.f90:950
This module contains all string manipulation and transformation routines.
Definition: strings.f90:45
Contains information for a field defined on a region.
Definition: types.f90:1346
subroutine, public tentusscher06_integrate(cells, materials, t0, t1, err, error,)
subroutine, public exits(NAME)
Records the exit out of the named procedure.
This module contains all type definitions in order to avoid cyclic module references.
Definition: types.f90:70
This module contains all the low-level base routines e.g., all debug, control, and low-level communic...
Contains information for a field variable defined on a field.
Definition: types.f90:1289
subroutine tentusscher06_dydt(t, y, dydt, activ)
subroutine, public bueno_orovio_initialize(field, err, error,)
This module contains all kind definitions.
Definition: kinds.f90:45
subroutine, public tentusscher06_initialize(field, err, error,)