Mathematical models of electric activity in cardiac tissue are becoming increasingly powerful tools in the study of cardiac arrhythmias. Considered here are mathematical models based on ordinary differential equations (ODEs) and partial differential equations (PDEs) that describe the behaviour of this electrical activity. Generating an efficient numerical solution of these models is a challenging task, and in fact the physiological accuracy of tissue-scale models is often limited by the efficiency of the numerical solution process. In this thesis, we discuss two sets of experiments that test ideas for making the numerical solution process more efficient. In the first set of experiments, we examine the numerical solution of four single cell cardiac electrophysiological models, which consist solely of ODEs. We study the efficiency of using implicit-explicit Runge-Kutta (IMEX-RK) splitting methods to solve these models. We find that variable step-size implementations of IMEX-RK methods (ARK3 and ARK5) that take advantage of Jacobian structure clearly outperform most methods commonly used in practice for two of the models, and they outperform all methods commonly used in practice for the remaining models. In the second set of experiments, we examine the solution of the bidomain model, a model consisting of both ODEs and PDEs that are typically solved separately. We focus these experiments on numerical methods for the solution of the two PDEs in the bidomain model. The most popular method for this task, the Crank-Nicolson method, produces unphysical oscillations; we propose a method based on a second-order L-stable singly diagonally implicit Runge-Kutta (SDIRK) method to eliminate these oscillations.<p> We find that although the SDIRK method is able to eliminate these unphysical oscillations, it is only more efficient for crude error tolerances.

Mathematical models of electric activity in cardiac tissue are becoming increasingly powerful tools in the study of cardiac arrhythmias. Considered here are mathematical models based on ordinary differential equations (ODEs) and partial differential equations (PDEs) that describe the behaviour of this electrical activity. Generating an efficient numerical solution of these models is a challenging task, and in fact the physiological accuracy of tissue-scale models is often limited by the efficiency of the numerical solution process. In this thesis, we discuss two sets of experiments that test ideas for making the numerical solution process more efficient. In the first set of experiments, we examine the numerical solution of four single cell cardiac electrophysiological models, which consist solely of ODEs. We study the efficiency of using implicit-explicit Runge-Kutta (IMEX-RK) splitting methods to solve these models. We find that variable step-size implementations of IMEX-RK methods (ARK3 and ARK5) that take advantage of Jacobian structure clearly outperform most methods commonly used in practice for two of the models, and they outperform all methods commonly used in practice for the remaining models. In the second set of experiments, we examine the solution of the bidomain model, a model consisting of both ODEs and PDEs that are typically solved separately. We focus these experiments on numerical methods for the solution of the two PDEs in the bidomain model. The most popular method for this task, the Crank-Nicolson method, produces unphysical oscillations; we propose a method based on a second-order L-stable singly diagonally implicit Runge-Kutta (SDIRK) method to eliminate these oscillations.<p> We find that although the SDIRK method is able to eliminate these unphysical oscillations, it is only more efficient for crude error tolerances.

<p>Cardiac arrhythmias triggered by both reentrant and focal sources are closely correlated with regions of tissue characterized by significant structural heterogeneity. Experimental and modeling studies of electrical activity in the heart have shown that local microscopic heterogeneities which average out at the macroscale in healthy tissue play a much more important role in diseased and aging cardiac tissue which have low levels of coupling and abnormal or reduced membrane excitability. However, it is still largely unknown how various combinations of microheterogeneity in the intracellular and interstitial spaces affect wavefront propagation in these critical regimes. </p> <p>This thesis uses biophysically realistic 1-D and 2-D computer models to investigate how heterogeneity in the interstitial and intracellular spaces influence both the initiation of ectopic beats and the escape of multiple ectopic beats from a poorly coupled region of tissue into surrounding well-coupled tissue. An approximate discrete monodomain model that incorporates local heterogeneity in both the interstitial and intracellular spaces was developed to represent the tissue domain. </p> <p>The results showed that increasing the effective interstitial resistivity in poorly coupled fibers alters the distribution of electrical load at the microscale and causes propagation to become more like that observed in continuous fibers. In poorly coupled domains, this nearly continuous state is modulated by cell length and is characterized by decreased gap junction delay, sustained conduction velocity, increased sodium current, reduced maximum upstroke velocity, and increased safety factor. In inhomogeneous fibers with adjacent well-coupled and poorly coupled regions, locally increasing the effective interstitial resistivity in the poorly coupled region reduces the size of the focal source needed to generate an ectopic beat, reduces dispersion of repolarization, and delays the onset of conduction block that is caused by source-load mismatch at the boundary between well-coupled and poorly-coupled regions. In 2-D tissue models, local increases in effective interstitial resistivity as well as microstructural variations in cell arrangement at the boundary between poorly coupled and well-coupled regions of tissue modulate the distribution of maximum sodium current which facilitates the unidirectional escape of focal beats. Variations in the distribution of sodium current as a function of cell length and width lead to directional differences in the response to increased effective interstitial resistivity. Propagation in critical regimes such as the ectopic substrate is very sensitive to source-load interactions and local increases in maximum sodium current caused by microheterogeneity in both intracellular and interstitial structure.</p> / Dissertation

