Solving the cardiac bidomain equations for discontinuous conductivities

Graph-based homogenisation for modelling cardiac fibrosis

Fibrosis, the excess of extracellular matrix, can affect, and even block, propagation of action potential in cardiac tissue. This can result in deleterious effects on heart function, but the nature and severity of these effects depend strongly on the localisation of fibrosis and its by-products in cardiac tissue, such as collagen scar formation. Computer simulation is an important means of understanding the complex effects of fibrosis on activation patterns in the heart, but concerns of computational cost place restrictions on the spatial resolution of these simulations. In this work, we present a novel numerical homogenisation technique that uses both Eikonal and graph approaches to allow fine-scale heterogeneities in conductivity to be incorporated into a coarser mesh. Homogenisation achieves this by deriving effective conductivity tensors so that a coarser mesh can then be used for numerical simulation. By taking a graph-based approach, our homogenisation technique functions naturally on irregular grids and does not rely upon any assumptions of periodicity, even implicitly. We present results of action potential propagation through fibrotic tissue in two dimensions that show the graph-based homogenisation technique is an accurate and effective way to capture fine-scale domain information on coarser meshes in the context of sharp-fronted travelling waves of activation. As test problems, we consider excitation propagation in tissue with diffuse fibrosis and through a tunnel-like structure designed to test homogenisation, interaction of an excitation wave with a scar region, and functional re-entry.

Cardiac Conduction Velocity, Remodeling and Arrhythmogenesis

Cardiac electrophysiological disorders, in particular arrhythmias, are a key cause of morbidity and mortality throughout the world. There are two basic requirements for arrhythmogenesis: an underlying substrate and a trigger. Altered conduction velocity (CV) provides a key substrate for arrhythmogenesis, with slowed CV increasing the probability of re-entrant arrhythmias by reducing the length scale over which re-entry can occur. In this review, we examine methods to measure cardiac CV in vivo and ex vivo, discuss underlying determinants of CV, and address how pathological variations alter CV, potentially increasing arrhythmogenic risk. Finally, we will highlight future directions both for methodologies to measure CV and for possible treatments to restore normal CV.

Approaches for determining cardiac bidomain conductivity values: progress and challenges

Modelling the electrical activity of the heart is an important tool for understanding electrical function in various diseases and conduction disorders. Clearly, for model results to be useful, it is necessary to have accurate inputs for the models, in particular the commonly used bidomain model. However, there are only three sets of four experimentally determined conductivity values for cardiac ventricular tissue and these are inconsistent, were measured around 40 years ago, often produce different results in simulations and do not fully represent the three-dimensional anisotropic nature of cardiac tissue. Despite efforts in the intervening years, difficulties associated with making the measurements and also determining the conductivities from the experimental data have not yet been overcome. In this review, we summarise what is known about the conductivity values, as well as progress to date in meeting the challenges associated with both the mathematical modelling and the experimental techniques. Graphical abstract Epicardial potential distributions, arising from a subendocardial ischaemic region, modelled using conductivity data from the indicated studies.

Scalable and Accurate ECG Simulation for Reaction-Diffusion Models of the Human Heart

Realistic electrocardiogram (ECG) simulation with numerical models is important for research linking cellular and molecular physiology to clinically observable signals, and crucial for patient tailoring of numerical heart models. However, ECG simulation with a realistic torso model is computationally much harder than simulation of cardiac activity itself, so that many studies with sophisticated heart models have resorted to crude approximations of the ECG. This paper shows how the classical concept of electrocardiographic lead fields can be used for an ECG simulation method that matches the realism of modern heart models. The accuracy and resource requirements were compared to those of a full-torso solution for the potential and scaling was tested up to 14,336 cores with a heart model consisting of 11 million nodes. Reference ECGs were computed on a 3.3 billion-node heart-torso mesh at 0.2 mm resolution. The results show that the lead-field method is more efficient than a full-torso solution when the number of simulated samples is larger than the number of computed ECG leads. While the initial computation of the lead fields remains a hard and poorly scalable problem, the ECG computation itself scales almost perfectly and, even for several hundreds of ECG leads, takes much less time than the underlying simulation of cardiac activity.

Stability analysis of the POD reduced order method for solving the bidomain model in cardiac electrophysiology

