lifex-ep: a robust and efficient software for cardiac electrophysiology simulations

Background Simulating the cardiac function requires the numerical solution of multi-physics and multi-scale mathematical models. This underscores the need for streamlined, accurate, and high-performance computational tools. Despite the dedicated endeavors of various research teams, comprehensive and user-friendly software programs for cardiac simulations, capable of accurately replicating both normal and pathological conditions, are still in the process of achieving full maturity within the scientific community. Results This work introduces \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\texttt {life}^{\text{x}}$$\end{document}lifex-ep, a publicly available software for numerical simulations of the electrophysiology activity of the cardiac muscle, under both normal and pathological conditions. \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\texttt {life}^{\text{x}}$$\end{document}lifex-ep employs the monodomain equation to model the heart’s electrical activity. It incorporates both phenomenological and second-generation ionic models. These models are discretized using the Finite Element method on tetrahedral or hexahedral meshes. Additionally, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\texttt {life}^{\text{x}}$$\end{document}lifex-ep integrates the generation of myocardial fibers based on Laplace–Dirichlet Rule-Based Methods, previously released in Africa et al., 2023, within \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\texttt {life}^{\text{x}}$$\end{document}lifex-fiber. As an alternative, users can also choose to import myofibers from a file. This paper provides a concise overview of the mathematical models and numerical methods underlying \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\texttt {life}^{\text{x}}$$\end{document}lifex-ep, along with comprehensive implementation details and instructions for users. \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\texttt {life}^{\text{x}}$$\end{document}lifex-ep features exceptional parallel speedup, scaling efficiently when using up to thousands of cores, and its implementation has been verified against an established benchmark problem for computational electrophysiology. We showcase the key features of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\texttt {life}^{\text{x}}$$\end{document}lifex-ep through various idealized and realistic simulations conducted in both normal and pathological scenarios. Furthermore, the software offers a user-friendly and flexible interface, simplifying the setup of simulations using self-documenting parameter files. Conclusions \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\texttt {life}^{\text{x}}$$\end{document}lifex-ep provides easy access to cardiac electrophysiology simulations for a wide user community. It offers a computational tool that integrates models and accurate methods for simulating cardiac electrophysiology within a high-performance framework, while maintaining a user-friendly interface. \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\texttt {life}^{\text{x}}$$\end{document}lifex-ep represents a valuable tool for conducting in silico patient-specific simulations.


Background
Cardiac electrophysiology focuses on the heart conduction system from both the normal and pathological perspectives: e.g., it involves the study, diagnosis and treatment planning of cardiac arrhythmias [1,2].
Nowadays, several clinical tools are widely employed to address these rhythm disorders.The Electrocardiogram (ECG) provides a recording of the electrical activity of the heart.Various deviations from sinus rhythm, such as Atrial Fibrillation (AF), Ventricular Tachycardia (VT), Left Bundle Branch Block (LBBB), can be monitored and identified based on the distinctive shape and morphology of the ECG on an individualized patient-specific level [3].The ECG may be combined with other imaging data, such as Magnetic Resonance Imaging (MRI) and Computed Tomography (CT), or electroanatomical maps that are directly recorded on the internal and external surfaces of the heart.Thanks to the aforementioned tools, clinicians can successfully reconstruct both the heart anatomy of a patient and multiple electrophysiology properties of the cardiac tissue.These information lay the foundations for traditional decision-making in cardiology.Indeed, medications and surgeries, such as the implantation of pacemakers or cardioverter defibrillators, are planned accordingly [4].
In recent years, the advent of computer models and in-silico simulations in cardiology has enabled the integration of novel scientific tools in standard clinical practice [5,6].Physics-based mathematical models and data-driven methods are combined to generate digital replicas of human hearts containing a detailed electrophysiology description, both at the cellular level and at the organ scale [7].Patient-specific clinical data are embedded inside numerical simulations of the cardiac function for precision medicine [8].For instance, personalized electrophysiology simulations are employed to assess risk stratification of arrhythmias, to define optimal catheter-based ablation targets, and to perform Cardiac Resynchronization Therapy (CRT) [9][10][11][12].
In this work, we present life x -ep, a publicly available software specifically designed for conducting numerical simulations of cardiac electrophysiology, encompassing both normal and pathological conditions.life x -ep is built upon the foun- dation of life x [13], an open-source, high-performance C++ Finite Element (FE) numerical solver capable of tackling multi-physics, multi-scale, and multi-domain differential problems.Leveraging the deal.II [14] FE core, life x was conceived as part of the iHEART project (refer to Section "Funding") with the primary aim of providing the scientific community with a cutting-edge FE solver for cardiac modeling.
Overall, developing a personalized computer model of the human heart electrophysiology that accurately represents all the underlying biological aspects has garnered significant attention.Crucial to achieving this goal is the development of efficient parallel numerical algorithms that can handle computationally intensive tasks with high accuracy.In this context, life x -ep proposes a software that offers a user-friendly yet very detailed and highly customizable interface, a variety of modeling and numerical options to choose from, while also ensuring accuracy and computational efficiency, representing a valuable tool for conducting in silico patient-specific simulations.
The present release focuses on the modeling of cardiac electrophysiology.In the following paragraphs, we provide a concise overview of cardiac electrophysiology and the prevailing mathematical methods employed to model it.

Cardiac electrophysiology: physiology and modeling
The heart wall consists of three distinct layers: the internal thin endocardium, the external thin epicardium and the thick muscular cardiac tissue known as the myocardium.The latter is predominantly composed of cardiomyocytes, which are specialized, striated excitable muscle cells responsible for the essential cardiac function.When these cardiomyocytes are stimulated by an electrical impulse, a change in the electro-chemical balance of the cell membrane results in a series of biochemical reactions that determine a large variation of the transmembrane potential, namely the voltage differential between the intra and extracellular spaces of the cell.The rapid depolarization and subsequent slow repolarization mechanism is known as the action potential.This process is triggered and controlled by the opening and closing of voltage-gated ion channels that make the cell membrane permeable to specific ionic species, like sodium, potassium, and calcium.The transmembrane potential changes as a result of the ionic fluxes, which, in turn, are driven by the voltage difference itself.
The cells of the heart tissue are connected to each other through gap junctions, i.e. intercellular low-resistance ionic channels, which allow the electric signal to travel from cell to cell across the whole cardiac muscle (known as myocardium).The propagation of the electric signal within the myocardium is highly anisotropic, with preferential directions for conduction determined by the presence of muscle fibers [17].These fibers are organized in sheets, that determine a second preferential direction of propagation.For Fig. 1 The life x library encompasses essential features and a framework for numerically solving Finite Element problems.life x -ep is a publicly released package designed for cardiac electrophysiology simulations based on life x .Left picture icons sourced from https:// fonta wesome.com/ licen se further details on the biophysical mechanisms of cardiac electrophysiology, we refer the interested reader to [1,2,41].A visual representation of the different spatial scales involved in this complex chain of physical processes is given in Fig. 2.
The mathematical description of these processes is made of two building blocks: a ionic model, describing the chemical processes taking place at the cellular scale, and an action potential propagation model, representing the spatial propagation of the action potential wavefront at the tissue scale.