A comparative study of fully implicit staggered and monolithic solution algorithms for the coupled bidomain equations of cardiac electrophysiology was discussed in Part I (JCAM 407:114021, 2022). In this follow-up study, we further employ the similar approach and present an extensive comparative analysis between the staggered and monolithic solution techniques for the coupled partial differential equations (PDEs) of cardiac electromechanics both in the monodomain and bidomain setting. For the monolithic solution technique, we adopt the work of Göktepe et al. (CM 45: 227-243, 2010) and Dal et al. (CMAME 253: 353-336, 2013), where the governing PDEs of the excitation–contraction problem are solved simultaneously and all rate equations are treated with an implicit time integration scheme. For the staggered scheme, however, we simply suggest a decoupled solution of the governing PDEs while keeping other aspects of the monolithic algorithm the same, such as utilization of the implicit time integration scheme and coupling of ordinary differential equations, which describe the state variables of cardiac electrophysiology, to the parabolic PDE of the electrophysiology. Both solution algorithms are applied to several problems, where we simulate regular planar wave propagations and scroll waves for different spatial and temporal resolution and material parameters. The comparison between the solution schemes is performed in terms of accuracy, efficiency and stability. We reveal that for regular wave propagation the results of different solution algorithms are identical. The only deviation is observed during scroll wave propagation, if the parameters controlling the mutual interaction between the mechanical and electrical field are relatively large and at the same time if the time increment is big, e.g., 𝛥𝑡 = 2.0 ms. Nevertheless, the suggested staggered solution scheme yields almost identical results even for high coupling if smaller time increments are used. Besides, the staggered scheme provides an enormous gain in computational efficiency, particularly, as the number of degrees of freedom is increased and we do not encounter any stability issue arising specifically from the decoupled solution of the governing PDEs. Similar to the conclusion that was made for cardiac electrophysiology in Part I, we report that the suggested fully implicit staggered solution algorithm has a tremendous application potential for computer simulations of excitation–contraction coupling of the heart both in the monodomain and bidomain setting.

The work developed in this thesis focusses on the electrical activity of the heart, from the modeling of the action potential originating from cardiac cells and propagating through the heart, as well as its electrical manifestation at the body surface. The study is divided in two main parts: modeling the action potential, and numerical simulations. For modeling the action potential a dimensional and asymptotic analysis is done. The key advance in this part of the work is that this analysis gives the steps to reliably control the action potential. It allows predicting the time/space scales and speed of any action potential that is to say the shape of the action potential and its propagation. This can be done as the explicit relations on all the physiological constants are defined precisely. This method facilitates the integrative modeling of a complete human heart with tissue-specific ionic models. It even proves that using a single model for the cardiac action potential is enough in many situations. For efficient numerical simulations, a numerical method for solving the heart-torso coupling problem is explored according to a level set description of the domains. This is done in the perspective of using directly medical images for building computational domains. A finite element method is then developed to manage meshes not adapted to internal interfaces. Finally, an anisotropic adaptive remeshing methods for unstructured finite element meshes is used to efficiently capture propagating action potentials within complex, realistic two dimensional geometries.

