Julien Sigüenza
email: j.siguenza@sim-and-cure.com

Desiree Pott

Simon Mendez

Simon J Sonntag

S Kaufmann

Ulrich Steinseifer

Franck Nicoud

T A S Kaufmann

Fluid-structure interaction of a pulsatile flow with an aortic valve model: A combined experimental and numerical study

Keywords: fluid-structure interaction, aortic valve, flexible leaflets, particle image velocimetry, immersed boundary method, unstructured fluid solver, multi-cycle analysis, turbulence, large eddy simulation

ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

INTRODUCTION

The aortic heart valve separates the left ventricle from the aorta. It is composed of three thin flexible leaflets that open and close passively during the cardiac cycle, preventing blood from flowing back into the left ventricle, and thus ensuring an unidirectional blood flow through the cardiovascular system. Aortic valves may degenerate and lead to either insufficiency or stenosis, which can cause the death of the patient if not treated. Usually, medication is not a sufficient treatment option, and the diseased aortic valve needs to be replaced by a prosthesis.

During the cardiac cycle, blood exerts a continuous loading on the aortic valve, which experiences compression, stretching, and bending stresses [START_REF] Sacks | On the biomechanics of heart valve function[END_REF]. Conversely, the interaction of the pulsatile blood flow with the valve flexible leaflets gives rise to complex flow structures in the vicinity of the valve, which could eventually lead to transition to turbulence [START_REF] Sotiropoulos | Fluid mechanics of heart valves and their replacements[END_REF]. Assessing these complex structural and flow features could ultimately allow to improve the understanding of the degeneration process leading to aortic valve calcification, but also drastically help the design of aortic valve prostheses. Indeed, many ex vivo studies suggest that mechanical factors such as pressure, leaflets tension, and fluid shear stresses play an important role in the activation of the inflammatory pathways leading to the aortic valve calcification [START_REF] Xing | Cyclic pressure affects the biological properties of porcine aortic valve leaflets in a magnitude and frequency dependent manner[END_REF][START_REF] Butcher | Transcriptional profiles of valvular and vascular endothelial cells reveal phenotypic differences[END_REF][START_REF] Balachandran | Elevated cyclic stretch alters matrix remodeling in aortic valve cusps: implications for degenerative aortic valve disease[END_REF][START_REF] Balachandran | Elevated cyclic stretch induces aortic valve calcification in a bone morphogenic protein-dependent manner[END_REF][START_REF] Sucosky | Elevated cyclic stretch alters matrix remodeling in aortic valve cusps: implications for degenerative aortic valve disease[END_REF][START_REF] Smith | Cyclic strain inhibits acute pro-inflammatory gene expression in aortic valve interstitial cells[END_REF][START_REF] Sun | Single-cell microfluidic impedance cytometry: a review[END_REF][START_REF] Sun | Ex vivo evidence for the contribution of hemodynamic shear stress abnormalities to the early pathogenesis of calcific bicuspid aortic valve disease[END_REF][START_REF] Miller | Calcific aortic valve stenosis: methods, models, and mechanisms[END_REF][START_REF] Fisher | Calcific nodule morphogenesis by heart valve interstitial cells is strain dependent[END_REF]. Likewise, mechanical features such as valve durability, transvalvular pressure gradient, valve regurgitation, flow turbulence, fluid shear stresses, blood stagnation and flow separation regions have been pointed out as being relevant indicators to evaluate the performances of aortic valve prostheses [START_REF] Yoganathan | Fluid mechanics of heart valves[END_REF][START_REF] Claessens | Image-Based Computational Modeling of the Human Circulatory and Pulmonary Systems, chap. Mechanical Valve Fluid Dynamics and Thrombus Initiation[END_REF].

Assessing such detailed mechanical features is hardly feasible by experimental means [START_REF] Yap | Experimental measurement of dynamic fluid shear stress on the aortic surface of the aortic valve leaflet[END_REF][START_REF] Moore | Spatiotemporal complexity of the aortic sinus vortex[END_REF], whereas numerical simulation constitutes a powerful tool which provides a full access to the details of the flow and mechanical quantities of interest, such as the turbulent characteristics of the flow or the stresses undergone by the valve leaflets. As a result, important research efforts have been made on the fluid-structure interaction (FSI) modeling of flexible native and prosthetic aortic valves, either by developing dedicated FSI models [START_REF] Hart | A two-dimensional fluid-structure interaction model of the aortic value[END_REF][START_REF] Hart | A three-dimensional computational analysis of fluid-structure interaction in the aortic valve[END_REF][START_REF] Griffith | Immersed boundary model of aortic heart valve dynamics with physiological driving and loading conditions[END_REF][START_REF] Borazjani | Fluid-structure interaction, immersed boundary-finite element method simulations of bio-prosthetic heart valves[END_REF][START_REF] Kamensky | An immersogeometric variational framework for fluid-structure interaction: Application to bioprosthetic heart valves[END_REF][START_REF] Hsu | Fluid-structure interaction analysis of bioprosthetic heart valves: significance of arterial wall deformation[END_REF][START_REF] Hsu | Dynamic and fluid-structure interaction simulations of bioprosthetic heart valves using parametric design with T-splines and Fung-type material models[END_REF][START_REF] Gilmanov | A numerical approach for simulating fluid structure interaction of flexible thin shells undergoing arbitrarily large deformations in complex domains[END_REF] or by using commercial FSI software packages [START_REF] Ranga | Computational simulations of the aortic valve validated by imaging data: evaluation of valve-sparing techniques[END_REF][START_REF] Weinberg | Transient, three-dimensional, multiscale simulations of the human aortic valve[END_REF][START_REF] Weinberg | A multiscale computational comparison of the bicuspid and tricuspid aortic valves in relation to calcific aortic stenosis[END_REF][START_REF] Sundaram | Aortic valve dynamics using a fluid structure interaction model -The physiology of opening and closing[END_REF][START_REF] Luraghi | Evaluation of an aortic valve prosthesis: Fluid-structure interaction or structural simulation[END_REF]. Such models however need to be supported by experimental data in order to establish the validity and reliability of the numerical results. Nonetheless, as pointed out by Kheradvar et al. [START_REF] Kheradvar | Emerging trends in heart valve engineering: Part iv. computational modeling and experimental studies[END_REF], there are still very few experimental validations of FSI valve models regarding flexible native or prosthetic aortic valves. Two main limitations are highlighted by Kheradvar et al. [START_REF] Kheradvar | Emerging trends in heart valve engineering: Part iv. computational modeling and experimental studies[END_REF]: the ability of FSI valve models to perform multiple cardiac cycles, probably due to the high computational cost of such FSI computations, and the difficulty to simulate the contact between the thin deformable leaflets when the valve is fully closed.

