Numerical Simulation of Calcium Waves in Heart Cells

Matthias K. Gobbert and Bradford E. Peercy

Department of Mathematics and Statistics
University of Maryland, Baltimore County


This page can be reached via the homepage http://www.math.umbc.edu/~gobbert.

Overview

Calcium ions play a crucial role in driving the heart beat. Calcium is released into the heart cell at a lattice of discrete positions throughout the cell, known as calcium release units (CRUs). The probability that calcium will be released from a CRU depends on the calcium concentration at that CRU. When a CRU releases calcium (it 'fires'), the local concentration of calcium ions increases sharply and the calcium that diffuses raises the probability for release at neighboring sites. As a sequence of CRUs begins to release calcium throughout the cell, the release can self-organize into a wave of increasing concentration. If the wave is triggered amongst the normal physiological signaling (i.e., the cardiac action potential) then it can lead to irregular heart beats and possibly a life threatening ventricular fibrillation.

A mathematical model for this effect has been developed by Leighton T. Izu and co-workers [Biophysical Journal, pp. 88-102 and pp. 103-120, 2001]. It consists of a system of three time-dependent partial differential equations, so-called reaction-diffusion equations, which are coupled through non-linear reaction terms. The most important component in the driving terms of the equation for the calcium concentration is the superposition of Dirac delta distributions that models the CRUs as point sources.

Realistic simulations can only be obtained using mathematical models that discretize the the three-dimensional cell domain with a high-resolution numerical mesh. Therefore, a special-purpose computer code is needed that pays great attention to the efficient use of computer memory and that is able to pool the memory from several nodes in a parallel computing cluster. Moreover, large final times are required to reach laboratory time scales and to simulate for long enough to allow several waves to propagate through the cell. The numerical solution of this problem is made additionally challenging by the discontinuities caused by the Dirac delta distributions, for which the classical finite element theory does not apply. This means that standard numerical analysis cannot guarantee the convergence of the numerical method, which is necessary to trust the computed solutions.

Initial work by graduate student Alexander Hanhart [M.S. thesis, 2002], funded by a UMBC DRIF grant, demonstrated that a parallel implementation can achieve the required high resolution of the numerical mesh. Moreover, it confirmed computationally that a finite element method with tri-linear nodal basis functions converges for this problem [Journal of Computational and Applied Mathematics, pp. 431-458, 2004].

Work involving several students, including the undergraduate student Kevin Allen and the high school senior Greg McGlynn, created a full special-purpose code for the application that uses state-of-the-art time-stepping and matrix-free linear solvers. Its confirmed numerical convergence, scalable parallel performance, ability to achieve spontaneous wave-initiation, and correct representation of vital physical effects such as mass conservation make it suitable for long-time simulations on the scales needed for the application [Gobbert, SIAM Journal on Scientific Computing, pp. 2922-2947, 2008]. Such simulations are made possible by the computing cluster tara with state-of-the-art Intel Nehalem processors and extended memory of 24 GB per node that is part of the UMBC High Performance Computing Facility (www.umbc.edu/hpcf), funded in part by the National Science Foundation through MRI and SCREMS grants with additional significant support from UMBC.

Very recent work in the department has allowed us to solve two key issues: The convergence of the finite element method is now rigorously established together with departmental colleague Thomas I. Seidman [Seidman, Gobbert, Trott, and Kruzik, Finite Element Approximation for Time-Dependent Diffusion with Measure-Valued Source, submitted (2011)], in addition to previous computational evidence.

However, the original values of the experimental parameters for the code resulted in unrealistic physiological behavior. This issue was now resolved by graduate student Zana Coulibaly, funded by a National Science Foundation and a UMBC DRIF grant, and his advisor Bradford E. Peercy. A directed search of parameter space indicated by a lower order model has established that the 3-D model and its implementation can yield waves with recovery [Coulibaly, Peercy, and Gobbert, 3-D Cardiac Cell Based on Analysis of a 1-D Deterministic Model, in preparation (2011)].

Moreover, further studies produced spontaneous spiral waves [Coulibaly, Peercy, and Gobbert, Spontaneous Spiral Wave Initiation in a 3-D Cardiac Cell, in preparation (2011)], which are an experimentally observed phenomenon. The first plot shows a simulated confocal-like image of the calcium concentration with a spiral shape that is rotating clockwise. This image of a simulation result is designed to mirror the appearance of an experimental result.


Simulated confocal-like image of the calcium concentration with a spiral wave.
The second plot shows an iso surface that brings out the three-dimensional shape of the region of highly elevated calcium concentration and which also shows the bending spiral shape in the right half of the cell.

Three-dimensional iso surface plot of elevated calcium concentration.

