Search

Deep Neural Networks are widely used for solving complex problems in several scientific areas, such as speech recognition, machine translation, image analysis. The strategies employed to investigate their theoretical properties mainly rely on Euclidean geometry, but in the last years new approaches based on Riemannian geometry have been developed. Motivated by some open problems, we study a particular sequence of maps between manifolds, with the last manifold of the sequence equipped with a Riemannian metric. We investigate the structures induced through pullbacks on the other manifolds of the sequence and on some related quotients. In particular, we show that the pullbacks of the final Riemannian metric to any manifolds of the sequence is a degenerate Riemannian metric inducing a structure of pseudometric space. We prove that the Kolmogorov quotient of this pseudometric space yields a smooth manifold, which is the base space of a particular vertical bundle. We investigate the theoretical properties of the maps of such sequence, eventually we focus on the case of maps between manifolds implementing neural networks of practical interest and we present some applications of the geometric framework we introduced in the first part of the paper.

In this paper we consider the spatial semi-discretization of conservative PDEs. Such finite dimensional approximations of infinite dimensional dynamical systems can be described as flows in suitable matrix spaces, which in turn leads to the need to solve polynomial matrix equations, a classical and important topic both in theoretical and in applied mathematics. Solving numerically these equations is challenging due to the presence of several conservation laws which our finite models incorporate and which must be retained while integrating the equations of motion. In the last thirty years, the theory of geometric integration has provided a variety of techniques to tackle this problem. These numerical methods require solving both direct and inverse problems in matrix spaces. We present three algorithms to solve a cubic matrix equation arising in the geometric integration of isospectral flows. This type of ODEs includes finite models of ideal hydrodynamics, plasma dynamics, and spin particles, which we use as test problems for our algorithms., 14 pages, 1 figure

This work deals with optimal control problems as a strategy to drive bifurcating solution of nonlinear parametrized partial differential equations towards a desired branch. Indeed, for these governing equations, multiple solution configurations can arise from the same parametric instance. We thus aim at describing how optimal control allows to change the solution profile and the stability of state solution branches. First of all, a general framework for nonlinear optimal control problem is presented in order to reconstruct each branch of optimal solutions, discussing in detail the stability properties of the obtained controlled solutions. Then, we apply the proposed framework to several optimal control problems governed by bifurcating Navier–Stokes equations in a sudden-expansion channel, describing the qualitative and quantitative effect of the control over a pitchfork bifurcation, and commenting in detail the stability eigenvalue analysis of the controlled state. Finally, we propose reduced order modeling as a tool to efficiently and reliably solve parametric stability analysis of such optimal control systems, which can be challenging to perform with standard discretization techniques such as Finite Element Method.

Background Modeling the whole cardiac function involves the solution of several complex multi-physics and multi-scale models that are highly computationally demanding, which call for simpler yet accurate, high-performance computational tools. Despite the efforts made by several research groups, no software for whole-heart fully-coupled cardiac simulations in the scientific community has reached full maturity yet. Results In this work we present $$\texttt {life}^{\texttt {x}}$$ life x -fiber, an innovative tool for the generation of myocardial fibers based on Laplace-Dirichlet Rule-Based Methods, which are the essential building blocks for modeling the electrophysiological, mechanical and electromechanical cardiac function, from single-chamber to whole-heart simulations. $$\texttt {life}^{\texttt {x}}$$ life x -fiber is the first publicly released module for cardiac simulations based on $$\texttt {life}^{\texttt {x}}$$ life x , an open-source, high-performance Finite Element solver for multi-physics, multi-scale and multi-domain problems developed in the framework of the iHEART project, which aims at making in silico experiments easily reproducible and accessible to a wide community of users, including those with a background in medicine or bio-engineering. Conclusions The tool presented in this document is intended to provide the scientific community with a computational tool that incorporates general state of the art models and solvers for simulating the cardiac function within a high-performance framework that exposes a user- and developer-friendly interface. This report comes with an extensive technical and mathematical documentation to welcome new users to the core structure of $$\texttt {life}^{\texttt {x}}$$ life x -fiber and to provide them with a possible approach to include the generated cardiac fibers into more sophisticated computational pipelines. In the near future, more modules will be successively published either as pre-compiled binaries for systems or as open source software.