One of the main challenges when attempting to establish a proper FSI valve model is thus to be able to effectively handle the complex interaction of the flexible valve leaflets with the pulsatile blood flow. The present paper aims at demonstrating that such a complex interaction can be handled using the immersed thick boundary method (ITBM), previously developed by Sigüenza et al. [START_REF] Sigüenza | Validation of an immersed thick boundary method for simulating fluid-structure interactions of deformable membranes[END_REF]. In this previous work, the ITBM was described and extensively validated on a wide variety of FSI benchmarks, and was finally used to compute the interaction of a pulsatile flow with a flexible trileaflet aortic valve at a moderate Reynolds number (Re = 4500). Such complex FSI computations however require further validation, and thus need to be supported by experimental data, as pointed out before. One difficulty regarding the validation of solvers dedicated to aortic valve applications is the diversity of challenges met: 1) the geometry and the flow are complex; 2) there are fluid-structure interactions between the valve leaflets and the blood; 3) the Reynolds number is sufficiently high to yield turbulence transition or at least intermittency; 4) the aorta is deformable; 5) blood and valves (either native or bioprosthetic valves) are biological materials whose characterization is far from trivial [START_REF] Robertson | Cardiovascular Mathematics. Modeling and simulation of the circulatory system[END_REF] [START_REF] Kheradvar | Emerging trends in heart valve engineering: Part iv. computational modeling and experimental studies[END_REF] In the present study, experimental data are generated thanks to a bio-inspired experimental test rig which attempts to mimic the interaction of the blood flow with the aortic valve, namely a pulsatile flow driving a model of aortic valve to alternatively open and close at a moderate Reynolds number. It must however be emphasized that such a configuration does not attempt to perfectly reproduce the physiological flow conditions which the aortic valve is subjected to, neither to faithfully model the anisotropic and nonlinear material properties of the aortic valve tissue or the non-Newtonian behavior of the blood. In particular, the valve is not made of biological material, but is made of polyurethane. The objective is to set up a well-controlled experimental configuration easy to model and simulate, in order to perform a meaningful experiment-vs-simulation comparison, and thus effectively validate FSI numerical methods. This comparison between the simulation and the experiment includes: 1) the global and detailed dynamics of the valve, whose motion is experimentaly recorded using a high-speed camera; 2) the flow field in the vicinity of the valve, which is experimentally reconstructed by means of particle image velocimetry (PIV) measurements;

FLUID-STRUCTURE INTERACTION OF A PULSATILE FLOW WITH AN AORTIC VALVE MODEL 3
3) the turbulence of the flow downstream of the valve, which is characterized by performing a multicycle analysis (over 24 cycles).

MATERIALS AND METHODS

Numerical method

The Fluid-Structure Interaction model used in the present study is based on the immersed thick boundary method (ITBM), which is presented and extensively validated in [START_REF] Sigüenza | Validation of an immersed thick boundary method for simulating fluid-structure interactions of deformable membranes[END_REF]. The ITBM was adapted from the original immersed boundary method (IBM) introduced by Peskin [START_REF] Peskin | The immersed boundary method[END_REF] to deal with 3D membranes having a finite thickness. The computational domain is decomposed into the solid part (i.e. the valve) and the fluid part, which are both discretized with independent meshes. The mechanics of the valve is solved on a Lagrangian mesh allowed to move to follow the valve displacements, whereas the fluid motion is solved on a fixed Eulerian unstructured mesh. The different steps of the ITBM are the following:

1. The mechanical force -→ F resulting from the valve deformation is calculated on the Lagrangian mesh.

The volumetric force

-→ f , which represents the force exerted by the valve on the fluid, is calculated by regularizing the mechanical force -→ F on the Eulerian mesh, such as

-→ f ( -→ x , t) = Ωs -→ F -→ X , t δ -→ x - -→ X dX,
where -→ x and -→ X respectively denote the coordinates vectors of the Eulerian nodes and Lagrangian nodes, Ω s denotes the solid domain defining the valve, and δ is the well known Dirac function. 3. The Navier-Stokes equations (forced by the regularized mechanical forces) are solved on the Eulerian mesh, yielding the velocity of the fluid -→ v .

The velocity of the valve

-→ V is calculated by interpolating the fluid velocity -→ v on the Lagrangian mesh, such as

-→ V -→ X , t = Ω f -→ v ( -→ x , t) δ -→ x - -→ X dx,
where Ω f denotes the fluid domain.

The Dirac function δ used in the procedures of regularization and interpolation of steps ( 2) and ( 4) is numerically represented by a smooth dicrete Dirac function [START_REF] Peskin | The immersed boundary method[END_REF], which is adapted to unstructured meshes using the Reproducing Kernel Particle Method (RKPM) [START_REF] Sigüenza | Validation of an immersed thick boundary method for simulating fluid-structure interactions of deformable membranes[END_REF][START_REF] Pinelli | Immersed-boundary methods for general finitedifference and finite-volume Navier-Stokes solvers[END_REF][START_REF] Mendez | An unstructured solver for simulations of deformable particles in flows at arbitrary Reynolds numbers[END_REF]. Imposing the exact impermeability of the valve leaflets treated through the IBM framework is challenging. Indeed, interpolation procedures induce small leakage of the flow through the leaflets. Moreover, there is no particular treatment performed to solve the contact between the flexible leaflets during the valve closure. The present ITBM intrinsically prevents any penetration of the leaflets, the minimum gap between the leaflets depending on the size of the discrete Dirac function, which is typically chosen as being 4 times the fluid mesh resolution h.

The computation of the mechanical force -→ F is performed by the in-house LMGC90 solid mechanics solver [START_REF]Discrete Numerical Modeling of Granular Materials[END_REF], using the classical finite-element method. In the present work, a quasiincompressible Neo-Hookean material is used to model the valve, defined by its strain energy function:

W = G 2 Ī1 -3 + K 2 (ln J) 2 , (1)
where G and K are the shear modulus and the bulk modulus, respectively. Ī1 is the first invariant of the isochoric Right Cauchy-Green deformation tensor C, and J is the Jacobian of the transformation. Provided that fluid and solid densities are of the same order of magnitude, and that the thickness of the valve is small enough, the inertia of the valve can be neglected.

The in-house YALES2BIO flow solver is used to solve the forced Navier-Stokes equations over the Eulerian unstructured mesh using a projection method [START_REF] Chorin | Numerical solution of the Navier-Stokes equations[END_REF]. The momentum conservation equations reads:

∂ -→ v ∂t + -→ ∇. ( -→ v ⊗ -→ v ) = - -→ ∇p ρ + ν∆ -→ v + -→ f ρ , (2)
where -→ v and p are the velocity vector and pressure, ρ the density and ν the kinematic viscosity. For an incompressible flow, the mass conservation equation becomes:

-→ ∇. -→ v = 0.
(

) 3
The fluid velocity is advanced using a 4th-order centred scheme in space and a 4th-order Runge-Kutta scheme in time. A divergence-free velocity field is obtained at the end of the time-step by solving a Poisson equation for pressure to correct the predicted velocity. A Deflated Preconditioned Conjugate Gradient algorithm is used to solve this Poisson equation. The employed numerical method is described in more details in appendix, and additional information can be found in [START_REF] Moureau | Design of a massively parallel CFD code for complex geometries[END_REF][START_REF] Malandain | Optimization of the deflated conjugate gradient algorithm for the solving of elliptic equations on massively parallel machines[END_REF][START_REF] Chnafa | Image-based large-eddy simulation in a realistic left heart[END_REF][START_REF] Sigüenza | Fluid-structure interaction problems involving deformable membranes: application to blood flows at macroscopic and microscopic scales[END_REF]. Large Eddy Simulation (LES) is used to account for possible turbulence in the flow. LES has the ability to predict intermittent flows [START_REF] Chnafa | Image-based large-eddy simulation in a realistic left heart[END_REF][START_REF] Chnafa | Image-based simulations show important flow fluctuations in a normal left ventricle: What could be the implications?[END_REF]. In particular, the so-called Sigma subgrid model [START_REF] Nicoud | Using singular values to build a subgrid-scale model for large eddy simulations[END_REF][START_REF] Toda | Assessment of subgrid-scale models with a large-eddy simulation-dedicated experimental database: The pulsatile impinging jet in turbulent cross-flow[END_REF] is well adapted to configurations with weak turbulence, as it has been designed not to add subgrid viscosity in a number of laminar canonical cases. The present fluid solver was validated on the benchmark of idealized cardiovascular device proposed by the U.S. Food and Drug Administration (FDA) [START_REF] Zmijanovic | About the numerical robustness of biomedical benchmark cases: Interlaboratory FDA's idealized medical device[END_REF], as well as on a variety of test cases [START_REF] Mendez | An unstructured solver for simulations of deformable particles in flows at arbitrary Reynolds numbers[END_REF][START_REF] Martins Afonso | On the damped oscillations of an elastic quasicircular membrane in a two-dimensional incompressible fluid[END_REF].

Experimental setup

A self-made heart valve model with a rigid frame made with PEEK material and leaflets manufactured as thin polyurethane foil was used. The cylindrical shaped leaflet was heat treated and mechanically closed, enabling to obtain an almost closed design of the valve (see Fig. 1 (a)). This valve model was placed inside an aortic anatomy, cast inside a rigid silicone model. The aortic anatomy is based on physiological data, including the three sinuses of Valsalva after the calculations by Reul et al. [START_REF] Reul | The geometry of the aortic root in health, at valve disease and after valve replacement[END_REF]. The aortic root diameter is 25 mm, chosen equal to the outer diameter of the valve model. The ascending aorta is designed with a diameter of 31 mm.

The silicone model with the inserted valve model was connected to a pulse duplicator which consists of a rotatory pump (dp3, Medos Medizintechnik AG, Stolberg, Germany), an adjustable compliance and resistance, tubing and a reservoir with a heating unit. A sketch is shown in Fig. 1. Data acquisition and controlling of the pulse duplicator was achieved by in-house written programs (LabVIEW, Austin, TX, USA). The rotatory pump was controlled by a program which delivers a square-wave signal, leading to an oscillatory rotation of the pump which enables to generate a pulsatile flow allowing the full opening and closure of the valve. As shown in Fig. 1, the flow measured by an ultrasonic flow meter downstream of the aortic valve is quasi-periodic with reduced cycle-to-cycle variations (less than 5% of variations). One typical cycle is displayed in Fig. 2 (a). This study focuses on the most challenging phase of the cycle, the systole, when the aortic valve opens and closes [START_REF] Yoganathan | Fluid mechanics of heart valves[END_REF]. The pulsatile flow rate thus imposes the periodic opening and closure of the valve, with a small negative flow at the end of the cycle.

In large vessels, the Newtonian assumption for blood rheology is generally considered to be acceptable [START_REF] Robertson | Cardiovascular Mathematics. Modeling and simulation of the circulatory system[END_REF]. As blood analogue fluid, a transparent solution of 56.4% glycerol and 43.6% water by mass was used inside the flow circuit as working fluid. At 45 • C, it matches the refractive index of the silicone model to eliminate optical distortions and it possesses a constant viscosity of 3.6 mPa.s. Both the flow and the valve can be characterized by a set of parameters summarized in Tab. I.

Determination of the material parameters of the leaflet was carried out by means of experimental tensile tests and parameter identification using the LMGC90 solid mechanics solver. First, an uniaxial tensile test was performed on a stick of polyurethane, giving the stress-strain relationship of the material (see Fig. 2 (b)). Then, the experimental tensile test was simulated with the LMGC90 solid mechanics solver. A dedicated algorithm enabled to determine the material parameters which provide the best fit between the experiment and the simulation, as shown in Fig. 2 (b). The quasiincompressible Neo-Hookean strain energy function introduced in Eq. ( 1) enables a very good representation of the material mechanical behaviour over the range of 0-30% deformation.

