Virtual Heart Simulator Revolutionizes Cardiovascular Science and Healthcare

According to the World Heart Federation, cardiovascular diseases remain a leading cause of global mortality and are responsible for more than a third of deaths worldwide [1]. Given the severity of this situation, we developed a comprehensive mathematical and computational model called the iHEART simulator that employs innovative mathematical techniques to advance our understanding of cardiac function, aid in clinical decision-making, reduce computational costs, and impact far-reaching medical applications. Developed at Politecnico di Milano’s MOX laboratory with funding from the European Union, the iHEART simulator offers a significant opportunity for personalized medicine to enhance the diagnosis and treatment of cardiovascular pathologies.

Revolutionizing Heart Modeling

The complexity of the human heart presents a unique and intricate healthcare challenge. This complexity, combined with the variability of individual patients, underscores the need to better understand the heart’s functionality. The computational model that underlies the iHEART simulator is a powerful tool that illuminates the heart's intricate mechanisms. This model finds extensive applications in the medical field, allowing healthcare professionals to make informed judgments and devise strategies that address a myriad of heart-related issues. By comprehensively modeling cardiac functions—including electromechanics, fluid dynamics, coronary perfusion, and blood circulation—the iHEART simulator provides a holistic view that can revolutionize clinical decision-making and treatment.

Advancing Mathematical Techniques in Heart Research

To develop this virtual heart, our team analyzed various mathematical and physical aspects, built core cardiac models, and consolidated them into a single comprehensive model. The core models encompass electrophysiology, active and passive mechanics, fluid dynamics, valve kinematics, blood circulation, myocardial perfusion, and the propagation of electric signals throughout the torso. 

Cardiac electrophysiology is modeled via monodomain or bidomain partial differential equations (PDEs), which describe the propagation of action potentials and the excitation of cardiac cells throughout the tissue. These equations are interconnected with appropriate ionic models that offer insights into ionic currents and the subcellular mechanisms that are responsible for excitation:

\[\begin{equation}
\begin{aligned}
&\left\{\begin{aligned}
& \begin{array}{ll}
J \chi C_{\mathrm{m}} \frac{\partial v}{\partial t}-\nabla \cdot\left(J F^{-1} D_{\mathrm{i}} F^{-T} \nabla\left(v+v_{\mathrm{e}}\right)\right) \\
\qquad \qquad  +J \chi I_{\text {ion }}(v, \mathbf{w}, \mathbf{z})=J \chi I_{\text {app }}(\mathbf{x}, t) 
\end{array} & \text { in } \Omega
\\
&-\nabla \cdot\left(J F^{-1} D_{\mathrm{i}} F^{-T} \nabla v\right)-\nabla \cdot\left(J F^{-1}\left(D_{\mathrm{i}}+D_{\mathrm{e}}\right) \nabla v_{\mathrm{e}}\right)=0 & \text { in } \Omega
\end{aligned}\right. \\
& \begin{cases}\frac{\partial \mathbf{w}}{\partial t}=\mathbf{F}_{\text {ion }}^{\mathbf{w}}(v, \mathbf{w}) & \text { in } \Omega \\
\frac{\partial \mathbf{z}}{\partial t}=\mathbf{F}_{\mathrm{ion}}^{\mathbf{z}}(v, \mathbf{w}, \mathbf{z}) & \text { in } \Omega\end{cases}
\end{aligned}    
\end{equation}\]

An important variable in the electrophysiology model is the intracellular calcium concentration: a pivotal factor in the excitation-contraction coupling that causes muscle cells to generate contractile force. We represent this process with a new, highly detailed force generation model that explicitly describes subcellular mechanisms. At the microscale level, these mechanisms are determined by the local calcium concentration and the sarcomere elongation, the latter of which is a function of the local tissue strain along the cardiac fibers:

\[\begin{equation}
\begin{aligned}
& \begin{cases}\frac{\partial \mathbf{s}}{\partial t}=\mathbf{F}_{\mathrm{act}}\left(\mathbf{s},\left[\mathrm{Ca}^{2+}\right]_{\mathrm{i}}, \mathrm{SL}, \frac{\partial \mathrm{SL}}{\partial t}\right) & \text { in } \Omega \\
\mathrm{SL}=\mathrm{SL}_0 \sqrt{I_{4 \mathrm{f}}} & \text { in } \Omega\end{cases} \\
& P_{\mathrm{act}}(\mathbf{d}, \mathbf{s})=T_{\mathbf{act}}\left(n_{\mathrm{f}} \frac{F \mathbf{f}_0 \cdot \mathbf{f}_0}{\sqrt{I_{4 \mathrm{f}}}}+n_{\mathrm{s}} \frac{F \mathbf{s}_0 \cdot \mathbf{s}_0}{\sqrt{I_{4 \mathrm{s}}}}+n_{\mathbf{n}} \frac{F \mathbf{n}_0 \cdot \mathbf{n}_0}{\sqrt{I_{4 \mathrm{n}}}}\right) \\
& T_{\text {act }}=T_{\text {act }}(\mathbf{s}, \mathrm{SL}), \quad I_{4 \mathrm{i}}=F \mathbf{i}_0 \cdot F \mathbf{i}_0 \quad \mathrm{i} \in\{\mathbf{f}, \mathbf{s}, \mathbf{n}\} 
\end{aligned}
\end{equation}\]

The law of momentum balance governs the mechanical strain within the myocardium. The corresponding equations account for myocardial displacement due to interactions with blood, active forces that are triggered by electrical impulses, and the heart muscle’s passive response:

\[\begin{equation}
\begin{gathered}
\rho_{\mathrm{s}} \frac{\partial^2 \mathbf{d}}{\partial t^2}-\nabla \cdot P_{\mathrm{s}}(\mathbf{d}, \mathbf{s})=\mathbf{0} \quad \text { in } \Omega \\
P_{\mathrm{s}}(\mathbf{d}, \mathbf{s})=P_{\mathrm{pas}}(\mathbf{d})+P_{\mathrm{act}}(\mathbf{d}, \mathbf{s}) \\
P_{\mathrm{pas}}(\mathbf{d})=\frac{\partial \mathcal{W}}{\partial F} \quad F=I+\nabla \mathbf{d}
\end{gathered}
\end{equation}\]

The heartbeat arises from the coordinated interaction of the four chambers, which occurs in two primary stages. The first is ventricular systole, including an isochoric phase and ensuing ventricular contraction phase. This stage causes the aortic and pulmonary valves to open when the ventricular pressure surpasses aortic and pulmonary pressures, thereby facilitating blood ejection into the systemic circulation. The second stage is the relaxation phase (i.e., ventricular diastole), which includes a second isochoric phase due to the closure of all heart valves and the subsequent opening of mitral and tricuspid valves, along with ventricular expansion and gradual filling. This phase also encompasses atrial systole, which contributes to ventricular filling. The model uses the incompressible Navier-Stokes equations in a moving domain to describe these processes. 

We utilize the Resistive Immersed Implicit Surface method to incorporate cardiac valves: 

\[\begin{equation}
\begin{aligned}
& \left\{\begin{array}{ll}
-\nabla \cdot P_{\mathrm{ALE}}\left(\mathbf{d}_{\mathrm{ALE}}\right)=\mathbf{0} & \text {in } \hat{\Omega} \\
\mathbf{d}_{\mathrm{ALE}}=\mathbf{d} & \text {on } \hat{\Sigma}
\end{array} \quad \mathbf{u}_{\mathrm{ALE}}=\frac{\partial \mathbf{d}_{\mathrm{ALE}}}{\partial t}\right. \\
& \rho_{\mathrm{f}}\left[\frac{\partial \mathbf{u}}{\partial t}+\left(\left(\mathbf{u}-\mathbf{u}_{\mathrm{ALE}}\right) \cdot \nabla \right) \mathbf{u}\right]-\nabla \cdot \sigma_{\mathrm{f}}(\mathbf{u}, p)+\boldsymbol{\mathcal{R}}\left(\mathbf{u}, \mathbf{u}_{\mathrm{ALE}}\right)=\mathbf{0} \quad \text { in } \Omega \\
& \boldsymbol{\mathcal{R}}\left(\mathbf{u}, \mathbf{u}_{\mathrm{ALE}}\right)=\sum_{\mathrm{k} \in \mathcal{V}} \frac{R_{\mathrm{k}}}{\varepsilon_{\mathrm{k}}} \delta_{\varepsilon_k}\left(\varphi_{\mathrm{k}}^t(\mathbf{x})\right)\left(\mathbf{u}-\mathbf{u}_{\mathrm{ALE}}-\mathbf{u}_{\Gamma_k}\right)
\end{aligned}
\end{equation}\]

We also model the cardiac perfusion process, which represents blood and oxygen supply through blood flow in the coronary tree. This process extends from epicardial coronary arteries—comparable to the organ scale and modeled by three-dimensional (3D) Navier-Stokes equations—to capillaries, which are comparable to the cellular scale and modeled by a multicompartment Darcy model for the intramyocardial coronaries and microcirculation:

\[\begin{equation}
\begin{cases}
\mathrm{NS}(\mathbf{u}, p)=\mathbf{0} & \text{in }\Omega_{\mathrm{cor}} 
\\
\mathbf{K}_{i}^{-1} \mathbf{u}_{\mathrm{myo},i}+\boldsymbol{\nabla} p_{\mathrm{myo},i}=\mathbf{0}, \quad i=1,2,3 & \text{in } \Omega_{\mathrm{myo}} 
\\
\boldsymbol{\nabla} \cdot \mathbf{u}_{\mathrm{myo},1} = \sum_{j=1}^{J} \frac{\chi_{\Omega_{\mathrm{myo}}^{j}}}{\mid \Omega_{\mathrm{myo}}^{j} \mid} \int_{\Gamma_{j}^{\mathrm{cor}}} \mathbf{u} \cdot \textbf{n}-\beta_{1,2} \left( p_{\mathrm{myo},1}-p_{\mathrm{myo},2} \right) & \text{in } \Omega_{\mathrm{myo}} 
\\
\boldsymbol{\nabla} \cdot \mathbf{u}_{\mathrm{myo},2} =-\beta_{2,1}\left( p_{\mathrm{myo},2} -  p_{\mathrm{myo},1} \right)-\beta_{2,3}\left( p_{\mathrm{myo},2} - p_{\mathrm{myo},3} \right) & \text{in } \Omega_{\mathrm{myo}} 
\\
\boldsymbol{\nabla} \cdot \mathbf{u}_{\mathrm{myo},3} =-\gamma\left(  p_{\mathrm{myo},3} -  p_{\mathrm{veins}} \right) - \beta_{3,2}\left( p_{\mathrm{myo},3} - p_{\mathrm{myo},2} \right) & \text{in } \Omega_{\mathrm{myo}}
\end{cases}
\end{equation}\]

Furthermore, we created a coupled model that considers myocardial displacement in electrocardiogram generation by combining a cardiac electromechanical (EM) model with a simulation of the electric potential in the torso. This final step yields the fully coupled 3D cardiac EM model that includes blood hemodynamics and the circulatory system. 

From a numerical standpoint, we made original contributions to improve the underlying mathematical models by developing new multiscale and variable-order numerical approximations, employing operator splitting and staggering techniques for time advancement, and incorporating ad hoc domain decomposition and multigrid parallel preconditioners for high-performance computing. Additionally, we combined physics-based algorithms with novel machine learning algorithms to address data and model uncertainty and achieve reduced order models (ROMs). 

Given the diverse range of core mathematical models and their heterogeneity within PDEs, our cardiovascular problem required an ad hoc approach for numerical approximation. We used a baseline Galerkin finite element method for spatial approximation that employs both tetrahedra and hexahedra. The local polynomial degree is nonuniform and utilizes high-degree polynomials when necessary for compatibility conditions and dispersion error control in fast propagating wave problems. 

Ultimately, we captured the complex variability of cardiac function by developing parameterized models that align with clinical data. The efficiency of ROMs facilitated sensitivity analysis and uncertainty quantification, which revealed underlying mechanisms and parameter interplays. Additionally, we addressed parameter estimation and inverse problems that could potentially contribute to personalized medicine. Our iHEART simulator harmoniously encapsulates all of these facets into a biophysically detailed model of cardiac function.

The iHEART Simulator

The iHEART simulator seamlessly replicates intricate cardiac processes within a unified platform and simulates heart function and heart-related ailments with an unparalleled level of biophysical precision. This groundbreaking advancement in the field of computational cardiology offers future researchers a comprehensive means of simulating the heart.

In conjunction with the simulator, we created a freely available parallel C++ library for high-performance finite element simulations called life^x. This library encompasses complex cardiac functions and is accessible to a wide community of users, including those with backgrounds in medicine and bioengineering; it can capably handle complex cardiac simulations as well as multiphysics, multiscale, and multidomain problems. The versatility, efficiency, and user-friendliness of life^x make it a valuable computational tool that can optimally scale up to thousands of computer cores. The interplay between physics-based numerical models and data-driven machine learning algorithms characterizes the scientific machine learning approach that empowers the iHEART simulator.

The roadmap that guides our simulator begins with a preprocessing phase that uses clinical images for tasks such as image segmentation, surface and volume meshing, and fiber reconstruction. We then draw on multiphysics and multiscale techniques to approximate fundamental cardiac models, after which we employ solution algorithms for time advancement, problem linearization, and the effective parallel implementation of scientific machine learning strategies through the life^x library. Afterwards, the phases of software verification, numerical verification, and validation ensure the results’ quality and reliability. See Figure 1 for an example of a numerical simulation of a beating heart from the iHEART simulator.

Beyond Science: Impacting Lives

The significance of this endeavor transcends the realm of research laboratories. The iHEART project envisions a future where digital replicas of individual hearts help clinicians diagnose conditions, plan surgeries, and predict outcomes with unprecedented accuracy. And the journey does not stop with scientific discovery; it extends into tangible medical applications. Collaborations with clinical divisions around the world have already begun to yield promising results. We are actively partnering with leading healthcare institutions to study cardiac arrhythmias, optimize treatments like cardiac resynchronization therapy (CRT) and transcatheter aortic valve implantation (TAVI), and address conditions like atrial fibrillation, cardiac myocardial perfusion, and obstructions in hypertrophic cardiomyopathy.

Collaborative research in partnership with the cardiac arrhythmology and electrophysiology units at IRCCS Ospedale San Raffaele and Humanitas Research Hospital, both in Italy, has revealed valuable quantitative insights into the factors that contribute to the onset and persistence of arrhythmias — specifically post-ischemic ventricular tachycardia. Our simulator integrates data from patients' electrical measurements to enable the visualization and in-depth analysis of reentry mechanisms and their progression, which are usually difficult to assess in clinical practice. Our efforts highlight the significance of variations in conduction and refractoriness properties of the atrial myocardium for sustained arrhythmic events. 

Furthermore—and in cooperation with the cardiology and radiology departments at Santa Maria del Carmine di Rovereto in Italy—we are working to enhance CRT: a cardiac treatment that implants a device to restore proper heart rhythm synchronization after a disruption from conduction issues or scarring. Our research has generated a mathematical tool that streamlines the mapping process and reduces both the number of invasive procedures and the time that is required to calculate the latest activation point. 

We are also partnering with the Monzino Cardiology Center in Italy to investigate the functionality of TAVI, which has become an increasingly popular treatment for aortic stenosis. We have identified predictive risk indicators prior to surgery to help physicians make informed decisions about surgical procedures and post-implantation monitoring. In addition, we are likewise actively engaged in cardiac perfusion modeling with the Monzino Cardiology Center. Our iHEART simulator—which enables us to gauge the impact of factors such as coronary artery stenosis on blood supply in different myocardial regions—has the potential to complement or even replace exercise testing with an alternative way to evaluate perfusion under stress conditions.

A Visionary Future

The effective preoperative use of the iHEART simulator can model the cardiac functions of specific patients to characterize pathophysiological mechanisms and assist in surgical decision-making. Numerical simulations allow us to virtually explore different therapeutic scenarios—like cardiac ablation, valve implantation, or myectomy surgery—and evaluate the potential efficacy and risks of various interventional strategies. The simulator also facilitates the estimation of important quantities—such as tissue stresses, blood pressure within chambers, blood particle residence time, and reentrant drivers—that are difficult to directly measure without invasive acquisition methods. High-resolution reconstructions of these parameters contribute to more comprehensive preoperative and postoperative evaluations.

Ultimately, our iHEART simulator embodies a synergy between scientific innovation and healthcare needs. It foresees a future where digital simulations guide medical interventions, heralding a new era of personalized heart healthcare. We aim to empower clinicians with tools that transform patient care, one heartbeat at a time.

Alfio Quarteroni received the 2023 Lagrange Prize from the International Council for Industrial and Applied Mathematics for his groundbreaking work in finite element and spectral methods, domain decomposition methods, discontinuous Galerkin methods, and numerical solution of incompressible Navier-Stokes equations. He delivered a corresponding prize lecture about the iHEART simulator at the 10th International Congress on Industrial and Applied Mathematics, which took place in Tokyo, Japan, last year.

References 
[1] Di Cesare, M., Bixby, H., Gaziano, T., Hadeed, L., Kabudula, C., McGhie, D.V., … Pinto, F. (2023). World heart report 2023: Confronting the world’s number one killer. Geneva, Switzerland: World Heart Federation.
[2] Fedele, M., Piersanti, R., Regazzoni, F., Salvador, M., Africa, P.C., Bucelli, M., … Quarteroni, A. (2023). A comprehensive and biophysically detailed computational model of the whole human heart electromechanics. Comput. Methods Appl. Mech. Eng., 410, 115983.
[3] Zingaro, A., Bucelli, M., Piersanti, R., Regazzoni, F., Dede', L., & Quarteroni, A. (2024). An electromechanics-driven fluid dynamics model for the simulation of the whole human heart. J. Comput. Phys., 504, 112885.

Further Reading
• Bucelli, M., Zingaro, A., Africa, P.C., Fumagalli, I., Dedè, L., & Quarteroni, A. (2023). A mathematical model that integrates cardiac electrophysiology, mechanics and fluid dynamics: Application to the human left heart. Int. J. Numer. Method. Biomed. Eng., 39(3), e3678. 
• Piersanti, R., Africa, P.C., Fedele, M., Vergara, C., Dedè, L., Corno, A.F., & Quarteroni, A. (2021). Modeling cardiac muscle fibers in ventricular and atrial electrophysiology simulations. Comput. Methods Appl. Mech. Eng., 373, 113468.
• Regazzoni, F., Dedè, L., & Quarteroni, A. (2020). Machine learning of multiscale active force generation models for the efficient simulation of cardiac electromechanics. Comput. Methods Appl. Mech. Eng., 370, 113268.
• Zingaro, A., Vergara, C., Dedè, L., Regazzoni, F., & Quarteroni, A. (2023). A comprehensive mathematical model for cardiac perfusion. Sci. Rep., 13, 14220.