NUMERICAL CHARACTERIZATION OF HEMODYNAMICS CONDITIONS NEAR AORTIC VALVE AFTER IMPLANTATION OF LEFT VENTRICULAR ASSIST DEVICE

Left Ventricular Assist Devices (LVADs) are implantable mechanical pumps that temporarily aid the function of the left ventricle. The use of LVADs has been associated with thrombus formation next to the aortic valve and close to the anastomosis region, especially in patients in which the native cardiac function is negligible and the aortic valve remains closed. Stagnation points and recirculation zones have been implicated as the main fluid dynamics factors contributing to thrombus formation. The purpose of the present study was to develop and use computer simulations based on a fluid-structure interaction (FSI) solver to study flow conditions corresponding to different strategies in LVAD ascending aortic anastomosis providing a scenario with the lowest likelihood of thrombus formation. A novel FSI algorithm was developed to deal with the presence of multiple structures corresponding to different elastic properties of the native aorta and of the LVAD cannula. A sensitivity analysis of different variables was performed to assess their impact of flow conditions potentially leading to thrombus formation. It was found that the location of the anastomosis closest to the aortic valve (within 4 cm away from the valve) and at the angle of 30 minimizes the likelihood of thrombus formation. Furthermore, it was shown that the rigidity of the dacron anastomosis cannula plays almost no role in generating pathological conditions downstream from the anastomosis. Additionally, the flow analysis presented in this manuscript indicates that compliance of the cardiovascular tissue acts as a natural inhibitor of pathological flow conditions conducive to thrombus formation and should not be neglected in computer simulations.


