Effect of linear and Mooney–Rivlin material model on carotid artery hemodynamics

Atherosclerosis is one of the most common cardiovascular diseases leading to high morbidity. The study of arterial dynamics using fluid–structure interaction (FSI) technique by taking into account the physiology of flow, the critical hemodynamic parameters can be determined which plays a crucial role in predictive medicine. Due to advances in the computational facilities, coupled field analysis such as FSI can facilitate understanding of the mechanics of stenosis progression and its early diagnosis. In this study a two-way FSI analysis is carried out using modified Navier–Stokes equations as the governing equations of blood flow for determining hemodynamic parameters. The arterial wall has been described at different linear elastic modulus and compared with hyperelastic Mooney–Rivlin model to evaluate the effect of different arterial stiffness on hemodynamics. The Mooney–Rivlin model predicts flow reduction with the severe backflow at arterial bifurcation resulting in decreased shear stress and oscillatory behavior. Furthermore, these findings may be used in understanding the advantages and disadvantages of using hyperelastic artery model in numerical simulations to better understand and predict the variable that causes cardiovascular diseases and as a diagnostic tool. In the present study, variation due to change in arterial wall properties such as linear elastic and Mooney Rivlin hyperelastic and its influence on hemodynamics are investigated.


Introduction
Atherosclerosis is a cardiovascular disease, which is the main cause of death in many advanced countries [1]. It is generally characterized by the accumulation of fatty particles, cholesterol, and other particles in the inner layer of the arterial wall. This leads to blockage of the lumen termed as arterial stenosis followed by ischemic stroke [2]. The formation of plaque is a localized phenomenon that generally occurs at locations such as bifurcations and curvatures observed in carotid arteries [3,4]. These highly vulnerable sites are conditioned by hemodynamic parameters such as low wall shear stress (WSS), oscillatory flow, or turbulent flow. Many studies are conducted to characterize the local hemodynamic nature of atherosclerosis development. Regions downstream of stenosis, post-stenosis region, experience different hemodynamic environment than normal regions due to flow reversals and sometimes turbulent flow. This may lead to high WSS at the neck of the stenosis due to high acceleration due to a reduction in the area which may lead to plaque rupture and also damage the endothelial cells [5,6]. In addition to this, a higher spatial WSS gradient associated with turbulence enhances mass transport into the artery wall, which may weaken the plaque and make it vulnerable to rupture [7].
Presently, many researchers considered low WSS is directly associated with arterial disease [8][9][10]. Yonis et al. [11] observed that WSS magnitude less than 1.5 Pa stimulates deposition of plaque generally observed in atherosclerotic prone areas. The location of the carotid bulb at the base of the internal carotid artery (ICA) is the key characteristic of carotid bifurcation. The carotid sinus along with the resistance at the outlets of the branches of the carotid artery contributes to a complex flow mechanism. Flow recirculation is observed along the outer wall of the carotid bulb at the opposite wall of the flow divider. This causes the velocity offsetting toward the inner wall of the ICA. These behaviors lead to low WSS regions along the outer wall of the carotid artery. The cause flow reversals at the decelerating phase of the cardiac cycle. The maximum WSS magnitude is at the apex of the carotid bifurcation as the flow gets divided [12].
Fluid-structure interaction (FSI) is often used to study atherosclerotic plaque and its effect on local hemodynamics, initiation of plaque, growth, and rupture. Chhai et al. [13] investigated the factors contributing to plaque rupture using FSI by examining asymmetric plaque models. The effect of maximum cap stress, WSS, was studied and concluded that the risk of plaque rupture was proportional to plaque asymmetry. The plaque inflammation was studied by Tang et al. [14] and its effect on rupture risk using FSI in anatomically realistic models. They simulated plaque weakening due to inflammation in normal and hypertensive circumstances. They observed that the inflammation of the cap causes increased deformation when the cap is thin and the patient is hypertensive. He et al. [15] studied the effects of initial development of atherosclerosis on neighboring arteries and found WSS is low and strain is high in arteries near atherosclerotic arteries. This indicates that atherosclerosis does not harm the adjacent arteries. Kafi et al. [16] studied the stress profile in a rigid artery and compared the non-Newtonian and Newtonian blood viscosity models. It was observed that the arterial rigidity has a significant effect on von Mises and wall shear stress. Kallekar et al. [17] compared the effect of rigid, linear elastic, neo-Hookean, Mooney-Rivlin, and Holzapfel material models on FSI results under steady-state and pulsatile flow conditions. He found that the results are highly dependent on the wall model used showed the importance of FSI in simulation.
Blood comprises different elements such as platelets, red blood cells, and white blood cells due to which the consideration of non-Newtonian model to mimic the shear thinning behavior is important [18]. However, researchers have interpreted blood as a Newtonian fluid because the shear rates are greater than 100 s −1 in large vessels (lumen diameter ranging from 1 to 3 mm), which makes the non-Newtonian behavior negligible [19][20][21][22].
The arterial dilation is regulated by the linear momentum equations that are affected by the anatomical model used and the design of the geometry [23]. In cardiovascular biomechanics, structural artery modeling is explored in depth [24,25]. A simpler way to mimic an artery is to consider it as a linearly elastic material. The advantage of this assumption is that it is inexpensive to compute and the results obtained by this method agree well with the experimental results within the physiological pressure range [26]. However, Holzapfel et al. [27] found that the arterial wall shows highly hyperelastic behavior. The use of Mooney-Rivlin model showed more reliable indicators.
In this study, a numerical study of artery and blood in an anatomically realistic phantom carotid artery bifurcation is studied using fluid-structure interaction. Hemodynamic parameters such as velocity, pressure, WSS, vortex core, and oscillatory shear index (OSI) are compared. The two-way fully coupled FSI analysis is performed by considering four different arterial stiffness models, viz., E 0.3, E 0.5, and E 0.9 with an elastic modulus of 0.3 MPa, 0.5 MPa, and 0.9 MPa, respectively, and two constant Mooney-Rivlin hyperelastic model to simulate the effect of different arterial stiffness on hemodynamics [28,29].

