Numerical Study of the Blood Flow in a Deformable Human Aorta

Featured Application: This article developed a computer simulation based on the finite element method to investigate the interaction between blood flow and arterial wall deformation. The numerical results show that the blood flow has an unsteady and complex behavior. The velocity field, pressure, and displacement distributions were investigated during a cardiac cycle. The method described in this study will be useful for investigating cardiovascular diseases and performing hemodynamics analyses. Abstract: In this work, we present a numerical investigation of blood ﬂow in a portion of the human vascular system. More precisely, the present work analyzed the blood ﬂow in the upper portion of the aorta. The aorta and its ramiﬁed blood vessels are surrounded by the cardiac muscle. The blood ﬂow generates pressure on the internal surfaces of the artery and its ramiﬁcations, thereby causing deformation of the cardiac muscle. The numerical analysis used the Navier–Stokes equations as the governing equations of blood ﬂow for the calculation of the velocity ﬁeld and pressure distribution in the blood. The neo-Hookean hyperelastic model was used for the description of the behavior of the vessel walls. The velocity and pressure distributions were analyzed. The deformation of the vessel was also investigated. The numerical results could be used to better understand and predict the factors that trigger cardiovascular diseases and distortions of the aorta and as a diagnostic tool in clinical applications.


Introduction
The cardiovascular system consists of the heart and vessels (arteries and veins). It has the function to distribute to the organs oxygen and nutrients essential for their life through the blood. The heart and its various organs can be affected by different pathologies, called cardiovascular diseases, namely, stenosis, aneurysm, occlusion, atherosclerosis, etc. According to recent statistics, cardiovascular diseases are the third leading cause of death in the western world [1]. The understanding of the vascular function requires knowledge of the behavior and the material properties of blood vessels [2,3]. In addition, the deep characterization of blood flow is crucial for the understanding of the behavior of the whole human system and for the improvement of human health [4,5]. Even though a few experimental works are available on the histological structure of arteries, they are still limited by the difficulty of making rapid analyses on the heart [6,7].
Numerical simulation is a widely used tool for understanding many phenomena and especially, it was confirmed to be useful in studying the blood flow in the cardiovascular system, to diagnose cardiovascular diseases [8][9][10][11][12], evaluate the hemodynamics of healthy and diseased blood vessels [13][14][15] and design and optimize vascular medical devices [16,17]. Over the last decades, various computational techniques have been used extensively by a large number of researchers seeking to simulate and analyze blood flow and to understand the relationship between vascular diseases and hemodynamics [18,19]. In addition, numerical simulation and computer models provide a fast, simple, safe way to study phenomena involved in vascular diseases.
Recently, there have been numerous studies of blood flow related to the aorta that use computer simulation to better investigate blood flow variables. Panta et al. [20] analyzed, through their numerical simulation, blood velocity as a function of the shape of human arteries, namely straight, bent, T-shaped, and the main arterial network models. Wood et al. [21] combined magnetic resonance imaging (MRI) and computational fluid dynamic (CFD) simulation of flow to analyze the human descending aorta. Their objective was to investigate the potential of combining CFD simulations with in vivo MRI scans. Tokuda et al. [22] simulated blood flow numerically using the finite element method to cognize the mechanism of stroke during cardio-pulmonary bypass (CPB). The location of turbulence, wall pressure and flow distribution within the aortic arch were investigated. They succeeded in simulating and visualizing the blood flow during CPB using CFD techniques. Torii et al. [23] developed a computer modeling technique based on the deforming-spatial-domain/stabilized space-time method to study fluid-structure interactions between flowing blood and arterial walls. They demonstrated that the wall shear stress is a very important factor in cardiovascular diseases, and their results proved that numerical modeling is an operative tool in hemodynamic studies of cardiovascular diseases.
Bazilevs et al. [24] presented a computational framework for the simulation of patient-specific cerebral aneurysm configurations and also described a new approach for the computation of blood vessel tissue prestress. They studied some variables, namely, hemodynamics, wall shear stress, and wall tension, in order to check the relevance of fluid-structure interactions when adapting the rigidity of the arterial wall and the relevance of a flexible wall. Crosetto et al. [25] performed a numerical algorithm based on finite elements to evaluate the mean velocity and temporal evolution of pressure and the distribution of wall shear stress on diverse sections of the aorta. More specifically, they reported the results of 3D simulations using a physiological geometry over several heartbeats. They changed the inlet condition of flux by pressure which caused a backward diastolic flow through the section and thereafter led to the closure of the aortic valve and the circulation in the coronary arteries. However, these observations were only predicted numerically. Garje et al. [26] studied numerically the changes in blood flow vessels and the growth of a blockage, which is an occlusion of the coronary arteries typical of coronary artery disease (CAD), due to an accumulation of cholesterol. Their work was focused on the size and position of the blockage. Priyadharshini and Ponalagusamy [27] studied numerically the effects of magneto-hydrodynamics on blood flow parameters in a stenotic artery under an applied external magnetic field and body acceleration.
In this work, we intended to simulate the interaction of blood flow with the artery wall. We simulated the pressure along the surface of the artery, the displacements and stresses of the walls of part of a child's aorta and some of its ramifications embedded in the heart muscle when subjected to the pressure caused by heart beating.
In Section 2, the geometry of the aorta and its ramifications are presented. Then, the Navier-Stokes equations and the equation governing the displacement of the blood vessels are introduced. The initial and boundary conditions, as well as the numerical method, are also presented. The results are discussed in Section 3. Figure 1 shows the simulated geometry. The domain was divided into three parts: blood, aorta and its ramifications (branching vessels), and cardiac muscle. The cardiac muscle is very important because it protects the artery from the deformation due to the pressure. The inlet and outlet orifices are also shown in this figure.