This result demonstrates the significance of the accomplishments in that our simulator produces physiologically correct behavior in full three-dimensional simulations and can provide a better understanding of physiological processes than available from other simulators to date. With this ability, we hope to provide input to help guide in what ranges of parameters to focus experiments, thus making research significantly more effective.

Work on improving the numerical efficiency of the code continues, which will help make runs faster and thus allow for more rapid parameter studies as would be necessary for the support of experiments. In collaboration with graduate student Yu Wang, supported by UMBC as HPCF RA, and Marc Olano from the Department of Computer Science and Engineering at UMBC, we are working on leveraging cutting-edge GPGPUs on each node in the parallel computing cluster to improve the computational kernel of the simulator.


People Involved

Students Formerly Involved


Spontaneous Calcium Waves

For accurate representation of the experiments, it is vital that computationally created waves are both self-initiated and sustained, while also returning the calicum concentrations to basal levels in between waves. This is borne out by the waves with recovery below.
Spontaneous spiral waves are phenomenon that can be seen experimentally. It is extremely exciting that our simulation tool has produced this type of wave below!
Both results by Zana Coulibaly and Dr. Peercy (publications in preparation (2011)).

Long-Time Simulations on High Resolution Meshes to Model Calcium Waves

Long-time simulation of calcium flow in a cell of size 12.8-by-12.8-by-64.0 um with a 15-by-15-by-31 lattice of CRUs using a 64-by-64-by-256 element mesh from time 0 ms to 1000 ms without any initial stimulation. The simulations use the parameters and boundary conditions of the three-dimensional model as specified in Gobbert (2008).

Parallel Performance Studies for Scalar Test Problem

Scalability studies for a time-dependent linear scalar test problem with realistic diffusion coefficients on the domain and on time scales representative of the application problem, studies up to 32 processors using one CPU per node on the distributed-memory cluster kali with progressively finer meshes up to final time 1000 ms. The test problem and more details are contained in Gobbert (CiSE, 2005).

Three-Dimensional Results Using the Finite Element Method


Publications Resulting From This Research

In reverse chronological order