The work developed in this thesis focusses on the electrical activity of the heart, from the modeling of the action potential originating from cardiac cells and propagating through the heart, as well as its electrical manifestation at the body surface. The study is divided in two main parts: modeling the action potential, and numerical simulations. For modeling the action potential a dimensional and asymptotic analysis is done. The key advance in this part of the work is that this analysis gives the steps to reliably control the action potential. It allows predicting the time/space scales and speed of any action potential that is to say the shape of the action potential and its propagation. This can be done as the explicit relations on all the physiological constants are defined precisely. This method facilitates the integrative modeling of a complete human heart with tissue-specific ionic models. It even proves that using a single model for the cardiac action potential is enough in many situations. For efficient numerical simulations, a numerical method for solving the heart-torso coupling problem is explored according to a level set description of the domains. This is done in the perspective of using directly medical images for building computational domains. A finite element method is then developed to manage meshes not adapted to internal interfaces. Finally, an anisotropic adaptive remeshing methods for unstructured finite element meshes is used to efficiently capture propagating action potentials within complex, realistic two dimensional geometries.

The work developed in this thesis focusses on the electrical activity of the heart, from the modeling of the action potential originating from cardiac cells and propagating through the heart, as well as its electrical manifestation at the body surface. The study is divided in two main parts: modeling the action potential, and numerical simulations. For modeling the action potential a dimensional and asymptotic analysis is done. The key advance in this part of the work is that this analysis gives the steps to reliably control the action potential. It allows predicting the time/space scales and speed of any action potential that is to say the shape of the action potential and its propagation. This can be done as the explicit relations on all the physiological constants are defined precisely. This method facilitates the integrative modeling of a complete human heart with tissue-specific ionic models. It even proves that using a single model for the cardiac action potential is enough in many situations. For efficient numerical simulations, a numerical method for solving the heart-torso coupling problem is explored according to a level set description of the domains. This is done in the perspective of using directly medical images for building computational domains. A finite element method is then developed to manage meshes not adapted to internal interfaces. Finally, an anisotropic adaptive remeshing methods for unstructured finite element meshes is used to efficiently capture propagating action potentials within complex, realistic two dimensional geometries.

The work developed in this thesis focusses on the electrical activity of the heart, from the modeling of the action potential originating from cardiac cells and propagating through the heart, as well as its electrical manifestation at the body surface. The study is divided in two main parts: modeling the action potential, and numerical simulations. For modeling the action potential a dimensional and asymptotic analysis is done. The key advance in this part of the work is that this analysis gives the steps to reliably control the action potential. It allows predicting the time/space scales and speed of any action potential that is to say the shape of the action potential and its propagation. This can be done as the explicit relations on all the physiological constants are defined precisely. This method facilitates the integrative modeling of a complete human heart with tissue-specific ionic models. It even proves that using a single model for the cardiac action potential is enough in many situations. For efficient numerical simulations, a numerical method for solving the heart-torso coupling problem is explored according to a level set description of the domains. This is done in the perspective of using directly medical images for building computational domains. A finite element method is then developed to manage meshes not adapted to internal interfaces. Finally, an anisotropic adaptive remeshing methods for unstructured finite element meshes is used to efficiently capture propagating action potentials within complex, realistic two dimensional geometries.