Simulated Geometry
Appl. Sci. 2019, 9, x FOR PEER REVIEW 3 of 12 Figure 1 shows the simulated geometry. The domain was divided into three parts: blood, aorta and its ramifications (branching vessels), and cardiac muscle. The cardiac muscle is very important because it protects the artery from the deformation due to the pressure. The inlet and outlet orifices are also shown in this figure.

Governing Equations
In this study, we aimed to compute the flow of blood in the aorta and to investigate the deformation of the cardiac muscle due to blood pressure and the behavior of the vessel walls. The numerical study consisted of two parts. Navier-Stokes equations: The blood flow fluid dynamics was described using the timedependent Navier-Stokes equations and the equations for the conservation of mass. The fluid is assumed to be incompressible and homogeneous, with a Newtonian behavior. The blood was assumed to be Newtonian in this study, as its viscosity is almost constant in arteries with relatively large diameters (5 mm); therefore, non-Newtonian effects were neglected [28].
The Navier-Stokes equations may be written as follows: where is the blood density, is the blood velocity, and is the stress tensor [29]. The relationship describing the stress tensor and the strain rate tensor is defined as: where is the blood pressure, is the dynamic viscosity of the blood, is the identity tensor, and denotes the transpose tensor.
The fluid properties such as the density and the dynamic viscosity were assumed to be constant. Their values were, respectively = 1060 / and = 0.005 . .

Governing Equations
In this study, we aimed to compute the flow of blood in the aorta and to investigate the deformation of the cardiac muscle due to blood pressure and the behavior of the vessel walls. The numerical study consisted of two parts.

Part 1:
Fluid dynamics analysis which required to solve the Navier-Stokes equations. It included the calculation of the velocity field and the distribution of blood pressure (variable over time and in space). Part 2: Mechanical analysis of the deformation of the cardiac tissue and artery.

Navier-Stokes equations:
The blood flow fluid dynamics was described using the time-dependent Navier-Stokes equations and the equations for the conservation of mass. The fluid is assumed to be incompressible and homogeneous, with a Newtonian behavior. The blood was assumed to be Newtonian in this study, as its viscosity is almost constant in arteries with relatively large diameters (5 mm); therefore, non-Newtonian effects were neglected [28].
The Navier-Stokes equations may be written as follows: where ρ is the blood density, u is the blood velocity, and σ is the stress tensor [29]. The relationship describing the stress tensor and the strain rate tensor is defined as: where p is the blood pressure, µ is the dynamic viscosity of the blood, I is the identity tensor, and T denotes the transpose tensor.
The fluid properties such as the density and the dynamic viscosity were assumed to be constant. Their values were, respectively ρ blood = 1060 kg/m 3 and µ blood = 0.005 Pa.s.