The aim of electromagnetic (EM) sounding methods in geophysics is to obtain information about the subsurface of the earth by recorded measurements taken at the surface. In particular, the goal is to determine variations in the electrical conductivity of the earth with depth by employing an inversion procedure. In this work we focus on one technique, that consists of placing a magnetic dipole above the surface, composed of a transmitter coil and different couples of adjacent receiver coils. The receiver couples are placed at different distances (offsets) from the transmitter coil. In this setting, the electromagnetic induction effect, encoded in the first-order linear Maxwell’s differential equations, produce eddy alterning currents in the soil which induce response electromagnetic fields, that can be used to determine the conductivity profile of the ground by meaning of an inversion algorithm. A typical inversion strategy consists in an iterative procedure involving the computation of the EM response of a layered model (forward modelling) and the solution of the inverse problem. Then, the algorithm attempts to minimize the mismatch between the measured data and the predicted data, by updating the model parameters at each iteration. By assuming that the local subsurface structures are composed by horizontal and homogeneous layers, general integral solutions of Maxwell equations (i.e., the EM fields) for vertical and horizontal magnetic dipoles, can be derived and represented as Hankel transforms, which contain the subsurface model parameters, i.e., the conductivity and the thickness of each layer. By a mathematical point of view, in general, these Hankel transforms are not analytically computable and therefore it is necessary to employ a numerical scheme. Anyway, the slowly decay of the oscillations determined by the Bessel function makes the problem very difficult to handle, because traditional quadrature rules typically fail to converge. In this work we consider two different approaches. The first one is based on the decomposition of the integrand function in a first function for which the corresponding Hankel transform is known exactly, and an oscillating function decays exponentially. For realistic sets of parameters, the oscillations are quite rapidly damped, and the corresponding integral can be accurately computed by a classical quadrature rule on finite intervals. The second approach consists in the application of a Gaussian quadrature formula. We develop a Gaussian rule for weight functions involving fractional powers, exponentials and Bessel functions of the first kind. Moreover, we derive an analytical approximation of these integrals that has a general validity and allows to overcome the limits of common methods based on the modelling of apparent conductivity in the low induction number (LIN) approximation. Having at disposal a reliable method for evaluating the Hankel transforms, by assuming as forward model a homogeneous layered earth, here we also consider the inverse problem of computing the model parameters (i.e., conductivity and thickness of the layers) from a set of measured field values at different offsets. We focus on the specific case of the DUALEM system. We employ two optimization algorithms. The first one is based on the BFGS line-search method and, in order to reduce as much as possible the number of integral evaluations, the analytic approximation of these integrals is used to have a first estimate of the solution. For the second approach we employ the damped Gauss-Newton method. To avoid the dependence on the initial guess of the iterative procedure, we consider a set of different initial models, and we use each one to solve the optimization problem. The numerical experiments, carried out for the study of river-levees integrity, are obtained by employing a virtual machine equipped with the NVIDIA A100 Tensor Core GPU. The aim of electromagnetic (EM) sounding methods in geophysics is to obtain information about the subsurface of the earth by recorded measurements taken at the surface. In particular, the goal is to determine variations in the electrical conductivity of the earth with depth by employing an inversion procedure. In this work we focus on one technique, that consists of placing a magnetic dipole above the surface, composed of a transmitter coil and different couples of adjacent receiver coils. The receiver couples are placed at different distances (offsets) from the transmitter coil. In this setting, the electromagnetic induction effect, encoded in the first-order linear Maxwell’s differential equations, produce eddy alterning currents in the soil which induce response electromagnetic fields, that can be used to determine the conductivity profile of the ground by meaning of an inversion algorithm. A typical inversion strategy consists in an iterative procedure involving the computation of the EM response of a layered model (forward modelling) and the solution of the inverse problem. Then, the algorithm attempts to minimize the mismatch between the measured data and the predicted data, by updating the model parameters at each iteration. By assuming that the local subsurface structures are composed by horizontal and homogeneous layers, general integral solutions of Maxwell equations (i.e., the EM fields) for vertical and horizontal magnetic dipoles, can be derived and represented as Hankel transforms, which contain the subsurface model parameters, i.e., the conductivity and the thickness of each layer. By a mathematical point of view, in general, these Hankel transforms are not analytically computable and therefore it is necessary to employ a numerical scheme. Anyway, the slowly decay of the oscillations determined by the Bessel function makes the problem very difficult to handle, because traditional quadrature rules typically fail to converge. In this work we consider two different approaches. The first one is based on the decomposition of the integrand function in a first function for which the corresponding Hankel transform is known exactly, and an oscillating function decays exponentially. For realistic sets of parameters, the oscillations are quite rapidly damped, and the corresponding integral can be accurately computed by a classical quadrature rule on finite intervals. The second approach consists in the application of a Gaussian quadrature formula. We develop a Gaussian rule for weight functions involving fractional powers, exponentials and Bessel functions of the first kind. Moreover, we derive an analytical approximation of these integrals that has a general validity and allows to overcome the limits of common methods based on the modelling of apparent conductivity in the low induction number (LIN) approximation. Having at disposal a reliable method for evaluating the Hankel transforms, by assuming as forward model a homogeneous layered earth, here we also consider the inverse problem of computing the model parameters (i.e., conductivity and thickness of the layers) from a set of measured field values at different offsets. We focus on the specific case of the DUALEM system. We employ two optimization algorithms. The first one is based on the BFGS line-search method and, in order to reduce as much as possible the number of integral evaluations, the analytic approximation of these integrals is used to have a first estimate of the solution. For the second approach we employ the damped Gauss-Newton method. To avoid the dependence on the initial guess of the iterative procedure, we consider a set of different initial models, and we use each one to solve the optimization problem. The numerical experiments, carried out for the study of river-levees integrity, are obtained by employing a virtual machine equipped with the NVIDIA A100 Tensor Core GPU.