Modeling and computational setup
In this study, blood is assumed to be incompressible, homogenous, Newtonian fluid, and the flow is described by the Navier-Stokes equation [30]. The laminar fluid flow is governed by the Navier-Stokes equations for incompressible flows given as.
for fluid domain where f and v represent density and velocity fields, respectively.
The flow is defined by pressure and velocity fields p f and v f , respectively. The relation between stress tensors is given as: where represents dynamic viscosity, I represents the identity matrix, and the Lagrange multiplier p f corresponds to the incompressibility constraint in (1). ϵ v f is the strain-rate tensor given as: The governing equation for the structural domain (arterial vessel), denoted as s, is given as: where s represents arterial vessel density, g represents the externally applied body force on the structure, and s is the Cauchy stress tensor [12].
To describe the structural configuration, a common Lagrangian description with respect to fixed reference state is used, given as where v s is the velocity field and u s is the displacement.
In a FSI problem, if the fluid flow through a solid domain is conformal, the displacement experienced by the fluid volume must be completely transferred to the solid, and vice versa. Since fluid flow is described using Eulerian motion, while solid using Lagrangian, the interaction at the fluid-solid interface is described using arbitrary Lagrange-Eulerian (ALE) framework. The fluid domain deforms arbitrarily and the boundaries follow the deformation of the structure. Mass and linear momentum are conserved at this interface, satisfied by the compatibility of displacement and traction equilibrium. The interaction at the interface is given as: where def s,Γ represents structural domain deformation, def f ,Γ represents fluid domain deformation, and F s,Γ and F f ,Γ , respectively, represent the interface forces of structural and fluid domains.

Material properties
The dynamic viscosity of blood is assumed to be 0.0035 Pa. s and density as 1060 kg/m 3 [31]. The non-Newtonian property was found to always be accounted for only if the lumen diameter is less than 0.1 mm or if the shear rate is less than 100 s −1 [32]. The wall was assumed to be isotropic, homogenous, incompressible, and linearly. The Young's modulus or the arterial wall was taken as 0.3 MPa, 0.5 MPa, and 0.9 MPa and Poisson's ratio as 0.45 [29]. The structural equations of motion for elastic materials are satisfied using Hooke's law in case of small deformations known as Lame's equations. It is derived from the continuum equation of motion by substituting the equations of stress and small strain tensor [33]. To study the hyperelastic behavior of the arterial wall, Mooney-Rivlin nonlinear material model is used. For a neo-Hookean material, strain energy density function is given as where C 1 ,C 2 , and d are material parameters,I 1 = trC , , and I 3 = detC is the first, second, and third invariant of modified right Cauchy-Green tensor 4 Boundary condition Figure 1 shows the fluid and structural model of the carotid artery used in this study. A generously shared carotid artery model by the GrabCAD community is used in this study [35]. The standard intima-media carotid range is in the 0.6-0.8 mm range. The thickness of the arterial wall in this study was believed to be 0.7 mm [36][37][38]. The arterial model is generated using the fluid model with Materialise 3 Matic software, by offsetting the fluid domain by 0.7 mm to represent the arterial thickness. At the inlet, pulsatile velocity derived from Doppler ultrasound was applied. At the outlet, the impedance boundary condition was applied. Blood flow velocity was acquired after the data were aligned with the cardiac cycle and the velocity waveform was integrated over the lumen inlet by calculating the velocity components of the typical carotid artery. Figure 2 illustrates the velocity profile of the selected typical carotid artery [12]. A pulsatile pressure waveform was used at the outlets, as shown in Fig. 3 [39]. On the walls of the fluid model, a no-slip boundary condition is defined. For outlets, an impedance method developed was implied by modeling the vascular system as a fractal network [40]. In this method, the vessels beyond the imaging domain are approximated with the help of 1D flow approach [41].