Structural mechanics:
The motion of the vessel wall can be mathematically described by the following equation [30]: where ρ w is the wall density, f s represents the external body forces acting on the solid, and σ s is the Cauchy stress tensor. The wall can be the artery or the cardiac muscle. The value of the wall density is ρ artery = 960 kg/m 3 for the artery and ρ muscle = 1200 kg/m 3 for the cardiac muscle. It should be noted that u is a velocity (unit m/s); however, u solid is a displacement (unit meter). The artery is assumed to be subjected to large strains and undergo large deformations. Its behavior is non-linear. More precisely, it behaves as a hyperelastic neo-Hookean material. The constitutive law of this material involves the material shear modulus µ s and the bulk modulus κ s . The values of these parameters are provided in Table 1. In the present study, the material is assumed to be isothermal, and no effect of the temperature is considered. For a three-dimensional, compressible, neo-Hookean material such as the artery, the strain energy density is given by [24]: where J = detF, F is the deformation gradient tensor, and C = F T F is the Cauchy-Green deformation gradient tensor. From the above constitutive law, the expressions of various tensors can be deduced by means of more or less arduous calculations.

Initial and Boundary Conditions
As mentioned above, the model uses three domains, namely, the blood, the artery and the cardiac muscle. The numerical simulation considers the solution of the 3D Navier-Stokes equations and the fluid-structure interactions. The initial and boundary conditions required to solve the equations of the model are given hereafter.
The human heart beat is cyclically time-varying. Thus, the velocity at the inlet of the aorta is not set to be constant. It has a periodic profile as a function of time. During a cycle, the velocity varies between a minimal and a maximal value. The duration of each period is 0.8 s. We introduced the following function [31]: The inlet velocity for the artery in this simulation was fixed to 0.3 [m/s]*f(t). The temporal evolution of the inlet velocity is represented in Figure 2. The outlet condition pressure was considered to be 125 [mmHg]*f(t). These conditions are similar to those adopted by Choudhari and Panse [31]. This function has two peaks which have different amplitudes and different periods.
On the aortic walls, a non-slip boundary condition is applied, and the wall is assumed to be impermeable since the fluid is supposed not to penetrate through the walls as it flows along the vessel. On the aortic walls, a non-slip boundary condition is applied, and the wall is assumed to be impermeable since the fluid is supposed not to penetrate through the walls as it flows along the vessel.

Numerical Method
In this study, the finite element method was used to solve the set equations describing the problem [32]. The simulated geometry was meshed in tetrahedral discretization as shown in Figure  3a. Figure 3b depicts the radial discretization of the blood aorta. The central part stands for the blood. The artery and its thickness appear clearly in this figure. However, the muscle was considered as having no thickness, and only a boundary condition was applied to represent this surface. Indeed, the muscle was considered a free surface, with no external force applied on it.
The model used to solve the governing equations is two-way coupled. Actually, any deformation of blood vessel walls has an effect on the fluid flow and vice versa [33,34]. Mesh sizes are summarized in Table 2. The total number of elements were 195,050, and the total number of nodes was 33,506.

Numerical Method
In this study, the finite element method was used to solve the set equations describing the problem [32]. The simulated geometry was meshed in tetrahedral discretization as shown in Figure 3a. Figure 3b depicts the radial discretization of the blood aorta. The central part stands for the blood. The artery and its thickness appear clearly in this figure. However, the muscle was considered as having no thickness, and only a boundary condition was applied to represent this surface. Indeed, the muscle was considered a free surface, with no external force applied on it.
The model used to solve the governing equations is two-way coupled. Actually, any deformation of blood vessel walls has an effect on the fluid flow and vice versa [33,34]. On the aortic walls, a non-slip boundary condition is applied, and the wall is assumed to be impermeable since the fluid is supposed not to penetrate through the walls as it flows along the vessel.