The movement of the valve leaflets was recorded from the top of the valve by a high-speed camera (CMOS, 1280x1024, 10 bit). Furthermore, two-dimensional PIV measurements were performed to determine the experimental flow field in the centerplane (through the commissure of the valve) of the aortic geometry behind the valve model. A two-cavity Nd:Ylf laser (Pegasus, New Wave For data acquisition and PIV post-processing, the software Dynamic Studio 2.3 (Dantec Dynamics A/S, Skovlunde, Denmark) was used. 24 cycles were recorded with a time interval ∆t between the image pairs of 950 µs to obtain an optimal pixel displacement for further analysis. The data acquisition of the particle images was set to 100 Hz, i.e. 100 pair images were recorded per cycle. As calibration method, a target with a millimeter scale aligned with the centerplane, was recorded. With knowledge of the displacement of the particles and the time delay between the image pairs, the velocity can be calculated. For the evaluation, the recorded pair images were divided into small subareas called interrogation areas. By means of statistical methods, i.e., adaptive correlation, the local velocity vectors for the pair images were determined. The final interrogation area size is 32×32 pixels with an overlap of 50%. The resolution of the resulting vector field is 0.0425×0.0425 mm 2 /px 2 . The resulting mesh of the PIV data is 63×79 points (radial/axial).

Computational setup

From the computational point of view, the geometry of interest is separated into the solid and the fluid domains. The solid domain is composed of the polyurethane leaflet whose geometry after manufacturing is cylindrical of radius R and thickness e. After performing the heat treatment on the valve, the geometry of the leaflet is closed. To reproduce this heat treatment in the numerical simulation, several preprocessing steps were performed: first a CAD model of the open leaflet before heat treatment were designed; then a structural computation performed with the LMGC90 solver enabled to close the leaflet by imposing a normal pressure onto the outer surface of the leaflet, thus mimicking the mechanical closure performed in the experiment (see Fig. 1 (b)). The obtained geometry was then used as the initial and unconstrained leaflet geometry for the FSI computation.

Both the aortic and the stentframe geometry are given by the CAD models of the experimental setup, and were used to design the geometry of the fluid domain. The silicone aorta and the stentframe are defined as rigid bodies whose borders are no-slip wall boundary conditions for the fluid domain. The flow is generated using inlet and outlet conditions, imposed far from the zone of interest (displayed in Fig. 3). One typical cycle of the flow rate signal measured during the experiment is imposed as a periodic inlet condition (see Fig. 2 (a)). Given the relatively high value of the Womersley number (W 0 = 17), the assumption of a flat velocity profile is made to impose the inlet flow rate signal. The inlet velocity is imposed as a Dirichlet condition normal to the inlet boundary:

--→ V in (t) = - Q(t) S -→ n 0 , (4)
where Q(t) corresponds to the flow rate displayed in Fig. 2 (a), S is the surface area of the inlet boundary, and -→ n 0 is the outward normal at the inlet boundary. The outlet is defined as a convective outlet boundary condition:

δ -→ v δt + V conv δ -→ v δ -→ n = 0, (5)
where -→ n is the outward normal at the outlet boundary, and V conv is the convective velocity, imposed in such a way to meet the global mass conservation. In order to avoid incoming flow at the outlet boundary, the inlet and outlet boundaries are switched depending if the flow rate is positive or negative. The forward flow is thus imposed upstream of the valve, whereas the backward flow is imposed downstream. The computational zone of interest is displayed in Fig. 3. Note that the computational fluid domain is extended both upstream and downstream of this zone of interest, in order to avoid any undesirable border effects. Two different levels of mesh were considered in the present study. Mesh M 1 is the coarse mesh: the fluid mesh is composed of 470 791 tetrahedral elements with a uniform mesh resolution of 1 mm, and the solid mesh is composed of 384 quadratic hexahedral elements with a mesh resolution of 1 mm. Mesh M 2 is the fine mesh: the fluid mesh is composed of 3 381 583 tetrahedral elements with a uniform mesh resolution of 0.5 mm, and the solid mesh is composed of 1536 quadratic hexahedral elements with a mesh resolution of 0.5 mm. Note that similar mesh resolutions were used for the fluid and the solid mesh, which is recommended when using the ITBM. We can also evaluate the Pope criterion [START_REF] Pope | Turbulent flows[END_REF][START_REF] Pope | Ten questions concerning the large-eddy simulation of turbulent flows[END_REF] which is used to estimate the subgrid dissipation caused by the LES filtering, expressed as:

k SGS k = 3 2 C ∆ πL 2/3 , (6)
where C = 1.5, ∆ = 0.5 mm being the mesh resolution, and L = 12.5 mm being the characteristic length of the largest structures. For the present LES, k SGS /k = 12% which remains below the 15 -20% threshold usually used to evaluate if a LES is sufficiently resolved [START_REF] Pope | Ten questions concerning the large-eddy simulation of turbulent flows[END_REF]. 3. RESULTS

Valve dynamics

Figure 4 (a) shows several snapshots of the valve movement recorded from the top by the highspeed camera, extracted over different instants of one cardiac cycle (depicted with red dots on the flow rate curve). Corresponding snapshots of the valve deformation extracted from a typical cycle of the simulation are displayed in Fig. 4 (b). The 3 leaflets are labeled to ease results presentation. The evolution of the opening area of the valve along time is also provided in Fig. 4 (c), as well as the evolution of the tips locations (midpoints of leaflets edges) of the 3 leaflets in Fig. 4 (d).

In both the experiment and the simulation, the valve shows a strong resistance to flexural deformation. The openings of the experimental and simulated valve are comparable. As seen in Fig. 4 (d), the 3 leaflets open simultaneously at the same opening times (between t 0 and t 2 ). After complete valve opening (after t 2 ), the 3 leaflets of the experimental valve behave slightly differently. Leaflets 1 and 3 start closing directly after the peak flow rate (which is reached at t 3 ), with an earlier closure for leaflet 3. Conversely, leaflet 2 reaches an equilibrium state indicating a buckling phenomenon. Closure of the leaflet suddenly occurs between t 6 and t 7 , due to the unbuckling of the leaflet resulting from the backflow. In the simulation, the valve opening is followed by a fluttering of the leaflets tips (between t 2 and t 4 ) and all the 3 leaflets then reach an equilibrium buckled state. For the cardiac cycle considered, leaflet 1 is the first one to close, followed by leaflet 3, and finally leaflet 2, which closes at the same time as in the experiment (see also snapshots of t 7 in Fig. 4 (a)-(b)). Although the detailed dynamics of the experimental and simulated valves significantly differ regarding the evolution of the tips locations of the 3 leaflets, a good quantitative agreement is obtained regarding the global dynamics given by the evolution of the opening area of the valve (Fig. 4 (c)), even when the coarse mesh M 1 is considered. The opening area first quickly increases during valve opening, due to the simultaneous opening of the leaflets. The opening section of the valve then gradually decreases, due to the progressive closure of leaflets 1 and 3 in the experiment, and due to a reduction of the opening area near the commissures of the valve in the simulation (see the difference between instants t 4 and t 5 in Fig. 4 (b)). While a strong fluttering of the leaflets tips is observed between t 2 and t 4 (in Fig. 4 (d)), a fluttering phenomenon on the opening area evolution is still observed between t 4 and t 6 , induced by the commissures of the valve which fluctuate (not shown).

The dynamics of the leaflet edges during valve opening, valve open position and valve closure is further analyzed in Fig. 5. The shapes of the 3 leaflets are first compared for a given cardiac cycle in Fig. 5 (a). In the experiment, the leaflets open with 2 curvature modes with a clockwise propagating curvature. In the simulation however, the deformation of the leaflets is symmetric and composed of 3 modes of curvature. Note that such 3-mode openings were already observed in previous experimental studies performed on bioprosthetic heart valves [START_REF] Gao | Bioprosthetic heart valve leaflet motion monitored by dual camera stereo photogrammetry[END_REF][START_REF] Watton | Dynamic modelling of prosthetic chorded mitral valves using the immersed boundary method[END_REF]. When the valve is fully open, the curvature of the valve is completely reversed with one unique mode of deformation. A small asymmetry of the experimental valve is observed whereas the simulated valve is perfectly symmetric. During the experimental valve closure, the leaflets recover their closed configuration again with 2 curvature modes, this time with a counterclockwise propagating curvature. In the simulation the leaflets also close with 2 curvature modes but the curvature propagates either toward one commissure or the other.

In order to describe the reproductibility of this dynamics over several cardiac cycles, the shape of leaflet 1 is compared for 3 consecutive cycles. It is seen that the dynamics of the experimental valve is virtually reproducible, but the simulated valve exhibits some cycle-to-cycle variations, especially during the valve closure. There is actually no privileged direction, from one cycle to the other, for the propagation of the curvature when the leaflets unbuckle. It is also observed that the order of closure of the leaflets (visible in Fig. 4 (d)) can change from one cycle to the other (not shown). These cycle-to-cycle variations come from the buckling of the leaflets which is a very sensitive and unstable phenomenon. In the experiment, the dynamics of the valve is more reproducible due to the asymmetry of the valve which enables one leaflet only to properly buckle (as seen in Fig. 4 (d)).

Figure 6 shows the Von Mises stress distribution in the leaflets obtained from the simulation when the valve is opening, fully open and closing. Stress patterns are mostly located in the region where the leaflets experience high flexural deformation. Maximum stresses are observed along the junction between leaflets and the rigid frame, particularly in the basal region. In this sense, the 2mode curvature observed during the opening of the experimental valve (in Fig. 5) would minimize the stresses undergone by the valve leaflets compared to the 3-mode curvature observed on the simulated valve.

Flow dynamics

The two-dimensional velocity field (u, v) downstream of the valve is reconstructed from the PIV measurements in the centerplane of the aortic geometry, where u is the streamwise velocity component and v the transverse velocity component. Instantaneous streamwise velocity fields obtained from the experiment and the simulation are displayed in Fig. 7. Four characteristic instants of the cardiac cycle are shown: t 2 = 0.20 s, just before the flow rate reaches its maximum value, referred to as Early Systole (ES) for convenience; t 3 = 0.26 s, when the flow rate is maximum, referred to as Peak Systole (PS); t 4 = 0.35 s, just after the peak flow rate, referred to as Mid-Systole (MS); t 5 = 0.55 s, when the flow rate is decreasing, referred to as Late Systole (LS). The color scale was chosen such that red regions correspond to positive velocity, whereas negative velocity regions appear in blue. The velocity fluctuations downstream of the valve are also investigated by calculating the experimental and numerical streamwise root-mean square (RMS) velocity fields on 24 successive cardiac cycles, as shown in Fig. 8. The local extrema of the velocity RMS fields are indicated by different labels (black labels for the experiment, and white labels for the simulation), enabling to identify different flow features. The same labels are reported in the velocity fields of Fig. 7 in order to characterize the source of these fluctuations.

At ES, a jet starts to emerge from the valve. In the experiment, the jet is widely spread compared to the simulation, and seems to be slightly tilted towards the sinus side of the aorta. Fluctuations are observed along the contour of the jet (see labels 1a, 2a, 3a). In the simulation, most of the fluctuations are observed along the leaflet (see label 1b).

At PS, the jet propagates downstream of the valve and grows broader. Flow recirculation is observed along both the sinus and commissure sides of the aorta, as indicated by the negative velocities. Similar fluctuations are observed at different critical locations of the jet structure: At the tip of the leaflet (see labels 5a vs. 2b), near to the shear layers generated by the recirculation of the jet (see labels 6a vs. 3b and 7a vs. 4b), and at the apex of the jet towards the sinus side of the aorta (see labels 8a vs. 5b). Regarding the structure of the jet, the flow seems to remain laminar, and the fluctuations of velocity are likely related to cycle-to-cycle variations in the shear layer positions and in the dynamics of coherent structures. At MS, the jet progresses further downstream of the valve, and the recirculation zones get stronger near the aorta wall, both at the commissure and the sinus sides. The flow seems to transition from laminar to turbulent along the shear layers: larger fluctuations are indeed observed at the commissure side (see labels 11a, 12a, 13a vs. 6b, 7b) and at the sinus side of the aorta (see labels 14a, 15a vs. 8b, 9b), in the wake of the valve. These are indicators of a production of turbulence downstream of the valve. Both for the experiment and the simulation, the maximum fluctuations are observed at the commissure side of the aorta. At LS, when the flow rate decreases, the jet decelerates and the recirculation velocities decrease. The two areas where fluctuations were observed at MS are still identifiable at the commissure side (see labels 16a, 17a, 18a vs. 10b, 11b) and sinus side of the aorta (see labels 19a, 20a, 21a vs. 12b). These areas are however more spread and start to dissipate downstream of the valve.

The experiment and the simulation results are thus qualitatively in good agreement. The main flow features are similar (both for the coarse mesh M 1 and the fine mesh M 2 ). However, the simulation seems to underestimate the backflow induced by the recirculation of the jet compared to the experiment. Note also that the magnitude of the fluctuations observed in the experiment is much higher than in the simulation (with maximum fluctuations of 0.2 m/s in the experiment against 0.04 m/s in the simulation, both for meshes M 1 and M 2 ). Moreover, the fluctuations are much more concentrated in the simulation, and are conversely more spread in the experiment. This difference is further commented in the discussion section. Note that additional comparisons for a plane with a 6 mm offset from the centerplane are provided in the thesis of Sigüenza [START_REF] Sigüenza | Fluid-structure interaction problems involving deformable membranes: application to blood flows at macroscopic and microscopic scales[END_REF].

In order to quantitatively compare the experimental flow and the numerical flow, phase-averaged streamwise velocity profiles are displayed in Fig. 9. Three successive velocity profiles downstream of the valve are compared, referred to as profile 1 (P1) which is located directly behind the valve, profile 2 (P2) and profile 3 (P3) which are located further downstream of the valve (see Fig. 7 for the exact locations). In spite of the high complexity of the flow, the experimental and numerical velocity profiles are generally in a very satisfactory agreement, even when the coarse mesh M 1 is used. The fine mesh M 2 however seems to provide a better agreement with the experiment. However, few discrepancies are observed. A more pinched numerical velocity profile is observed on P1 when the jet emerges from the valve at ES (see P1-ES). In addition, the simulation does not show the strong backflow observed along the commissure side of the aorta in the experiment (see P2-PS and P1-MS). Note that these differences are robust to the change of mesh, so that they are not the result of an insufficient mesh resolution. These discrepancies were expected given that the detailed dynamics (notably the buckling dynamics) of the three leaflets differ between the experiment and the simulation (see Fig. 4 (d)). Transverse velocity profiles for the same locations are provided in the thesis of Sigüenza [START_REF] Sigüenza | Fluid-structure interaction problems involving deformable membranes: application to blood flows at macroscopic and microscopic scales[END_REF].

Figure 10 (a) shows iso-surfaces of the instantaneous Q-criterion [START_REF] Hunt | Eddies, streams, and convergence zones in turbulent flows[END_REF], to visualize instantaneous vortical structures. At ES, vortical structures are shed downstream of the valve commissures, at the tips of the rigid frame. These structures propagate between ES and MS over an annular region in the wake of the valve. At LS, the vortical structures are smaller and spread over the lumen of the aorta.

Figure 10 (b) shows iso-surfaces of the fluctuating kinetic energy, defined as k = ρ 2 u 2 rms + v 2 rms + w 2 rms . This enables to better visualize the structure of the velocity fluctuations, partially presented in Fig. 8 over cutting planes. At ES, velocity fluctuations are related to the fluctuations of position of the shear layers along the leaflets. On the contrary, at PS, high k values are related to strong vortices shed downstream of the valve commissures. At PS, the iso-surface of k demonstrates fluctuations are much higher downstream of the commissure than downstream of the leaflet tips: remarkably, the fixed part of the valve generates more turbulence than the leaflets. Between PS and MS, the transition from laminar to turbulent is clearly suspected by the strong increase and presence of high k levels in the whole wake of the valve. The iso-surface of k is rather toroidal at MS, showing that the turbulence region is localized in the valve wake. At LS, the k iso-surface spreads over a wider volume of the aorta, as anticipated from the Q-criterion data.

DISCUSSION

In the present work, a complex fluid-structure interaction problem relevant to the flow of blood passing through the aortic valves has been investigated, both experimentally and numerically. An experimental test rig was set up, consisting of a pulse duplicator which generates a pulsatile flow through a realistic aortic root geometry, in which a model of aortic valve was placed. This constitutes a simplified and numerically reproducible configuration which was simulated with the immersed thick boundary method previously developed by Sigüenza et al. [START_REF] Sigüenza | Validation of an immersed thick boundary method for simulating fluid-structure interactions of deformable membranes[END_REF]. This configuration is thus numerically tractable, while taking into account the complex geometries of both the aorta and the tri-leaflet aortic valve. Experimental and numerical results can thus be meaningfully compared as regards to the fluid-structure interaction between the flow and the flexible leaflets of the valve.

Investigation of the valve dynamics shows similar global dynamics between the experimental and simulated valves. Discrepancies are nonetheless observed when analyzing the individual dynamics of the leaflets, with a 2-mode opening observed in the experiment against a 3-mode opening observed in the simulation. This symmetric 3-mode opening is however expected given the symmetry of the geometry and flow conditions of the computational setup. The non-symmetric swirl motion observed in the experiment can be explained by the dissymmetry observed in the experimental valve, probably induced by the manual heat treatment performed on the leaflets to obtain the closed design of the valve.

Investigation of the flow downstream of the valve shows similar features between the experiment and the simulation. A flow transition to turbulence seems to occur after peak systole, when the flow starts to decelerate. Spots of fluctuating kinetic energy are observed in the wake of the valve, and then spread in the lumen of the aorta. The authors originally expected a maximum level of fluctuations at the sinus side of the aorta, in the wake of the leaflets tips. However, the present study shows that fluctuations are maximum at the commissure side of the aorta. This should be induced by the significant thickness of the rigid frame maintaining the leaflets, which is seen as a backward-facing step known to generate turbulence [START_REF] Kasagi | Three-dimensional particle-tracking velocimetry measurement of turbulence statistics and energy budget in a backward-facing step flow[END_REF]. Note that experimental and numerical phase-averaged velocity profiles are in good agreement. The simulation however underestimates the levels of velocity fluctuations measured in the experiment, even if it provides a good qualitative agreement with the experiment in terms of velocity fluctuations structure. In the experiment, higher levels of fluctuations are observed upstream of the valve, which indicates a possible production of turbulence of the rotatory pump which is used to generate the flow. In addition, the flow rate measured in the experiment exhibits some cycle-to-cycle variations which are not taken into account in the simulation (the same flow rate is imposed for all the cycles simulated). Another possible explanation for this discrepancy is the LES filtering which may cut off a part of the turbulent kinetic energy. However, the evaluation of the Pope criterion [START_REF] Pope | Turbulent flows[END_REF][START_REF] Pope | Ten questions concerning the large-eddy simulation of turbulent flows[END_REF] has revealed that the LES was sufficiently resolved, which cuts off only 12% of the total turbulent kinetic energy.

The few discrepancies observed in the present comparison between the experiment and the simulation point out the complexity of carrying out such a combined experimental and numerical study, especially in such a very sensitive configuration involving fluid-structure interaction and buckling. In spite of these discrepancies, the results obtained in this study are very promising, demonstrating the ability of the ITBM [START_REF] Sigüenza | Validation of an immersed thick boundary method for simulating fluid-structure interactions of deformable membranes[END_REF] to handle such a complex FSI problem. It is important to note that despite several FSI modeling of aortic valve exist in the literature [START_REF] Hart | A two-dimensional fluid-structure interaction model of the aortic value[END_REF][START_REF] Hart | A three-dimensional computational analysis of fluid-structure interaction in the aortic valve[END_REF][START_REF] Griffith | Immersed boundary model of aortic heart valve dynamics with physiological driving and loading conditions[END_REF][START_REF] Borazjani | Fluid-structure interaction, immersed boundary-finite element method simulations of bio-prosthetic heart valves[END_REF][START_REF] Kamensky | An immersogeometric variational framework for fluid-structure interaction: Application to bioprosthetic heart valves[END_REF][START_REF] Hsu | Fluid-structure interaction analysis of bioprosthetic heart valves: significance of arterial wall deformation[END_REF][START_REF] Hsu | Dynamic and fluid-structure interaction simulations of bioprosthetic heart valves using parametric design with T-splines and Fung-type material models[END_REF][START_REF] Gilmanov | A numerical approach for simulating fluid structure interaction of flexible thin shells undergoing arbitrarily large deformations in complex domains[END_REF][START_REF] Ranga | Computational simulations of the aortic valve validated by imaging data: evaluation of valve-sparing techniques[END_REF][START_REF] Weinberg | Transient, three-dimensional, multiscale simulations of the human aortic valve[END_REF][START_REF] Weinberg | A multiscale computational comparison of the bicuspid and tricuspid aortic valves in relation to calcific aortic stenosis[END_REF][START_REF] Sundaram | Aortic valve dynamics using a fluid structure interaction model -The physiology of opening and closing[END_REF][START_REF] Luraghi | Evaluation of an aortic valve prosthesis: Fluid-structure interaction or structural simulation[END_REF], none of them have been subjected to such an extended experimental validation, including both the dynamics of the valve and the flow field in the vicinity of the valve. The most advanced experimental validation of aortic valve modeling found in the literature is the study of Luraghi et al. [START_REF] Luraghi | Evaluation of an aortic valve prosthesis: Fluid-structure interaction or structural simulation[END_REF], which was however only focused on the validation of the valve dynamics. One of the main added value of the present work remains the characterization of the turbulent features of the flow downstream of the valve. Indeed, most of the papers presenting aortic valve computations only consider a few number of cycles, because of the high computational challenge that represents such computations. This makes the investigation of the turbulent flow features quite challenging. This computational challenge has nonetheless been taken up for the first time by the present numerical method, which has been able to simulate 24 cardiac cycles. Still, the present study has some limitations discussed in the following section.

LIMITATIONS

One of the main limitation of the present work is the way the experimental flow conditions are imposed in the simulation. Indeed, the computations have been performed imposing the experimental flow rate signal during the whole cycle. However, as pointed out by Sundaram et al. [START_REF] Sundaram | Aortic valve dynamics using a fluid structure interaction model -The physiology of opening and closing[END_REF], the use of such boundary conditions introduces a bias in the analysis of the valve dynamics, since only its response to the prescribed flow rate can be analyzed. Studying the valve in its natural environment would require to impose the impedance of the pumping device and of the downstream line, and was not considered by lack of acurate data on the one hand, and because of numerical stability issues on the other hand. Furthermore, the immersed thick boundary method presently used to handle the FSI coupling suffers from some leakage problems, characterized by a large under-estimation of the transvalvular pressure gradient (TPG) compared to the experiment during the backward flow phase (see Tab. II, TPG min ). There are two possible sources for these leakage problems: 1) the poor closure of the valve observed in the simulation, meaning that the valve leaflets never directly get into contact (see Fig. 11 (a)); 2) the small permeability of the valve leaflets, which is inherent to the IBM formalism.