(Communicated by the associate editor name)
Abstract.Left Ventricular Assist Devices (LVADs) are implantable mechanical pumps that temporarily aid the function of the left ventricle.The use of LVADs has been associated with thrombus formation next to the aortic valve and close to the anastomosis region, especially in patients in which the native cardiac function is negligible and the aortic valve remains closed.Stagnation points and recirculation zones have been implicated as the main fluid dynamics factors contributing to thrombus formation.The purpose of the present study was to develop and use computer simulations based on a fluid-structure interaction (FSI) solver to study flow conditions corresponding to different strategies in LVAD ascending aortic anastomosis providing a scenario with the lowest likelihood of thrombus formation.A novel FSI algorithm was developed to deal with the presence of multiple structures corresponding to different elastic properties of the native aorta and of the LVAD cannula.A sensitivity analysis of different variables was performed to assess their impact of flow conditions potentially leading to thrombus formation.It was found that the location of the anastomosis closest to the aortic valve (within 4 cm away from the valve) and at the angle of 30 • minimizes the likelihood of thrombus formation.Furthermore, it was shown that the rigidity of the dacron anastomosis cannula plays almost no role in generating pathological conditions downstream from the anastomosis.Additionally, the flow analysis presented in this manuscript indicates that compliance of the cardiovascular tissue acts as a natural inhibitor of pathological flow conditions conducive to thrombus formation and should not be neglected in computer simulations.
1. Introduction.Left Ventricular Assist Devices (LVADs) are implantable mechanical pumps that temporarily aid the function of the left ventricle as a bridge to recovery or transplant.In some patients they are used as destination therapy.The use of LVADs has been associated with thrombus formation next to the aortic valve and close to the anastomosis region, especially in patients in which the native cardiac function is negligible and the aortic valve remains closed [9,10,14,16,20,23].Stagnation points and recirculation zones have been implicated as the main fluid dynamics factors contributing to thrombus formation [7,12,13,22].The purpose of our study was to develop and use computer simulations based on a fluid-structure interaction solver to investigate formation of stagnation points and recirculation zones near the anastomosis region and the aortic valve for several locations and angles of the anastomisis cannula.In contrast with other works which include only a study of the region near the anastomosis, ignore the aortic valve location, and almost exclusively ignore the compliance of the arterial wall, the present study shows that all these factors are important in the optimal LVAD cannula placement analysis.
A fluid-structure interaction (FSI) solver, originally developed by Badia, Quaini and Quarteroni in [2,4], was extended in the present work to include structures with different elastic wall properties modeling different compliant behaviors of the ascending aorta and of the anastomosis cannula.In this manuscript a semi-implicit monolithic approach is proposed to handle all the interactions at the same time.At all the interfaces (fluid-structure and structure-structure) the coupling conditions are treated in the same way: weak treatment of the continuity of stresses and strong treatment of the continuity of velocity.The monolithic fluid-structure interaction system, obtained after space-time discretization and linearization in a fixed point approach, is preconditioned with a diagonal scaling and an ILUT preconditioner [2].The system is then solved by a Krylov method.This approach is stable and convergent even in the case when the structure and the fluid have comparable densities, which is the case in the blood flow application.
The analysis of flow conditions between the anastomosis and the aortic leaflets suggests that the optimal placement of the anastomosis cannula is the one which is closest to the aortic valve (lower that the usual placement of 5 cm away from the valve region) and sutured at the angle of 30 • to the ascending aorta.This strategy provides the lowest likelihood of stagnant flow formation near the (closed) aortic valve, and the smallest deviation from laminar flow downstream from the anastomosis.Furthermore, it was shown that the rigidity of the dacron anastomosis cannula plays almost no role in generating pathological conditions downstream from the anastomosis.Additionally, the flow analysis presented in this manuscript indicates that compliance of the cardiovascular tissue acts as a natural inhibitor of pathological flow conditions conducive to thrombus formation and should not be neglected in the computer simulations.
Details are presented in the following sections.The numerical method is described in Section 2. Section 3 discusses the LVAD ascending aortic anastomosis problem.Sections 4 and 5 discuss numerical findings regarding the influence of the placement of LVAD cannula on the blood flow conditions quantified by the introduction of the vorticity index and the stagnation index.Section 6 contains the conclusions.
2. The Numerical Method.In this work a semi-implicit monolithic method is proposed for the solution of fluid-structure interaction problems in the presence of structures with different material properties.The basic features of the proposed approach concern the mathematical formulation of fluid and structure problems in terms of velocities, the use of a single finite element partition for the whole domain and the use of the same finite element space for fluid and structure subproblems.At both fluid-structure and structure-structure interfaces the subproblems are coupled by two transmission conditions: continuity of velocity (kinematic condition) and continuity of stress (dynamic condition).At all the interfaces those conditions are imposed in the same way: the dynamic coupling condition is satisfied in the weak sense, and the kinematic one in the strong sense.
Nonlinearities in the resulting coupled problem are dealt via a fixed point approach in an explicit way.This means that only one iteration of the linearization technique is performed at every time step.It was shown in [17] that nonlinear iterations can be disregarded without endangering stability thereby allowing important computational time savings.A FSI algorithm in which the fluid-structure coupling is treated implicitly and the nonlinearities are treated explicitly is said to be semi-implicit.
The solution of the resulting linearized monolithic system is preconditioned in two steps.First a diagonal scaling of the FSI matrix is performed, and then the GMRES method is used to solve the resulting system which is preconditioned by an incomplete LU factorization of the diagonally scaled matrix (the ILUT preconditioner).It was shown in [2] that the ILUT-GMRES method performs well even in the case when the density of the fluid is close to that of the structure, which is the case in the blood flow applications.
A detailed description of the method is presented next.1).The initial configuration Ω 0 at t = 0 is considered as the reference one.We are interested in the flow of an incompressible, viscous, Newtonian fluid in a domain bounded by linearly elastic structures.
The fluid problem is governed by the incompressible Navier-Stokes equations , where u is the fluid velocity, σ f the Cauchy stress tensor and f f the body force.For Newtonian fluids σ f has the following expression where p is the pressure, µ is the fluid viscosity, and is the strain rate tensor, with ∇ denoting the spatial gradient operator.The mechanical behavior of the i-th structure is governed by the elastodynamics equation where η i is the structure displacement, f i the body force, and D t denotes the material derivative.Equation ( 1) must be supplemented with a constitutive law that relates Cauchy stress tensor σ i with strain ǫ.As a simple example, in our numerical simulations we adopt the linear Saint-Venant Kirchhoff elastic model (isotropic, hyperelastic material), for which the consitutive law reads as follows: where µ i ℓ and λ i ℓ are the Lamé constants of the i-th structure.We employ the linear strain-displacement relationship given by ǫ(η i ) = (∇η i + (∇η i ) T )/2 which then leads to the Cauchy stress-displacement relationship of the following form It will be assumed in this manuscript that the mechanical properties of the structures differ only by the values of the coefficients ρ i , µ i ℓ , and λ i ℓ .Of course, other structure models can be chosen according to the specific problem under consideration.
At Σ f i t , the fluid and structure subproblems are coupled via the kinematic and dynamic coupling conditions: the continuity of velocity (kinematic coupling condition) u = ∂ t η i on Σ f i t × (0, T ), i = 1, ..., N s and balance of stresses (dynamic coupling condition): Here, n f denotes the outward normal to Ω f t , and n i denotes the outward normal to the i-th structure domain boundary.The same contact conditions hold at the interface between two structures: To capture the time-evolution of the whole domain Ω t an Arbitrary Lagrangian-Eulerian (ALE) approach is considered.More precisely, define the following two families of mappings: and The map L i t = L i (•, t) tracks the i-th structure domain in time while A t = A(•, t) tracks the fluid domain.The two must agree on the interfaces: See Figure 1.The mapping L i t is defined (in Lagrangian framework) as: where the "hat" denotes a variable evaluated at the reference configuration.The mapping A t can be defined as an appropriate extension operator to Ω f t of its values on the interfaces: A classical choice for the extension operator is to consider a harmonic extension in the reference domain.A t is called the Arbitrary Lagrangian-Eulerian (ALE) mapping [11].The ALE derivative of the fluid velocity is then defined as: where the domain velocity w is calculated using the following expression: The fluid-structure interaction problem with multiple structures in strong form is defined as follows: (i) Geometry: find the fluid domain displacement so that (ii) Fluid-structure itneraction: find velocity u, pressure p and displacement η i , i = 1, ..., N s , such that Problem ( 6) needs to be supplemented by initial and boundary condition.Since the particular form of the initial and boundary conditions does not influence the rest of the discussion in this section, it will be omitted.The initial and boundary conditions for the specific problem discussed in this manuscript will be specified in Section 3.