Numerical Method
In this study, the finite element method was used to solve the set equations describing the problem [32]. The simulated geometry was meshed in tetrahedral discretization as shown in Figure  3a. Figure 3b depicts the radial discretization of the blood aorta. The central part stands for the blood. The artery and its thickness appear clearly in this figure. However, the muscle was considered as having no thickness, and only a boundary condition was applied to represent this surface. Indeed, the muscle was considered a free surface, with no external force applied on it.
The model used to solve the governing equations is two-way coupled. Actually, any deformation of blood vessel walls has an effect on the fluid flow and vice versa [33,34]. Mesh sizes are summarized in Table 2. The total number of elements were 195,050, and the total number of nodes was 33,506.  Mesh sizes are summarized in Table 2. The total number of elements were 195,050, and the total number of nodes was 33,506.

Results
The results concern mainly the temporal evolution of the velocity field, the pressure, and the wall displacement. Figure 4 shows the blood speed distribution inside the vessel at different times. To discuss this figure, one should keep in mind the profile of the inlet velocity depicted in Figure 2. At t = 0.1 s, the fluid was mainly directed to the first inlet. The section between the first and second outlets started to move. However, the remaining part of the aorta was almost at rest. At t = 0.2 s, the fluid also reached the second and the third exits. Progressively, the blood flowed into the vessel accelerating or decelerating. The un-steady flow remained laminar.

