Abstract
The paper presents a numerical investigation of non-Newtonian modeling effects on unsteady periodic flows in a two-dimensional (2D) constricted channel with moving wall using finite volume method. The governing Navier-Stokes equations have been modified using the Cartesian curvilinear coordinates to handle complex geometries, such as, arterial stenosis. The physiological pulsatile flow has been used at the inlet position as an inlet velocity. The flow is characterized by the Reynolds numbers 300, 500, and 750 that are appropriate for large arteries. The investigations have been carried out to characterize four different non-Newtonian constitutive equations of blood, namely, the (i) Carreau, (ii) Cross, (iii) Modified-Casson, and (iv) Quemada. In these four models, blood viscosity is a nonlinear function of shear rates. The Newtonian model has been investigated to study the physics of fluid and the results are compared with the non-Newtonian viscosity models. The numerical results are presented in terms of streamwise velocity, wall shear stress, pressure distribution as well as the vorticity, streamlines, and vector plots indicating recirculation zones at the poststenotic region. Comparison has also been illustrated in terms of wall pressure and wall shear stress for the Cross model considering different amplitudes of wall oscillation.
1. Introduction
Atherosclerosis is known as a major arterial disease. In atherosclerosis, localized deposits and accumulation of cholesterol and lipid compounds as well as proliferation of connective tissues originate a partial decline in the arterial cross-sectional area which, in particular, is called stenosis. Atherosclerotic lesions mainly occur in arterial segments with high curvature or bifurcations and junctions initiating notable alterations in flow structure and fluid loading on vessel walls [1]. Diabetes, smoking, inflammation, ischemia, and so forth are the main risk factors of stenosis development. Although the precise mechanisms responsible for the initiation of this phenomenon are not apparently known, it has been established that once a mild stenosis is developed, the resulting flow disorder plays a significant function in the further development of the disease that eventually changes the regional blood rheology as well [2, 3].
Non-Newtonian viscosity can be employed to characterize the rheological behavior of blood. The rheology and the fluid dynamical properties of blood flow can play a significant role in the fundamental understanding, diagnosis, and treatment of many cardiovascular and arterial diseases [4]. Therefore, many researchers have paid their attention to studying the hemodynamics for various viscous conditions. For example, Tu and Deville [5] implemented Galerkin finite-element method simulations for physiological pulsatile flow through a severe stenosis. They treated blood as non-Newtonian fluid employing a Herschel-Bulkley model which roughly behaves like blood. They illustrated results for steady and pulsatile flow conditions in terms of velocity profile, formation of vortex in separate regions, pressure drop across the stenosis, wall shear stress, and the vorticity contours. Non-Newtonian behavior of the blood is a crucial factor affecting the primary and secondary flow patterns near the junction between the bypass graft and the stenosed artery. Chen et al. [6] explored the non-Newtonian fluid flow in a stenosed coronary bypass numerically employing the Carreau-Yasuda model where they revealed significant differences in axial velocity profiles, secondary flow streamlines, and wall shear stress (WSS) between the non-Newtonian and Newtonian fluid flows. Moreover, the effects of both non-Newtonian behavior and the pulsation of blood flow on the distributions of luminal surface low-density lipoprotein (LDL) concentration and oxygen flux along the wall of the human aorta were numerically analyzed by Liu et al. [7] whereas Razavi et al. [8] performed a numerical analysis to study the viscous effects of blood. The power-law model demonstrated higher deviations in terms of velocity and wall shear stress compared to Newtonian and six other non-Newtonian viscosity models in their investigation. They also found that increasing stenosis intensity causes more disturbed flow patterns in the downstream of the stenosis and WSS develops remarkably at the stenosis throat. Sriram et al. [9] reveal the importance of incorporating non-Newtonian blood properties into estimates of WSS in microvessels.
Taking the physiological inlet condition of the blood flow into consideration, Long et al. [10] numerically simulated pulsatile blood flow in straight tube stenosis models with area reduction of 25%, 50%, and 75%. A measured human common carotid artery blood flow waveform was used as the upstream flow condition which has a mean Reynolds number of 300. The results indicated WSS oscillations (between positive and negative values) at various downstream locations in some models. On the other hand, comparison between numerical solution for simple pulsatile and physiologically pulsatile flow through a 61% stenosed artery was worked out by Zendehbudi and Moayeri [11]. They considered laminar flow, Newtonian and axisymmetric rigid tube for flow field computation.
Effects of compliance on diagnostic parameters have also been studied by Moayeri and Zendehbudi [2] where they numerically investigated pulsatile blood flow through stenotic arteries considering physiological flow in a dog femoral artery assuming Newtonian fluid. An isotropic elastic and incompressible material was assumed for the wall at each axial section but a nonuniform distribution of the shear modulus in axial direction to model the high stiffness of the wall at the stenotic location. Their results indicated that deformability of the wall causes an increase in the time average of pressure drop, but a decrease in the maximum wall shear stress. Furthermore, flow field and stress field for different degrees of stenoses under physiological conditions were investigated by Li et al. [12] which suggests that severe stenoses cause considerably large pressure drop across the throat that inhibits wall motion resulting in higher blood velocities and higher peak wall shear stress and localization of hoop stress.
A different kind of enquiry was carried out to study the influence of arterial wall-stenosis compliance on the coronary diagnostic parameters by Konala et al. [13]. Three parameters were defined as fractional flow reserve (FFR), pressure drop coefficient (CDP), and lesion flow coefficient (LFC) to be assessed for varying degrees of epicardial stenoses. The paper concluded that the differences in diagnostic parameters with compliance at intermediate stenosis (78.7–82.7% area blockage) could lead to misinterpretation of the stenosis severity. Again, the critical rates of blood flow acceleration and deceleration at sites of artificially induced stenosis (vessel side-wall compression or ligation) are a function of tissue elasticity; this was found by Tovar-Lopez et al. [14] when they investigated the relationship between the local hydrodynamic strain-rates and the severity of arteriolar stenosis in the small bowel mesenteric vessels of mice. Vahidi and Fatouraee [15] presented a computational model using fluid structure interactions (FSI) to investigate the physical motion of a blood clot inside the human common carotid artery. They simulated transportation of a buoyant embolus in an unsteady flow within a finite length tube having stenosis. The maximum magnitude of arterial wall shear stress during embolism occurred at a short distance proximal to the throat of the stenosis.
Paul et al. [16, 17] and Molla et al. [18] investigated the simple sinusoidal pulsatile flow in a planer channel with a cosine shape single stenosis for maximum using the LES technique. They [19–21] also have studied the physiological pulsatile flows in a channel with a single stenosis for the Newtonian and non-Newtonian fluids using the large-eddy simulation technique. Again, rheological properties of blood have been studied for transition-to-turbulent condition using LES technique in [22].
From the physiological point of view, wall pressure and wall shear stress caused by physiological pulsatile flow in a stenosed artery play significant roles in hemodynamics. Fry [23] revealed that high wall shear stress initiated by atherosclerosis is a strong parameter for endothelial or inner side damage in an artery. It can also overstimulate platelet thrombosis causing blockage [24]. Therefore, it is important to study the hemodynamic factors to realize the basic scenario behind the physiology of arterial diseases. Moreover, from the above discussion it is apparent that many researchers performed various simulations to study the fluid flow pattern both physically and physiologically. In the present study, we also simulated the pulsatile flow through an axisymmetric artery considering Newtonian, Carreau [25], Cross [26], Modified Casson [27], and Quemada [28] models as molecular shear thinning viscosity models and studied different interesting flow phenomena varying the fluid viscosity, wall condition, flow velocity, and geometry. The objective of this study is to analyze the modeling consequences of various non-Newtonian fluids in terms of streamlines, wall pressure, and wall shear stress along with the Newtonian one. Efforts have been made to obtain a very good flow insight in a stenotic artery considering physiological inlet and moving arterial wall with a very small degree of oscillation.
2. Hypothesis
The wall motion is generated using (A.8) that is given in the appendix. In this study, the wall moves sinusoidally along the height of the channel. Apparently, the wall oscillation is prescribed by the authors. As a result, it influences the flow pattern but the wall motion is not induced by the fluid flow. This is how our study differs from fluid structure interaction where both fluid and structure influence each other.
3. Governing Equations
Physiological pulsatile flow is simulated for Reynolds number 300. The arterial wall is assumed moving and blood is modeled as both Newtonian and non-Newtonian fluids for the flow field computation. The governing momentum equations for non-Newtonian 2D flows are
The blood viscosity, , depends on the shear rate , and its magnitude is defined as . When blood is treated as a Newtonian fluid, its viscosity tends to become constant value which is denoted by Pa·s. Moreover, constitutive relations used for the apparent viscosity of the blood are presented in Section 6 for non-Newtonian models. The above governing equations have been modified using the general Cartesian curvilinear coordinate system which is described in the following section.
4. Computational Geometry
The geometry is a two-dimensional (2D) channel with a cosine-shaped blockage or stenosis. Owing to the presence of the stenosis, the cross-sectional area of the channel, , is a variable in the streamwise direction [i.e., ]. Away from the stenosis, the height of the channel is constant and is represented here using (i.e., in the region either upstream or downstream of the stenosis). The stenosis is centered at downstream of the channel inlet (i.e., the inlet location is ) and upstream from the channel outlet. The stenosis is centered at and length of the stenosis is . The mathematical form of the stenosis chosen for this study iswhere is a parameter that controls the height of the stenosis. In the present study, is fixed to 1/2 that results in a 50% reduction of the cross-sectional area at the center of the stenosis. However, a schematic view of the model has been depicted in Figure 1(a) along with a subsequent illustration of the grid in Figure 1(b).
(a) Schematic view
(b) Grid arrangement
5. Physiological Flow
The physiological pulsatile velocity profile is obtained from the solution of a one-dimensional momentum equation where the pressure gradient is the Fourier series of time [29]. Womersley first calculated the physiological velocity profile for a tube by using the pressure gradient.
The steady part of the solution of the velocity field is obtained as
The oscillatory part of the solution takes the following form:
Using the definition of Womersley number, , the full solution including the steady and oscillatory part for harmonics can be written asThe real part of this solution is used to generate physiological velocity profile at the inlet of the channel which is shown below:where
In the above equation, constants and correspond to the steady and oscillatory parts of the pressure gradient; and represent coefficients and the phase angle; is the number of harmonics of the flow set to 4 (considering the first four harmonics of the pressure pulse that are necessary for modeling realistic arterial blood waveforms); is the frequency of the pulsations, is the time period of a pulsation cycle, and is the unit imaginary number, where bulk velocity depends on the Reynolds number. In present study, for streamwise velocity no slip condition has been used at the wall but the normal velocity component changes as . The zero gradient condition is applied at the outlet of artery, where velocity gradients and are zero along the streamwise direction.
6. Non-Newtonian Viscosity Models
The shear rate varies from 1 to 1200 s−1 over a cardiac cycle in human arteries [12]. Blood behaves as non-Newtonian fluid in some parts of the cardiac cycle as well. Hence, it is necessary to consider the non-Newtonian behavior of blood while studying hemodynamics. There are four non-Newtonian viscosity models that are incorporated for the numerical simulation. The models are summarized in Table 1. Newtonian fluid viscosity is always constant against the shear rate. On the other hand, non-Newtonian fluid viscosity changes depending on the shear rate. However, Boger [30] showed in his experiment that non-Newtonian fluids should exhibit Newtonian behavior or their viscosity should remain constant in both lower and higher shear rate.
Figure 2 represents the relationship between blood viscosity and the shear rate for Newtonian and four different non-Newtonian models. In these models, there is a limitation in the lower shear rate limits; they are not showing constant viscosity or Newtonian behavior according to Boger [30].
The relationship between the shear rates and viscosity for the non-Newtonian blood viscosity models, that is, Carreau, Cross, Modified Casson, and Quemada models along with the Newtonian viscosity model, is presented in Figure 2. Blood viscosity is constant in the Newtonian model shown by the solid line. On the other hand, the viscosity of blood produced by non-Newtonian models for low shear rates (less than 100 s−1) is higher than that of the Newtonian model. Viscosity in the Carreau and Modified Casson models tends to asymptotic constant viscosity, , at the shear rate s−1. The Quemada and the Cross models exhibit the non-Newtonian properties of blood at shear rates s−1. Particularly, the Cross model asymptotically matches the constant viscosity at the shear rates s−1 whereas the Quemada model shows the asymptotic nature below the constant viscosity, .
7. Numerical Procedures
The FORTRAN code employed in our computation uses a finite volume approach where the governing equations are integrated over the mesh control surface to obtain a system of algebraic equations. To facilitate calculation, the governing equations are transformed into curvilinear coordinates.
To discretize the governing equations (A.9) a second-order central difference is used for the spatial derivatives. Time derivatives are discretized by a three-point backward difference scheme with a constant timestep of .
Using this pressure correction algorithm, the pressure and velocity components are stored at the center of a control volume according to the collocated grid arrangement. The Poisson like pressure correction equation is discretized by using the pressure smoothing approach, which prevents the even odd node uncoupling in the pressure and velocity fields. A BI-CGSTAB [31] solver is used for solving the matrix of velocity vectors, while for the Poisson like pressure correction equation a ICCG [32] solver is applied due to its symmetric and positive definite nature.
8. Code Validation
In the present study, a 2D simulation is performed using the same code and for the code validation purpose, the 2D code is compared with the experimental results of Ahmed and Giddens [33] and numerical results of Damodaran et al. [34] in Figure 3. Moreover, for this specific study, two grid independence tests have also been performed as a part of code validation.
(a)
(b)
Here, the grid independence test has been carried out in Figure 4 to establish a suitable combination of grid configuration to accurately predict the physiological flow behavior for different viscous fluids in a model arterial stenosis. Different combination of grids has been chosen here as different cases that are Case 1: ; Case 2: ; and Case 3: . The number of streamwise grid points upstream of the stenosis is always fixed at 50 but the rest of the grid points are distributed nonuniformly within and downstream of the stenosis. Figure 1 shows the grid system and here it is clearly seen that the grid is significantly refined to accurately resolve the wall shear stress in the near wall region and after the center of the stenosis.
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
(j)
(k)
(l)
In this section, the grid independence test is done for the streamwise velocity at different locations for above three cases considering . The results are identical for Cases 1 and 2 in streamwise velocity up to the stenosed portion. On the other hand, a bit variation is observed between the three cases in some poststenotic regions and this percentage of error is quite acceptable since the high Reynolds number and also the physiological inlet condition drive the flow to be transitional to turbulent. The effect of wall oscillation (the amplitude of oscillation, ) also contributes to this slight variation in magnitudes. Therefore, it can be concluded that the present grid is capable of providing convergent solution independent of different grid sizes. Moreover, it is convenient to investigate the characteristics of the hemodynamic factors using any of the cases as they provide appropriate and convergent solutions; hence, Case 2 (330 × 110) has been chosen in this investigation.
9. Results and Discussion
9.1. Laminar Flow Behavior
In the following parts of the discussion, we have particularly studied the hemodynamic parameters, that is, streamlines, wall pressure and wall shear stress, vorticity, and so forth, for two different Reynolds numbers; that is, and considering various rheological models. However, the bulk velocity is computed here considering to study the laminar flow behavior.
9.1.1. Effect of Wall Oscillation
The wall of the model geometry has been considered oscillating in this study. Hence, different results considering four different amplitudes of oscillation, that is, (rigid wall condition), , , and , are presented here.
Figures 5(a) and 5(b) represent the wall pressure and wall shear stress (WSS), respectively, for different amplitudes of wall oscillation. Non-Newtonian physiological pulsatile flow is considered as an inlet flow for Modified Casson model. From the figure it is clearly observed that pressure drop is significantly higher at the stenosis throat in the rigid wall condition . As the amplitude of wall oscillation increased, the pressure drop at the center of stenosis is also optimized. Moreover, the oscillating wall conditions create highly oscillating pressure distribution throughout the arterial segment. On the other hand, the wall shear stress is the maximum at the stenosis center for highly oscillating wall condition that is illustrated in Figure 5(b). Flow with different oscillation amplitudes causes more or less closer range of shear stress in the poststenosis region whereas the rigid wall condition contributes to less fluctuating and less negative shear stress in the downstream of the stenosis.
(a)
(b)
The centerline velocity pattern in the streamwise direction considering different wall conditions is presented in Figure 6 for . All the flows cause increasing the velocity at the stenosis throat. The rigid wall model then decreases the velocity in the downstream regions that fluctuates in a small range between . On the contrary, as the amplitude is increased, the fluctuation increases tremendously in the downstream regions. Moreover, he highest oscillating wall model causes the velocity to oscillate between at to .
The effect of wall oscillation is also presented here in terms of streamlines considering Modified Casson model in Figure 7. Figure 7(a) depicts the streamlines for the rigid wall model that corresponds to the presence of a pair of wide recirculation regions in the downstream of the stenosis. However, the streamtraces are straight near the upper and lower walls of the geometry in this model. Conversely, as we increase the amplitude of wall oscillation, the streamtraces start bending in the near wall regions and the sinusoidal shaped streamtraces are apparent in the oscillating wall models.
(a)
(b)
(c)
(d)
9.1.2. Comparison between Re = 300 and Re = 500 Flow Characteristics
Since studying the flow intensity in a stenosis can help to predict various biomedical issues as blockage of the arterial segment, magnitude of velocity at every point in the geometry is an important indicator as well. In the following part, the results of streamwise velocity caused by physiological pulsatile flow with a critical stenosis are going to be discussed.
The nondimensionalized streamwise velocity () recorded at different axial locations is presented in Figures 8(a)–8(l) for the Newtonian and different non-Newtonian models considering . The investigation has been started with the inlet position and ended with outlet but between these two positions several points were also taken to show how velocity develops gradually in the artery in these points. The -axis represents the velocity profile in some streamwise locations where the -axis shows different locations on the artery. A fully developed physiological pulsatile profile is observed at the inlet that is shown in Figure 8(a). The velocity increases and attains the highest value at the neck of the stenosis in Figure 8(b) where the Newtonian and non-Newtonian models show almost identical patterns; the non-Newtonian models cause the velocity profiles to be flattened. Velocity of the flow starts decreasing gradually throughout the downstream of the stenosis and reaches the lowest peak in Figure 8(l). Significant variations in velocity between the different viscous models are observed in the poststenotic regions through Figures 8(c)–8(h) where Newtonian model causes the highest and Cross model causes the lowest velocity. Negative values of velocity in the downstream locations of the stenosis near the upper and lower wall correspond to the presence of permanent recirculation zones shown in Figures 8(c)–8(e).
(a)
(b)
(c)
(d)
(e)
(f)
(g)
(h)
(i)
(j)
(k)
(l)
In the presence of severe stenosis, resistance, or impedance to flow, pressure drop becomes highly significant and, in detail, can finally restrict the flow to human body parts [35]. Thus, knowledge of the relationship between the flow and the pressure drop across a stenosis is a prerequisite to comprehending the effect of stenotic obstructions on the distribution of blood flow to peripheral vascular beds. However, thrombus and atherosclerotic plaque form a stenosis and modify the local hemodynamics including shear which is well known to affect thrombosis characteristics. Therefore, thrombus formation and embolization are shear-dependent. Likewise, atherothrombosis can induce acute myocardial infarction and stroke by progressive stenosis of a blood vessel lumen to full occlusion [36]. As a result, it is also important to study the effect of wall shear stress inside an idealized constriction caused by physiological pulsatile inlet flow.
According to Ku [37], the cyclic nature of the heart pump creates pulsatile conditions in all arteries. Moreover, heart always maintains two cyclic phases. Blood is ejected and filled in the heart in alternating cycles called systole and diastole. Blood is pumped out of the heart during systole and rests during diastole when no blood is ejected. Wall pressure through an oscillating wall arterial segment with stenosis for considering different phases of a cardiac cycle has been estimated here which is depicted in Figure 10. These diagrams are illustrated based on a physiological cardiac cycle which is also shown in Figure 9. Six distinct phases of a cardiac cycle are plotted here for wall pressure which is again shown inset.
(a)
(b)
(c)
(d)
(e)
(f)
Wall pressure exhibits a gradual decreasing manner coupled with slight drop in the constriction in the early systolic phase that is shown in Figure 10(a). It drops maximum for the Carreau model at the center of the stenosis near the peak systole that is shown in Figure 10(b). Oscillating nature of the wall pressure pattern is prominent in these phases. In the late systole, wall pressure follows again a decreasing manner accompanied by less oscillations throughout the arterial segment. Wall pressure for the Newtonian model is remarkably visible in the early diastole (Figure 10(d)) since it is pretty higher in magnitude than the other models. Moreover, wall pressure caused by all the models shows sharp zigzags in the prestenotic zones that drop in the bulge and then oscillate almost in sinusoidal manner in the poststenotic regions. Wall pressure becomes nearly steady and oscillates slowly in the middiastolic phase that is shown in Figure 10(e). The Modified Casson model shows relatively high oscillating magnitudes whereas the Cross and the Carreau models cause lower magnitudes of wall pressure. Oscillation becomes again high for all the rheological models near the late diastolic phase shown in Figure 10(f).
Wall pressure variations for various viscous models considering are presented in Figure 11. Different phases of a physiological pulsatile cycle are considered here to visualize the exact scenario of fluid pressure created upon the arterial walls. In the beginning of the pulse, the wall pressure varies in a range of () for all the models while the highest wall pressure occurs at the proximal end of the stenosis and the lowest pressure occurs at the stenosis throat by the viscosity varying models. On the contrary, multiple peak values are observed to occur by the Newtonian model at different axial locations; that is, , , and . However, Carreau model maintains a lower wall pressure and Newtonian model maintains a higher pressure throughout the arterial locations. All the models cause a lower range of wall pressure during the peak systolic phase at when the wall pressure becomes most negative at the stenosis throat. However, it recovers the magnitude during the end systolic phase when the highest wall pressure occurs at by the Newtonian model. The wall pressure again falls to a lower range during the cardiac phase, where the Carreau model shows a comparatively higher magnitude and the Newtonian model maintains the lower range of wall pressure. The wall pressure recovers during the middiastolic phase where the stenosis throat experiences the lowest pressure like the other phases and the Carreau model shows the lower range of pressure. Near the end diastolic phase, the wall pressure follows the previous pattern creating the lowest pressure at the center of the stenosis (−4 in magnitude) where the Newtonian model causes the highest and the Carreau contributes to the lowest wall pressure.
(a)
(b)
(c)
(d)
(e)
(f)
Similar results for wall shear stress considering different phases of a cardiac cycle have been evaluated in this work which is depicted in Figure 12. Six distinct phases of a cardiac cycle, that is, early systole, peak systole, late systole, mid diastole, and late diastole, are plotted here considering for WSS that is also shown inset.
(a)
(b)
(c)
(d)
(e)
(f)
In the early systolic phase, it is clearly seen that all five viscous fluids follow a oscillating pattern and the peak WSS occurs at the throat of the stenosis. In the peak systole, the similar pattern is maintained by all the rheological models as well as the Newtonian one. The peak value occurs sharply at the center of the stenosis followed by oscillating zigzag pattern during peak systole that is also the highest WSS (with a magnitude 0.8) among all other phases. All the viscous models show roughly closer values of WSS throughout the arterial segment during these two phases. However, the Carreau and the Cross models cause relatively higher range of shear stress while the Newtonian model shows lower range of WSS. In the late systole and middiastolic phases that are shown in Figures 12(c) and 12(d), the WSS follows the same pattern as the previous phases but the peak value decreases gradually here and drops to a magnitude around 0.5 at the early diastolic phase. Moreover, the viscous models exhibit some variations in magnitude in the arterial locations; that is, here each curve for the models can easily be identified separately. WSS for another diastolic phase and the late diastole is depicted in Figures 12(e) and 12(f). The peak value occurs at the throat location. The oscillation also ceases at the middiastole but it starts increasing near the late diastolic phase and tries to catch the WSS pattern of the early systolic one. In simple words, shear stress is very high in the systolic phases that is again erratically sharp but it becomes small in magnitude and also slower in frequency in the diastolic phases especially in the middiastole.
Wall shear stress distribution for various rheological models considering has been illustrated in Figure 13 where the scenario is presented for different phases of a cardiac cycle. The WSS follows almost similar pattern throughout the arterial positions during different phases of the cardiac pulse. The peak stress generally occurs at the stenosis throat and an oscillating WSS is observed at different poststenotic arterial locations. However, the highest WSS is observed during peak systolic and end diastolic phases; that is, and with a magnitude near 0.3. Moreover, all the models cause the highest peak at the center position of the stenosis and the Cross and the Carreau models behave in a similar manner when .
(a)
(b)
(c)
(d)
(e)
(f)
Figure 14(a) illustrates the streamlines caused by physiological pulsatile flow in a stenosed arterial segment for . Figure 14(a)(i) represents the streamline for Newtonian case and Figure 14(a)(ii) to Figure 14(a)(v) represent the streamlines for four different non-Newtonian viscosity models (Carreau, Cross, Modified Casson, and Quemada, resp.). In the case of Newtonian fluid, the flow is fully developed in the constriction. The magnitude of velocity field ranges from 1.60 to 1.80 from the proximal region to at least downstream region from the stenosis neck. For the Carreau and Cross models, the length of high velocity field is significantly shorter than the Newtonian model. It ends at position and no flow recirculation is observed in these regions. Unlike the previous two viscosity models, the Modified Casson model creates relatively longer high velocity field in the poststenotic region. The Quemada model acts differently than the Modified Casson model; the velocity field is very weak in the constriction. In fact, recirculation region is absent inside the stenosed arterial segment when the . The difference in the length of high velocity field for each model can be explained from the fact that the behavior of the Quemada model is the most viscous followed in order of decreasing viscous behavior by the Cross, the Carreau, the Modified Casson, and finally the Newtonian model.
(a) Streamlines for
(b) Streamlines for
Vorticity highly depends on viscosity. More insight into the flow separation seen in Figure 15 is given through the streamwise vorticity contours, . Figure 15(i) to Figure 15(v) represent five different viscosity models, that is, the Newtonian, Carreau, Cross, Modified Casson, and Quemada models, respectively. It is noted that a total of 15 different contour levels are plotted between their maximum and minimum values that can be viewed easily through the legend color bar. The vortex units rotated in the clockwise and anticlockwise direction give positive and negative values of , respectively. The clockwise rotations are represented by the solid lines where the anticlockwise rotations are represented by dashed lines. Vorticity is very high in low viscosity models. As a result, more vortices are generated from the nose of the stenosis where the flow separation begins in the Newtonian fluid shown in Figure 15(a)(i). Multiple vortical structures are observed in the downstream regions between ; one acts in the anticlockwise direction and the other acts in the clockwise direction near the upper and lower walls. Both of them interact with each other and then roll up to downstream region. According to Berger and Jou [1], multiple vortices propagate at higher speed than that of a single vortex that happens in the Newtonian fluid. Furthermore, small clockwise vortices are seen in both the upper and lower wall region even throughout the distant downstream regions; vortex multiplication and flexibility of the wall might contribute to this phenomenon. In case of the non-Newtonian models, the Modified Casson model creates comparatively larger area of vortex propagation from that is mostly contributed by the near wall small vortex structures. Other rheological models show small vortical structures in the downstream regions. The strength of the vortices is also weaker in these models than the Newtonian fluid due to more viscosity. Among all of the viscous models, Modified Casson model shows more and the Cross model shows fewer vortex structures.
(a) Vorticity for
(b) Vorticity for
Vorticity for various viscous models considering is illustrated in Figure 15(b) that presents the flow separation clearly. A separate region of comparatively high velocity is present along the central axis of the artery of the Newtonian model and multiple negative velocities with anticlockwise movement of the fluid particle is found near the upper and lower walls of the geometry that covers a larger region than that of the Newtonian fluid. Clockwise movements of the fluid particles are dominant in almost every regions of the arterial segment for the Carreau, Quemada, and Cross models except in the near wall regions. A small area with negative vorticity exits there and the viscous behavior follows a descending order in these three models. However, the Modified Casson model corresponds to the least viscous behavior among the four rheological models; as a result, both the high and low velocity regions are relatively larger in this model.
9.2. Transitional Flow Behavior
High Reynolds numbers are incorporated here to characterize the transitional behavior of the flow studying the hemodynamic parameters, that is, wall shear stress, streamtraces, and vorticity for different heights of the model stenosis and different Reynolds numbers. However, all the simulations have been done here considering Modified Casson model and has been used to compute the bulk velocity, .
Wall shear stress distribution for Modified Casson model considering and different percentages of the stenosis height is shown in Figure 16. It is apparently seen in the figure that as the stenosis height increases the WSS distribution varies more in different parts of the geometry and the 60% stenosed model generally causes the highest peaks at the throat location during different phases of the cardiac cycle. A peak is found at the center of the stenosis by all the models during the beginning phase of the cycle which is associated by a drop at about by the 60% stenosed model. During the peak systole, two sharp peaks are observed at (a positive rise) and (a negative drop) by the mostly stenosed model. At the end of the systole, the negative WSS is absent in the downstream region. WSS upswings the most by all the stenosed models during the middiastolic phase shown in Figure 16(e). However, the stress variation is evident for the 60% stenosed model throughout the artery while the two other models show relatively close values in different parts of the artery. The end diastolic phase shows nearly similar patterned WSS for all the models as the beginning phase of the cardiac cycle.
(a)
(b)
(c)
(d)
(e)
(f)
The effect of increasing the percentage of stenosis in terms of streamlines for is depicted in Figure 17(a) where Figures 17(a)(i), 17(a)(ii), and 17(a)(iii) show 40%, 50%, and 60% stenosed models, respectively. The streamtraces are straight near the central axis of the artery as a result of the absence of flow recirculation in these regions in the 40% stenosed model. Velocity field covers a region of higher magnitude extending from up to the downstream region and a pair of small vortices are found at the distal end of the stenosis. When the percentage of stenosis is increased to 50%, the density of flow recirculation increases more than the previous case. Total four pairs of moderate sized vortices form in different parts of the downstream regions where the velocity range is also very high (1.70–1.90). Very large recirculation regions are observed near the upper wall of the arterial segment extending from . to . due to the effect of 60% stenosed model. Velocity also becomes very high in these regions and the vortex distribution becomes asymmetric about the central axis. As a consequence, the flow becomes transitional in this stenosis height for .
(a) Streamlines
(b) Vorticity
Vorticities for the Modified Casson model considering and different heights of the stenosis, that is, 40%, 50%, and 60%, are depicted in Figure 17(b). Positive velocity field exists near the central axis of the artery throughout the artery while negative vorticity is found in the near wall regions at about . The scenario completely changes when the stenosis height is increased to 50% that is shown in Figure 17(b)(ii). Comparatively high velocity region is observed near the central axis under the cover of the constriction while a relatively large negative velocity field is observed in the downstream region. Both the high and low velocity regions extend; as a consequence, both the clockwise and anticlockwise rotations increase. If the percentage of stenosis is increased more, to 60%, the flow loses its symmetric behavior and becomes transitional. Negative vorticity appears near the upper wall and the distal end of the constriction covering a large region whereas, multiple small regions of anticlock wise fluid movement are observed in the lower wall regions.
Wall shear stress distribution considering a range of Reynolds number for the 50% stenosed Modified Casson model is illustrated in Figure 18. WSS maintains zero magnitude in the prestenosis regions for all the Reynolds number fluids. It increases as a sharp peak at the throat location and then drops to a range between (). WSS maintains an oscillating pattern in the downstream region up to the outlet. The notable observation is that as we increase the Reynolds number, the WSS magnitude also increases at the stenosis throat and the highest peak is found during the middiastolic phase; that is, . Moreover, multiple negative peaks are found at different downstream locations of the geometry that are also contributed by the highest Reynolds number.
(a)
(b)
(c)
(d)
(e)
(f)
The effect of increasing the Reynolds number in terms of streamlines considering Modified Casson model is depicted in Figure 19(a) where Figures 19(a)(i), Figure 19(a)(ii), and Figure 19(a)(iii) show , , and models, respectively. The streamtraces are straight near the central axis of the artery due to the absence of flow recirculation in these regions in the model. Velocity field covers a region of higher magnitude extending from up to the . When the Reynolds number is increased to 500, the flow recirculation increases. Total two pairs of vortices form in the distal end region where the velocity range is (0.1–0.3). The size and length of the recirculation region increase near the upper wall of the arterial segment covering from . to . owing to the effect of Reynolds number, . Velocity becomes very high (up to 3.20) in these regions and the vortex distribution becomes asymmetric about the central axis as a number of small vortices are found near the lower wall of the geometry. As a consequence, the flow becomes transitional.
(a) Streamlines
(b) Vorticity
The effects of flow separation caused by a range of Reynolds numbers are depicted in Figure 19(b) in that the presence of very low velocity in different parts of the arterial segment is visible due to the effect of . Negative velocity and thus anticlockwise rotation of the fluid particles are observed near the upper and lower walls of the geometry. A small region of high velocity is observed at the stenosis throat when the Reynolds number is increased to 500. However, negative vorticity is found near the upper and lower walls of the geometry in this model. Multiple negative vortical structures form near the upper wall and in downstream region that are absent near the lower wall. The asymmetric distribution of the vortex generation and propagation corresponds to transitional flow characteristics when the .
10. Conclusion
Finite volume numerical simulations of unsteady, incompressible, homogeneous blood flow in two-dimensional stenosed model with sinusoidally oscillating wall condition have been presented in this paper. Flow field, flow induced wall pressure, and wall shear stress have been evaluated for Newtonian and non-Newtonian models (Carreau, Cross, Modified Casson, and Quemada) under physiological pulsatile condition. The findings of this investigation can be summarized in following ways.(i)Pressure drop is insignificant in the model predicting less risk of limiting the blood supply to body parts. Moreover, the Newtonian model maintains relatively higher wall pressure during the systolic phases throughout the arterial positions.(ii)The wall shear stress reaches its highest peak at the throat of the stenosis in the peak systole leaving the atheromatous plaques prone to rupture. The phenomenon of plaque disruption and thrombosis is not ended by luminal occlusion only and it may lead to extensive surface ulcerations comprising large areas of the aortic wall. The Modified Casson and the Newtonian models consistently follow a lower range of shear stress and the Carreau and the Cross models show relatively higher range of WSS.(iii)Streamlines demonstrate the presence of recirculation zones in the flow field that is found in the largest scale by the Newtonian fluid where the flow is fully developed in the constricted region. This again increases the possibility of thrombosis. However, the smallest recirculation region is caused by the Carreau model.
It can be stated in essence the effect of wall oscillation is evident in both the Newtonian and the non-Newtonian models. Both the flow parameter and the wall condition, being time variant, remarkably affect the wall pressure, wall shear stress, and flow recirculation. As a result, the Newtonian fluid is more likely to cause heart attack or blockage due to its characteristics of low pressure and largest recirculation region at the throat locations of the stenosis than the non-Newtonian models. On the other hand, the Carreau and Cross models exhibit higher risk of atherothrombosis owing to the high wall shear stress in the arterial segment.
Appendices
A. Coordinate Transformation
Thompson et al. [38] introduced an approach where the finite difference equations are formulated in a transformed curvilinear coordinate system that coincides with the boundaries of the fluid domain. In this approach, flow domain in physical space is mapped onto a rectangular domain in computational space, as shown in Figure 20, where a two-dimensional case is represented for simplicity. Mapping where represents the elements of the Jacobian matrix, , of the transformation, The determinant of the Jacobian matrix, is denoted by ; that is,where are the elements of the cofactor matrix, , of the Jacobian, defined asThe derivatives can be expressed in transformed time and space in the following way:where is a generic variable. Here a moving wall condition is used in the radial direction and the streamwise coordinate () does not depend on time. So (A.4) can be written asThe radial variable is defined as the function of time and space asSo the time derivative of the radial coordinate iswhere is the amplitude of oscillation and is the diameter of pipe. The governing equations (1) becomewhere and are used to represent coordinates along and directions, respectively.
(a)
(b)
B. Physiological Flow
B.1. Steady Part of the Solution
Let us consider a one-dimensional fluid flow with density and dynamic viscosity between two parallel plates separated by a distance . The governing equation of motion of the fluid is given by the following form of Navier-Stokes equation:whereHere is the constant corresponding to the steady pressure gradient.
If the flow is in steady state,
Now from (B.1),
The solution of the velocity field is obtained as (3).
B.2. Oscillatory Part of the Solution
For the oscillatory part,where is the constant corresponding to the oscillatory pressure gradient. Let us assume that the solution of (B.1) isFrom (B.1) and (B.6),
Solving the auxiliary equation corresponding to (B.7), , the complementary function can be written as
For the particular integral, the corresponding equation can be written as
Now combining the auxiliary equation and complementary function
The corresponding boundary conditions are at and at . Applying these two boundary conditions into (B.10), we get and . Based on these, the solution can be found in (4).
Using the definition of Womersley number, , the full solution including the steady and oscillatory part for harmonics can be written as (5).
The real part of this solution is used to generate physiological velocity profile at the inlet of the channel.
B.3. Real Part of the Solution
In order to separate the real part of the solution from (B.10), De Moivre’s theorem of complex numbers is needed together with some trigonometric formulas. De Moivre’s theorem givesIf , thenUsing De Moivre’s theorem, we obtain
From (5), (B.14), and (B.15), we obtain
Separating from the above equation, the real part of the solution is expressed in (6).
Nomenclature
| : | Amplitude of the wall oscillation (m) |
| : | Elements of the cofactor matrix |
| : | Diameter of artery (m) |
| : | Jacobian of transform coordinates |
| : | Pressure (Pa) |
| : | Radius of the pipe (m) |
| : | Reynolds number |
| : | Time (s) |
| : | Bulk velocity (m·s−1) |
| : | Velocity along the streamwise direction (m·s−1) |
| : | Velocity along the radial direction (m·s−1). |
| : | Viscosity of blood (kg·m−1 s−1) |
| : | Height of stenosis (m) |
| : | Kinematic viscosity (m2·s−1) |
| : | Vorticity (s−1) |
| : | Density of blood (kg·m−3) |
| : | Wall shear stress (kg·m2·s−1) |
| : | Coordinate along the streamwise direction (m) |
| : | Coordinate along the radial direction (m). |
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgment
The authors wish to thank the North South University, Bangladesh, for the financial support during this research.