CEPS
APPbilayerAtria

Bisurfacic monodomain model for atria

This model is suited for atrial geometries, where the tissue is considered as a single 2D surface. However, to take into account the fact that fibers are not oriented in the same way on the endo and epicardial sides, two layers located at the same place are considered. Those layers are coupled on several regions of the mesh to represent the conduction paths between the layers.

There are then two unknowns $v_1$ and $v_2$ that are the transmembrane potential on the endo- and epicardial layers respectively. The usual parameters can be associated to each layer indepedently: fiber orientations (of course!), conductivities, ionic models, volume fraction, etc.

Equations

Let us denote each layer by $\Omega_k,\ k=1,2$, with their associated parameters and ionic state variables $\sigma_k$ and $w_k$,

\[ \left\{\begin{array}{ll} A_\mathrm{m}(C_\mathrm{m} \partial_t v_1 + I_\mathrm{ion}(v_1,w_1) + I_{\mathrm{app},1}) - \dfrac{1}{1+\lambda}\nabla\cdot(\sigma_{i,1}\nabla v_1) + \alpha\chi(x)(v_1-v_2) = 0 & \text{on }\Omega_1,\\ A_\mathrm{m}(C_\mathrm{m} \partial_t v_2 + I_\mathrm{ion}(v_2,w_2) + I_{\mathrm{app},2}) - \dfrac{1}{1+\lambda}\nabla\cdot(\sigma_{i,2}\nabla v_2) + \alpha\chi(x)(v_2-v_1) = 0 & \text{on }\Omega_2,\\ \partial_t w_k + g(v_k,w_k) = 0 &\text{on }\Omega_k,\\ \sigma_{i,k}\partial_{\bm{n}} v_k = 0&\text{on }\partial\Omega_k,\\ + \text{ initial conditions.} \end{array}\right. \]

where $\alpha$ is a coupling strength coefficient (strictly positive), and $\chi(x)$ is the indicator function of the coupling region(s). The parameter $\alpha$ is given in the input file, and the region is defined through a mesh attribute.

How to use

The parameters are almost the same as for the monodomain equation. Here are the additional or modified entries required for the Bilayer model:

intracellular conductivity layer 1 : <options>
extracellular conductivity layer 1 : <options>
intracellular conductivity layer 2 : <options>
extracellular conductivity layer 2 : <options>

The conductivities defined on each layer. See the iftissue inputs page for details of options.

coupling strength : <real>
coupling attribute: <integer>

The mesh attribute locating the coupling region, and the associated strength of the coupling.

Here is an example of input file

# =========================
# Bilayer monodomain model:
# component index 0 represents the transmembrane potential on first layer
# component index 1 represents the transmembrane potential on second layer
# ===============
# Output controls
# ===============
# Output name
output file name : ./results/bilayerAB2
# Output each 0.5ms
snapshot time : 0.5
# Uncomment to plot currents
#save ionic current : true
#save stimulation current : true
# Data used to generate the activation time map. Values are
# - AP detection threshold
# - Minimum voltage used for peak detection
# - AP duration percentage
activation time data : 0.2 0.2 0.75 0.2 0.2 0.75
# Probe points: values of the unknowns at those points*
# along time will be written in a separate file
# (*at closest point in given mesh to be exact)
probe points : 2.68655 -15.8401 20.9223 , 1.34332 -16.1237 19.8474
# ===============
# Simulation time
# ===============
# Times are given in milliseconds
PDE start time : 0.0
PDE end time : 25
# ========
# Geometry
# ========
# Meshes
2d mesh : ./meshes/atriaTwoLayers.msh
# The given mesh has units of mm, so we scale by 0.1
# to obtain centimetres which is CEPS unit of length
geometry scale : 0.1
# Partitioning method can be
# - ParmetisElem
# - ParmetisElemGeom
# - PTScotchNode
partitioning method : PTScotchNode
# ===============
# Layers coupling
# ===============
# All regions except those with attribute 5 will be coupled: strength*(u0-u1) + monodomainEq(u0,u1) = 0
coupling attribute : NOT 5
coupling strength : 200.
# ======
# Tissue
# ======
# Surface to volume ratio, in cm-1
Am : CONSTANT 1.0
# Surfacic membrane capacitance, in muF per cm2
Cm : CONSTANT 1.0
# Conductivity, in mS per cm
# Format: regionAttribute g_l g_t g_n conductivityType
# * -1 for the region attribute means everywhere
# * g_l g_t g_n are the terms of the diagonalized conductivity tensor
# * conductivity types: 0 for first layer, 1 for second layer
conductivity : -1 2.4 0.8 0.8 0
conductivity : -1 2.4 0.8 0.8 1
# ===========
# Ionic model
# ===========
# Format: <ionicModelName> <regionAttribute> <component> <odeSolver>
# Here we use the Mitchell Schaeffer model everywhere (-1)
# for transmembrane potential (component 0) and solve the ODEs
# using the exponential Adams-Bashforth scheme of order 2.
ionic model : MS -1 0 EAB 2
# ===========
# Stimulation
# ===========
# Format: <attr> <amp> <tstart> <duration> <period> <component>
stimulus : -1 15. 1.1 1 -1 0 stimFunc1
# definition of the spatial variation of stimulation
# format: <label> <C4|C5|CINF|CONSTANT> <center_x> <center_y> <center_z> <diameter> <amplitude>
# here we use a ~5mm wide bell-shaped C-infty profile, centerend around the
# given coordinates
stimulation profile : stimFunc1 CINF -1.48443 -19.2075 18.6753 1.
# ==================
# PDE solver options
# ==================
finite elements : P1
PDE solver : SBDF 2
PDE time step : 0.05
# Linear solver options
linear solver absolute tolerance : 1.E-12
linear solver relative tolerance : 1.E-12
linear solver solver type : GMRES
@ CONSTANT
constant field, aka real_t or tensor [future]
@ GMRES
well, gmres

Outputs

The same outputs as the monodomain model are generated, but for both unknowns.

References

An asymptotic two-layer monodomain model of cardiac electrophysiology in the atria: derivation and convergence Yves Coudière, Jacques Henry, Simon Labarthe SIAM Journal on Applied Mathematics, Society for Industrial and Applied Mathematics, 2017, 77 (2), pp.409 – 429. ⟨10.1137/15M1016886⟩ Read it on HAL