In this paper we show the numerical stability of the Proper Orthogonal Decomposition (POD) reduced order method used in cardiac electrophysiology applications. The difficulty of proving the stability comes from the fact that we are interested in the bidomain model, which is a system of degenerate parabolic equations coupled to a system of ODEs representing the cell membrane electrical activity. The proof of the stability of this method is based on a priori estimates controlling the gap between the reduced order solution and the Galerkin finite element one. We present some numerical simulations confirming the theoretical results. We also combine the POD method with a time splitting scheme allowing a faster solving of the bidomain problem and show numerical results. Finally, we conduct numerical simulation in 2D illustrating the stability of the POD method in its sensitivity to the ionic model parameters. We also perform 3D simulation using a massively parallel code. We show the computational gain using the POD reduced order model. We also show that this method has a better scalability than the full finite element method.

Predicting electromyographic signals under realistic conditions using a multiscale chemo-electro-mechanical finite element model

This paper presents a novel multiscale finite element-based framework for modelling electromyographic (EMG) signals. The framework combines (i) a biophysical description of the excitation-contraction coupling at the half-sarcomere level, (ii) a model of the action potential (AP) propagation along muscle fibres, (iii) a continuum-mechanical formulation of force generation and deformation of the muscle, and (iv) a model for predicting the intramuscular and surface EMG. Owing to the biophysical description of the half-sarcomere, the model inherently accounts for physiological properties of skeletal muscle. To demonstrate this, the influence of membrane fatigue on the EMG signal during sustained contractions is investigated. During a stimulation period of 500 ms at 100 Hz, the predicted EMG amplitude decreases by 40% and the AP propagation velocity decreases by 15%. Further, the model can take into account contraction-induced deformations of the muscle. This is demonstrated by simulating fixed-length contractions of an idealized geometry and a model of the human tibialis anterior muscle (TA). The model of the TA furthermore demonstrates that the proposed finite element model is capable of simulating realistic geometries, complex fibre architectures, and can include different types of heterogeneities. In addition, the TA model accounts for a distributed innervation zone, different fibre types and appeals to motor unit discharge times that are based on a biophysical description of the α motor neurons.

Cardiac response to low-energy field pacing challenges the standard theory of defibrillation

BACKGROUND

The electric response of myocardial tissue to periodic field stimuli has attracted significant attention as the basis for low-energy antifibrillation pacing, potentially more effective than traditional single high-energy shocks. In conventional models, an electric field produces a highly nonuniform response of the myocardial wall, with discrete excitations, or hot spots (HS), occurring at cathodal tissue surfaces or large coronary vessels. We test this prediction using novel 3-dimensional tomographic optical imaging.

METHODS AND RESULTS

Experiments were performed in isolated coronary perfused pig ventricular wall preparations stained with near-infrared voltage-sensitive fluorescent dye DI-4-ANBDQBS. The 3-dimensional coordinates of HS were determined using alternating transillumination. To relate HS formation with myocardial structures, we used ultradeep confocal imaging (interrogation depths, >4 mm). The peak HS distribution is located deep inside the heart wall, and the depth is not significantly affected by field polarity. We did not observe the strong colocalization of HS with major coronary vessels anticipated from theory. Yet, we observed considerable lateral displacement of HS with field polarity reversal. Models that de-emphasized lateral intracellular coupling and accounted for resistive heterogeneity in the extracellular space showed similar HS distributions to the experimental observations.

CONCLUSIONS

The HS distributions within the myocardial wall and the significant lateral displacements with field polarity reversal are inconsistent with standard theories of defibrillation. Extended theories based on enhanced descriptions of cellular scale electric mechanisms may be necessary. The considerable lateral displacement of HS with field polarity reversal supports the hypothesis of biphasic stimuli in low-energy antifibrillation pacing being advantageous.

Numerical quadrature and operator splitting in finite element methods for cardiac electrophysiology

We study the numerical accuracy and computational efficiency of alternative formulations of the finite element solution procedure for the monodomain equations of cardiac electrophysiology, focusing on the interaction of spatial quadrature implementations with operator splitting and examining both nodal and Gauss quadrature methods and implementations that mix nodal storage of state variables with Gauss quadrature. We evaluate the performance of all possible combinations of 'lumped' approximations of consistent capacitance and mass matrices. Most generally, we find that quadrature schemes and lumped approximations that produce decoupled nodal ionic equations allow for the greatest computational efficiency, this being afforded through the use of asynchronous adaptive time-stepping of the ionic state variable ODEs. We identify two lumped approximation schemes that exhibit superior accuracy, rivaling that of the most expensive variationally consistent implementations. Finally, we illustrate some of the physiological consequences of discretization error in electrophysiological simulation relevant to cardiac arrhythmia and fibrillation. These results suggest caution with the use of semi-automated free-form tetrahedral and hexahedral meshing algorithms available in most commercially available meshing software, which produce nonuniform meshes having a large distribution of element sizes.

Modeling defibrillation of the heart: approaches and insights