Results
The results concern mainly the temporal evolution of the velocity field, the pressure, and the wall displacement.   From the results of Figure 4, one can also compute the flow rate in each output as well as the stroke volume. Figure 5 exhibits the temporal evolution of the maximal velocity. The maximum value reached varied from 0.1 to 0.32 m/s depending on time. The relative curve is coarsely similar to that of the input velocity. It is constituted by two peaks that have different amplitudes and different frequencies.
A least-square fit allowed to obtain the amplitude of these peaks. We obtained 0.35 m/s for the first peak and 0.26 m/s for the second one. The discrepancy between the computation results and the theoretical results obtained by the least-square fit was in a tolerable limit. The mean uncertainty was about 10%.
Appl. Sci. 2019, 9, x FOR PEER REVIEW 7 of 12 fluid was mainly directed to the first inlet. The section between the first and second outlets started to move. However, the remaining part of the aorta was almost at rest. At t = 0.2 s, the fluid also reached the second and the third exits. Progressively, the blood flowed into the vessel accelerating or decelerating. The un-steady flow remained laminar. From the results of Figure 4, one can also compute the flow rate in each output as well as the stroke volume. Figure 5 exhibits the temporal evolution of the maximal velocity. The maximum value reached varied from 0.1 to 0.32 m/s depending on time. The relative curve is coarsely similar to that of the input velocity. It is constituted by two peaks that have different amplitudes and different frequencies.
A least-square fit allowed to obtain the amplitude of these peaks. We obtained 0.35 m/s for the first peak and 0.26 m/s for the second one. The discrepancy between the computation results and the theoretical results obtained by the least-square fit was in a tolerable limit. The mean uncertainty was about 10%. It was also interesting to analyze the pressure distribution versus time. Figure 6 exhibits the pressure distribution at t = 0.1 s and t = 0.7 s. It is clear that the pressure has a complex distribution. Furthermore, this distribution depends strongly on time, since the input conditions are unsteady.
The values of the wall shear stress applied to the blood vessel walls during a cardiac cycle vary strongly in time and position. This has an important effect on cardiovascular diseases.  It was also interesting to analyze the pressure distribution versus time. Figure 6 exhibits the pressure distribution at t = 0.1 s and t = 0.7 s. It is clear that the pressure has a complex distribution. Furthermore, this distribution depends strongly on time, since the input conditions are unsteady.  Figure 7 shows the temporal evolution of the maximal pressure. This figure confirms the large variation in pressure during a cardiac cycle. Evidently, this has an important effect on the wall shear stress and consequently on cardiovascular diseases. Like the maximal velocity, the temporal evolution of the maximal pressure consisted of two pulses. The amplitudes of the pulses were 12.65 kPa and 9.62 kPa. The mean uncertainty between the computational results and the data obtained from the least-square fit was lesser than 8%. The values of the wall shear stress applied to the blood vessel walls during a cardiac cycle vary strongly in time and position. This has an important effect on cardiovascular diseases. Figure 7 shows the temporal evolution of the maximal pressure. This figure confirms the large variation in pressure during a cardiac cycle. Evidently, this has an important effect on the wall shear stress and consequently on cardiovascular diseases. Like the maximal velocity, the temporal evolution of the maximal pressure consisted of two pulses. The amplitudes of the pulses were 12.65 kPa and 9.62 kPa. The mean uncertainty between the computational results and the data obtained from the least-square fit was lesser than 8%.  Figure 7 shows the temporal evolution of the maximal pressure. This figure confirms the large variation in pressure during a cardiac cycle. Evidently, this has an important effect on the wall shear stress and consequently on cardiovascular diseases. Like the maximal velocity, the temporal evolution of the maximal pressure consisted of two pulses. The amplitudes of the pulses were 12.65 kPa and 9.62 kPa. The mean uncertainty between the computational results and the data obtained from the least-square fit was lesser than 8%.    The temporal evolution of maximal displacement during a cardiac cycle was investigated. Figure  9 presents the related results. Similarly to the maximal velocity and the maximal pressure, the temporal evolution of the maximal pressure consisted of two pulses. Using the least-square fit, we found that the amplitudes were 8.33 μm and 6.26 μm. The mean uncertainty due to the fit was about 3%. The temporal evolution of maximal displacement during a cardiac cycle was investigated. Figure 9 presents the related results. Similarly to the maximal velocity and the maximal pressure, the temporal evolution of the maximal pressure consisted of two pulses. Using the least-square fit, we found that the amplitudes were 8.33 µm and 6.26 µm. The mean uncertainty due to the fit was about 3%. The temporal evolution of maximal displacement during a cardiac cycle was investigated. Figure  9 presents the related results. Similarly to the maximal velocity and the maximal pressure, the temporal evolution of the maximal pressure consisted of two pulses. Using the least-square fit, we found that the amplitudes were 8.33 μm and 6.26 μm. The mean uncertainty due to the fit was about 3%. It appears clearly that the maximal displacement depended on time and could reach 8 μm. These values are in accordance with the important shear stress on the vessel. This has an impact on cardiovascular diseases such as stenosis.
Our results have the same order of magnitude as previous works [31,[35][36][37]. For instance, Choudhari and Panse obtained a radial deformation ranging from 1.5 to 2.2 μm and a maximal velocity about 0.6 m/s. In this study, the pressure, the velocity, and the deformation were predicted with respect to a given blood flow velocity. It could be interesting, in future work, to analyze the inverse problem. A numerical method such as the iterative gradient search method is expected to be useful to perform this analysis [38,39]. It appears clearly that the maximal displacement depended on time and could reach 8 µm. These values are in accordance with the important shear stress on the vessel. This has an impact on cardiovascular diseases such as stenosis.
Our results have the same order of magnitude as previous works [31,[35][36][37]. For instance, Choudhari and Panse obtained a radial deformation ranging from 1.5 to 2.2 µm and a maximal velocity about 0.6 m/s. In this study, the pressure, the velocity, and the deformation were predicted with respect to a given blood flow velocity. It could be interesting, in future work, to analyze the inverse problem. A numerical method such as the iterative gradient search method is expected to be useful to perform this analysis [38,39].

Conclusions
In this work, we developed a computer simulation based on the finite element method to investigate the interaction between blood flow and arterial wall deformation. This simulation tool was applied to the upper part of a child's aorta. The numerical results showed that the blood flow had an unsteady and complex behavior. The velocity field, pressure, and displacement distributions were investigated during a cardiac cycle. Further, the present study reveals that modeling the blood flow is helpful to the investigation of cardiovascular diseases and it can be an effective tool in hemodynamics studies of cardiovascular diseases.
Experimental validation on real human aorta and a silicon rubber model would be very useful. We intend to perform this validation in a future work.