Results
To choose the optimal grid size for the analysis, a grid sensitivity test was conducted. The variance of average deformation and average velocity for various grid sizes is shown in Table 1. The analysis selected an optimal grid size of 105,963 solid elements and 323,514 fluid volumes. The total cardiac cycle time is 0.8 s and it is divided into 100 times steps of 0.008 s.
In this analysis, pulsatile blood flow is taken into consideration and hemodynamic parameters such as wall deformation, wall shear stress, and flow velocity are discussed. Figure 3 shows the temporal variation in mass flow at the outlets of the carotid artery for one cardiac cycle. Maximum flow is through the ICA with 71% of the total flow through the CCA and 34% through ECA as shown in Fig. 4. This flow division has been validated by experimental studies [42]. Also, the arterial stiffness has no significant impact on the flow division between the artery branches. Further, it can be noted that the downstream peripheral resistance directly influences the flow separation through the arterial outlets, with preferential flow through the vessel with the least resistance. The difference in the flow rate observed here is mainly due to the change in the cross-sectional area of the branches due to different elasticity of the artery. Figure 5 shows the velocity contours at the mid-sagittal plane of the carotid artery model. Higher velocities are observed at the distal CCA, and the maximum velocity magnitude decreases with an increase in arterial stiffness and minimum for Mooney-Rivlin artery model. However, at the carotid sinus region, higher recirculation is observed with E 0.3, inversely varying with arterial stiffness, and minimum for Mooney-Rivlin model. At the flow divider, due to reduced cross-sectional area, blood flows at a higher velocity near the inner wall of the ECA. The velocity increases with an increase in stiffness for the linear elastic model, due to obvious reduced arterial deformation with an increase in the elasticity.
The oscillatory shear index (OSI) highlights the oscillatory flow intensity within a cardiac cycle [43], which is given as where T is the time period of the cardiac cycle. OSI measures the time-averaged variations in the shear stress. Figure 6 shows the oscillatory shear index for different arterial stiffness models. The oscillatory flow is more prominent in the E 0.3 model and reduces with decreased elasticity. The higher flow recirculation is concentrated at the carotid sinus region for all the models. As seen in the contours, the OSI follows a uniform trend among the linear elastic models. The low oscillation region is minimum at the outer wall of the ICA near the base of the carotid sinus where the flow separation takes place and increases with an increase in stiffness toward the outer wall of the ECA. However, the Mooney-Rivlin model shows a low OSI region at the base of the bifurcation, and higher OSI is concentrated at the outer wall of the carotid sinus. Figure 7 shows the vorticity, indicating the curl of the velocity vector. The higher magnitude regions seen in the  Higher OSI observed at the carotid sinus due to higher deformation in this region and also due to higher flow recirculation contour plot in Fig. 7 shows the flow separation zones after which the flow recirculated creating low WSS regions. These are the regions highly susceptible to plaque buildup. In these regions, low-velocity flow and gradients and low wall shear stress occur. These facts are of prime importance as the thrombus deposition and detection of early atherosclerotic lesions. WSS is related to atherosclerotic development and progression. Figure 8 shows the time-averaged WSS contours on the carotid artery models of different arterial stiffness. The maximum WSS is observed at the inner wall near the  root of the ECA. The magnitude increases with an increase in elastic modulus for the linear elastic artery. However, the Mooney-Rivlin model shows a higher magnitude as observed in the linear elastic model as well as near the carotid sinus where the flow separation is observed. Low WSS is at the carotid sinus in all the cases; however, the low WSS zone is reduced with increased stiffness due to lower deformation. However, our results refer to a factor that may be important in the genesis of atherosclerosis at the carotid bifurcation, namely the duration of the stay of blood and its components within the area of flow separation.