The synthesis of suboptimal feedback laws for controlling nonlinear dynamics arising from semi-discretized PDEs is studied. An approach based on the State-dependent Riccati Equation (SDRE) is presented for 2 and ∞ control problems. Depending on the nonlinearity and the dimension of the resulting problem, offline, online, and hybrid offline-online alternatives to the SDRE synthesis are proposed. The hybrid offline-online SDRE method reduces to the sequential solution of Lyapunov equations, effectively enabling the computation of suboptimal feedback controls for two-dimensional PDEs. Numerical tests for the Sine-Gordon, degenerate Zeldovich, and viscous Burgers’ PDEs are presented, providing a thorough experimental assessment of the proposed methodology.

We propose a test for multisample comparison studies that can be applied without strict assumptions, especially when the underlying population distributions are far from normal. The new test can detect differences not only in location or scale but also in shape parameters among parent population distributions. We are motivated by numerous medical studies, where the variables are not normally distributed and may present in the various groups more complex differences than simple differences in a particular aspect of underlying distributions, such as location or scale. In these situations, traditional ANOVA and Kruskal-Wallis tests are unreliable since the underlying assumptions are not valid. The proposed procedure also allows the researcher to determine which aspects are more responsible for a significant result. This is an important practical advantage over procedures that test for general differences among the distribution functions but cannot identify which aspects lead to significant results. The asymptotic distribution of the test statistic is analyzed along with its small sample behavior against several competing tests. The practical advantages of the proposed procedure are illustrated with a multisample comparison study of a biomarker for liver damage in patients with hepatitis C.