Cardiac defibrillation, as accomplished nowadays by automatic, implantable devices (ICDs), constitutes the most important means of combating sudden cardiac death. While ICD therapy has proved to be efficient and reliable, defibrillation is a traumatic experience. Thus, research on defibrillation mechanisms, particularly aimed at lowering defibrillation voltage, remains an important topic. Advancing our understanding towards a full appreciation of the mechanisms by which a shock interacts with the heart is the most promising approach to achieve this goal. The aim of this paper is to assess the current state-of-the-art in ventricular defibrillation modeling, focusing on both numerical modeling approaches and major insights that have been obtained using defibrillation models, primarily those of realistic ventricular geometry. The paper showcases the contributions that modeling and simulation have made to our understanding of the defibrillation process. The review thus provides an example of biophysically based computational modeling of the heart (i.e., cardiac defibrillation) that has advanced the understanding of cardiac electrophysiological interaction at the organ level and has the potential to contribute to the betterment of the clinical practice of defibrillation.

Computational modelling of cardiac electrophysiology: explanation of the variability of results from different numerical solvers

A recent verification study compared 11 large-scale cardiac electrophysiology solvers on an unambiguously defined common problem. An unexpected amount of variation was observed between the codes, including significant error in conduction velocity in the majority of the codes at certain spatial resolutions. In particular, the results of the six finite element codes varied considerably despite each using the same order of interpolation. In this present study, we compare various algorithms for cardiac electrophysiological simulation, which allows us to fully explain the differences between the solvers. We identify the use of mass lumping as the fundamental cause of the largest variations, specifically the combination of the commonly used techniques of mass lumping and operator splitting, which results in a slightly different form of mass lumping to that supported by theory and leads to increased numerical error. Other variations are explained through the manner in which the ionic current is interpolated. We also investigate the effect of different forms of mass lumping in various types of simulation.

Simulating a dual-array electrode configuration to investigate the influence of skeletal muscle fatigue following functional electrical stimulation

A novel, anatomically-accurate model of a tibialis anterior muscle is used to investigate the electro-physiological properties of denervated muscles following functional electrical stimulation. The model includes a state-of-the-art description of cell electro-physiology. The main objective of this work is to develop a computational framework capable of predicting the effects of different stimulation trains and electrode configurations on the excitability and fatigue of skeletal muscle tissue. Utilizing a reduced but computationally amenable model, the effects of different electrode sizes and inter-electrode distances on the number of activated muscle fibers are investigated and qualitatively compared to existing literature. To analyze muscle fatigue, the sodium current, specifically the K+ ion concentrations within the t-tubule and the calcium release from the sarcoplasmic reticulum, is used to quantify membrane and metabolic fatigue. The simulations demonstrate that lower stimulation frequencies and biphasic pulse waveforms cause less fatigue than higher stimulation frequencies and monophasic pulses. A comparison between single and dual electrode configurations (with the same overall stimulation surface) is presented to locally investigate the differences in muscle fatigue. The dual electrode configuration causes the muscle tissue to fatigue quicker.

Reduced-order modeling for cardiac electrophysiology. Application to parameter identification

A reduced-order model based on proper orthogonal decomposition (POD) is proposed for the bidomain equations of cardiac electrophysiology. Its accuracy is assessed through electrocardiograms in various configurations, including myocardium infarctions and long-time simulations. We show in particular that a restitution curve can efficiently be approximated by this approach. The reduced-order model is then used in an inverse problem solved by an evolutionary algorithm. Some attempts are presented to identify ionic parameters and infarction locations from synthetic electrocardiograms.

Construction of 3D MR image-based computer models of pathologic hearts, augmented with histology and optical fluorescence imaging to characterize action potential propagation