Weak formulation.
To define the variational formulation of the fluid-structure interaction problem (5)-( 6), denote with L 2 (Ω) the space of square integrable functions in Ω and with H 1 (Ω) the space of functions in L 2 (Ω) with first derivatives in L 2 (Ω).Furthermore, denote by (•, •) Ω and •, • Ω the L 2 -product and a duality pairing in L 2 , respectively.For simplicity, let the number of different structures equal to 2, N s = 2 (see Figure 1).The following spaces will be used bellow: Introduce the following notation: The weak formulation of the fluid-multi structure interaction problem is then given by: find (u, p, η1 , η2 ) for all (v f 0 , q, v1 , v2 The continuity of velocities at both fluid-structure and structure-structure interfaces has been enforced in a strong way by (7d) and (7e), respectively.On the other hand, the continuity of stresses at the interfaces is imposed in a weak way via (7a)-(7c).

The fully discrete problem: space and time discretization. Let
and V 2 0 at the reference configuration, respectively.With a slight abuse of notation denote by ) the finite element spaces for a given time step t n defined via the domain mappings (2)- (3).
It is well known that the velocity-pressure finite element space (Q h , V f h ) is required to satisfy a discrete inf-sup condition in order to obtain pressure stability, see [6].Thus, either inf-sup stable finite elements need to be used or a stabilization is required.In this work, a stabilized formulation, called orthogonal subgrid scales technique [8], is used to overcome both pressure instabilities and instabilities arising when the convective term is dominant.The stabilized version of the fluid problem is then obtained by replacing N in (7a) with where the perturbation term introduced by the orthogonal subgrid scales technique (in its quasi-static form) reads Here, Π ⊥ (•) is the L 2 orthogonal projection onto the finite element space, i. e.: , where Π(•) is the L 2 projection onto the finite element space and I(•) the identity operator.See [8] for the computation of the stabilization parameters τ 1 and τ 2 .The time discretization is obtained using the backward Euler scheme for both the fluid and the structure equations.By defining the backward Euler operator δ t as δ t f n+1 = (f n+1 − f n )/δt and denoting by Ext h (•) a discretized version of the extension operator Ext(•), at each time level t n+1 , the fully discretized fluidmulti structure interaction problem reads: (i) Geometry: find the fluid domain displacement so that (ii) Fluid-structure interaction: find (u (10a) 2.4.Solution of the monolithic system.To deal with the nonlinearities in the problem which are due to domain motion and to the convective term in (10a), a fixed point algorithm is employed.The linearization of the fluid-structure interaction problem ( 9)-( 10) by the fixed point algorithm uses η n+1 1,h and η n+1 2,h (predictions for η n+1 1,h and η n+1 2,h ) to compute Ω f t n+1 and u n+1 h (predictions for u n+1 h ) as convective velocity.More details can be found in [3].The monolithic linearized system which results at every fixed point iteration is then solved by an extension of the ILUT-GMRES method proposed in [2].
We rely on a single finite element partition of the overall domain, which implies matching grids on the fluid-structure and structure-structure interfaces.The structure equations are reformulated in terms of velocity instead of displacement in order to define the fluid-structure interaction problem on the entire domain in terms of velocity.The same finite element spaces are then used for fluid velocity and structure velocity.Moreover, since a stabilization technique is used to circumvent the discrete inf-sup condition, the same finite element interpolation space can be used for the pressure unknowns as well.In particular, for the numerical simulation of the problem discussed in Section 3, we used P 1 − P 1 finite elements [19] for the fluid and P 1 finite elements for the structures.
In this framework, the fluid-structure and structure-structure coupling conditions (the transmission conditions) are easily implemented.The kinematic coupling condition describing continuity of velocities on each interface is implicitly enforced by the finite element space interpolation used over the whole domain.In particular, we adopt continuous finite element spaces, so the continuity of velocities at the interface is enforced automatically by the fact that there is only one value of the velocity at the interface nodes.The dynamic coupling condition describing continuity of stresses is imposed weakly.The weak transmission of stresses simply arises from the fact that shape functions on the interface nodes have support on the fluid and structure subdomains (fluid-structure interface) or on both structure subdomains (structure-structure interface).This way the final system has the simple form similar to the one described in [3,2].
Standard preconditioned iterative methods were then used to solve the resulting monolithic system.In particular, we used the GMRES algorithm [21] preconditioned in two steps: first a diagonal scaling is applied and then the resulting system is preconditioned by an incomplete LU factorization (ILUT preconditioner).See [2] for details.
The numerical method presented in this section will be used in the next section to analyze the flow conditions associated with the ascending aortic anastomosis of LVAD cannula.

Numerical Characterization of Flow Conditions in LVAD Cannula
Anastomosis.Left Ventricle Assist Devices (LVADs) are implantable mechanical pumps that temporarily aid the function of the left ventricle.An LVAD is typically connected to the left ventricle apex and then anastomosed to the ascending aorta, as shown in Figure 2. Implantation of LVADs has been associated with thrombus formation close to the anastomosis region in the ascending aorta and next to the aortic valve, especially if native cardiac function is negligible and the aortic valve remains closed throughout cardiac cycle [10].Stagnation points and recirculation zones have been implicated as the main fluid dynamics factors contributing to thrombus formation.In the present work we employ the fluid-structure interaction algorithm described in Section 2.4 to simulate the fluid dynamics conditions in the anastomosis cannula and in the ascending aorta between the aortic valve and the anastomosis region.Clinicians have been trying to understand the optimal location of the cannula suture with respect to the aortic valve location that would minimize complications associated with thrombus formation, especially near the aortic valve leaflets.
In this work we investigate the stagnation and recirculation zones quantified by the stagnation index and vorticity index defined below.Of particular interest are the flow conditions generated by the cannula placement at 4, 5 and 6 cm from the aortic valve (clinical placement is typically 5 cm from the aortic valve), and with the cannula sutured to the ascending aorta at the angles of 30 • and 90 • (clinical suture is typically performed at 90 • ).Additionally, we have investigated how inclusion of the compliance of ascending aorta and a different compliant behavior of the anastomosis cannula influence the fluid flow simulations, stagnation zones and vorticity.
3.1.Problem set up.A realistic three-dimensional simulation of the fluid dynamics at the LVAD ascending aortic anastomosis with compliant aorta, LVAD cannula, and functional aortic leaflets, is exceedingly complicated.Often times various simplifying assumptions are made in order to obtain an insight into the corresponding fluid dynamics in a reasonable time frame.Although these simplifying assumptions bring a certain degree of idealization to the problem, the flow analysis of the simplified problem still provides important information about the main qualitative features of the underlying problem.
In this work we performed two-dimensional fluid flow simulations for the geometry corresponding to a section of the anastomosis cannula and the region of the ascending aorta from the closed aortic valve to an area downstream from the anastomosis, as shown in Figure 3(d).Three different combinations of structures compliant properties were considered: The 2d domain assumes an idealization of the aortic geometry and of the aortic valve, ignoring the tricuspid geometry of the aortic valve leaflets.The 2d Navier-Stokes equations for an incompressible, viscous Newtonian fluid were coupled with the following equation of linear elasticity to describe the motion of the compliant boundaries This structure model can be obtained from the equations of 3D linear elasticity assuming cylindrical geometry.The term which accounts for the effects of cylindrical geometry is the reaction term a i η in (11), where 2 , E i being the Young modulus, ν i the Poisson ratio of the i-th structure, and R the distance of the structure from the axis of symmetry.This term plays an important role in recovering the physiologically reasonable pressures when elastic structures are considered since it accounts for the circumferential stress coming from the 3d geometry and should therefore not be neglected in the 2d simulations.
In this work the case of severe cardiac failure is considered, where native cardiac function may be negligible, and the aortic valve remains closed as all systemic perfusion is provided by the LVAD.In these conditions thrombus formation is most likely.This is one reason why compliance of the aortic leaflets was not taken into account.Another reason why compliance of the aortic leaflets was not taken into account has to do with dimension reduction from 3d to 2d.As mentioned above, the reaction term in the structure model ( 11) contains a coefficient which becomes singular as R → 0 at the valve commissure where the leaflets meet.Even though there are ways to get around this difficulty in 2d, this was not done in the present paper as we are currently extending our investigation to involve 3d simulations, in which case the singularity in the structure model at the valve commissures is no longer present.Other factors, however, related to 3d FSI simulations of the aortic valve motion make this problem exceedingly complicated and will not be addressed in the present manuscript.The fluid and the structures physical parameters used in this manuscript are listed in Table 1.The LVAD cannula is made of polyvinyl chloride (dacron), whose Young modulus in much higher than that of the aortic wall.Thus, we expect that accounting for the cannula compliance will not make much difference in the fluid dynamics.However, problem (c) is useful to test how the properties of the cannula material affect the flow conditions.The diameter of the LVAD cannula is taken to be 1.2 cm, while the diameter of the aorta is assumed to be 2.5 cm.
At the inlet of the anastomosis cannula the following Dirichlet condition for the velocity is imposed describing time-periodic inflow parallel to the walls of the cannula.Here n and t are the unit vectors normal and tangent to the inflow boundary, respectively (see Figure 3(d)).The magnitude of the velocity has been chosen to recover a flow rate of 5.45 − 5.56 l/min, which is associated with a GYRO flow pump working at 2000 rpm [5].Moreover, with this inlet boundary condition a Reynolds number (based on the diameter of the aorta) of 4200 was recovered, which is reasonable for blood flow in the ascending aorta.The algorithm presented in Section 2 with the appropriate choice of discretization performs well for such Reynolds numbers as was shown in [18] where a classical test problem of a wall-driven cavity flow at Reynolds numbers up to 10 4 [1] was studied.Notice that condition ( 12) is periodic with period 1 s (we are assuming a pulsatile pump), which is approximately the duration of a cardiac cycle.On the outflow section we impose a non-homogeneous Neumann condition: this value being physiological.Notice that the use of this outlet condition does not imply constant pressure in the aorta, as the pressure is just one component of the normal stress.The pressure itself will be time-dependent as a consequence of boundary condition (12).In addition, since the valve is not displacing, a homogeneous Dirichlet condition is imposed on it.Whenever a compliant structure is present, its displacement is kept fixed at the inlet and at the outlet (proximal and distal ends), while a homogeneous Neumann condition is imposed on the rest of the physical boundary.The condition of fixed structure ends has been chosen for simplicity, even though it is unphysical and leads to small oscillations due to a boundary layer formed in a small neighborhood of the end points.More realistic structure boundary conditions will be studied in the future.
For the space discretization of the fluid we use stabilized P 1 − P 1 finite elements.Whenever one or more structures are present, matching meshes were implemented at the interfaces and P 1 finite elements were used for space discretization.The mesh of the fluid domain is the same for problems (a), (b), and (c), and consists of 6215 nodes and 11237 triangles.For the structure triangulations of problems (b) and (c) a diameter of 0.025 is taken.A time step δt = 0.02 s is set for all the simulations.
The enormous amount of information given by the numerical simulations needs to be synthesized in indexes, which aim at measuring some quantity of interest.For the problems under consideration we propose a downstream zone vorticity index Table 1.Fluid and structure physical properties for the numerical tests Fluid Density: ρ f = 1.0 g/cm 3  Viscosity: µ = 0.035 poise Aortic wall Density: ρ s = 1.2 g/cm 3  Thickness: h s = 0.1 cm Young modulus: E = 6.5 • 10 7 dyne/cm 2 Poisson ratio: ν = 0.49 LVAD cannula Density: ρ s = 1.37 g/cm 3  Thickness: h s = 0.1 cm Young modulus: E = 2.2 • 10 10 dyne/cm 2 Poisson ratio: ν = 0.43 to measure the mean vorticity downstream from the anastomosis in the region indicated by Ω d in Figure 3(d), and a valve zone stagnation index to quantify the region of stagnant flow close to the valve (denoted by Ω v in Figure 3(d)).Index I s represents the mean velocity over Ω v .Index (13) needs to be minimized in order to have a smooth laminar in the ascending aorta.On the other hand, index (14) has to be as high as possible to avoid flow stagnation, which is associated with clot formation.

Significance of wall compliance.
In the first set of simulations we have investigated the importance of wall compliance on downstream vorticity and on flow stagnation in the aortic valve region.We assume that the anastomosis cannula is inserted into the aorta at a distance of 5 cm from the proximal end, and at an angle of 30 • (see Figure 3(d)).It was shown in [15] that such an angle reduces local flow recirculation and high shear stress on the aortic wall, typical of high angle of insertion.
We have generated the flow conditions corresponding to the 3 different structure configurations: (a), (b), and (c), stated in the previous section.Fig. 4 shows the downstream zone vorticity index I ω and the valve zone stagnation index I s as functions of time, for problems (a), (b), and (c).Fig. 4(a) shows that I ω oscillates around a mean value, with different amplitudes and phases if the aortic wall and cannula are modeled as rigid or if their elasticity is taken into consideration.The mean value of I ω is 75.3 for problem (a), 73.7 for problem (b) and 71.9 for problem (c).Moreover, I ω for problem (a) displays oscillations of larger amplitude.This means that the compliance of the structures surrounding the fluid naturally tends to reduce the magnitude of the vorticity.
Similarly, the values of I s , shown in Figure 4(b), change for the three problems, depending on whether or not the compliance of the walls has been taken into account.The results shown in Figure 4(b) show that the compliance of the aortic wall and cannula increase the mean fluid velocity near the aortic valve, i.e., over Ω v , for the whole duration of a cardiac cycle.
This shows that compliance of the cardiovascular tissue acts as a natural inhibitor of pathological flow conditions conducive to thrombus formation and should not be neglected in the computer simulations.In the rest of the paper all the simulations will account for the compliance of both the aortic wall and the cannula.4. Sensitivity of the valve stagnation index to the position of the anastomosis, material of the cannula and insertion angle.Sensitivity of the stagnation index I s to the variations in the position of the anastomosis, material of the cannula, and angle of insertion is considered.Three positions of the anastomosis: at a distance of 4, 5, and 6 cm from the aortic valve; two materials for the cannula: dacron and a material with the same properties as the aortic wall (see Table 1); and two angles of insertion: 30 • and 90 • , were tested.
Figure 5 shows the graphs of the stagnation index I s for the three different positions and for the two different materials.The angle of insertion is kept at 30 • .The graphs in Figure 5 show that the lower position of the cannula anastomosis guarantees a value of I s above 1 cm/s for the entire cardiac cycle.This is due to the fact that the vortex generated by the blood stream coming from the cannula generates secondary vortices in the aortic sinuses (see Fig. 6 (a)).If the anastomosis is placed farther from the valve, the secondary vortices move progressively out of the sinuses (see Fig. 6 (b)).In the case of anastomosis located 6 cm away from the valve (upper anastomosis), shown in Fig. 6 (c), only the tertiary vortices enter the sinuses.The graphs in Figure 5 show that the insertion of the anastomosis cannula closest to the aortic valve provides the lowest likelihood for the formation of stagnant flow near the closed aortic valve, as expected.
Next, the sensitivity of I s to the cannula material was studied.It was found that the cannula made of a material with the same properties (density, Youngs modulus and Poisson ratio) as the aortic wall, instead of a dacron cannula, does not significantly reduce the likelihood of stagnant flow formation near the aortic valve.Fig. 5 shows that the values of I s associated with the softer cannula are comparable with those obtained for the dacron cannula for the upper and middle positions, while whey are smaller for the lower position.Therefore, the use of softer  cannula instead of the dacron one does not improve the behavior of the stagnation index I s near a closed aortic valve.Finally, the sensitivity of I s on the cannula insertion angle was studied.Fig. 7 shows I s for different angles of insertions and different materials for the cannula for all three anastomosis positions (4, 5, and 6 cm from the aortic valve).Fig. 7 shows that the angle of insertion of 90 • significantly lowers the values of I s , regardless of the cannula material.It can be seen in Fig. 8(a) that this configuration distorts the primary vortex at the anastomosis which then, in turn, distorts the secondary vortices.The secondary vortex close to the anastomosis is "pulled out" of the aortic sinus, leaving a large region of stagnant flow behind.The distortion of the vortices is clear when Fig. 8(a) is compared with Fig. 6(a).Fig. 7(b) also shows that for the middle anastomosis the angle of insertion of 90 • causes a reduction of I s , both for dacron and soft cannula.Notice that the curves of I s for the 90 • cannula in the middle position are slightly higher than the corresponding ones in the lower position.This is because the vortices get less and less distorted as the cannula, placed at 90 • , is moved up to the middle location in the ascending aorta (see Fig. 8).For the upper position of the anastomosis the values of I s , reported in Figure 7(c), are slightly higher for the whole simulated cardiac cycle when the angle of insertion is 90 • , independently of the cannula material.
Based on the results from Fig. 7 the optimal location for the anastomosis cannula that lowers the likelihood of stagnant flow near the aortic valve is the one closest to the aortic valve, inserted at the angle of 30 • with respect to the native aorta.The graphs in Fig. 7 also indicate that the stagnation index I s is not significantly influenced by the properties of the cannula material.5. Sensitivity of the downstream zone vorticity index to the position of the anastomosis, material of the cannula and insertion angle.A sensitivity analysis to the variations in the position of the anastomosis, material of the cannula, and angle of insertion was performed for the vorticity index I ω .We tested the same three positions of the anastomosis, two materials for the cannula, and two angles of insertion considered in Sec. 4 for the sensitivity index.
Figures 9(a) and (b) show the graphs of the vorticity index I ω for the three different positions with dacron cannula and soft cannula, respectively.The angle of insertion is kept at 30 • .The mean value of I ω for every configuration considered in Fig. 9 is reported in Table 2. Fig. 9 and Table 2 suggest that by moving the anastomosis farther from the valve, the magnitude of the vorticity is reduced and a flow closer to laminar is obtained.However, notice that the vorticity index (Eq.( 13)) is defined in the region downstream of the anastomosis all the way to the outlet flow boundary (see Fig. 3(d)).This means that in the case when the anastomosis is further away from the outlet boundary, the flow has more room to develop, showing secondary and tertiary vertices downstream from the anastomosis, and less influence of the outlet boundary on the flow conditions in the vicinity of the anastomosis.See Figure 10.A larger downstream flow region should be included for a thorough investigation of the vorticity influenced by the different positions of the cannula, which is beyond the scope of this manuscript.The sensitivity of the vorticity index on the cannula material and on the angle of insertion are much less sensitive to the outlet flow boundary, and are presented next.Concerning the sensitivity of I ω to the cannula material, it was found that the cannula made of a material with the same properties (density, Youngs modulus and Poisson ratio) as the aortic wall, instead of a dacron cannula, does not significantly reduce vorticity.Table 2 shows that the mean values of I ω associated with the softer cannula are comparable with those obtained for the dacron cannula for all the positions of the anastomosis.Moreover, from Fig. 9 we see that I ω displays oscillations of similar amplitude and only a small phase lag (oscillations are retarded by the softer material) when the position of the anastomosis is the same and the cannula material is different.The phase lag gets smaller as the anastomosis approaches the valve.We can conclude that the use of a softer cannula instead of the dacron one does not improve the quality of the flow in the aorta.
Finally, the sensitivity of I ω on the cannula insertion angle was studied.Fig. 11 shows I ω for different angles of insertions and different cannula materials for all three anastomosis positions (4, 5, and 6 cm from the aortic valve).Fig. 11 shows that the angle of insertion of 90 • significantly increases the values of I ω , regardless of the cannula material.The mean value of I ω for each configuration with angle of insertion of 90 • is reported in Table 3.Every value in Table 3 is between 35%  and 55% higher than the correspondent value in Table 2, proving that increasing the angle of insertion from 30 • to 90 • worsen the quality of the flow in the aorta.Notice from Table 3 that the mean value of I ω is almost independent of the position of the anastomosis, suggesting that the vorticity distribution does not change much with the distance from the anastomosis.The vorticity filed at a given time level for the three configurations at 90 • is shown in Fig. 12.
Based on the results from the sensitivity analysis of the vorticity index, the optimal configuration found in Sec. 4 which minimizes the risk of thrombus formation  based on the stagnation index values (lower anastomosis, 30 • angle of insertion) is found to also provide the lowest vorticity formation in the flow downstream from the anastomosis.In addition, our analysis shows that the vorticity index I ω , like the stagnation index, is not significantly influenced by the properties of the cannula material.
6. Conclusions.We have implemented a novel fluid-structure interaction algorithm to study flow conditions associated with LVAD ascending aortic anastomosis.The proposed numerical algorithm is capable of capturing different elastic properties of the native aortic wall and of the LVAD cannula.The numerical method is based on a monolithic approach where fluid-structure and structure-structure interaction are treated in the same way: continuity of stresses is imposed in a weak sense while the continuity of velocity is imposed in the strong sense.Flow conditions between the anastomosis location and the closed aortic valve exhibiting no function were studied.Two indexes were proposed to quantify the deviation from laminar flow downstream from anastomosis (the vorticity index) and to quantify the severity of stagnant flow near the aortic leaflets (the stagnation The study presented in this manuscript found that the cannula placement closest to the aortic valve (within 4 cm away from the valve) and at the angle of 30 • with respect to the native aorta close to the aortic valve, gives rise to the flow conditions which are least likely to lead to thrombus formation.Rigidity of the dacron cannula was found to play almost no role in generating pathological conditions downstream from the anastomosis.Compliance of the native aortic wall was found to be important in naturally reducing the magnitude of vorticity, thereby lowering the probability of thrombus formation.
In contrast with other works in this area which include a study of flow conditions near anastomosis, ignore the aortic valve location, and almost excusively ignore the compliance of the aortic wall, the present study shows that all these factors are important in the analysis of optimal LVAD cannula placement.

Figure 2 .
Figure 2. Left Ventricular Assist Device with anastomosis to ascending aorta.
(a) rigid aortic wall and cannula; (b) compliant aortic wall and rigid cannula; (c) compliant aortic wall and cannula.The corresponding domains are shown in Figure 3(a)-(c).

Figure 3 .
Figure 3. Computational domain for the simulations with (a) rigid boundaries, (b) deformable aortic wall and rigid cannula, (c) deformable aortic wall and cannula.(d) Explanation of the geometry.
(a) downstream zone vorticity index (b) Valve zone stagnation index

Figure 4 .
Figure 4. (a) downstream zone vorticity index I ω and (b) valve zone stagnation index I s for the problem with rigid boundaries (problem (a)), the problem with compliant aortic wall and rigid cannula (problem (b)), and the problem with compliant aortic wall and cannula (problem (c)).

Figure 5 .
Figure 5.Comparison of the valve zone stagnation index I s for three different positions of the anastomosis and two different materials for the cannula.The angle of insertion is 30 • in all the cases.

Figure 6 .
Figure 6.Streamlines at the end of one cardiac cycle for (a) lower anastomosis, (b) middle anastomosis, and (c) upper anastomosis.The angle of insertion is 30 • and the cannula is made of dacron.The streamlines are colored based on the velocity magnitude (in cm/s).

Figure 7 .
Figure 7.Comparison of the valve zone stagnation index I s for two different angles of insertion and two different materials for the cannula for (a) lower anastomosis, (b) middle anastomosis, and (c) upper anastomosis.The legend in (c) is common to the three subfigures.

Figure 8 .
Figure 8. Streamlines at the end of one cardiac cycle for (a) lower anastomosis, (b) middle anastomosis, and (c) upper anastomosis.The angle of insertion is 90 • and the cannula is made of dacron.The streamlines are colored based on the velocity magnitude (in cm/s).

Figure 9 .
Figure 9.Comparison of the downstream zone vorticity index I ω for three different positions of the anastomosis with (a) dacron cannula and (b) soft cannula.The angle of insertion is 30 • in all the cases.

Figure 10 .
Figure 10.Vorticity magnitude at the end of one cardiac cycle for (a) lower anastomosis, (b) middle anastomosis, and (c) upper anastomosis.The angle of insertion is 30 • and the cannula is made of dacron.

Figure 11 .
Figure 11.Comparison of the dowstream zone vorticity index I ω for two different angles of insertion and two different materials for the cannula for (a) lower anastomosis, (b) middle anastomosis, and (c) upper anastomosis.The legend in (c) is common to the three subfigures.

Figure 12 .
Figure 12.Vorticity magnitudes at the end of one cardiac cycle for (a) lower anastomosis, (b) middle anastomosis, and (c) upper anastomosis.The angle of insertion is 90 • and the cannula is made of dacron.

Table 2 .
Mean values of I ω for three different positions of the anastomosis and two cannula material with angle of insertion of 30 • .

Table 3 .
Mean values of I ω for three different positions of the anastomosis and two cannula material with angle of insertion of 90 • .