Figure 11 shows that the first option is hardly conceivable. When the valve is supposed to be closed at t = 0.95 s (see Fig. 11 (a)), a small gap between the leaflets still remains. This however does not imply that the flow goes through this gap, and that the valve is leaking. As a matter of fact, the fluid sees the valve through the volumetric force -→ f which comes from the regularization of the valve mechanical force -→ F (see section 2.1). This mechanical force -→ F is regularized on a patch of fluid elements located at the neighboring of the valve (see Fig. 11 (b)), via the smooth discrete Dirac function. One thus easily understands that when there is an overlapping of two patches of fluid elements, corresponding to two different leaflets, this induces a local repulsion of the valve leaflets which are not able to get closer. As seen in Fig. 11 (b), when the valve is closed at t = 0.95 s, the overlapped patches of fluid elements describe a fully closed valve without any gap.

It is thus likely that porous leakage through the leaflets is the origin of the underestimation of the TPG min . Indeed, the valve may experience continuous leakage during the whole cardiac cycle, because of the errors generated by interpolation procedures of the IBM. The valve leakage is however difficult to quantify, and is expected to increase when the valve highly resists to the flow. Given that the valve slightly resists to forward flow and highly resists to backward flow, this would explain the large underestimation of the TPG min and the slight underestimation of the TPG max (see Tab. II). Despite of these issues which should be addressed in the future, the results obtained in the present work remain valid. They demonstrate the ability of the method the properly represent the valve dynamics and associated TPG during the forward flow phase.