Cardiac computer models can help us understand and predict the propagation of excitation waves (i.e., action potential, AP) in healthy and pathologic hearts. Our broad aim is to develop accurate 3D MR image-based computer models of electrophysiology in large hearts (translatable to clinical applications) and to validate them experimentally. The specific goals of this paper were to match models with maps of the propagation of optical AP on the epicardial surface using large porcine hearts with scars, estimating several parameters relevant to macroscopic reaction-diffusion electrophysiological models. We used voltage-sensitive dyes to image AP in large porcine hearts with scars (three specimens had chronic myocardial infarct, and three had radiofrequency RF acute scars). We first analyzed the main AP waves' characteristics: duration (APD) and propagation under controlled pacing locations and frequencies as recorded from 2D optical images. We further built 3D MR image-based computer models that have information derived from the optical measures, as well as morphologic MRI data (i.e., myocardial anatomy, fiber directions and scar definition). The scar morphology from MR images was validated against corresponding whole-mount histology. We also compared the measured 3D isochronal maps of depolarization to simulated isochrones (the latter replicating precisely the experimental conditions), performing model customization and 3D volumetric adjustments of the local conductivity. Our results demonstrated that mean APD in the border zone (BZ) of the infarct scars was reduced by ~13% (compared to ~318 ms measured in normal zone, NZ), but APD did not change significantly in the thin BZ of the ablation scars. A generic value for velocity ratio (1:2.7) in healthy myocardial tissue was derived from measured values of transverse and longitudinal conduction velocities relative to fibers direction (22 cm/s and 60 cm/s, respectively). The model customization and 3D volumetric adjustment reduced the differences between measurements and simulations; for example, from one pacing location, the adjustment reduced the absolute error in local depolarization times by a factor of 5 (i.e., from 58 ms to 11 ms) in the infarcted heart, and by a factor of 6 (i.e., from 60 ms to 9 ms) in the heart with the RF scar. Moreover, the sensitivity of adjusted conductivity maps to different pacing locations was tested, and the errors in activation times were found to be of approximately 10-12 ms independent of pacing location used to adjust model parameters, suggesting that any location can be used for model predictions.

Tissue-specific mathematical models of slow wave entrainment in wild-type and 5-HT(2B) knockout mice with altered interstitial cells of Cajal networks

Gastrointestinal slow waves are generated within networks of interstitial cells of Cajal (ICCs). In the intact tissue, slow waves are entrained to neighboring ICCs with higher intrinsic frequencies, leading to active propagation of slow waves. Degradation of ICC networks in humans is associated with motility disorders; however, the pathophysiological mechanisms of this relationship are uncertain. A recently developed biophysically based mathematical model of ICC was adopted and updated to simulate entrainment of slow waves. Simulated slow wave propagation was successfully entrained in a one-dimensional model, which contained a gradient of intrinsic frequencies. Slow wave propagation was then simulated in tissue models which contained a realistic two-dimensional microstructure of the myenteric ICC networks translated from wild-type (WT) and 5-HT(2B) knockout (degraded) mouse jejunum. The results showed that the peak current density in the WT model was 0.49 muA mm(-2) higher than the 5-HT(2B) knockout model, and the intracellular Ca(2+) density after 400 ms was 0.26 mM mm(-2) higher in the WT model. In conclusion, tissue-specific models of slow waves are presented, and simulations quantitatively demonstrated physiological differences between WT and 5-HT(2B) knockout models. This study provides a framework for evaluating how ICC network degradation may impair slow wave propagation and ultimately motility and transit.

Models of cardiac tissue electrophysiology: progress, challenges and open questions

Models of cardiac tissue electrophysiology are an important component of the Cardiac Physiome Project, which is an international effort to build biophysically based multi-scale mathematical models of the heart. Models of tissue electrophysiology can provide a bridge between electrophysiological cell models at smaller scales, and tissue mechanics, metabolism and blood flow at larger scales. This paper is a critical review of cardiac tissue electrophysiology models, focussing on the micro-structure of cardiac tissue, generic behaviours of action potential propagation, different models of cardiac tissue electrophysiology, the choice of parameter values and tissue geometry, emergent properties in tissue models, numerical techniques and computational issues. We propose a tentative list of information that could be included in published descriptions of tissue electrophysiology models, and used to support interpretation and evaluation of simulation results. We conclude with a discussion of challenges and open questions.

A numerical guide to the solution of the bi-domain equations of cardiac electrophysiology

Simulation of cardiac electrical activity using the bi-domain equations can be a massively computationally demanding problem. This study provides a comprehensive guide to numerical bi-domain modelling. Each component of bi-domain simulations--discretization, ODE-solution, linear system solution, and parallelization--is discussed, and previously-used methods are reviewed, new methods are proposed, and issues which cause particular difficulty are highlighted. Particular attention is paid to the choice of stimulus currents, compatibility conditions for the equations, the solution of singular linear systems, and convergence of the numerical scheme.

Structural heterogeneity alone is a sufficient substrate for dynamic instability and altered restitution

BACKGROUND

Marked changes in ventricular APD restitution and associated alternans rhythm have been demonstrated in structural heart disease (SHD). However, whether this is due to structural heterogeneity or regional variation in cellular properties remains uncertain. In this study, we address the hypothesis that the structural heterogeneity associated with SHD is sufficient to alter dynamic restitution and increase the probability of electric instability.

METHODS AND RESULTS