The heart is one the most important organs in the human body and many other live creatures. The electrical activity in the heart controls the heart function, and many heart diseases are linked to the abnormalities in the electrical activity in the heart. Mathematical equations and computer simulation can be used to model the electrical activity in the heart. The heart models are challenging to solve because of the complexity of the models and the huge size of the problems. Several cell models have been proposed to model the electrical activity in a single heart cell. These models must be coupled with a heart model to model the electrical activity in the entire heart. The bidomain model is a popular model to simulate the propagation of electricity in myocardial tissue. It is a continuum-based model consisting of non-linear ordinary differential equations (ODEs) describing the electrical activity at the cellular scale and a system of partial differential equations (PDEs) describing propagation of electricity at the tissue scale. Because of this multi-scale, ODE/PDE structure of the model, splitting methods that treat the ODEs and PDEs in separate steps are natural candidates as numerical methods. First, we need to solve the problem at the cellular scale using ODE solvers. One of the most popular methods to solve the ODEs is known as the Rush-Larsen (RL) method. Its popularity stems from its improved stability over integrators such as the forward Euler (FE) method along with its easy implementation. The RL method partitions the ODEs into two sets: one for the gating variables, which are treated by an exponential integrator, and another for the remaining equations, which are treated by the FE method. The success of the RL method can be understood in terms of its relatively good stability when treating the gating variables. However, this feature would not be expected to be of benefit on cell models for which the stiffness is not captured by the gating equations. We demonstrate that this is indeed the case on a number of stiff cell models. We further propose a new partitioned method based on the combination of a first-order generalization of the RL method with the FE method. This new method leads to simulations of stiff cell models that are often one or two orders of magnitude faster than the original RL method. After solving the ODEs, we need to use bidomain solvers to solve the bidomain model. Two well-known, first-order time-integration methods for solving the bidomain model are the semi-implicit method and the Godunov operator-splitting method. Both methods decouple the numerical procedure at the cellular scale from that at the tissue scale but in slightly different ways. The methods are analyzed in terms of their accuracy, and their relative performance is compared on one-, two-, and three-dimensional test cases. As suggested by the analysis, the test cases show that the Godunov method is significantly faster than the semi-implicit method for the same level of accuracy, specifically, between 5 and 15 times in the cases presented. Second-order bidomain solvers can generally be expected to be more effective than first-order bidomain solvers under normal accuracy requirements. However, the simplest and the most commonly applied second-order method for the PDE step, the Crank-Nicolson (CN) method, may generate unphysical oscillations. We investigate the performance of a two-stage, L-stable singly diagonally implicit Runge-Kutta method for solving the PDEs of the bidomain model and present a stability analysis. Numerical experiments show that the enhanced stability property of this method leads to more physically realistic numerical simulations compared to both the CN and Backward Euler (BE) methods.

La modélisation du vivant, en particulier la modélisation de l'activité cardiaque, est devenue un défi scientifique majeur. Le but de cette thématique est de mieux comprendre les phénomènes physiologiques et donc d'apporter des solutions à des problèmes cliniques. Nous nous intéressons dans cette thèse à la modélisation et à l'étude numérique de l'activité électrique du cœur, en particulier l'étude des électrocardiogrammes (ECGs). L'onde électrique dans le cœur est gouvernée par un système d'équations de réaction-diffusion appelé modèle bidomaine ce système est couplé à une EDO représentant l'activité cellulaire. Afin simuler des ECGs, nous tenons en compte la propagation de l'onde électrique dans le thorax qui est décrite par une équation de diffusion. Nous commençons par une démonstrer l'existence d'une solution faible du système couplé cœur-thorax pour une classe de modèles ioniques phénoménologiques. Nous prouvons ensuite l'unicité de cette solution sous certaines conditions. Le plus grand apport de cette thèse est l'étude et la simulation numérique du couplage électrique cœur-thorax. Les résultats de simulations sont représentés à l'aide des ECGs. Dans une première partie, nous produisons des simulations pour un cas normal et pour des cas pathologiques (blocs de branche gauche et droit et des arhythmies). Nous étudions également l'impact de certaines hypothèses de modélisation sur les ECGs (couplage faible, utilisation du modèle monodomaine, isotropie, homogénéité cellulaire, comportement résistance-condensateur du péricarde,. . . ). Nous étudions à la fin de cette partie la sensibilité des ECGs par apport aux paramètres du modèle. En deuxième partie, nous effectuons l'analyse numérique de schémas du premier ordre en temps découplant les calculs du potentiel d'action et du potentiel extérieur. Puis, nous combinons ces schémas en temps avec un traîtement explicite du type Robin-Robin des conditions de couplage entre le cœur et le thorax. Nous proposons une analyse de stabilité de ces schémas et nous illustrons les résultats avec des simulations numériques d'ECGs. La dernière partie est consacrée à trois applications. Nous commençons par l'estimation de certains paramètres du modèle (conductivité du thorax et paramètres ioniques). Dans la deuxième application, qui est d'originie industrielle, nous utilisons des méthodes d'apprentissage statistique pour reconstruire des ECGs à partir de mesures ('électrogrammes). Enfin, nous présentons des simulations électro-mécaniques du coeur sur une géométrie réelle dans diverses situations physiologiques et pathologiques. Les indicateurs cliniques, électriques et mécaniques, calculés à partir de ces simulations sont très similaires à ceux observés en réalité.