Newtonian fluid, the Navier-Stokes equations and the mass conservation equation read:

         ∂ -→ v ∂t + -→ ∇. ( -→ v ⊗ -→ v ) = - -→ ∇p ρ + ν∆ -→ v + -→ f ρ -→ ∇. -→ v = 0 (7)
where -→ v is the velocity vector, p the pressure, and -→ f a possible volumetric force. ρ and ν are the properties, respectively the density and the kinematic viscosity.

The projection method

A projection method [START_REF] Chorin | Numerical solution of the Navier-Stokes equations[END_REF] is used to numerically solve Eq. [START_REF] Sucosky | Elevated cyclic stretch alters matrix remodeling in aortic valve cusps: implications for degenerative aortic valve disease[END_REF]. One can write the following time semi-discrete form of Eq. ( 7), considering an Euler explicit scheme:

         -→ v n+1 --→ v n ∆t = - -→ ∇. ( -→ v n ⊗ -→ v n ) - -→ ∇p n+1 ρ + ν∆ -→ v n + -→ f ρ -→ ∇. -→ v n+1 = 0 (8)
To advance both the velocity -→ v n+1 and the pressure p n+1 , the velocity -→ v * is first advanced considering the pressure p n at the previous time step:

-→ v * --→ v n ∆t = - -→ ∇. ( -→ v n ⊗ -→ v n ) - -→ ∇p n ρ + ν∆ -→ v n + -→ f ρ (9)
This first estimation of the velocity -→ v * is then substracted to the final velocity -→ v n+1 , which gives:

-→ v n+1 --→ v * = -∆t -→ ∇ p n+1 -p n ρ . (10)
Finally, taking the divergence of Eq. ( 10) and applying the mass conservation constraint (∇. -→ v n+1 = 0) leads to:

∆ p n+1 -p n ρ = -→ ∇. -→ v * ∆t . (11)
The pressure p n+1 is thus advanced by solving the Poisson equation of Eq. ( 11), using a Deflated Preconditioned Conjugate Gradient (DPCG) algorithm [START_REF] Moureau | Design of a massively parallel CFD code for complex geometries[END_REF][START_REF] Malandain | Optimization of the deflated conjugate gradient algorithm for the solving of elliptic equations on massively parallel machines[END_REF]. The first estimation of the velocity -→ v * is then corrected by the calculated pressure p n+1 (using Eq. ( 10)), yielding the final velocity -→ v n+1 :

-→ v n+1 = -→ v * -∆t -→ ∇ p n+1 -p n ρ . ( 12
)

Time advancement

The time advancement of the velocity -→ v * is performed using a 4th-order Runge-Kutta scheme in time (RK4) [START_REF] Williamson | Low-storage Runge-Kutta schemes[END_REF]. Rewritting Eq. ( 9) as:

-→ v * --→ v n ∆t = f ( -→ v n , p n ) , (13)

Figure 1 .

 1 Figure 1. In-vitro pulse duplicator designed to reproduce physiological flow conditions.

Figure 2 .

 2 Figure 2. (a) One typical cycle of the flow rate signal measured during the experiment, used as inlet boundary condition of the simulation. (b) Stress-strain relationship extracted from the uniaxial tensile test.The material parameters G and K are determinated by fitting the simulation with the experiment.

Flow parameters Density ρ = 1100 kg/m 3

 3 Dynamic viscosity µ = 3.6 10 -3 Pa.s Heart rate n bpm = 60 Mean cardiac output Q mean = 3.48 L/min Reynolds number R e = 1388 Womersley number W 0 = 17 Valve parameters Density ρ s = 1000 kg/m 3 Shear modulus G = 2.4 MPa Bulk modulus K = 1.6 MPa Thickness of the leaflet e = 0.15 mm Radius of the leaflet R = 12.5 mm Thickness of the frame e f = 1.45 mm Research Inc.) operated as light source. The laser light beam is expanded into a thin light sheet of approximately 1 mm thickness by an arrangement of lenses. The laser light illuminates the tracer particles (diameter of 10.5 µm, Intelligent Laser Applications, Jülich, Germany), seeded into the circuit and carried by the fluid, twice within a defined time interval ∆t to obtain pair images. A Nanosense MKIII 1280x1024 CMOS High-Speed camera (IDT, Redlake) with a 105 mm Nikon Makro Nikkor 105mm F2.8D lens recorded the particle images. The camera was positioned at an angle of 90 • to the illuminated laser light plane.

Figure 3 .

 3 Figure 3. Computational zone of interest composed of the valve solid mesh and the aortic fluid mesh. Two levels of mesh are considered, the coarse mesh M 1 and the fine mesh M 2 .