Activation was simulated in a 14x14 mm(2) domain in the presence and absence (control) of a central region containing nonuniform discontinuities resembling patchy fibrosis. A modified LR1 cardiac activation model was used in a bidomain formulation with isotropic conductivities. Bipolar stimulation was imposed above the central region with coupling intervals decreasing progressively from 500 ms and then maintained at 105 ms. Structural discontinuities had little effect on electric activation at low stimulus rates, but activation time and APD distributions became highly nonuniform within and adjacent to the discontinuous region at high rates. Discordant APD alternans occurred in both "fibrosis" and control, but at lower stimulus rates and with markedly greater extent in the former. Tortuous conduction through the discontinuous region resulted in large fluctuations of diastolic intervals giving rise to regional electric instability, which modulates dynamic conduction velocity and APD restitution. This led to heterogeneous conduction block and reentry not observed in control.

CONCLUSIONS

We show that structural discontinuities can amplify discordant alternans and provide a rate-dependent substrate for reentry. This work provides new insights into the mechanisms by which fibrosis may contribute to arrhythmogenesis.

Shock induced electrical activation in structurally detailed models of pig left-ventricular tissue

Detailed models of sample specific structures in pig left-ventricular tissue have been constructed. These models include epicardial and endocardial surfaces, fiber and sheet orientations, vessels and cleavage planes with significant dimensions. This work shows that it is possible to extract from 3D tissue images reduced dimension descriptions of cleavage planes in the heart wall. These descriptions are used to analyze the response of tissue to electrical shocks of varying strengths. The presence of explicit discontinuities in the heart significantly reduces the time required for transmural activation and provides a basis for understanding successful defibrillation.

CHASTE: incorporating a novel multi-scale spatial and temporal algorithm into a large-scale open source library

Recent work has described the software engineering and computational infrastructure that has been set up as part of the Cancer, Heart and Soft Tissue Environment (CHASTE) project. CHASTE is an open source software package that currently has heart and cancer modelling functionality. This software has been written using a programming paradigm imported from the commercial sector and has resulted in a code that has been subject to a far more rigorous testing procedure than that is usual in this field. In this paper, we explain how new functionality may be incorporated into CHASTE. Whiteley has developed a numerical algorithm for solving the bidomain equations that uses the multi-scale (MS) nature of the physiology modelled to enhance computational efficiency. Using a simple geometry in two dimensions and a purpose-built code, this algorithm was reported to give an increase in computational efficiency of more than two orders of magnitude. In this paper, we begin by reviewing numerical methods currently in use for solving the bidomain equations, explaining how these methods may be developed to use the MS algorithm discussed above. We then demonstrate the use of this algorithm within the CHASTE framework for solving the monodomain and bidomain equations in a three-dimensional realistic heart geometry. Finally, we discuss how CHASTE may be developed to include new physiological functionality--such as modelling a beating heart and fluid flow in the heart--and how new algorithms aimed at increasing the efficiency of the code may be incorporated.

Numerical solution of the bidomain equations

Knowledge of cardiac electrophysiology is efficiently formulated in terms of mathematical models. However, most of these models are very complex and thus defeat direct mathematical reasoning founded on classical and analytical considerations. This is particularly so for the celebrated bidomain model that was developed almost 40 years ago for the concurrent analysis of extra- and intracellular electrical activity. Numerical simulations based on this model represent an indispensable tool for studying electrophysiology. However, complex mathematical models, steep gradients in the solutions and complicated geometries lead to extremely challenging computational problems. The greatest achievement in scientific computing over the past 50 years has been to enable the solving of linear systems of algebraic equations that arise from discretizations of partial differential equations in an optimal manner, i.e. such that the central processing unit (CPU) effort increases linearly with the number of computational nodes. Over the past decade, such optimal methods have been introduced in the simulation of electrophysiology. This development, together with the development of affordable parallel computers, has enabled the solution of the bidomain model combined with accurate cellular models, on geometries resembling a human heart. However, in spite of recent progress, the full potential of modern computational methods has yet to be exploited for the solution of the bidomain model. This paper reviews the development of numerical methods for solving the bidomain model. However, the field is huge and we thus restrict our focus to developments that have been made since the year 2000.

Structure specific models of electrical function in the right atrial appendage

