The bidomain model and its monodomain simplification are a valuable tool to compute excitation conduction in cardiac tissue. The reaction-diffusion system is commonly solved using the finite differences or finite element (FE) methods. It is known that mesh resolution has a strong effect on the conduction velocity (CV) resulting from these simulations. Additionally, and particularly at coarse resolutions, CV is impacted by the choice of numerical schemes, such as mass lumping. In this work, we show that the (ir)regularity of node locations has a comparably strong influence.

Initial simulations were performed on 3-d tetrahedral grids conforming to the size and resolutions of the benchmark by Niederer et al. One grid was created with Gmsh, another with regularly arrayed nodes of distances $h\in\{0.5, 0.2, 0.1\}$ mm in $x$-, $y$-, and $z$-directions. Further grids were created by shifting nodes from this grid within an interval of $\pm \frac{q h}{2}$ in those directions, defined by a perturbation factor $q\in[0,1)$. The electrophysiological model by ten Tusscher & Panfilov was used for ionic currents with our PETSc-based simulation framework to solve the system posed by the monodomain model. All parameters were chosen as in the benchmark. Activation times (AT) were recorded in the node diagonally opposite of the stimulation corner.

Further simulations were performed on 2-d rectangular grids measuring 6 cm x 3.5 cm. A block of 3.5 cm x 2 cm to the left was meshed with linear triangular elements of resolution $h=2^{-3}~$mm. The remaining 4 cm of width were filled with linear triangles of sizes $h=1~\text{mm},~2^{-1}~\text{mm},~\ldots,~2^{-6}~\text{mm}~$ from top to bottom, with $h$ constant within successive 5 mm intervals. One of the grids was created using Gmsh, one grid was created by placing nodes in regular intervals according to $h$, and variations of that grid were created by randomly shifting all nodes like above. Fiber orientation was defined along the long axis and stimulation occurred from the left side.

Simulations on the 3-d benchmark geometry showed that the random shifting of nodes has a considerable impact at coarse resolutions. For $h=0.5$ mm, recorded AT could be reduced from 54 ms ($q=0$) to 41 ms with $q=0.5$ (reference solution from the benchmark: 42.82 ms). Higher values for $q$ did, however, not necessarily improve the CV (by bringing it closer to the reference solution) but instead generally increased CV.

The stabilization scheme proposed by Pezzuto et al. seems to be a good solution to the problem of CV depending on grid regularity. Simulations on grids with the same local resolution $h$ but a different node regularity do not require different values for the parameters $\theta_{lhs}$ and $\delta$ to have equally stable CVs. For instance, with $\theta_{lhs}=\theta_{rhs}=0.5$, $\delta=105$, and $p=2$, we observed that AT in the 2-d test geometries were within $\pm 2$ ms for all corresponding nodes of the regular ($q=0$) and the perturbed grid ($q=0.5$). CVs were almost equal for all regions with $h\le2^{-2}$ mm. With $ \delta=0$, the difference was 8 ms, and 5 ms with $\delta=\theta_{lhs}=0$. In those cases, additionally, the CVs differed significantly between areas of different local $h$.

Node regularity can for instance be recognized and assessed using (multi-dimensional) Fourier-transform. Directions with regular node sequences will show up as peaks in the spectrum of node coordinates. Regularity can be somewhat assessed as the ratio of spatial frequencies that are above average in value to those below average.

It remains to be studied, how node regularity affects physiological simulations, for instance, if rotor properties or ECG patterns can be changed significantly by simply perturbing nodes.

This is joint work with Dr. Gunnar Seemann. This work was in parts supported by a Visiting Scholar grant from the Simula School of Research and Innovation and by the European Research Council Advanced Grant CardioNect (ERC-AdG-323099).