Ionic models
Several ionic models with different levels of biophysical detail and suitable for specific types of cells have been proposed in the literature [41].Most of them are written in the form of the following Ordinary Differential Equations (ODEs) where the unknowns are u = u(t) , the transmembrane potential, and the vector w = w(t) = (w 1 , . . ., w M ) , collecting M ionic variables.Ionic variables may include concentrations of different ionic species and gating variables (describing the opening probability of ionic channels).The term I ion describes the electric current generated by the flux of ionic species across the cell membrane, while I app represents an externally applied current.The dynamics of ionic variables is modeled by the function H . Different ionic models are characterized by the number of ionic variables M and the definition of the two functions I ion and H.

Action potential propagation models
One of the most popular action potential propagation models is the monodomain equation, a Partial Differential Equation (PDE) in which the spatial propagation of the electric signal through gap junctions is accounted for by a diffusion term [41].
(1) Referring to Fig. 3, we consider a computational domain ⊂ R 3 , representing the region of the myocardium of interest (e.g., a slab of cardiac tissue, atria, ventricles or whole heart geometries), see Fig. 3(a).In the whole domain, we consider a local orthonormal triplet of vectors f 0 , s 0 and n 0 , defining the fiber, the sheetlet and the sheet-normal directions see Fig. 3(b), respectively [17].
A notable feature of life x -ep is its capability to divide the domain into a generic number N of subdomains, thereby allowing distinct electrophysiology properties to be assigned to each individual region (e.g.different ionic models or different electrical conductivities).In practical implementation, this goal is realized by providing the simulation with a computational mesh, wherein each subdomain is distinctly labeled.More formally, we introduce a partition of into N disjoint subdomains, namely � 1 , � 2 , . . ., � N (more precisely, we assume = ∪ N i=1 i and i ∩ j = ∅ for i = j ).We denote by Ŵ i = ∂� ∩ ∂� i the external boundary of the i-th subdomain, and by Ŵ ij = ∂� i ∩ ∂� j the interface between the i-th and the j-th subdomains, see Fig. 3(a).
For each subdomain i , we consider a ionic model, featuring M i ionic variables, and characterized by the functions H i and I i ion and initial conditions u i 0 and w i 0 .Moreover, we introduce the unknown w i : i × [0, T ] → R M i .The transmembrane potential, namely u : × [0, T ] → R , is instead defined in the whole computational domain.
The monodomain model reads as follows [5,41] where i, j = 1, . . ., N , and the diffusion tensor D i is defined as (2) where σ i l , σ i t , σ i n being positive coefficients denoting the longitudinal, transversal and normal conductivities, respectively [17].The operator ⊗ denotes the outer product acting on the vectors f 0 , s 0 and n 0 .Given two vectors v and w , their outer product is a tensor, component-wise defined as (v ⊗ w) ij = v i w j , for i, j = 1, 2, 3 .Homogeneous Neu- mann boundary conditions are prescribed on the whole boundary ∂� to impose the con- dition of electrically isolated domain, with n i denoting the outward normal unit vector to the boundary.On the interface between subdomains, we impose continuity of flux conditions.
The action potential is triggered by an external applied current I app (x, t) , that mim- ics the presence of a natural or artificial pacemaker, or the role of specialized conduction systems, such as the Purkinje network [1,2] or other relevant bundles.Since these anatomical entities are not explicitly modeled in this release, life x -ep provides great flexibility in the definition of I app (x, t) , allowing for the definition of multiple stimuli of various shape and with customized duration and delay.The applied current is assumed to be in the following form: wherein N stim ∈ N is the number of applied stimuli, g t i : (0, T ) → {0, 1} is an indica- tor function that only activates the i-th stimulus during a prescribed time interval, and g x i : � → [0, ∞) is a space-dependent function that describes the shape of the stimulus, defined either as the indicator function of a given subset of (for cubic, spherical or planar impulses) or as a Gaussian function.Location, shape and timing of each stimulus can be customized by the user.
The possibility of prescribing different properties to different regions of the domain enables several use cases.For instance, it allows multi-chamber simulations, in which different ionic models are assigned to atria and ventricles [24].Moreover, it also gives the possibility to simulate pathological scenarios, such as the presence of scars or fibrosis [18,25,46,47].Non-conductive zones can be introduced in the model, by labelling a given subregion of the domain as non-conductive.In this case, the corresponding portion of the domain will be excluded from the simulation, thus creating a block in the conduction [24].Furthermore, grey zones can be modeled as well, by suitably adjusting the tissue conductivities and the ionic model properties in grey zone and fibrotic subregions [18,25,48,49].
We remark that the monodomain equation is often expressed as where C m denotes the membrane cell capacitance and χ m is the membrane surface-to- volume ratio.This formulation can clearly be recast to the one of (2) by the rescaling (3)