Atrial fibrillation is the most common cardiac arrhythmia. It is prevalent in the elderly and contributes to mortality in congestive heart failure. Computer models of atrial electrical activation that incorporate realistic structure provide a means of investigating the mechanisms that initiate and maintain reentrant atrial arrhythmia. We are extending computational and experimental techniques that are well established in our laboratory to develop a detailed structure-based model of atrial electrical function. The 3D geometry of the atria and veno-atrial junctions is reconstructed from magnetic resonance images and detailed structure from specific regions of the atria is acquired using a semi-automated extended-volume imaging system. As an example of our approach, we present a reconstruction of the pig right atrial appendage (RAA), including the pectinate muscles (PM) and crista terminalis (CT). The RAA was embedded in wax and the block surface was serially etched, stained and imaged, then removed using an ultramiller to produce a uniformly-spaced image stack. Tissue was segmented and connected voxels were selected using a 3D region-growing algorithm to construct RAA geometry. Electrical activity has been modeled on this structure using the Courtemanche atrial cell activation model. A bidomain formulation was used employing a grid-based finite element solver. The RAA was activated by applying a stimulus (150 microA/mm3, 5 ms) to the 27 grid points at the top of the CT. Despite the complex structure of the PM, RAA activation was relatively uniform.

A reliability analysis of cardiac repolarization time markers

Only a limited number of studies have addressed the reliability of extracellular markers of cardiac repolarization time, such as the classical marker RT(eg) defined as the time of maximum upslope of the electrogram T wave. This work presents an extensive three-dimensional simulation study of cardiac repolarization time, extending the previous one-dimensional simulation study of a myocardial strand by Steinhaus [B.M. Steinhaus, Estimating cardiac transmembrane activation and recovery times from unipolar and bipolar extracellular electrograms: a simulation study, Circ. Res. 64 (3) (1989) 449]. The simulations are based on the bidomain - Luo-Rudy phase I system with rotational fiber anisotropy and homogeneous or heterogeneous transmural intrinsic membrane properties. The classical extracellular marker RT(eg) is compared with the gold standard of fastest repolarization time RT(tap), defined as the time of minimum derivative during the downstroke of the transmembrane action potential (TAP). Additionally, a new extracellular marker RT90(eg) is compared with the gold standard of late repolarization time RT90(tap), defined as the time when the TAP reaches 90% of its resting value. The results show a good global match between the extracellular and transmembrane repolarization markers, with small relative mean discrepancy (<or=1.6%) and high correlation coefficients (>or=0.92), ensuring a reasonably good global match between the associated repolarization sequences. However, large local discrepancies of the extracellular versus transmembrane markers may ensue in regions where the curvature of the repolarization front changes abruptly (e.g. near front collisions) or is negligible (e.g. where repolarization proceeds almost uniformly across fiber). As a consequence, the spatial distribution of activation-recovery intervals (ARI) may provide an inaccurate estimate of (and weakly correlated with) the spatial distribution of action potential durations (APD).

A tissue-specific model of reentry in the right atrial appendage

INTRODUCTION

Atrial fibrillation is prevalent in the elderly and contributes to mortality in congestive heart failure. Development of computer models of atrial electrical activation that incorporate realistic structures provides a means of investigating the mechanisms that initiate and maintain reentrant atrial arrhythmia. As a step toward this, we have developed a model of the right atrial appendage (RAA) including detailed geometry of the pectinate muscles (PM) and crista terminalis (CT) with high spatial resolution, as well as complete fiber architecture.

METHODS AND RESULTS

Detailed structural images of a pig RAA were acquired using a semiautomated extended-volume imaging system. The generally accepted anisotropic ratio of 10:1 was adopted in the computer model. To deal with the regional action potential duration heterogeneity in the RAA, a Courtemanche cell model and a Luo-Rudy cell model were used for the CT and PM, respectively. Activation through the CT and PM network was adequately reproduced with acceptable accuracy using reduced-order computer models. Using a train of reducing cycle length stimuli applied to a CT/PM junction, we observed functional block both parallel with and perpendicular to the axis of the CT.

CONCLUSION

With stimulation from the CT at the junction of a PM, we conclude: (a) that conduction block within the CT is due to a reduced safety factor; and (b) that unidirectional block and reentry within the CT is due to its high anisotropy. Regional differences in effective refractive period do not explain the observed conduction block.

Experiment-specific models of ventricular electrical activation: construction and application

Experimental intramural recordings of electrical activity at high resolution have been made in the in-vivo pig LV free wall. To analyze features of these recordings experiment-specific 3D computer models of tissue structures and electrical behavior around the recording sites were constructed. The construction of the models used novel tissue image registration, correction and feature extraction methods. Appropriate model conductivity parameters were deduced from measurements and used to replicate features of experimental recordings.

From mitochondrial ion channels to arrhythmias in the heart: computational techniques to bridge the spatio-temporal scales