Numerical stabilization is often used to eliminate (alleviate) the spurious oscillations generally produced by full order models (FOMs) in under-resolved or marginally-resolved simulations of convection-dominated flows. In this article, we investigate the role of numerical stabilization in reduced order models (ROMs) of marginally-resolved, convection-dominated incompressible flows. Specifically, we investigate the FOM-ROM consistency, that is, whether the numerical stabilization is beneficial both at the FOM and the ROM level. As a numerical stabilization strategy, we focus on the evolve-filter-relax (EFR) regularization algorithm, which centers around spatial filtering. To investigate the FOM-ROM consistency, we consider two ROM strategies: (i) the EFR-noEFR, in which the EFR stabilization is used at the FOM level, but not at the ROM level; and (ii) the EFR-EFR, in which the EFR stabilization is used both at the FOM and at the ROM level. We compare the EFR-noEFR with the EFR-EFR in the numerical simulation of a 2D incompressible flow past a circular cylinder in the convection-dominated, marginally-resolved regime. We also perform model reduction with respect to both time and Reynolds number. Our numerical investigation shows that the EFR-EFR is more accurate than the EFR-noEFR, which suggests that FOM-ROM consistency is beneficial in convection-dominated, marginally-resolved flows.

We propose a mathematical and numerical model for the simulation of the heart function that couples cardiac electrophysiology, active and passive mechanics and hemodynamics, and includes reduced models for cardiac valves and the circulatory system. Our model accounts for the major feedback effects among the different processes that characterize the heart function, including electro-mechanical and mechano-electrical feedback as well as force-strain and force-velocity relationships. Moreover, it provides a three-dimensional representation of both the cardiac muscle and the hemodynamics, coupled in a fluid-structure interaction (FSI) model. By leveraging the multiphysics nature of the problem, we discretize it in time with a segregated electrophysiology-force generation-FSI approach, allowing for efficiency and flexibility in the numerical solution. We employ a monolithic approach for the numerical discretization of the FSI problem. We use finite elements for the spatial discretization of partial differential equations. We carry out a numerical simulation on a realistic human left heart model, obtaining results that are qualitatively and quantitatively in agreement with physiological ranges and medical images. This article is protected by copyright. All rights reserved.

While ventricular electromechanics is extensively studied in both physiological and pathological conditions, four-chamber heart models have only been addressed recently; most of these works however neglect atrial contraction. Indeed, as atria are characterized by a complex anatomy and a physiology that is strongly influenced by the ventricular function, developing computational models able to capture the physiological atrial function and atrioventricular interaction is very challenging. In this paper, we propose a biophysically detailed electromechanical model of the whole human heart that considers both atrial and ventricular contraction. Our model includes: (i) an anatomically accurate whole-heart geometry; (ii) a comprehensive myocardial fiber architecture; (iii) a biophysically detailed microscale model for the active force generation; (iv) a 0D closed-loop model of the circulatory system, fully-coupled with the mechanical model of the heart; (v) the fundamental interactions among the different core models, such as the mechano-electric feedback or the fibers-stretch and fibers-stretch-rate feedbacks; (vi) specific constitutive laws and model parameters for each cardiac region. Concerning the numerical discretization, we propose an efficient segregated-intergrid-staggered scheme that includes a computationally efficient strategy to handle the non-conductive regions. We also propose extending recent stabilization techniques - regarding the circulation and the fibers-stretch-rate feedback - to the whole heart, demonstrating their cruciality for obtaining a stable formulation in a four-chamber scenario. We are able to reproduce the healthy cardiac function for all the heart chambers, in terms of pressure-volume loops, time evolution of pressures, volumes and fluxes, and three-dimensional cardiac deformation, with volumetric indexes within reference ranges for cardiovascular magnetic resonance. We also show the importance of considering atrial contraction, fibers-stretch-rate feedback and the proposed stabilization techniques, by comparing the results obtained with and without these features in the model. In particular, we show that the fibers-stretch-rate feedback, often neglected due to the numerical challenges that it entails, plays a fundamental role in the regulation of the blood flux ejected by ventricles. The proposed model represents the state-of-the-art electromechanical model of the iHEART ERC project - an Integrated Heart Model for the Simulation of the Cardiac Function - and is a fundamental step toward the building of physics-based digital twins of the human heart.(c) 2023 The Author(s). Published by Elsevier B.V. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/).