The following list is a partial list with a focus on student involvement in the work; please go to the publication list on my webpage for additional papers.
  1. Zana A. CoulibalyG, Bradford E. Peercy, and Matthias K. Gobbert. Spontaneous Spiral Wave Initiation in a 3-D Cardiac Cell. In preparation.

  2. Zana A. CoulibalyG, Bradford E. Peercy, and Matthias K. Gobbert. Insight Into Spontaneous Recurrent Calcium Waves in a 3-D Cardiac Cell Based on Analysis of a 1-D Deterministic Model. In preparation.

  3. Jonas SchäferG, Xuan HuangG, Stefan Kopecz, Philipp Birken, Matthias K. Gobbert, and Andreas Meister. A Memory-Efficient Finite Volume Method for Advection-Diffusion-Reaction Systems with Non-Smooth Sources. Submitted. Preprint in PDF-format.
  4. Thomas I. Seidman, Matthias K. Gobbert, David W. TrottG, and Martin Kružík. Finite Element Approximation for Time-Dependent Diffusion with Measure-Valued Source. Numerische Mathematik, vol. 122, no. 4, pp. 709-723, 2012. Link to this article at SpringerLink and preprint in PDF-format.
  5. David W. TrottG and Matthias K. Gobbert. Finite Element Convergence for Time-Dependent PDEs with a Point Source in COMSOL 4. In: Yeswanth Rao, editor, Proceedings of the COMSOL Conference 2011, Boston, MA, 6 pages, 2011. Reprint in PDF-format. David W. TrottG and Matthias K. Gobbert. Parallel Performance Studies for a Three-Species Application Problem on the Cluster tara, Technical Report number HPCF-2010-11, UMBC High Performance Computing Facility, University of Maryland, Baltimore County, 2010. Reprint in PDF-format.

  6. David W. TrottG and Matthias K. Gobbert. Finite Element Convergence Studies using COMSOL 4.0a and LiveLink for MATLAB. Technical Report number HPCF-2010-8, UMBC High Performance Computing Facility, University of Maryland, Baltimore County, 2010. Reprint in PDF-format.

  7. Michael MuscedereG, Andrew M. RaimG, and Matthias K. Gobbert. Parallel Performance Studies for a Parabolic Test Problem on the Cluster tara. Technical Report number HPCF-2010-4, UMBC High Performance Computing Facility, University of Maryland, Baltimore County, 2010. Reprint in PDF-format.

  8. Andrew M. RaimG and Matthias K. Gobbert. Parallel Performance Studies for an Elliptic Test Problem on the Cluster tara. Technical Report number HPCF-2010-2, UMBC High Performance Computing Facility, University of Maryland, Baltimore County, 2010. Reprint in PDF-format.

  9. Michael MuscedereG and Matthias K. Gobbert. Parameter Study of a Reaction-Diffusion System Near the Reactant Coefficient Asymptotic Limit. Dynamics of Continuous, Discrete and Impulsive Systems Series A Supplement, pp. 29-36, 2009. Preprint in PDF-format.

  10. Zana CoulibalyUG, Michael MuscedereG, Matthias K. Gobbert, and Bradford E. Peercy. Long-Time Simulation of Calcium Waves in a Heart Cell to Study the Effects of Calcium Release Flux Density and of Coefficients in the Pump and Leak Mechanisms on Self-Organizing Wave Behavior. Technical Report number HPCF-2009-6, UMBC High Performance Computing Facility, University of Maryland, Baltimore County, 2009. Reprint in PDF-format.

  11. Matthias K. Gobbert and Shiming YangG. Numerical Demonstration of Finite Element Convergence for Lagrange Elements in COMSOL Multiphysics. In: Vineet Dravid, editor, Proceedings of the COMSOL Conference 2008, Boston, MA, 6 pages, 2008. Reprint in PDF-format; see also the COMSOL area of my homepage.

  12. Shiming YangG and Matthias K. Gobbert. Convergence Order Studies for Elliptic Test Problems with COMSOL Multiphysics. Technical Report number HPCF-2008-4, UMBC High Performance Computing Facility, University of Maryland, Baltimore County, 2008. Reprint in PDF-format.

  13. Michael MuscedereG and Matthias K. Gobbert. Parallel Performance Studies for a Parabolic Test Problem. Technical Report number HPCF-2008-2, UMBC High Performance Computing Facility, University of Maryland, Baltimore County, 2008. Reprint in PDF-format.

  14. Matthias K. Gobbert. Parallel Performance Studies for an Elliptic Test Problem. Technical Report number HPCF-2008-1, UMBC High Performance Computing Facility, University of Maryland, Baltimore County, 2008. Reprint in PDF-format.

  15. Matthias K. Gobbert. Long-Time Simulations on High Resolution Meshes to Model Calcium Waves in a Heart Cell. SIAM Journal on Scientific Computing, vol. 30, no. 6, pp. 2922-2947, 2008. Direct link to the abstract and reprint in PDF-format.

  16. Matthias K. Gobbert. Configuration and Performance of a Beowulf Cluster for Large-Scale Scientific Simulations. Computing in Science and Engineering (CiSE), vol. 7, no. 2, pp. 14-26, March/April 2005. Link to IEEE's Xplore site for this article and reprint in PDF-format.

  17. Kevin P. AllenUG. Efficient Parallel Computing for Solving Linear Systems of Equations. UMBC Review: Journal of Undergraduate Research and Creative Works, vol. 5, pp. 8-17, 2004.

  18. Kevin P. AllenUG. A Parallel Matrix-Free Implementation of the Conjugate Gradient Method for the Poisson Equation. Senior Thesis, Department of Mathematics and Statistics, University of Maryland, Baltimore County, May 2003.

  19. Kevin P. AllenUG and Matthias K. Gobbert. A Matrix-Free Conjugate Gradient Method for Cluster Computing. Technical Report, University of Maryland, Baltimore County, 2003.

  20. Alexander L. HanhartG, Matthias K. Gobbert, and Leighton T. Izu. A Memory-Efficient Finite Element Method for Systems of Reaction-Diffusion Equations with Non-Smooth Forcing. Journal of Computational and Applied Mathematics, vol. 169, no. 2, pp. 431-458, 2004. Link to Elsevier's ScienceDirect Service for this article and reprint in PDF-format.

  21. Kevin P. AllenUG and Matthias K. Gobbert. Coarse-Grained Parallel Matrix-Free Solution of a Three-Dimensional Elliptic Prototype Problem. In: Vipin Kumar, Marina L. Gavrilova, Chih Jeng Kenneth Tan, and Pierre L'Ecuyer, editors, Computational Science and Its Applications - ICCSA 2003, Part II, Lecture Notes in Computer Science, vol. 2668, pp. 290-299, Springer-Verlag, 2003. Preprint in PDF-format.

  22. Alexander L. HanhartG. Coarse-Grained Parallel Solution of a Three-Dimensional Model for Calcium Concentration in Human Heart Cells. M.S. thesis, Department of Mathematics and Statistics, University of Maryland, Baltimore County, 2002.

  23. Alexander L. HanhartG, Matthias K. Gobbert, and Leighton T. Izu. Coarse-Grained Parallel Solution for a System of Reaction Diffusion Equations. Technical Report, University of Maryland, Baltimore County, 2002.


Disclaimer

The information provided on this page is related to the numerical solution of a system of partial differential equation that is designed to model physiological phenoma in human heart cells. This information does not constitute medical advice. Please see your physician or other medical professional, if you believe to have health problems.


Copyright © 2000-2013 by Matthias K. Gobbert. All Rights Reserved.
This page version 8.0, July 2013.