Computer simulations of electrical behaviour in the whole ventricles have become commonplace during the last few years. The goals of this article are (i) to review the techniques that are currently employed to model cardiac electrical activity in the heart, discussing the strengths and weaknesses of the various approaches, and (ii) to implement a novel modelling approach, based on physiological reasoning, that lifts some of the restrictions imposed by current state-of-the-art ionic models. To illustrate the latter approach, the present study uses a recently developed ionic model of the ventricular myocyte that incorporates an excitation-contraction coupling and mitochondrial energetics model. A paradigm to bridge the vastly disparate spatial and temporal scales, from subcellular processes to the entire organ, and from sub-microseconds to minutes, is presented. Achieving sufficient computational efficiency is the key to success in the quest to develop multiscale realistic models that are expected to lead to better understanding of the mechanisms of arrhythmia induction following failure at the organelle level, and ultimately to the development of novel therapeutic applications.

Construction and validation of a plunge electrode array for three-dimensional determination of conductivity in the heart

The heart's response to electrical shock, electrical propagation in sinus rhythm, and the spatiotemporal dynamics of ventricular fibrillation all depend critically on the electrical anisotropy of cardiac tissue. Analysis of the microstructure of the heart predicts that three unique intracellular electrical conductances can be defined at any point in the ventricular wall; however, to date, there has been no experimental confirmation of this concept. We report the design, fabrication, and validation of a novel plunge electrode array capable of addressing this issue. A new technique involving nylon coating of 24G hypodermic needles is performed to achieve nonconductive electrodes that can be combined to give moderate-density multisite intramural measurement of extracellular potential in the heart. Each needle houses 13 silver wires within a total diameter of 0.7 mm, and the combined electrode array gives 137 sites of recording. The ability of the electrode array to accurately assess conductances is validated by mapping the potential field induced by a point current source within baths of saline of varying concentration. A bidomain model of current injection in the heart is then used to test an approximate relationship between the monodomain conductivities measured by the array, and the full set of bidomain conductivities that describe cardiac tissue.

Multilevel homogenization applied to the cardiac bidomain equations

Accurate cardiac tissue-based modeling using the bidomain equations requires the incorporation of fine-scale structures observed at the 50-100 micron level. By including such features we can more easily observe how defibrillation shocks lead to total depolarization of the heart. Several modeling studies that have investigated the effect of fine scale structures on defibrillation success have been completed. Results have shown that such structures aid, through the creation of virtual electrodes, in total depolarization. An obstacle that occurs with this modeling style is the massive amount of data that must be incorporated into detailed tissue models for even a cubic millimeter sample of cardiac tissue. In this paper, we discuss our approach to generating upscaled, or homogenized, versions of these models that can be used to perform simulations at a more reasonable modeling scale. They have the advantage of incorporating fine scale structure into the model at a reduced modeling cost. We introduce and briefly explore the advantages of this upscaling method.

Models for mechanistic investigations of pacing arrthymogenesis and cardiac tissue structure

Two-dimensional in silico models have been constructed for investigating electrical dynamics arising due to pacing rates and tissue discontinuities. Discontinuities (structures) are introduced explicitly through flux discontinuities in the intracellular space. These models use isotropic electrophysiological properties with heterogeneity arising solely from tissue structure. Both simplified and tissue-specific structures can be used together with membrane models with varying restitution properties. Stimuli are delivered at increasing pacing rates, and solution restarts enable the investigation of repeat stimulation at any of the rates. Investigations using these models are providing new insights into mechanisms leading to conduction block and arrhythmogenesis in cardiac tissue.

Solvers for the cardiac bidomain equations

The bidomain equations are widely used for the simulation of electrical activity in cardiac tissue. They are especially important for accurately modeling extracellular stimulation, as evidenced by their prediction of virtual electrode polarization before experimental verification. However, solution of the equations is computationally expensive due to the fine spatial and temporal discretization needed. This limits the size and duration of the problem which can be modeled. Regardless of the specific form into which they are cast, the computational bottleneck becomes the repeated solution of a large, linear system. The purpose of this review is to give an overview of the equations and the methods by which they have been solved. Of particular note are recent developments in multigrid methods, which have proven to be the most efficient.

Myocardial segment-specific model generation for simulating the electrical action of the heart

BACKGROUND

Computer models of the electrical and mechanical actions of the heart, solved on geometrically realistic domains, are becoming an increasingly useful scientific tool. Construction of these models requires detailed measurement of the microstructural features which impact on the function of the heart. Currently a few generic cardiac models are in use for a wide range of simulation problems, and contributions to publicly accessible databases of cardiac structures, on which models can be solved, remain rare. This paper presents to-date the largest database of porcine left ventricular segment microstructural architecture, for use in both electrical and mechanical simulation.

METHODS