17 pages; International audience; Microfluidics proved to be a key technology in various applications, allowing to reproduce large-scale laboratory settings at a more sustainable small-scale. The current effort is focused on enhancing the mixing process of different passive species at the micro-scale, where a laminar flow regime damps turbulence effects. Chaotic advection is often used to improve mixing effects also at very low Reynolds numbers. In particular, we focus on passive micromixers, where chaotic advection is mainly achieved by properly selecting the geometry of microchannels. In such a context, reduced order modeling can play a role, especially in the design of new geometries. In this chapter, we verify the reliability and the computational benefits lead by a Hierarchical Model (HiMod) reduction when modeling the transport of a passive scalar in an S-shaped microchannel. Such a geometric configuration provides an ideal setting where to apply a HiMod approximation, which exploits the presence of a leading dynamics to commute the original three-dimensional model into a system of one-dimensional coupled problems. It can be proved that HiMod reduction guarantees a very good accuracy when compared with a high-fidelity model, despite a drastic reduction in terms of number of unknowns.

Convolutional Neural Network (CNN) is one of the most important architectures in deep learning. The fundamental building block of a CNN is a trainable filter, represented as a discrete grid, used to perform convolution on discrete input data. In this work, we propose a continuous version of a trainable convolutional filter able to work also with unstructured data. This new framework allows exploring CNNs beyond discrete domains, enlarging the usage of this important learning technique for many more complex problems. Our experiments show that the continuous filter can achieve a level of accuracy comparable to the state-of-the-art discrete filter, and that it can be used in current deep learning architectures as a building block to solve problems with unstructured domains as well.

Fast Field-Cycling Nuclear Magnetic Resonance relaxometry is a non-destructive technique to investigate molecular dynamics and structure of systems having a wide range of ap- plications such as environment, biology, and food. Besides a considerable amount of liter- ature about modeling and application of such technique in specific areas, an algorithmic approach to the related parameter identification problem is still lacking. We believe that a robust algorithmic approach will allow a unified treatment of different samples in several application areas. In this paper, we model the parameters identification problem as a con- strained L 1 -regularized non-linear least squares problem. Following the approach proposed in [ Analytical Chemistry 2021 93 (24)], the non-linear least squares term imposes data con- sistency by decomposing the acquired relaxation profiles into relaxation contributions as- sociated with 1 H −1 H and 1 H −14 N dipole-dipole interactions. The data fitting and the L 1 - based regularization terms are balanced by the so-called regularization parameter. For the parameters identification, we propose an algorithm that computes, at each iteration, both the regularization parameter and the model parameters. In particular, the regularization parameter value is updated according to a Balancing Principle and the model parameters values are obtained by solving the corresponding L 1 -regularized non-linear least squares problem by means of the non-linear Gauss-Seidel method. We analyse the convergence properties of the proposed algorithm and run extensive testing on synthetic and real data. A Matlab software, implementing the presented algorithm, is available upon request to the authors.