Discussion
In this study, hemodynamic parameters in a normal carotid artery are evaluated by varying the stiffness of the arterial material. It was observed that, given the results of Fig. 5, Fig. 6, and Fig. 8, it can be concluded that the carotid sinus is the most favorable area for the development of atherosclerosis. In this region, the flow of blood is reduced, leading to low WSS, which in turn results in higher OSI, which are the conditions that favor endothelial inflammation and consequently, plaque accumulation. Along with carotid sinus, plaque accumulation can also be observed at the inner wall of the ECA. The results in the present study play a significant role in predicting the hemodynamic mechanism of stenosis. The results refer to a factor that may be important in the genesis of atherosclerosis at the carotid bifurcation, namely the duration of the stay of blood and its components within the area of flow separation.
The percentage total area where the time-averaged WSS is less than or equal to 0.4 Pa was determined, with 7.97%, 7.33%, 7.04%, and 6.95% for E 0.3, E 0.5, E 0.9, and Mooney-Rivlin model, respectively, showing that Mooney-Rivlin model is almost similar to E 0.9. For simulating the behavior of lower stiffness of the carotid artery, the MR model is not suitable as it underestimates the low WSS region. Flow recirculation is maximum with E 0.3 model and decreases with increases in stiffness due to reduced arterial deformation with an increase in stiffness. This leads to variation inflow through the ECA satisfying the law of conservation of mass with low flow velocity at the inner wall of the ECA of E 0.3 and increasing with the increase in stiffness. Due to reduced deformation, the residence time of the blood is low at the carotid sinus and the flow is balanced with increased flow through the ECA. The maximum flow in the ECA is observed with the MR model, which suggests that the residence time is lowest with the MR model and hence low WSS surface area compared to the linear elastic models.
The low residence time also leads to low OSI at the carotid sinus as compared to the linear elastic models. At the bifurcation, the downstream flow leads to large contact areas leading to higher OSI. The increased rate of flow increases the contact intensity leading to the increased shear rate at the lumen surface making it a favorable region for plaque deposition. These observations encourage the more concentrated qualitative study of WSS-based indicators on the artery wall. WSS is dependent on flow patterns driven by the blood flow variation. The axial low leads to WSS to be higher in the axial direction compared to other directions. However, at increased flow rates, the oscillatory behavior is due to disturbed flow due to flow division and arterial deformation creating unstable flow domains which leading to increased flow in the direction perpendicular to the flow direction. Therefore, based on the WSS indicators, the likely location for plaque deposits can be evaluated. The WSS as observed in Fig. 5 shows higher WSS at the inner wall of the ECA near the bifurcation for all the models, with the lowest for E 0.3 and maximum for the MR model. Lower values are observed at the carotid sinus region for all the models.

Conclusion
In the present study, the pulsatile blood flow with the Newtonian assumption in the carotid bifurcation model under different arterial wall stiffness is investigated with two-way coupled FSI simulation. The results describe the important role in understanding the hemodynamic mechanism of stenosis. The obtained results demonstrated that linear elastic stiffness in comparison with the hyperelastic Mooney-Rivlin model showed significant variation in arterial wall stress distribution. However, with reduced arterial stiffness, a linear elastic model is more suitable as the Mooney-Rivlin model resembles an artery of the healthy adult with an average elastic modulus of 0.9 MPa. Increased velocities were found to be closer to the distal side of CCA, and it was observed that maximum velocity magnitude decreases with an increase in arterial stiffness, and minimum for Mooney-Rivlin artery model. The low OSI region is found to be toward the outer wall of the ICA and closer to the carotid bulb region which is prone to intense flow separation. This low OSI region increases with the increase in arterial wall stiffness and inclined toward the outer wall of ECA. However, the Mooney-Rivlin model shows that the low and high OSI region is found to be at the base and outer wall of the carotid bulb, respectively. The results of the present study demonstrate the potential of the numerical technique to understand the variation in hemodynamics under different arterial wall stiffness and its influence in understanding the atherosclerosis prognosis.
Funding Open access funding provided by Manipal Academy of Higher Education, Manipal.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.