Cryosectioning techniques were used to reconstruct the myofibre and myosheet orientations in tissue blocks of size ~15 x 15 x 15 mm, taken from the mid-anterior left ventricular freewall, of seven hearts. Tissue sections were gathered on orthogonal planes, and the angles of intersection of myofibres and myosheets with these planes determined automatically with a gradient intensity based algorithm. These angles were then combined to provide a description of myofibre and myosheet variation throughout the tissue, in a form able to be input to biophysically based computational models of the heart.

RESULTS

Several microstructural features were common across all hearts. Myofibres rotated through 141 +/- 18 degrees (mean +/- SD) from epicardium to endocardium, in near linear fashion. In the outer two-thirds of the wall sheet angles were predominantly negative, however, in the inner one-third an abrupt change in sheet angle, with reversal in sign, was seen in six of the seven hearts. Two distinct populations of sheets with orthogonal orientations often co-existed, usually with one population dominating. The utility of the tissue structures was demonstrated by simulating the passive and active electrical responses of two of the tissue blocks to current injection. Distinct patterns of electrical response were obtained in the two tissue blocks, illustrating the importance of testing model based predictions on a variety of tissue architectures.

CONCLUSION

This study significantly expands the set of geometries on which models of cardiac function can be solved.

Modeling cardiac electrical activity at the cell and tissue levels

Significant tissue structures exist in cardiac ventricular tissue, which are of supracellular dimension. It is hypothesized that these tissue structures contribute to the discontinuous spread of electrical activation, may contribute to arrhythmogenesis, and also provide a substrate for effective cardioversion. However, the influences of these mesoscale tissue structures in intact ventricular tissue are difficult to understand solely on the basis of experimental measurement. Current measurement technology is able to record at both the macroscale tissue level and the microscale cellular or subcellular level, but to date it has not been possible to obtain large volume, direct measurements at the mesoscales. To bridge this scale gap in experimental measurements, we use tissue-specific structure and mathematical modeling. Our models, which can incorporate ion channel models at the cell level into the reaction-diffusion equations at the tissue level, have enabled us to consider key hypotheses regarding discontinuous activation.

A comparison of multilevel solvers for the cardiac bidomain equations

Computing the extracellular potentials in a bidomain cardiac activation model is a computationally significant step in the solution process. Thus, using a fast solver can drastically reduce the overall time of simulation. Solving for the extracellular potentials involves inverting the matrix coming from the elliptic equation describing the extracellular-intracellular potential coupling. Elliptic equations are known to yield matrices that become progressively more ill-conditioned as the spatial resolution is increased. However, optimal multilevel solution methods are known to exist for these equations given enough effort is placed into developing the correct solution components. Two multilevel solvers that automatically perform much of this work are black box multigrid (BOXMG) and algebraic multigrid (AMG). In this paper, we compare the performance of BOXMG and AMG as solvers for the elliptic component of the bidomain equations. Our investigation is with respect to simulations of reentry in two-dimensional cardiac tissue.

A comparison of monodomain and bidomain reaction-diffusion models for action potential propagation in the human heart

A bidomain reaction-diffusion model of the human heart was developed, and potentials resulting from normal depolarization and repolarization were compared with results from a compatible monodomain model. Comparisons were made for an empty isolated heart and for a heart with fluid-filled ventricles. Both sinus rhythm and ectopic activation were simulated. The bidomain model took 2 days on 32 processors to simulate a complete cardiac cycle. Differences between monodomain and bidomain results were extremely small, even for the extracellular potentials, which in case of the monodomain model were computed with a high-resolution forward model. Propagation of activation was 2% faster in the bidomain model than in the monodomain model. Electrograms computed with monodomain and bidomain models were visually indistinguishable. We conclude that, in the absence of applied currents, propagating action potentials on the scale of a human heart can be studied with a monodomain model.

Cardiac electrophysiology and tissue structure: bridging the scale gap with a joint measurement and modelling paradigm

Significant tissue structures exist in cardiac ventricular tissue that are of supracellular dimension. It is hypothesized that these tissue structures contribute to the discontinuous spread of electrical activation, may contribute to arrhymogenesis and also provide a substrate for effective cardioversion. However, the influences of these mesoscale tissue structures in intact ventricular tissue are difficult to understand solely on the basis of experimental measurement. Current measurement technology is able to record at both the macroscale tissue level and the microscale cellular or subcellular level, but to date it has not been possible to obtain large volume, direct measurements at the mesoscales. To bridge this scale gap in experimental measurements, we use tissue-specific structure and mathematical modelling. Our models have enabled us to consider key hypotheses regarding discontinuous activation. We also consider the future developments of our intact tissue experimental programme.