Blood microcirculation is the site of control of tissue perfusion, blood-tissue exchange, and tissue blood volume. Despite the many irregularities, almost ubiquitously, one can recognize in microcirculation vessels a hierarchy of arterioles and venules, organized in tree-like structures, and capillary plexi, organized in net-like structures. Whilst for arterioles and venules it may be envisageable to obtain geometries needed for numerical simulations from imaging techniques, the size and numerosity of capillaries makes this task much more cumbersome. For this reason, it is interesting to study approaches to generate in silico-derived artifacts of capillary networks, even in view of machine-learning based approaches which require a large amount of samples for training. Artificial networks must correctly reflect proper metrics and topology, which in turn, will ensure with proper boundary conditions a physiological blood flux in the net and a sufficient nutrient distribution in the surrounding tissues. In this paper, we introduce the sequence of curves whose limit is the space filling Hilbert curve and we discuss its inherent properties and we obtain the backbone of the artificial capillary network from a suitable element of this sequence. The backbone represents a significant synthesis of basic metric features of the network and, in this context, its properties can be studied analytically. In this framework, the Hilbert curve is a malleable entity which allows to shape the backbone according to the physical indicators. In particular, two significant factors are shown to control the network topology and scaling: the iteration step of the Hilbert curve generation and the characteristic length of the REV, respectively. Based on the points we generate for a certain iteration step, we then obtain via spline interpolation a smoothed version of the curve, which fine–tunes the tortuosity. A volumetric construction is obtained building a tubular neighborhood of the backbone, whose metrics can be computed and tuned as well. Numerical simulations of the blood flow in the obtained geometry show the physical fields occurring in the artificial network.

In this paper we develop an adaptive algorithm for determining the optimal degree of regression in the constrained mock-Chebyshev least-squares interpolation of an analytic function to obtain quadrature formulas with high degree of exactness and accuracy from equispaced nodes. We numerically prove the effectiveness of the proposed algorithm by several examples.

In this manuscript a POD-Galerkin based Reduced Order Model for unsteady Fluid-Structure Interaction problems is presented. The model is based on a partitioned algorithm, with semi-implicit treatment of the coupling conditions. A Chorin–Temam projection scheme is applied to the incompressible Navier–Stokes problem, and a Robin coupling condition is used for the coupling between the fluid and the solid. The coupled problem is based on an Arbitrary Lagrangian Eulerian formulation, and the Proper Orthogonal Decomposition procedure is used for the generation of the reduced basis. We extend existing works on a segregated Reduced Order Model for Fluid-Structure Interaction to unsteady problems that couple an incompressible, Newtonian fluid with a linear elastic solid, in two spatial dimensions. We consider three test cases to assess the overall capabilities of the method: an unsteady, non-parametrized problem, a problem that presents a geometrical parametrization of the solid domain, and finally, a problem where a parametrization of the solid’s shear modulus is taken into account.

We extend the applicability of the popular interior penalty discontinuous Galerkin method discretizing advection-diffusion-reaction problems to meshes comprising extremely general, essentially arbitrarily-shaped element shapes. In particular, our analysis allows for curved element shapes, without the use of non-linear elemental maps. The feasibility of the method relies on the definition of a suitable choice of the discontinuity penalization, which turns out to be explicitly dependent on the particular element shape, but essentially independent on small shape variations. This is achieved upon proving extensions of classical trace and Markov-type inverse estimates to arbitrary element shapes. A further new H 1 − L 2 H^1-L_2 -type inverse estimate on essentially arbitrary element shapes enables the proof of inf-sup stability of the method in a streamline-diffusion-like norm. These inverse estimates may be of independent interest. A priori error bounds for the resulting method are given under very mild structural assumptions restricting the magnitude of the local curvature of element boundaries. Numerical experiments are also presented, indicating the practicality of the proposed approach.