Figure 4 .

 4 Figure 4. (a) Snapshots extracted from the valve movement recording at various instants of the cardiac cycle (depicted with red dots on the flow rate curve). (b) Corresponding snapshots of the valve deformation extracted from a typical cycle of the simulation. (c) Evolution of the opening area of the valve during the cardiac cycle. (d) Evolution of the tips locations of the 3 leaflets. Solid line corresponds to leaflet 1, dashed line to leaflet 2, and dotted line to leaflet 3. (red: experiment, blue: simulation for the fine mesh M 2 ).

Figure 5 .

 5 Figure 5. (a) Shape of the 3 leaflet edges observed during valve opening (t 1 = 0.15 s), valve open position (t 4 = 0.35 s) and valve closure (t 6 = 0.80 s). Solid line corresponds to leaflet edge 1, dashed line to leaflet edge 2, and dotted line to leaflet edge 3. (b) Shape of the leaflet edge 1 observed during valve opening, valve fully open and valve closure, for 3 consecutive cycles. (+: outer radial direction, -: inner radial direction) (red: experiment, blue: simulation for the fine mesh M 2 ).

6 .

 6 Leaflets Von Mises stress distribution during valve opening (t 1 = 0.15 s), at the fully open position (t 4 = 0.35 s) and during valve closure (t 6 = 0.80 s).

Figure 7 .

 7 Figure 7. (a) Fields of instantaneous streamwise velocity obtained from the experiment. (b) Fields of instantaneous streamwise velocity obtained from the simulation for the fine mesh M 2 . Four different instants of the cardiac cycle are depicted (ES, PS, MS and LS).

Figure 8 .

 8 Figure 8. (a) Fields of streamwise RMS velocity obtained from the experiment. (b) Fields of streamwise RMS velocity obtained from the simulation for the fine mesh M 2 . Four different instants of the cardiac cycle are depicted (ES, PS, MS and LS).

9 .

 9 Comparison of the experiment (red) and the simulations for meshes M 1 (black) and M 2 (blue), for 3 successive phase-averaged velocity profiles downstream of the valve. Four different instants of the cardiac cycle are depicted (ES, PS, MS and LS).

Figure 10 .

 10 Figure 10. (a) Iso-surfaces of the instantaneous Q-criterion for the iso-value of 5000 s -2 . (b) Iso-surfaces of the fluctuating kinetic energy for the iso-value of 0.08 J.m -3 . Four different instants of the cardiac cycle are depicted (ES, PS, MS and LS).

Figure 11 .

 11 Figure 11. (a) Deformation of the valve at two different instants of the cardiac cycle (t = 0.15 s when the valve is opening, and t = 0.95 s when the valve is fully closed). (b) Patch of fluid elements on which the mechanical force is regularized, illustrating how the valve is seen by the fluid.

Table I .

 I Parameters of the flow and the valve leaflet.

Table II .

 II Comparison of the experimental and numerical transvalvular pressure gradients (TPG) during the forward flow phase (TPGmax) and backward flow phase (TPG min ).

		Experiment Simulation
	TPG max (mmHg)	5.28	4.90
	TPG min (mmHg)	-93.35	-21.05

J. SIG ÜENZA ET AL.

Copyright c 2010 John Wiley & Sons, Ltd. Int. J. Numer. Meth. Biomed. Engng. (2010) Prepared using cnmauth.cls DOI: 10.1002/cnm

ACKNOWLEDGEMENT

J. Sigüenza and D. Pott acknowledge the European Society for Artificial Organs (ESAO) for funding the yESAO Exchange Program which enabled to achieve this collaborative work. V. Moureau and G. Lartigue from the CORIA lab, and the SUCCESS scientific group are acknowledged for providing the YALES2 solver which constitutes the basis of the YALES2BIO solver. D. Ambard, F. Dubois, F. Jourdan and R. Mozul from the LMGC lab are also acknowledged for providing the LMGC90 solid mechanics solver, and for their active collaboration in the coupling of YALES2BIO with LMGC90. The authors acknowledge supports from ANR (FORCE project ANR-11-JS09-0011), from BPIfrance (DAT@DIAG project N • I1112018W) and from the NUMEV Labex (ANR-10-LABX-20). This work was performed using HPC resources from GENCI-CINES (Grants 2015-c2014037194 and 2016-c2015037194). Desiree Pott is supported by the Deutsche Forschungsgemeinschaft (DFG) grant STE1680/5-1.

APPENDIX

In this section, the numerical method employed in the in-house YALES2BIO solver to solve the Navier-Stokes equations is described. In the case of an incompressible flow and assuming a the velocity -→ v * is advanced such as:

Finite-volume spatial discretization

The finite-volume method (FVM) is used for the spatial discretization, which is based on the integration of the equations on small polyhedral control volumes, called dual volumes. The definition of the dual volume V j associated to node j of an unstructured hybrid mesh is illustrated in Fig. 12. The integration of a quantity φ over the dual volume V j is expressed as:

Integrating the momentum equation over each dual volume V j gives:

where I is the identity matrix. Using the Green-Ostrogradski theorem, Eq. ( 16) can be written as:

with:

S j represents the surface bordering the dual volume V j and -→ n the vector normal to this surface.

The surface integral of the flux ψ represents a transport of momentum towards neighboring control volumes, and illustrates the conservative nature of the Navier-Stokes equations. The right-hand side of Eq. ( 17) is approximated as follow:

where -→ S jk represents the non-normalized vector normal to the surface patch of the dual volume V j attributed to the pair (j, k), and ψ jk the flux through this surface patch, which can be basically expressed as:

This expression of the flux ψ jk makes the scheme second-order accurate for regular meshes, and only first-order accurate for irregular meshes [START_REF] Vantieghem | Numerical simulations of quasi-static magnetohydrodynamics using an unstructured finite volume solver: development and applications[END_REF], potentially implying important diffusion and dispersion errors. High-order schemes enable to reduce these numerical errors, and are thus generally less dissipative.

A fourth-order and low-dissipative spatial scheme has thus been developed within the YALES2BIO flow solver. It is based on a deconvolution of the finite-volume integration [START_REF] Malandain | Simulation massivement parallèle des écoulement turbulents faible nombre de mach[END_REF]. Indeed, when using the FVM, the nodal values φ j associated to the quantity φ are not directly known, only their integral φ V j over the control volumes V j are known. The development of a highorder scheme however requires writing this nodal value φ j in an accurate way. The deconvolution thus enables to write the nodal value φ j as a function of the integrated quantity φ V j and its successive derivatives, using a Taylor-series development. The nodal values φ j and φ k are then used to express the value φ jk associated to the surface patch of the dual volume V j attributed to the pair (j, k). The high-order scheme necessitates to pre-calculate the derivatives of φ at the nodes, which increases the stencil of the spatial scheme.