Cardiac electrophysiology: numerical discretization
We partition the temporal domain (0, T] into N T subintervals with a time step t = t n+1 − t n , where t n = n t for n = 0, 1, . . ., N T .We denote with a subscript n the approximation of a variable at time t n (e.g., u n ≈ u(t n ) ).Time derivatives are approxi- mated by means of a Backward Differentiation Formula (BDF) scheme of order σ , with σ ∈ {1, 2, 3} .Given a generic time-dependent function f, its time derivative is thus approxi- mated by where f BDF,n is a linear combination of f n , f n−1 , . . ., and α BDF is a coefficient depend- ing on the order of the scheme [50].We also define explicit extrapolations of order σ , denoted by the subscript EXT (e.g., f EXT,n+1 ≈ f (t n+1 )).
The numerical discretization of the system (2) is obtained by means of the implicitexplicit time advancing scheme described in [26].Given the solution up to time step t n , to compute the solution at time t n+1 : 1 solve the time-discrete ionic model equations: for i = 1, . . ., N , We treat some of the ionic variables appearing in H i with an implicit formulation, and others with an explicit formulation, so that the resulting problem can be solved by directly inverting the equations [26]. 2 solve the time-discrete monodomain equation: for i, j = 1, . . ., N , Notice that the ionic current term is treated explicitly by using the extrapolated potential u EXT,n+1 , so that the resulting problem is linear in u n+1 .For the spatial discretization, we introduce a tetrahedral or hexahedral mesh over , and use the FE method [51] to approximate the solution variables u and w i as piecewise polyno- mials of order p .life x -ep supports polynomials of orders 1 and 2 on tetrahedral meshes, and polynomials of arbitrary degree on hexahedral meshes.
The ionic model ( 4) is solved at each support point of the degrees of freedom of the FE space.The ionic current I ion is also evaluated at every support point, and then interpo- lated onto quadrature nodes on the interior of the mesh elements, in the approach known as Ionic Current Interpolation (ICI) [52,53].The time-discrete monodomain equation ( 5) is discretized in space using the FE method, leading to the following algebraic linear system of equations ( where, denoting by ϕ j the basis functions of the FE space, M is the mass matrix of entries M jk = � ϕ k ϕ j dx , K is the stiffness matrix of entries K jk = N i=1 � i D i ∇ϕ k • ∇ϕ j dx and f n+1 is the vector arising from the applied current term I app , of entries (f n+1 ) i = � I app,n+1 ϕ i dx .The vector s n+1 arises from the ionic current terms I i ion .Exploiting the ICI formulation, it can be computed as where M i is the mass matrix for subdomain i , whose entries are M i jk = � i ϕ k ϕ j dx , and I i ion,n+1 is the vector of the evaluations of I i ion at the support points of the FE space.The integrals that arise from the FE discretization are numerically approximated using the Gauss-Legendre quadrature rule, with the minimum number of points to ensure exact integration of the mass matrix.

Overview of existing software
One of the first landmarks in computational cardiology has been defined by the pioneering work of Hodgkin and Huxley in the mid-1950s [54].The first cardiac action potential model has been developed by D. Noble in 1962 [55].Since then, a plethora of cardiac electrophysiology models has been proposed in the literature [6,41,56,57].These mathematical models feature a different degree of biophysical complexity and act either at the microscopic level, by describing the behavior of single cardiomyocytes, or at the organ scale, where an ensemble of many myocardial cells is considered.
The development of software to perform electrophysiology simulations is still mainly steered by academia and public institutions, despite this variety of mathematical models, the presence of increasingly elaborated numerical methods, and the evolution of computer hardware.Indeed, several research tools for multi-physics and multi-scale cardiac simulations have been proposed in the last two decades.Among them, an important role is played by openCARP [58], an open-source C++ simulation environment integrated with cellML, a public repository that encompasses many cell-based mathematical models [59], Chaste [60], an open-source C++ library mainly developed at the University of Oxford, Cardioid [61], a highly efficient and scalable tool for high resolution electrophysiology simulations mainly developed at the Lawrence Livermore National Laboratory, and Alya [62], a high-performance code from the Barcelona Supercomputing Center.Prominent examples of industrial projects that demonstrate translational research efforts for the heart are given by the Living Heart Project by Dassault Systémes [63,64] and the services and software provided by NumeriCor GmbH 1 [65].
The aforementioned tools allow to perform single-chamber, bi-atrial, bi-ventricular or four-chamber heart electrophysiology simulations by means of accurate, yet computationally expensive physiologically-based models, such as the bidomain or monodomain equation coupled with the TTP06 and CRN ionic models [41,44,66], or the (6) 1 https:// www.numer icor.at/ more efficient reaction-eikonal equation [67].Model parameters can be calibrated on a patient-specific basis to match ECG or Body Surface Potential Mapping (BSPM) [68].Different pathological scenarios involving AF [69,70], VT [9,12] and LBBB [11,71] have been addressed by employing these software tools.Moreover, heterogeneity in the tissue and cellular properties can be prescribed to incorporate the presence of scar, grey zones and fibrosis importing measurements from clinical data, such as Late Gadolinium Enhancement-Magnetic Resonance Imaging (LGE-MRI) [10], contrast enhanced CT [72], or the Imaging Itensity Ratio (IIR) [73].This allows to locally vary Conduction Velocity (CV) and Action Potential (AP) morphology.Another application is related to the in-silico assessment of drugs efficacy by means of numerical simulations [74,75], where model parameters can be tuned to replicate the effects of pharmacological therapies.Electrophysiology simulations are also used to evaluate gender differences in healthy and pathological conditions involving arrhythmias [75,76].
Compared to the aforementioned software and libraries, the life x -ep solver stands out with several distinctive features and advantages, mostly inherited from the life x core structure [13].It is designed to be user-friendly and easy to use, even for biomedical researchers without extensive experience in numerical methods.The solver is implemented in C++ using advanced programming paradigms and leverages Message Passing Interface (MPI) for distributed memory parallelism.Moreover, it supports the possibility to import arbitrary meshes with either hexahedral or tetrahedral elements, and incorporates advanced numerical solvers based on the Trilinos linear algebra backend, thus ensuring precise control over the numerical setting and accuracy.The solver also exhibits ideal scalability up to thousands of cores, as demonstrated in Section Strong scalability test, allowing to efficiently simulate large-scale scenarios.In addition to the numerical and programming features stemming from its foundation on life x , life x -ep offers two options for prescribing myocardial fibers, which can be either imported from a file or generated online taking advantage of the previous release life x -fiber [15], based on the LDRBMs presented in [17].Moreover, it supports spatial heterogeneity in the choice of both models and physical coefficients, easily configurable through a convenient parameter file, without the need to access and modify the source code.
In general, life x -ep stands out in its ease of use, performance, and compatibility with common I/O (input/output) file formats, as well as its comprehensive and self-contained infrastructure, achieved by combining sophisticated mathematical models with accurate numerical schemes.To the best of our knowledge, none of the packages mentioned above exhibits similar features altogether.

Implementation
In this section, we present the technical specifications of life x -ep and provide a com- prehensive documentation of its user interface.The aim is to guide users through the entire process, from downloading the software to successfully running a full cardiac electrophysiology simulation.life x -ep offers a numerical solver tailored for cardiac electrophysiology, leveraging the mathematical models and numerical algorithms discussed in the previous section.The linear algebra backend is provided by Trilinos [77], integrated into deal.II.This incorporates the implementation of various linear solvers (CG and GMRES) and flexible black-box preconditioners (AMG, additive Schwarz, block Jacobi) supported by life x .The code is inherently parallel and designed to run efficiently on a diverse range of architectures, spanning from personal laptops to High-Performance Computing (HPC) facilities and cloud platforms.To ensure reliability and performance, life x -ep has been thoroughly tested on multiple sys- tems, including a cluster node equipped with 192 cores based on Intel Xeon Gold 6238R (2.20 GHz) at MOX, Dipartimento di Matematica, Politecnico di Milano, as well as the GALI-LEO100 supercomputer available at CINECA (Intel CascadeLake 8260, 2.40GHz, technical specifications available at https:// wiki.u-gov.it/ confl uence/ displ ay/ SCAIUS/ UG3. 3 For further insights into the core functionalities of life x , we recommend referring to [13].Additionally, in the following sections, we provide a concise guide on how to quickly get started and run simulations using life x -ep.

Running simulations in life x -ep
The distribution and installation process of life x -ep is designed to be user-friendly and platform-independent.The software is conveniently provided as a binary AppImage 2 executable, which can be obtained from https:// doi.org/ 10. 5281/ zenodo.80852 66.Along with the executable, all the necessary input files required to reproduce the numerical results presented in the subsequent sections are included.
The adoption of the AppImage format ensures a universal package that is compatible with x86-64 Linux operating systems, eliminating the need for multiple distribution-specific versions.From the user's perspective, this translates to a seamless and straightforward download-then-run experience, without the hassle of manually managing system dependencies.
As an AppImage, life x -ep has been built on Debian Buster, 3 which corre- sponds to the current oldoldstable version.This follows the principle of "Build on old systems, run on newer systems". 4Consequently, the software is expected to function on virtually any recent x86-64 Linux distribution, provided that glibc 5 version 2.28 or higher is installed.
Once downloaded and extracted the life x -ep archive, the AppImage file needs to be made executable by typing the following command in a terminal: Finally, lifex_electrophysiology-1.5.0-x86_64.AppImage can be executed with: No root permissions are required for the commands mentioned above to run successfully.However, it is important to note that the AppImage relies on the userspace filesystem framework called FUSE. 6Please ensure that FUSE is installed on your system.If you encounter any errors, the following commands may be helpful in resolving the issue: Additionally, we recommend referring to the AppImage troubleshooting guide. 7he following command will provide an inline help that includes detailed information about all the available command line options and their purpose: The executable allows to run test cases with an arbitrary number of disjoint subdomains � i , i = 1, . . ., N , which are also referred to as volumes.The configuration of the simulation is supplied through a parameter file.The user can generate a template parameter file using the following command: To match different user needs, the level of detail in the parameter files can be adjusted using the optional minimal or full option after the -g flag.The minimal option reduces the level of detail, making it suitable for initial usage of life x -ep.On the other hand, the full option increases the level of detail, exposing advanced options, such as parameter choices and detailed options on linear algebra and preconditioning, among others.If the user does not specify either the minimal or full option, an intermediate verbosity level is selected by default.In all cases, parameters that are not present in the file will retain their default values.This flexibility allows users to customize the level of detail in the parameter files according to their specific requirements and familiarity with the software.
The parameters are written in a plain text file, organized as a list of key-value pairs grouped in subsections, which describe the configuration for the simulation to be run.Each parameter is accompanied by a brief documentation within the parameter file itself, explaining its meaning.
In the provided command, the optional argument -vol <volume labels> allows to specify a list of user-defined labels.These labels are used to differentiate each subdomain, enabling the selection of heterogeneous model options such as ionic model type, coefficients, and electrical conductivities for each subdomain.If not specified, a single subdomain characterized by its global Volumetric parameters is assumed to exist.
The following example illustrates a parameter file specifying three subdomains: Healthy, Fibrosis, and Scar: The input mesh is expected to include at least one volumetric tag corresponding to each of the subdomains that need to be differentiated.These subdomain tags are specified in the Material IDs list, located under each respective subdomain section.
The parameter file includes a dedicated section called Fiber generation, which serves the purpose of enabling the importing from a file of the myocardial fibers or the online generation of them on various geometry types, such as slabs, ventricles, and atria.To accomplish this, life x -ep incorporates the functionalities of its prede- cessor, life x -fiber [15], which utilizes the Laplace-Dirichlet Rule-Based Methods presented in [17].This unique feature of life x -ep sets it apart from other existing software alternatives.
Once the user has edited the parameter file, the simulation can be started using the following command: When executing the command, the volume labels provided must match the ones used for generating the parameter file.It is essential to use consistent volume labels throughout the process to ensure proper identification and configuration of the subdomains within the simulation.
A parallel simulation is started prepending the command with the mpirun or mpiexec wrapper (which may vary depending on the MPI implementation available), e.g.: where N represents the desired number of parallel processes.The binary package supports parallel execution using MPICH (https:// www.mpich.org/) version 4.0 or higher.
The parameter file also includes options that enable the serialization of the solution, allowing the simulation to be paused or stopped at any point and then resumed at a later time using the serialized data.This feature is particularly useful when dealing with longrunning simulations or when unexpected interruptions occur.

License and third-party software
This work is copyrighted by the life x -ep authors and licensed under the Creative Commons Attribution Non-Commercial No-Derivatives 4.0 International License. 8 It should be noted that life x -ep incorporates several third-party libraries, which are separately copyrighted by their respective authors and whose use is covered by various permissive licenses.
Third-party software bundled with (in binary form), required by, copied, modified, or explicitly used in life x -ep include: life x9 [13]: the open-source, high-performance software providing the core functionalities for the numerical solution of the Finite Element problems described in the previous section; deal.II 10  [14]: it provides support to mesh handling, assembling and solving Finite Element problems (compiled with enabled support to Trilinos 11 for linear algebra data structures and solvers) and to input/output functionalities; Boost 12 [78]: its modules Filesystem and Math are used for manipulating files/directories and for advanced mathematical functions and interpolators, respectively; VTK 13  [79]: it is used for importing external surface or volume input data and coefficients appearing in the mathematical formulation.
Some of the packages listed above, as stated by their respective authors, rely on additional third-party dependencies that may also be bundled (in binary form) with life x -ep, although not used directly.These dependencies include: ADOL-C, 14 ARPACK-NG, 15  BLACS, 16 Eigen, 17 FFTW, 18 GLPK, 19 HDF5, 20 HYPRE, 21 METIS, 22 MPICH, 23 MUMPS, 24  NetCDF, 25 OpenBLAS, 26 PETSc, 27 ParMETIS, 28 ScaLAPACK, 29 Scotch, 30 SuiteSparse, 31 SuperLU, 32 oneTBB, 33 p4est. 34 These libraries are free software and have relatively few restrictions on their use.However, please note that different terms may apply.For detailed information on the licenses and copyright statements for these packages, please refer to the content of the folder share/doc/licenses/ in the life x -ep archive.

Results and discussion
To highlight the versatility of life x -ep, this section presents the results obtained in various numerical examples conducted on a range of idealized and realistic geometries, encompassing both normal and pathological scenarios.To facilitate the reproduction of these test cases, this section provides a getting started guide, detailing the following steps: • Generating and importing the input data (e.g., computational meshes and fibers); • Configuring the parameter files specific to each test case and executing the corresponding simulation; • Performing post-processing of the results and visualizing the output.
Finally, a strong scalability test is presented in Section Strong scalability test.

Input data
All parameter files and meshes related to the numerical simulations described below can be downloaded from the life x -ep release archive https:// doi.org/ 10. 5281/ zenodo.80852 66.
We emphasize that the provided example meshes are solely used to illustrate the life x -ep features, as users have the flexibility to input any (idealized or patient-specific) meshes into life x -ep.The cardiac tissue slab, idealized left atrial, and left ventricular meshes are characterized by a single volumetric tag representing the entire myocardium, Fig. 4(a-c).On the other hand, the ventricular slab, realistic left atrial, and left ventricular meshes have multiple distinct volumetric tags, Fig. 4(d-i).The ventricular slab is divided into subendocardial, myocardial, and sub-epicardial layers, Fig. 4(d, g).The realistic left atrium and left ventricle include regions with scars, grey zones, and fibrosis, Fig. 4(e, f, h, i).The former meshes are used for single volume simulations, while the latter for multi-volume simulations.In all cases, volumetric tags must be defined during the mesh generation process [80].
The geometrical models for the idealized cardiac meshes were created using the builtin CAD engine of gmsh, 35 an open-source 3D FE mesh generator.The gmsh scripts used to generate these meshes are included in the life x -ep release archive (https:// doi.org/ 10. 5281/ zenodo.80852 66).For detailed information on the syntax of the scripts, we refer to the online documentation of gmsh.Tetrahedral mesh generation for the slab tissue, idealized left atrium, and ventricular slab are also performed using gmsh.On the other hand, the idealized and realistic left ventricular hexahedral meshes were instead generated using the cubit 36 mesh generation software toolkit.Finally, the realistic left atrial mesh was perfomed using the semi-automatic meshing tools developed in [80], based on the Vascular Modelling Toolkit (vmtk) software [81].
The realistic left atrium and left ventricle, also containing the scar and fibrotic regions, have been produced starting from the openly available meshes used in [82] (for the left atrium 37 ) and in [83] (for the left ventricle 38 ).For the latter, we used the (vmtk) software [81], in conjunction with the cubit mesh generator.
Mesh files can be specified in the life x -ep parameter file for electrophysiology sim- ulation within the section named Mesh and space discretization, by setting the proper element type (tetrahedra or hexahedra), the polynomial order to be used in the FE discretization and an appropriate mesh rescaling factor (e.g., if the coordinates in the input mesh file are given in millimeters, Scaling factor = 1e-3, since life x -ep internally handles all physical quantities in the International System of Units.
Regarding the prescription of the myocardial fiber architecture, an essential building block for cardiac electrophysiology simulations, life x -ep provides two options.Users can either import the myofibers from a file or generate them online using LDRBMs [17] by incorporating the life x -fiber release package, which was recently published in [15].
To generate the myocardial fibers using LDRBMs, users can select the appropriate Geometry type within the parameter file under the section named Fiber generation.Each Geometry type corresponds to the specific LDRBM applicable for different geometries, such as (ventricular and spherical) slabs, (cut at base and complete) left ventricular and left atrial geometries.The parameters related to the fiber generation for a particular Geometry Type are located within a subsection with the same name.We refer to [15] for a more comprehensive guide on configuring LDRBMs, and to [17] for an in-depth mathematical description.
As an alternative, users can choose to import myofibers from a file by setting Geometry type = Import from file.The myofibers are imported from a VTU file format with unstructured grid data.In this file, the three fiber directions ( f 0 , s 0 , and n 0 representing the fiber, sheet, and sheet-normal directions, respectively) must be normalized and embedded as point-data arrays.Additionally, users need to set an appropriate geometry rescaling factor (e.g., if the coordinates in the input fiber file are given in millimeters, Scaling factor = 1e-3).This ensures that the imported myofibers align properly with the geometry of the cardiac electrophysiology model.
The myocardial fibers, generated using LDRBM, for the electrophysiology simulations presented hereafter are illustrated in Fig. 5.

Simulations of normal cardiac electrophysiology
We present normal electrophysiology simulations applied to the set of idealized geometries, namely a rectangular slab of cardiac tissue, an idealized left atrium, an idealized left ventricle and a layered ventricular slab.

Slab benchmark
To perform software verification, we consider the N-version slab benchmark proposed in [84].This benchmark involves a rectangular slab domain of size (3 × 7 × 20) × 10 −3 m, as depicted in Fig. 4(a).The fiber directions within the domain are oriented along the long axis (0.02 m), as shown in Fig. 5(a).For comprehensive modeling and geometrical information regarding the benchmark definition, we refer to the original paper [84].
The simulation can be run by first generating the parameter file using then configuring the simulation by editing the lifex_electrophysiology.prm default parameter file, and finally running the simulation by typing The same can be obtained by directly using the already prepared parameter file par-ams_slab.prmuploaded in the release archive: To match with [84], precise specifications are provided for the ionic model employed (in this test case, the TTP06 model [44]), including initial conditions and conductivity values along the myofiber directions (longitudinal, transversal, and normal).
The applied stimulus current, delivered to a volume of (1.5 × 1.5 × 1.5) × 10 −3 m, situ- ated at one corner of the slab, is prescribed in the subsection Applied current.The stimulus has a duration of 2 × 10 −3 s and an amplitude of 35.714 V/s.In Fig. 6(a) and Fig. 7(a), we display a snapshot of the transmembrane potential and the total activation time computed as output of the numerical simulation, respectively.The computation of the activation time is evaluated in a given point in the cardiac muscle as the time when the transmembrane potential derivative ∂u ∂t reaches its maximum value.This can be enabled in the parameter file under the subsection named Activation time.
The problem was solved using both tetrahedral and hexahedral conforming meshes and the BDF2 scheme, using eight combinations of spatial resolutions ( dx = [0.5, 0.2, 0.1, 0.05] × 10 −3 m) and time steps ( t = [0.05,0.01, 0.005, 0.001] × 10 −3 s), as reported in Fig. 8.The results showcase the activation time at points along the diagonal line of the slab, as also shown in Fig. 9(a).For both the tetrahedral and hexahedral meshes, the numerical solutions converge towards the finer space-time discretization ( t = 0.001 × 10 −3 s, dx = 0.05 × 10 −3 m), yielding a latest activation time of 41.8 ms and 42.0 ms, respectively.These findings align with the values reported in the original N-version benchmark paper [84].To further evaluate the life x -ep results in comparison to the other eleven codes participating in the N-version benchmark [84], we report in Fig. 9 the activation time (for t = 0.005 × 10 −3 s, dx = 0.1 × 10 −3 m) at specific points along the slab diagonal (namely P1-P8-P9, see Fig. 9(a)) for all the codes, including life x -ep.The activa- tion times at points P1-P8-P9 for all the other eleven codes are available in the electronic supplementary material of [84].As shown in Fig. 9(b), the life x -ep results fall within the range of activation time values computed by the other codes.Moreover, at such refinement level ( dt = 0.005 × 10 −3 s, dx = 0.1 × 10 −3 m), the life x - ep numerical solutions (in both tetrahedral and hexahedral meshes) are in excellent agreement with the reference converged activation time values of 42-43 ms computed at point P8, where the accumulation of errors tends to occur as the wave propagates across the cuboid.

Idealized left atrium
We simulate the propagation of the electrophysiology wavefront in an idealized left atrial geometry (see Fig. 4(b)) using the APf ionic model [42].To run the simulation, the user can modify the generic parameter file lifex_electrophysiology.prm or use the predefined file params_ideal_la.prmavailable in the release archive.
We employ a second-order BDF temporal scheme with a time step t = 10 −4 s and a final time T = 0.15 s.Furthermore, we initialize the ionic model by running a 1000-cycle long single-cell simulation applying a stimulus with period 0.8 s.
Fig. 9 a Activation times evaluted along the cuboid diagonal and in the points P1, P9 and P8 in the N-version benchmark problem [84].b Comparison of the life x -ep numerical solutions (with Red line=Hexahedral mesh and Blue line=Tetrahedral mesh) with respect to the other codes partecipating to the benchmark problem [84] The simulation is initiated using a single spherical impulse with a radius of r = 3 × 10 −3 m.The parameters related to the impulse site, amplitude, duration, and initial time are specified in the Applied current/Spherical subsection of the parameter file.Notice that, when utilizing the APf model, it is essential to properly rescale the impulse amplitude and also conductivity values embedded within the diffusion tensor (3) in accordance with the model's formulation [42].Finally, the myofiber architecture is prescribed using the atrial LDRBM described in [17], as depicted in Fig. 5(b), by specifying Geometry type = Left atrium in the subsection Fiber generation.It is important to remark that in the APf model formulation, the transmembrane potential u is a dimensionless variable ranging from 0 to 1.However, for visualization purposes, the actual transmembrane potential u mV in millivolts is obtained by postprocessing the numerical results using the formula u mV = (100u − 80)[mV] [42].This conversion is applied to the data shown in Fig. 6(b).

Idealized left ventricle
We simulate the electrophysiology wavefront propagation in an idealized left ventricular hexahedral mesh, shown in Fig. 4(c), using the Bueno-Orovio (BO) ionic model [43].The simulation can be run either by modifying the generic parameter file lifex_electrophysiology.prm or by utilizing the dedicated parameter file params_ideal_lv.prm.This test case can be executed at various levels of hierarchical grid refinements by modifying the parameter value Number of refinements in the Mesh and space discretization subsection.We remark that this is possible for all life x simula- tions conducted with hexahedral meshes, but is not available for tetrahedral meshes, due to lack of support in deal.II.
For this simulation, we use a second-order BDF temporal scheme with a time step of t = 5 × 10 −5 s and a final time of T = 0.15 s.The ionic model is initialized using 1000- cycle long single-cell simulations with a cardiac period of 0.8 s.We employ a pacing protocol where three ventricular endocardial areas are activated using Gaussian impulses [85].The characteristics of these impulses, such as amplitude, duration, and initial time, can be set in the file parameter subsection Applied current/Gaussian.
The myofiber architecture is prescribed using the Rossi-Lassila (RL) ventricular LDRBM [17], see Fig.The simulation results are reported in Fig. 6(c) and Fig. 7(c), where a snapshot of the transmembrane potential and the total activation time are displayed, respectively.Note that the transmembrane potential u mV shown in Fig. 6(c) is obtained by postprocessing the numerical results using the formula u mV = (85.7u− 84) [mV ] [43].

Ventricular slab
This test case serves as an explanatory example for the multi-volume simulation framework of life x -ep.We simulate a portion of an idealized left ventricular geometry, also referred to as ventricular slab [15], where the computational domain is divided into three volumetric regions: Sub-endocardial, Myocardial, and Sub-epicardial layers, as shown in Fig. 4(d).
To perform multi-volume simulations, the parameter template is created by specifying the volume labels on the command line as follows: Doing so, the file lifex_electrophysiology.prm will contain three subsections named Sub Endocardium, Myocardium and Sub Epiucardium, where the volumetric tags (Material IDs) and all volume-specific parameters are prescribed.We use the TTP06 ionic model with a different Cell type for each layer (Endocardium, Myocardium and Epicardium) [44].
The simulation can then be run by providing on the command line the same volume labels used when generating the parameter file: A second-order BDF temporal scheme is employed in this simulation, with a time step of t = 5 × 10 −5 s and a final time of T = 0.12 s.The ionic model is initialized in each volume by conducting 1000-cycle long single-cell simulations with a cardiac period of 0.8 s.One single-cell simulation is run for every volume, so that the initialization is consistent with the spatially-varying parameters.The simulation onset is obtained with a single spherical impulse, of radius r = 2.5 × 10 −3 m.The myofiber architecture is prescribed in the Fiber generation subsection by utilizing the Slab LDRBM, as shown in Fig. 5(d).The fiber orientations are set to exhibit a linear transmural variation from 60 • to −60 • , passing from the endocardial to the epicardial surface.
The simulation results are shown in Fig. 6(d) and Fig. 7(d), which display a snapshot of the transmembrane potential and the total activation time, respectively.

Simulations of pathological cardiac electrophysiology
In the following section, we present the pathological electrophysiology simulations to realistic left atrial and left geometries.

Realistic pathological left atrium
We perform a simulation of reentrant drivers typical of AF in a realistic left atrial tetrahedral mesh [45,82], as shown in Fig. 4(e, h).
In this pathological scenario, four volumes have been introduced, named Healthy, Fibrosis Mild, Fibrosis and Scar.We use the CRN ionic model adopting different ionic conductances and conductivity values to model the pathophysiological behaviour in the fibrotic regions, that are labelled according to their IIR [82].Specifically, in the Fibrosis Mild volumetric region, we changed, with respect to the default CRN values, the parameters related to the transient outward current conductance gto, the L-type intracellular Ca 2+ current conductance gCaL, and the delayed rectifier current represented by gKur_fix and gKur_var, to model the effects of chronic AF.In the Fibrosis region, we set specific current conductance values for gK1, gCaL and gNa, which stand for the inward rectifier, the L-type calcium, potassium and sodium current conductances, respectively.Conductances and conductivities were adjusted according to [48].An additional region, labeled as Healthy, models healthy tissue, without modification to ionic parameters.Finally, we consider a region labelled as Scar, which we model as non-conductive, thus excluding that region from the domain of the monodomain and ionic models and creating a block in the conduction (see Section Action potential propagation models).
For the spatial discretization we employ second-order FEs by setting FE space degree = 2.We use a first-order BDF temporal scheme with a time step t = 5 × 10 −5 s and a final time T = 1.5 s.We initialize the ionic model by running a 24-cycle long single-cell simulation using an impulse period of 0.5 s [48].We use a pacing protocol with a sequence of four multiple spherical impulses delivered every 220 × 10 −3 s.This timing is set under the parameter set Impulse initial times in the Applied current/Spherical subsection of the parameter file.Finally, the myofiber architecture is prescribed using the atrial LDRBM [17], as depicted in Fig. 5(e), by specifying Geometry type = Left atrium in the subsection Fiber generation.
The simulation can be run using The simulation results are presented in Fig. 6(e) and Fig. 7(e), which display a snapshot of the transmembrane potential and the total activation time, respectively.In this video 39  (online version) we show the evolution of the transmembrane potential.

Realistic pathological left ventricle
We simulate a Ventricular Tachycardia (VT) macro-reentrant circuit in a realistic left ventricle [83], as shown in Fig. 4(f, i).
To account for the presence of grey zone and scar, we consider three volumes: Healthy, Border Zone and Scar.Different ionic conductances in the TTP06 ionic model and conductivity values in the monodomain equation are used to characterize the electrophysiology properties within each region.In particular, in the Border Zone, we reduce the conductances of the peak sodium, L-type calcium, and potassium currents GNa, GCaL, Gkr and Gks_myo, according to [49].The region has a normal configuration, while the Scar region is modeled as non-conductive (setting Disable conduction = true).
39 https:// polim i365-my.share point.com/: v:/g/ perso nal/ 10594 253_ polimi_ it/ EcRfI yP2ya 5EsD6 IixPi_ 4ABCa ZfJVi Twwge DpKBz Csquw?e= VUMMP9 We use second-order FEs for the spatial discretization, with one level of hierarchical grid refinement (Number of refinements = 1), and a second-order BDF temporal scheme with a time step of t = 5 × 10 −5 s and a final time of T = 2.5 s.The ionic model is initialized using 1000-cycle long single-cell simulations with a cardiac period of 0.8 s.We employ a pacing protocol with a sequence of multiple spherical impulses delivered in a specific ventricular endocardial area.The myofiber architecture is prescribed using the RL ventricular algorithm [17], see Fig.The simulation can be run using The simulation results are presented in Fig. 6(f ) and Fig. 7(f ), which display a snapshot of the transmembrane potential and the total activation time, respectively.In this video 40  (online version) we show the evolution of the transmembrane potential.

Strong scalability test
We consider the slab benchmark of [84], and discretize the domain with a structured hexahedral mesh of 47185920 elements and 47744577 nodes, corresponding to a mesh size of approximately dx = 3.6 × 10 −5 m .We set t = 1 × 10 −4 s and T = 5 × 10 −2 s , and we use linear FEs and the second-order BDF scheme for time discretization.
We run a strong scalability test, varying the number of parallel processes used in the computation and measuring the wall time necessary for the solution of the ionic models, the assembly of the monodomain system and the solution of the linear system.The test was run on the CINECA GALILEO100 supercomputer.
The wall times, plotted in Fig. 10 for the different steps of the solver, scales linearly up to over 1000 cores, confirming the results of [13] on the scalability properties of the life x core components.We report in Table 1 a breakdown of the computa- tional cost of the different sections for the simulation with 960 cores.Most of the computational time is spent for the solution of the linear system (6).The matrix of the linear system is assembled only once (since it is the same at every time step), and the right-hand-side is efficiently recomputed at every time step by means of matrixvector products, resulting in a very small computational cost for the system assembly phase.Moreover, since the ionic model is solved independently at each degree of freedom, the associated computational cost scales almost perfectly, and becomes very small if a sufficiently large number of parallel processes is employed.

Memory occupation
We consider the same setting used for the scalability test, with three differently refined meshes of 772497, 6038305 and 47744577 nodes (corresponding to a mesh size of approximately 1.5 × 10 −4 m , 7.3 × 10 −5 m and 3.6 × 10 −5 m , respectively).We measure the peak memory occupation of the simulation, varying the number of parallel processes.Results are reported in Fig. 11.We observe that the total memory occupation increases as the number of processes increases.This is expected due to distributed memory parallelization, that requires some duplicated information across processes.However, as shown in Fig. 11 (right), this overhead grows more slowly than the number of processes itself.Therefore, increasing the number of

Output and visualization
Two different types of output file formats are available in the life x -ep release: HDF5 and csv.Both of them can be enabled and configured in the Output subsections: The HDF5 output is available in the following subsections: • Electrophysiology/Output, • Electrophysiology/Activation time, • Fiber generation/Output.
This generates an XDMF file named output_filename.xdmf(which links to a corresponding HDF5 output file output_filename.h5).These files can be visualized using ParaView, 41 an open-source multi-platform data analysis and visualization application, see e.g., Fig. 5, 6 and 7.The HDF5 format ensures that the output can be easily postprocessed, not only for visualization purposes but also as input for more advanced computational pipelines.The Comma-Separated Values (csv) format consists of delimited text files where values are separated by commas, with each line representing a specific data record.The csv files can be found in different subsections: • Ionic model parameters/Output, • Ionic model parameters/0D Output, • Electrophysiology/Output.These csv output files can be conveniently used to plot electrophysiology variable (min, max and pointwise) values over time in the computed numerical simulation.

Conclusions
In this work, we introduced life x -ep, a robust and advanced software specifically designed for simulating the electrophysiology activity of the cardiac muscle.With the goal of addressing the computational challenges associated with cardiac simulations, life x -ep provides efficient numerical methods while maintaining precision and accuracy.life x -ep incorporates a numerical solver for the monodomain equa- tion coupled with both phenomenological and second-generation ionic models, namely Aliev-Panfilov [42], Bueno-Orovio [43], TTP06 [44], and CRN [45].These models are discretized in time using the BDF scheme and in space using the FE method of orders 1 and 2 on tetrahedral meshes, and of arbitrary degree on hexahedral meshes, thus providing a comprehensive framework for modeling the electrical activity of the heart under both normal and pathological conditions.
Leveraging the capabilities of life x , life x -ep provides users with a user-friendly and flexible interface, facilitated by self-documenting parameter files for easy simulation setup.For enhanced accessibility, life x -ep is distributed in an AppImage binary format, rendering it universally compatible with any recent x86-64 Linux system.Researchers from diverse backgrounds, such as medicine and bio-engineering, can readily access and utilize life x -ep for in-silico simulations.The underlying principles and structure of life x -ep can be readily understood thanks to the comprehensive techni- cal and mathematical documentation.
As unique and distinctive features, life x -ep provides two options for prescribing myocardial fibers.Users can either import them from a file or generate them online by exploiting the LDRBMs presented in [17] and implemented in the previous release life x -fiber [15].Moreover, it supports spatial heterogeneity in the choice of both models and physical coefficients, easily configurable through a convenient parameter file, without the need to access and modify the source code.life x -ep benefits from its high-performance computing capabilities, achieving ideal parallel speedup on thousands of cores.The accuracy and reliability of life x -ep have been verified through a benchmark for computational electrophysiology, ensuring the validity of its results.Furthermore, a range of idealized and realistic cardiac simulations in both normal and pathological settings highlights its capabilities and versatility in capturing complex cardiac dynamics and its potential for patient-specific studies.life x - ep offers the capability to facilitate the simulation of pathological scenarios, allowing the creation of scars, grey zones, and an arbitrary number of conduction "levels".This potential impact in the study of pathologies can prove to be highly valuable.
In conclusion, life x -ep provides to the scientific community a comprehensive, high- performance, and user-friendly software for conducting in-silico cardiac electrophysiology simulations.
In the future, efforts will be made to further improving the accuracy and efficiency of life x -ep.One possible approach is the adoption of lookup tables instead of repeat- edly evaluating the expensive functions in the ionic current term, which could provide a speed-up, although with a potential impact on accuracy [86].Inexact solvers, such as those based on domain decomposition methods, show promising optimality and scalability properties [87].Recent studies have also highlighted the advantages of employing high-order discretization schemes for accurately capturing the intricate electrical wavefront propagation observed in cardiac electrophysiology [19].These schemes not only offer improved accuracy but also enable the implementation of efficient matrix-free solvers that require minimal memory usage [19].This opens up the possibility of leveraging GPU architectures for accelerated computations [88].Additionally, the use of higher-order or adaptive time stepping schemes and hp-adaptive FE holds the potential to achieve greater accuracy in simulations while optimizing computational efficiency [89].These developments aim to further enhance the capabilities of life x -ep and expand its range of features.

Fig. 2
Fig. 2 Multiscale cardiac electrophysiology model.From the largest to the smallest spatial scale: a organ; b tissue; c cell; d membrane

Fig. 3 a
Fig. 3 a Computational domain ⊂ R 3 in an idealized slab of ventricular myocardium, partitioned in three subdomains ( = 1 ∪ 2 ∪ 3 ), where Ŵ 12 and Ŵ 23 denote the interfaces among the corresponding subdomains; (b) Representation of the orthonormal triplet for the myofiber orientations f 0 (fiber, in red), s 0 (sheet, in blu) and n 0 (sheet-normal, in green) in an idealized slab of ventricular tissue ion = I ion /C m , I app = I app /C m and D = D/(χ m C m ) .These equivalences are helpful when comparing different literature sources.

Fig. 4
Fig. 4 Domains and meshes used in the numerical examples.Domains in the top row a-c are composed of a single subdomain, while for domains in the middle row d-f, zoomed on the bottom g-i, colors are used to differentiate the subdomains

Fig. 5
Fig. 5 Fiber field computed using Laplace-Dirichlet Rule-Based Methods [15, 17] visualized as streamlines: a Slab tissue; b Idealized left atrium; c Idealized left ventricle; d Ventricular slab; e Realistic left atrium; f Realistic left ventricle.The Laplace solution φ is the transmural function where φ = 0 on the endocardium and φ = 1 on the epicardium

Fig. 6
Fig. 6 Snapshots of the transmembrane potential for all test cases: a Slab tissue; b Idealized left atrium; c Idealized left ventricle; d Ventricular slab; e Realistic left atrium; f Realistic left ventricle

Fig. 7 Fig. 8
Fig. 7 Activation maps computed for all the test cases: a Slab tissue; b Idealized left atrium; c Idealized left ventricle; d Ventricular slab; e Realistic left atrium; f Realistic left ventricle

5
(c), by specifying Geometry type = Left ventricle and setting Algorithm type = RL in the subsection Fiber generation.

5
(f ), by specifying Geometry type = Left ventricle and setting Algorithm type = RL in the subsection Fiber generation.

Table 1
Computational time (left) and parallel speedup (right) against the number of cores for the strong scalability test.Dashed lines indicate the ideal linear scaling Summary of the computational costs for the strong scalability test, using 960 parallel cores.For each section, we report the total wall time, the wall time for each time step and the wall time relative to the total.Sections are sorted in descending order of cost Total memory occupation (left) and memory occupation per core with three differently refined meshes computing cores, which generally implies increasing the number of available RAM, does not lead to out-of-memory issues.