Cardiac tissue at a microscopic size scale is a discrete structure.
From an electrical point of view the most important building block are myocytes,
in which the sarcolemmnal membrane separates the intracellular space inside the myocyte
from the interstitial space in between myocytes.
The intracellular spaces of adjacent myocytes are tightly coupled by gap junctions,
therefore cardiac tissue can be considered a functional syncytium.
Electrical behavior of cardiac tissue when observed from a more macroscopic point of view
can be considered a continuum.
Figure 1: Discrete nature of cardiac tissue and its representation as a continuum in a homogenized model
Tissue: Cardiac tissue consists of discrete myocytes which are interconnected via gap junctions.
Homogenization: Discrete structures in the intracellular and extracellular spaces within the myocardium are homogenized
to arrive at a continuum representation, the electrical properties at the discrete microscopic scale
and at the continuous macroscopic scale match.
Both homogenized domains are overlapped, separated by a membrane at each point in space (interpenetrating domains).
Spatial discretization: Discrete representation of both conductive domains, intra- and extracellular, and the membrane.
In terms of modeling cardiac electrophysiology at the tissue scale continuum representations prevail.
While discrete models representing tissue as interconnected myocytes are, in principle, feasible
and have been reported in the literature,
as the heart is composed of around 5 billion myocytes such models are computationally extremely
expensive,
limiting their application to smaller tissue specimens.
For larger tissue specimens up to the whole heart
the macroscopic cardiac bidomain equations are considered the most complete mathematical model
describing the spread of cardiac electrical activity.
The bidomain model can be derived based on geometric assumptions on cellular structures
through a homogenization procedure,
yielding a continuum formulation where tissue properties are taken into account in an averaged sense.
The domains of interest - intracellular and extracellular - are preserved in the formulation.
The membrane separates intracellular from interstitial space.
Both spaces - intracellular and extracellular - interpenetrate each other
and co-exist everywhere.
Current can flow from one space to the other by crossing the separating membrane.
Charge conservation in both domains yields then the following equations:
Cardiac Bidomain Equation
$$
\begin{align}
\nabla \cdot \sigma_{\rm i} \nabla \phi_{\rm i} &= \beta I_{\rm m} \\
\nabla \cdot \sigma_{\rm e} \nabla \phi_{\rm e} &= -\beta I_{\rm m} + I_{\rm }
\end{align}
$$
where
\(\beta\) is the surface area of membrane contained within a unit volume,
\(I_{\mathrm{e}}\) is an extracellular current density, specified on a per volume basis,
which serves as a stimulus to initiate propagation,
\(\phi_{\mathrm{X}}\) and \( \sigma_{\mathrm{X}} \) are the electrical potential and homogenized
conductivity in the domain \({\mathrm{X}}\) respectively, where \({\mathrm{X}}\) is \({\mathrm{i}}\) or \({\mathrm{e}}\) for intracellular or extracellular.
The conductivity tensors describe the orthotropic eigendirections of the tissue
at a given location, \(\mathbf{x}\), and are constructed as
$$
\sigma_{\mathrm{X}} = \sigma_{\mathrm f, \mathrm X} \, \mathbf{f} \otimes \mathbf{f} +
\sigma_{\mathrm s, \mathrm X} \mathbf{s} \otimes \mathbf{s} +
\sigma_{\mathrm n, \mathrm X} \mathbf{n} \otimes \mathbf{n},
$$
where \(\sigma_{\mathrm f, \mathrm X}\), \( \sigma_{\mathrm s, \mathrm X}\) and \( \sigma_{\mathrm n, \mathrm X} \)
are conductivities of the domain \( {\mathrm{X}} \) along the local eigenaxes, \( \mathbf{f} \),
along the prevailing myocyte orientation referred to as fiber orientation,
perpendicular to the fiber orientation but within a sheet, \(\mathbf{s}\), and the sheet normal direction, \( \mathbf{n} \).
The transmembrane voltage \(V_m\) is given as the difference between intra- and extracellular potential
$$
V_{\rm{m}} = \phi_{\rm{i}} - \phi_{\rm{e}}
$$
and the transmembrane current \(I_{\rm{m}}\) is given by
$$
\begin{align}
I_{\rm{m}} &= C_{\rm{m}} \frac{\partial V_{\rm{m}}}{\partial t} + I_{\rm{ion}}(V_{\rm{m}},\boldsymbol{\eta}) \\
\frac{\partial \boldsymbol{\eta}}{\partial t} &= f(V_{\rm m}, \boldsymbol{\eta})
\end{align}
$$
where \(C_{\rm{m}}\) is the membrane capacitance, and \(\boldsymbol{\eta}\) is a vector of state variables
describing the cellular dynamics governed by gating mechanisms and ion fluxes.
That is, the transmembrane current \(I_{\rm{m}}\) consists of a capacitive current \(I_{\rm{c}}\)
and a current due to ions crossing the membrane, \(I_{\rm{ion}}\).
The fast upstroke of the action potential lasting less than 1 ms translates into steep depolarization wave fronts in space
of a spatial extent of less than 1 mm. This spatio-temporal dynamics render solving the bidomain equations computationally expensive
since high spatio-temporal resolutions are needed to achieve sufficent accuracy.
Typically, spatio-temporal resolutions of 50 μs and 250 μm or finer are used
as shown in standard electrophysiology verification benchmarks.
The Eikonal Equation
In the early nineties, the eikonal model has been proposed as an efficient way
of computing arrival times of depolarization wavefronts in the myocardium,
and its potential limitations have been studied extensively
[Franzone and Rovida, 1990],
[Keener, 1991],
[Franzone et al., 1998].
Wavefront arrival times based on the eikonal model have also been used
to prescribe the spatio-temporal evolution of the transmembrane voltage
to predict extracellular potentials and electrograms
[Franzone et al., 1993],
[Franzone et al., 1998],
[Franzone et al., 2000].
Eikonal Equation
In the eikonal model, wavefront arrival times \(t_{\rm a} \) in the myocardium \(\Omega\) are described based on a spatially
heterogeneous orthotropic velocity function encoded as \(\boldsymbol{V}(\mathrm{x})\)
and initial activations \(t_0\) at locations \( \Gamma \).
The eikonal equation is of the form
$$
\left\{
\begin{array}{rcll}
\sqrt{\nabla t_{\rm a}^\top \, \boldsymbol{V} \, \nabla t_{\rm a}} & = & 1 \qquad & \text{in } \Omega \\
t_{\rm a} & = & t_0 & \text{in } \Gamma
\end{array}
\right. \label{eq:_eikonal}
$$
where \( \Omega \) is the myocardium and \( t_{\rm a}\) is a positive function
describing the wavefront arrival time at location \( \mathrm{x}\).
The symmetric positive definite \( 3 \times 3 \) tensor \( \boldsymbol{V}(\boldsymbol{x}) \) holds the squared velocities
\( \left( v_{\rm f}(\boldsymbol{x}), v_{\rm s}(\boldsymbol{x}), v_{\rm n}(\boldsymbol{x}) \right) \) associated to
\( \left( \boldsymbol{f}(\boldsymbol{x}), \boldsymbol{s}(\boldsymbol{x}), \boldsymbol{n}(\boldsymbol{x}) \right) \),
respectively the longitudinal, transversal and normal fiber directions:
$$
\boldsymbol{V} := v_{f}^2 \boldsymbol{f} \, \boldsymbol{f}^\top +
v_{s}^2 \boldsymbol{s} \, \boldsymbol{s}^\top +
v_{n}^2 \boldsymbol{n} \, \boldsymbol{n}^\top
$$
This definition is analogous to the that of the conductivity tensors \( \sigma_{\mathrm{X}}\) above in bidomain equation,
and since \( \boldsymbol{V}(\boldsymbol{x}) \) and \( \boldsymbol{\sigma_{\rm i}}(\boldsymbol{x}) \)
are based on the same orthotropic fiber and sheet arrangement
\( \left( \boldsymbol{f}, \boldsymbol{s}, \boldsymbol{n} \right) \)
it follows that the velocities \( \left( v_{\rm f}, v_{\rm s}, v_{\rm n} \right) \) and the conductivities
\(\sigma_{\mathrm f, \mathrm X}\), \( \sigma_{\mathrm s, \mathrm X}\) and \( \sigma_{\mathrm n, \mathrm X} \)
forming \( \boldsymbol{\sigma}_{\mathrm{X}} \) must be proportionally coupled.
While conductivities and velocities for a uniformly propagating planar depolarization wave front
are proportionally related by
$$
v_{\rm X} \propto \sqrt{\frac{\sigma_{\rm{m,X}}}{\beta}}
$$
with \(\sigma_{\rm{m,X}} \) being twice the harmonic mean conductivity tensor given as
$$
\sigma_{\rm{m,X}} = \frac{\sigma_{\rm{i,X}} \, \sigma_{\rm{e,X}}}{\sigma_{\rm{i,X}} +
\sigma_{\rm{e,X}}}
$$
fitting of conductivities to match a desired conduction velocity is a non-trivial task
since simulated velocities are also influenced by a number of technical factors
such as spatio-temporal discretization, chosen element type or integration rules.
Methods for automated fitting have been discussed elsewhere
[Mendonca Costa et al. 2013].
Video 1: Eikonal activation sequence
Activation starts at the sino-atrial node and spreads from the right atrium
to the left atrium. Excitation propagates slowly then through the atrioventricular
node into the His-Purkinje system which spreads the
impulse quickly over the ventricular endocardium. After a delay of about 10 ms -
the time needed for thin 1D Purkinje strands to excite the thick-walled 3D myocardium -
the depolarization wavefront is transduced to the ventricular myocardium.
Wavefronts traverse then the ventricular walls mostly in a transmural direction.