The thesis is about how perturbations in the initial value $y_{0}$ or in the coefficient matrix $A$ propagate along the solutions of $n$-dimensional linear ordinary differential equations (ODE) \begin{equation*} \left\{ \begin{array}{l} y^\prime(t) =Ay(t),\ t\geq 0,\\ y(0)=y_0, \end{array} \right. \end{equation*} where $A \in \mathbb{R}^{n\times n} $ and $y_0\in \mathbb{R}^n$ and $y(t)=e^{tA}y_0$ is the solution of the equation.\\ We begin by considering a perturbation analysis when the initial value $y_0$ is perturbed to $\tilde{y}_0$ with relative error \[ \varepsilon=\frac{\norm{\tilde{y}_0-y_0}}{\norm{y_0}}, \] where $\norm{\cdot}$ is a vector norm on $\mathbb{R}^n$. Due to perturbation in the initial value, the solution $y(t)=e^{tA}y_0$ is perturbed to $\tilde{y}(t)=e^{tA}\tilde{y}_0$ with relative error \[ \delta(t)=\frac{\left\| e^{tA}\tilde{y}_0-e^{tA}y_0\right\|}{\left\| e^{tA}y_0\right\|}. \] In other words, it is the (relative) conditioning of the problem \begin{equation*} y_0\mapsto e^{tA}y_0. \end{equation*} The relation between the error $\varepsilon$ and the error $\delta(t)$ is described by three condition numbers namely: the condition number with the direction of perturbation, the condition number independent of the direction of perturbation and the condition number not only independent of the specific direction of perturbation but also independent of the specific initial value. How these condition numbers behave over a long period of time is an important aspect of the study. In the thesis, we move towards perturbations in the matrix as well as componentwise relative errors, rather than normwise relative errors, for perturbations of the initial value. About the first topic of the thesis, we look over how perturbations propagate along the solution of the ODE, when it is the coefficient matrix $A$ rather than the initial value that perturbs. In other words, the interest is to study the conditioning of the problem $$ A\mapsto e^{tA}y_0. $$ In case when the matrix $A$ perturbs to $\tilde{A}$, the relative error is given by \[ \epsilon=\frac{\vertiii{\tilde{A}-A}}{\vertiii{A}} \] and the relative error in the solution of the ODE is given by \[ \xi(t)=\frac{\left\| e^{t\widetilde{A}}y_0-e^{tA}y_0\right\|}{\left\| e^{tA}y_0\right\|}. \] We introduce three condition numbers as before. The analysis of the condition numbers is done for a normal matrix $A$ and by making use of $2$-norm. We give very useful upper and lower bounds on these three condition numbers and we study their asymptotic behavior as time goes to infinity. There could be cases when someone is interested in the relative errors \[ \delta_l(t)=\frac{\vert \tilde{y}_l(t)-y_l(t) \vert}{\vert y_l(t) \vert}, \quad l=1,\dots,n, \] of the perturbed solution components. With the motivation that componentwise relative errors give more information than the normwise relative error, we make a componentwise relative error analysis, which is the other topic of this thesis. We consider perturbations in initial value $y_0$ with normwise relative error $\varepsilon$ and the relative error in the components of the solution of the equation given by $\delta_l(t)$. The interest is to study, for the $l$-th component, the conditioning of the problem $$ y_0\mapsto y_l(t)=e_l^Te^{tA}y_0, $$ where $e_l^T$ is the $l$-th vector of the canonical basis of $\mathbb{R}^n$. We make this analysis for a diagonalizable matrix $A$, diagonalizability being a generic situation for the matrix $A$. We give two condition numbers in this part of the thesis and study their asymptotic behavior as time goes to infinity. The thesis is about how perturbations in the initial value $y_{0}$ or in the coefficient matrix $A$ propagate along the solutions of $n$-dimensional linear ordinary differential equations (ODE) \begin{equation*} \left\{ \begin{array}{l} y^\prime(t) =Ay(t),\ t\geq 0,\\ y(0)=y_0, \end{array} \right. \end{equation*} where $A \in \mathbb{R}^{n\times n} $ and $y_0\in \mathbb{R}^n$ and $y(t)=e^{tA}y_0$ is the solution of the equation.\\ We begin by considering a perturbation analysis when the initial value $y_0$ is perturbed to $\tilde{y}_0$ with relative error \[ \varepsilon=\frac{\norm{\tilde{y}_0-y_0}}{\norm{y_0}}, \] where $\norm{\cdot}$ is a vector norm on $\mathbb{R}^n$. Due to perturbation in the initial value, the solution $y(t)=e^{tA}y_0$ is perturbed to $\tilde{y}(t)=e^{tA}\tilde{y}_0$ with relative error \[ \delta(t)=\frac{\left\| e^{tA}\tilde{y}_0-e^{tA}y_0\right\|}{\left\| e^{tA}y_0\right\|}. \] In other words, it is the (relative) conditioning of the problem \begin{equation*} y_0\mapsto e^{tA}y_0. \end{equation*} The relation between the error $\varepsilon$ and the error $\delta(t)$ is described by three condition numbers namely: the condition number with the direction of perturbation, the condition number independent of the direction of perturbation and the condition number not only independent of the specific direction of perturbation but also independent of the specific initial value. How these condition numbers behave over a long period of time is an important aspect of the study. In the thesis, we move towards perturbations in the matrix as well as componentwise relative errors, rather than normwise relative errors, for perturbations of the initial value. About the first topic of the thesis, we look over how perturbations propagate along the solution of the ODE, when it is the coefficient matrix $A$ rather than the initial value that perturbs. In other words, the interest is to study the conditioning of the problem $$ A\mapsto e^{tA}y_0. $$ In case when the matrix $A$ perturbs to $\tilde{A}$, the relative error is given by \[ \epsilon=\frac{\vertiii{\tilde{A}-A}}{\vertiii{A}} \] and the relative error in the solution of the ODE is given by \[ \xi(t)=\frac{\left\| e^{t\widetilde{A}}y_0-e^{tA}y_0\right\|}{\left\| e^{tA}y_0\right\|}. \] We introduce three condition numbers as before. The analysis of the condition numbers is done for a normal matrix $A$ and by making use of $2$-norm. We give very useful upper and lower bounds on these three condition numbers and we study their asymptotic behavior as time goes to infinity. There could be cases when someone is interested in the relative errors \[ \delta_l(t)=\frac{\vert \tilde{y}_l(t)-y_l(t) \vert}{\vert y_l(t) \vert}, \quad l=1,\dots,n, \] of the perturbed solution components. With the motivation that componentwise relative errors give more information than the normwise relative error, we make a componentwise relative error analysis, which is the other topic of this thesis. We consider perturbations in initial value $y_0$ with normwise relative error $\varepsilon$ and the relative error in the components of the solution of the equation given by $\delta_l(t)$. The interest is to study, for the $l$-th component, the conditioning of the problem $$ y_0\mapsto y_l(t)=e_l^Te^{tA}y_0, $$ where $e_l^T$ is the $l$-th vector of the canonical basis of $\mathbb{R}^n$. We make this analysis for a diagonalizable matrix $A$, diagonalizability being a generic situation for the matrix $A$. We give two condition numbers in this part of the thesis and study their asymptotic behavior as time goes to infinity.

We present a loosely coupled scheme for the numerical simulation of the cardiac electro-fluid-structure interaction problem, whose solution is typically computationally intensive due to the need to suitably treat the coupling of the different submodels. Our scheme relies on a segregated treatment of the subproblems, in particular on an explicit Robin-Neumann algorithm for the fluid-structure interaction, aiming at reducing the computational burden of numerical simulations. The results, both in an ideal and a realistic cardiac setting, show that the proposed scheme is stable at the regimes typical of cardiac simulations. From a comparison with a scheme with implicit fluid-structure interaction, it emerges that, while conservation properties are not fully preserved, computational times significantly benefit from the explicit scheme. Overall, the explicit discretization represents a good trade-off between accuracy and cost, and is a valuable alternative to implicit schemes for fast large-scale simulations.