Predicting Hemodynamic Alterations Due to Bifurcation Stenosis: A Mathematical Model and Numerical Analysis of Newtonian Blood Flow

Predicting Hemodynamic Alterations Due to Bifurcation Stenosis: A Mathematical Model and Numerical Analysis of Newtonian Blood Flow

Esam A. Alnussairy* Laheeb Muhsen Noman Ahmed Bakheet

Department of Biology, College of Science, Wasit University, Al Kut 52001, Iraq

Department of Mathematics, College of Science, Mustansiriyah University, Baghdad 10052, Iraq

Corresponding Author Email: 
eahmed@uowasit.edu.iq
Page: 
4122-4128
|
DOI: 
https://doi.org/10.18280/mmep.121202
Received: 
12 August 2025
|
Revised: 
17 September 2025
|
Accepted: 
22 September 2025
|
Available online: 
31 December 2025
| Citation

© 2025 The authors. This article is published by IIETA and is licensed under the CC BY 4.0 license (http://creativecommons.org/licenses/by/4.0/).

OPEN ACCESS

Abstract: 

Modeling blood flow through a stenosed bifurcated artery is crucial in the biomedical field, as it provides valuable insights into how these conditions affect the body and contribute to the diagnosis, treatment, and prevention of cardiovascular diseases. The Marker and Cell (MAC) method is used to numerically solve the equations, and the code was created using MATLAB programming to obtain and analyze the flow characteristics—including velocity profiles, flow patterns, and pressure distribution—that have not been thoroughly investigated in earlier studies in this field. The benefit of the MAC approach is that it eliminates the need for pressure boundary conditions at either the inlet or outflow. A graphic representation of these characteristics, along with relevant physical parameters, is provided in this study. The results indicate that the severity of the stenoses leads to an increase in the recirculation zones in the junction and the outer wall in the daughter artery, as well as the wall pressure rapidly drops at the throat of the stenosis. Also, as the blood passes through the stenosis, its flow is accelerated. Then, after the narrowest point of the stenosis, the artery returns to its normal, wider diameter. High-velocity flow stream, however, continues to move forward in a straight line, leaving a low-pressure region between the streaming blood and the artery wall. The blood in the low-pressure zone is pulled backward, and it starts to swirl in a circular motion to fill this space. This swirling motion generates the recirculation zones, which are widely recognized as a significant factor in the progression of cardiovascular disease. The current findings show good agreement with established results.

Keywords: 

numerical simulation, MAC method, bifurcated stenosed artery, Newtonian fluid

1. Introduction

Stenosis is a plaque that constricts the artery walls, resulting in significant modifications to the blood flow structure. Areas characterised by curvatures and bifurcations are the most encouraging locations for this plaque [1]. Thus, the investigation of blood circulation in the bifurcated artery is an important point of scientific interest regarding the diagnosis of atherosclerosis. This helps to gain a deeper comprehension of the pathogenic mechanism involved. Not only for researchers to gain a better understanding of how blood characteristics are altered through the bifurcated arteries, but it also assists physicians in making treatment predictions, which ultimately results in a more effective therapy for the patient. Moreover, in the narrow artery, the blood exhibits a non-Newtonian behaviour at low shear rates. When shear rates exceed 100 s-1 (reciprocal seconds), which is common in big arteries, the Newtonian model is proven to be a good approximation [2]. Consequently, the current study provides an examination of blood flow in a stenosed bifurcated artery. The artery was presumed to be symmetric along its axis, and blood was modeled as Newtonian fluid.

Assuming that there is a mild stenosis in the lumen of the parent aorta, Srinivasacharya and Rao [3] hypothesized that the stenosis had evolved in an axi-symmetric manner to study the flow phenomena in the bifurcated artery caused by a pulsatile pressure gradient. Considering the complex flow phenomena around the curvatures, junctions, and bifurcations of arteries, Zain and Ismail [4] proposed the flow phenomena in an overlapping stenosed bifurcated artery, which may have triggered the deposition of atherosclerotic plaque. In a bifurcated artery, the effect of various types of stenosis on blood flow is investigated, and results are obtained using COMSOL Multiphysics software for the generalised power law model [5] and for the Newtonian model [6, 7]. An investigation of the pulsatile flow occurring in a model of aortic bifurcation in the presence of a pressure gradient has been carried out by Chakravarty and Mandal [8]. The Carreau-Yasuda biorheological blood model was utilized, and a finite element computation on the magnetohydrodynamic flow with heat transfer in a bifurcated artery has been done [9]. Within permeable bifurcated arteries that are subjected to a magnetic field, the influence of heat transfers on magnetohydrodynamic blood flow, observed and addressed as a Newtonian rheological model, was studied in reference [10]. In the rigid model of the carotid bifurcation under pulsatile flow conditions, Meem et al. [11] explored the influence of the Biomagnetic field on pulsatile non-Newtonian features of blood. Specifically, they looked at how these features affected dynamics in the presence of an aneurysm. A study conducted by Shaw et al. [12] showed the impact of shear stress on the bifurcation angle of the Casson model to represent the non-Newtonian characteristics of blood viscosity under pulsatile flow conditions. An investigation into the flow of blood via a bifurcated artery with minor stenosis under variation of gravity acceleration was carried out by considering blood as a micropolar fluid [13] and the couple stress model [14]. Fan et al. [15] studied the pulsatile non-Newtonian flow in the carotid artery bifurcation; they reported that the flow behaviour produced by the Casson model did not differ from the flow characterizations provided by the Newtonian model.

Limited numerical works have been performed, which employed bifurcated arteries, and all these works used software programs or given equations of pressure gradient; these software programs are typically quite pricey and difficult to acquire in terms of availability [16-20]. Furthermore, the pressure distribution in the arterial segment is an unknown parameter that requires an appropriate method for calculating and to provide this, the current work investigates the effect of stenosis on the blood flow pattern in parent artery and main bronchi under laminar Newtonian fluid conditions using the Marker and Cell (MAC) method is advantageous because it eliminates the requirement for pressure boundary conditions at both the inlet and the outflow. The vector representing the velocity is determined, and the results of wall pressure and blood streamline are obtained with the requisite level of precision.

2. Governing Equations

A 2-dimensional, unsteady, incompressible, laminar, Newtonian model examines blood flow in the arterial lumen. In the system of cylindrical polar coordinates (r, x), the axisymmetric blood flow's governing momentum and continuity equations in dimensionless form are expressed as follows:

$\frac{1}{r} \frac{\partial}{\partial r}(r u)+\frac{\partial}{\partial x}(w)=0$                          (1)

$\frac{\partial w}{\partial t}+\frac{\partial(w u)}{\partial r}+\frac{\partial w^2}{\partial x}+\frac{w u}{r}=-\frac{\partial p}{\partial x}+\frac{1}{\operatorname{Re}}\left(\frac{1}{r} \frac{\partial w}{\partial r}\left(r \frac{\partial w}{\partial r}\right)+\frac{\partial^2 w}{\partial x^2}\right)$                         (2)

$\frac{\partial u}{\partial t}+\frac{\partial u^2}{\partial r}+\frac{\partial(u w)}{\partial x}+\frac{u^2}{r}=-\frac{\partial p}{\partial r}+\frac{1}{\operatorname{Re}}\left(\frac{1}{r} \frac{\partial u}{\partial r}\left(r \frac{\partial u}{\partial r}\right)-\frac{u}{r^2}+\frac{\partial^2 u}{\partial x^2}\right)$                 (3)

where, the Reynolds number $R e=\frac{\rho U_0 R_0}{\mu}$, and $t, \rho, \mu, p, U_0, R_0$ respectively are the time, density, viscosity of blood, pressure, the velocity of blood at the inlet, and the radius of the artery. The axial and radial equations of momentum (2), (3) are imposed [13]. The geometry model of this study is a segment of a bifurcated artery with mild stenosis occurred as proposed by studies [7, 8, 12, 13]. The functions R1(x) and R2(x), which stand for the artery's inner and outer walls, respectively (see Figure 1), are provided by:

$R_1(x)=\left\{\begin{array}{ccc}1 & ; & 0 \leq x \leq d ; d+L_0 \leq x \leq x_1 \\ 1-\frac{4 \delta}{L_0^2}\left\{(x-d) L_0-(x-d)^2\right\} & ; & d \leq x \leq d+L_0 \\ 1+r_1-\sqrt{r_1^2-\left(x-x_1\right)^2} & ; & x_1 \leq x<x_2 \\ 2 r_1^{\prime} \sec \beta+\left(x-x_2\right) \tan \beta & ; & x_2 \leq x \leq x_{\max }\end{array}\right.$                       (4)

$R_2(x)= \begin{cases}0, & 0 \leq x<x_3 \\ \sqrt{r_2^2-\left(x-x_3-r_2\right)^2}, & x_3 \leq x<x_3+r_2(1-\sin \beta) \\ r_2 \cos \beta+x_4, & x_3+r_2(1-\sin \beta) \leq x \leq x_{\max }\end{cases}$                         (5)

where,

$r_1=\frac{\left(1-2 r_1^{\prime} \sec \beta\right)}{(\cos \beta-1)} ; r_2=\frac{\left(x_3-x_2\right) \sin \beta}{1-\sin \beta}$

$x_2=x_1+\left(1-2 r_1^{\prime} \sec \beta\right) \frac{\sin \beta}{\cos \beta-1}$

$x_3=x_2+0.5 ; x_4=\left(x-\left(x_3+r_2(1-\sin \beta)\right)\right) \tan \beta$

Figure 1. Axisymmetric geometry for a bifurcated stenosed artery

2.1 Initial and boundary conditions

Due to the rigid wall, a no-slip condition is applied, meaning the blood flow velocity components are zero at the wall. The inlet is assumed to have a fully developed parabolic velocity profile, indicating maximum velocity at the center.

$w(r, x, t)=U_{\max }\left(1-\left(\frac{r}{R_0}\right)^2\right), u(r, x, t)=0$, at $x=0$                       (6)

The problem's boundary and initial conditions are established as follows:

$\begin{aligned} & w(r, x, t)=0=u(r, x, t)\ { on }\ r=R_1(x) \ { and }\ r=R_2(x), x_3 \leq x \leq x_{\max }\end{aligned}$                          (7)

$\frac{\partial w(r, x, t)}{\partial r}=0$ on $r=0,0 \leq x \leq x_3$                       (8)

$w(r, x, t)=u(r, x, t)=p(r, x, 0)=0$, for $x>0$                        (9)

2.2 Radical transformation

The radical transformation is introduced:

$\xi=\frac{r-R_2(x)}{R(x)}, R(x) \neq 0$                       (10)

where, the arterial wall is immobilized in the transformed coordinate by R1(x) and R2(x). Applying Eq. (10), consequently, Eqs. (1)-(3) have the following form:

$\begin{gathered}\left(\xi R+R_2\right) \frac{\partial w}{\partial x}-\frac{\xi R+R_2}{R}\left(\xi \frac{\partial R}{\partial x}+\frac{\partial R_2}{\partial x}\right) \frac{\partial w}{\partial \xi}+\frac{\partial\left(u\left(\xi R+R_2\right) / R\right)}{\partial \xi}=0\end{gathered}$                       (11)

$\frac{\partial w}{\partial t}=-\frac{\partial p}{\partial x}+\frac{1}{R}\left(\xi \frac{\partial R}{\partial x}+\frac{\partial R_2}{\partial x}\right) \frac{\partial p}{\partial \xi}+\operatorname{Con} w+\frac{1}{\operatorname{Re}} \operatorname{Diff} w$                      (12)

$\frac{\partial u}{\partial t}=-\frac{1}{R} \frac{\partial p}{\partial \xi}+\operatorname{Con} u+\frac{1}{\operatorname{Re}} \operatorname{Diff} u$                   (13)

where,

$\begin{aligned} \operatorname{Con} w & =-\frac{1}{R} \frac{\partial(u w)}{\partial \xi}-\frac{u w}{\xi R+R_2}-\frac{\partial w^2}{\partial x} +\frac{1}{R}\left(\xi \frac{\partial R}{\partial x}+\frac{\partial R_2}{\partial x}\right) \frac{\partial w^2}{\partial \xi}\end{aligned}$                      (14)

$\begin{aligned} \text { Diff } w & =\frac{\partial^2 w}{\partial x^2}-\frac{2}{R}\left(\xi \frac{\partial R}{\partial x}+\frac{\partial R_2}{\partial x}\right) \frac{\partial^2 w}{\partial x \xi}+\frac{1}{R^2} \left(1+\left(\xi \frac{\partial R}{\partial x}+\frac{\partial R_2}{\partial x}\right)^2\right) \frac{\partial^2 w}{\partial \xi^2} +\frac{1}{R}\left(\frac{1}{\xi R+R_2}+\frac{2}{R} \frac{\partial R}{\partial x}\left(\xi \frac{\partial R}{\partial x}+\frac{\partial R_2}{\partial x}\right)-\left(\xi \frac{\partial^2 R}{\partial x^2}+\frac{\partial^2 R_2}{\partial x^2}\right)\right) \frac{\partial w}{\partial \xi}\end{aligned}$                    (15)

$\begin{aligned} \operatorname{Con} u & =-\frac{1}{R} \frac{\partial u^2}{\partial \xi}-\frac{u^2}{\xi R+R_2}-\frac{\partial(w u)}{\partial x} +\frac{1}{R}\left(\xi \frac{\partial R}{\partial x}+\frac{\partial R_2}{\partial x}\right) \frac{\partial(w u)}{\partial \xi}\end{aligned}$                    (16)

$\begin{aligned} \operatorname{Diff} u= & \frac{\partial^2 u}{\partial x^2}-\frac{2}{R}\left(\xi \frac{\partial R}{\partial x}+\frac{\partial R_2}{\partial x}\right) \frac{\partial^2 u}{\partial z \xi}+\frac{1}{R^2}\left(1+\left(\xi \frac{\partial R}{\partial x}+\frac{\partial R_2}{\partial x}\right)^2\right) \frac{\partial^2 u}{\partial \xi^2} +\frac{1}{R}\left(\frac{1}{\xi R+R_2}+\frac{2}{R} \frac{\partial R}{\partial x}\left(\xi \frac{\partial R}{\partial x}+\frac{\partial R_2}{\partial x}\right)-\left(\xi \frac{\partial^2 R}{\partial x^2}+\frac{\partial^2 R_2}{\partial x^2}\right)\right) \frac{\partial u}{\partial \xi} -\frac{u}{\left(\xi R+R_2\right)^2}\end{aligned}$                 (17)

Similarly, Eq. (10) is used to transform the initial and boundary conditions (6)-(9) appropriately with $\xi \in[0,1]$.

$w(\xi, x, t)=U_{\max }\left(1-\xi^2\right), u(\xi, x, t)=0$ for $x=0$                  (18)

$\begin{aligned} w(\xi, x, t) & =0=u(\xi, x, t) \text { on } \xi=R_1(x) \ { and }\ \xi =R_2(x), x_3 \leq x \leq x_{\max }\end{aligned}$                   (19)

$\frac{\partial w(\zeta, x, t)}{\partial \zeta}=0$ on $\zeta=0,0 \leq x \leq x_3$                     (20)

$w(\zeta, x, t)=u(\zeta, x, t)=p(\zeta, x, 0)=0$, for $x>0$                     (21)

The MAC method is used to discretize the aforementioned unsteady governing Eqs. (11)-(13) [21]. At various points, the pressure and velocities are calculated. By defining $\zeta= j \Delta \zeta, x=i \Delta x, t=n \Delta t \quad$ and $\quad p(\zeta, x, t)=p(j \Delta \zeta, i \Delta x, n \Delta t)= p_{i, j}^n$, where, $t$ is the time, $\Delta t$ is the increment of $t$, and $\Delta \zeta, \Delta x$ represent the length and width of the control volume cell $(i, j)$ respectively.

2.3 Discretization of governing equations

At the (i, j) th cell, the continuity and momentum equations are discretized to yield:

$\begin{aligned} & \left(\xi_{l j} R_{l i}^n+R_{2 l i}^n\right)\left(\frac{w_{i, j}^n-w_{i-1, j}^n}{\Delta x}\right) -\frac{\xi_{l j} R_{l i}^n+R_{2 l i}^n}{R_{l i}^n}\left(\xi_{l j}\left(\frac{\partial R}{\partial x}\right)_{l i}^n+\left(\frac{\partial R_2}{\partial x}\right)_{l i}^n\right)\left(\frac{w_{a t}-w_{a b}}{\Delta \xi}\right) +\frac{\left(\xi_j R_i^n+R_{2 i}^n\right) u_{i, j}^n-\left(\xi_{j-1} R_i^n+R_{2 i}^n\right) u_{i, j-1}^n}{R_i^n \Delta \xi}=0\end{aligned}$                       (22)

where,

$w_{a t}=\frac{w_{i, j}^n+w_{i-1, j}^n+w_{i-1, j+1}^n+w_{i, j+1}^n}{4}$               (23)

$\begin{aligned} w_{a b} & =\frac{w_{i, j}^n+w_{i-1, j}^n+w_{i-1, j-1}^n+w_{i-1, j-1}^n}{4} \\ \zeta_{l j} & =\zeta_j-\frac{\Delta \zeta}{2}, R^n=R\left(x_{l j}\right), x_{l j}=x_i-\frac{\Delta x}{2}\end{aligned}$                       (24)

$\begin{aligned} \frac{w_{i, j}^{n+1}-w_{i, j}^n}{\Delta t} & =\left(\frac{p_{i, j}^n-p_{i+1, j}^n}{\Delta x}\right) +\frac{1}{R_i^n}\left(\xi_{l j}\left(\frac{\partial R}{\partial x}\right)_i^n+\left(\frac{\partial R_2}{\partial x}\right)_i^n\right)\left(\frac{p_t-p_b}{\Delta \xi}\right) +w c d_{i, j}^n\end{aligned}$                     (25)

where,

$w c d_{i, j}^n=\operatorname{Con} w_{i, j}^n+\frac{1}{\operatorname{Re}}\left(\operatorname{Diff} w_{i, j}^n\right)$                      (26)

$\frac{u_{i, j}^{n+1}-u_{i, j}^n}{\Delta t}=-\frac{1}{R_{l j}^n}\left(\frac{p_{i, j}^n-p_{i, j+1}^n}{\Delta \xi}\right)+u c d_{i, j}^n$                        (27)

where,

$u c d_{i, j}^n=\operatorname{Con} u_{i, j}^n+\frac{1}{\operatorname{Re}}\left(\operatorname{Diff} u_{i, j}^n\right)$                     (28)

with $u c d_{i, j}^n$ and $\operatorname{Diff} u_{i, j}^n$ stand for equation's convective and diffusive terms at the nth time level in the cell (i, j). The terms are differences in a similar manner to the w- equation of momentum. The entire numerical process and stability restriction have already been covered in references [22-24].

3. Results and Discussions

The numerical calculation has been carried out to quantitatively estimate the velocity profiles and wall pressure at the parent and daughter arteries for varying degrees of stenosis severity. For a comprehensive quantitative examination of how Newtonian blood rheology affects the phenomenon of stenotic bifurcated artery flow, the dimensional data for validation purpose have been taken from references [7, 8]: $\mathrm{R}_0=0.0075 \mathrm{~m}, L_0=0.015 \mathrm{~m}, d=0.005 \mathrm{~m}, \rho =105 \mathrm{kgm}^{-3}, \mu=0.0035 \mathrm{Pas}^{-1}, \beta=30^{\circ}$ (angle of branching between left branch and daughter artery), $r_1=0.51 \mathrm{R}_0$ (radius of bronchi), $x_1=0.092 \mathrm{~m}, \delta=0.0 \mathrm{R}_0, \delta=0.25 \mathrm{R}_0, \delta=0.37 \mathrm{R}_0$ (0.00%, 0.4375% and 0.6031% areal occlusions), $\Delta x$ = 0.1 and $\Delta \xi$ = 0.025. The Reynolds number is taken to be 155. Additionally, an arterial finite length xmax = 0.14 m is used to determine the streamlines pattern. The solutions are generated by using a staggered grid of size 1400 × 40 at constant, Re = 300 and U0 = 0.5 × 10⁻2 m·s⁻¹.

The MAC, a finite-difference method, is used to solve the unsteady governing PDE equations of motion. When the simulation reaches a steady state with dimensionless t = 60, the results are determined. After solving the momentum equations, the pressure is used to figure out the velocity profile.

3.1 Grid independence

The grid independence study simulations were run in order to investigate the error related to the grid and time step. Figure 2 shows how the profiles for the three different grid sizes nearly overlap with each other. In the current numerical simulation context, the grid independence study is therefore crucial to proving the accuracy of the outcomes.

Figure 2. Axial velocity profile variation at Re = 300 for various grid sizes

3.2 Numerical validation

For the purpose of numerical validation, it is done by comparing results with Jamali and Zuhail [7], and also Chakravarty and Mandal [8]. Figure 3 indicates the result of dimensionless axial velocity profile achieved using MAC method, which at the maximum constriction. The outcome shows that they are in good agreement after it was brought to the same platform as studies [7, 8]. When the Reynolds number was set to 155, the current method produced the desired result. All results returned parabolic curves at the z-axis, with the highest value roughly. Because of the various methods and inlet velocities used to solve the problem, there is an error between the results, which is very small (Table 1).

Figure 3. Validation of the dimensionless axial velocity [7, 8] at (Re = 155, $\delta$ = 0.4375)

Table 1. Comparison of numerical results of maximum velocity and their errors

Authors

Maximum Velocity

Error

Jamali and Zuhail [7]

0.07194

0.00276

Chakavarty and Mandal [8]

0.06918

-

Present study (MAC)

0.07053

0.00135

3.3 Pressure loss

For the pressure loss, Figure 4 illustrates how the wall pressure varies along the x-axis.

Figure 4. Wall pressure variation with axial position for various $\delta$

For stenosis severity ($\delta$ = 0 normal artery, $\delta$ = 0.4375 occlusion and $\delta$ = 0.6031 occlusion) at Re = 300. As the severity of the stenosis increases, the wall pressure decreases, in other words, due to flow acceleration (velocity increasing), at the location of the stenosis' throat, the wall pressure rapidly drops. Following this position, flow deceleration (velocity decreasing) tends to cause the pressure to recover, but it does not get as high as the starting points. The fact that the wall pressure stays negative on both sides of the apex is an important finding. This happens as a result of secondary flow close to the daughter artery's apex and backflow of blood at the lateral junction. Additionally, as the severity of stenosis increases, the narrowing of the stenotic region prevents blood from passing through the daughter artery, causing the velocity profile to increase up to its peak value. This means the impact of stenosis on pressure losses because there is not enough pressure to keep the lumen open, and the artery may collapse at the throat [25], resulting in a serious pathological state from a clinical standpoint. The heart may have to work harder to maintain a sufficient blood supply if the pressure drop is too large because there might not be enough pressure downstream to perfuse essential organs and tissues. Serious cardiovascular issues, such as heart muscle hypertrophy and ultimately heart failure, may result from this over time.

3.4 Streamlines

The understanding of blood flow patterns through the parent and daughter arteries is essential for understanding localized cardiac diagnostic disorders and aerosolized drug delivery mechanisms. Figures 5-7 show how blood flow streamlines behave based on their axial position for different levels of stenosis with Re = 300. The streamline patterns for blood flow through parent and daughter arteries in Figure 5 are based on a normal artery (δ = 0%). No recirculation zones appear for the straight artery (parent artery). Furthermore, it is observed that the recirculation zones increased in the junction and the outer wall of the daughter artery. This happens as a result of secondary flow in the daughter artery close to the apex and backflow (an increase in velocity) of the streaming blood at the lateral junction's onset.

Figures 6 and 7 show the effect of the severity of the stenosis ($\delta$ = 0.4375 and $\delta$ = 0.6031 areal reduction) on the streamlines behavior of blood flow at Re = 300. The velocity of blood increases as the severity of stenosis escalates, with no recirculation zones downstream of the stenosis for $\delta$ = 0.4375 areal occlusion in a straight artery. Moreover, the elevated values of vortices occurred in the parent and daughter artery with increasing stenosis for $\delta$ = 0.6031 areal occlusion, which aligns well with the experimental results obtained by Ahmed and Giddens [25]. The formation of the eddy depended on the severity of the stenosis. In other words, these recirculation zones are augmented in the outer wall of the artery as the severity of stenosis increases. It is important to pay attention to flow separation and recirculation zones that appear on downstream constrictions. The blood in the low-pressure zone is pulled backward, and it starts to swirl in a circular motion to fill this space. This swirling motion is what forms the recirculation zones. The more severe the stenosis, the faster the blood jet has to move to get through the narrow opening. This higher velocity causes a more pronounced separation from the wall and a larger, more intense low-pressure region. This, in turn, leads to the formation of a larger and more powerful recirculation zone.

These zones can cause the disease to get worse and cause new intimal thickening.

Figure 5. Streamlining of blood flow through the main and a single bifurcated artery for $\delta$ = 0.0

Figure 6. Streamlining of blood flow through the main and a single bifurcated artery for (Re = 300 and $\delta$ = 0.4375 occlusion)

Figure 7. Streamlining of blood flow through the main and a single bifurcated artery for (Re = 300 and $\delta$ = 0.6031 occlusion)

4. Conclusion

This paper presents a mathematical model to simulate Newtonian blood flow through a bifurcated rigid artery with varying degrees of stenosis. The MAC method is developed in MATLAB code to numerically solve the governing unsteady, nonlinear, two-dimensional equations, which were discretized on a uniform staggered grid. The model computed blood flow characteristics, including velocity profiles, pressure loss, and streamlines. The results showed good agreement with previous numerical studies for mild stenosis. A key finding was that pressure loss decreases with increasing stenosis severity, a phenomenon linked to secondary flow and backflow. As stenosis becomes more severe, the increased narrowing leads to a higher peak velocity and the amplification of recirculation zones. The formation of these backflow regions may have significant implications for cardiovascular diseases.

Standard diagnostic tools like angiography or ultrasound provide images, but mathematical models can take these images and calculate precise metrics that are invisible to the naked eye. This includes pressure drop, flow velocity, a streamline across a stenosis. A clinician can use this data to determine the severity of a blockage with greater accuracy than visual inspection alone. By simulating flow patterns, models can identify areas at risk of disease progression even before significant anatomical changes are visible. This can help clinicians intervene earlier to prevent a severe event like a heart attack or stroke.

Ultimately, this model can be extended to be non-Newtonian involves replacing a single, constant parameter (dynamic viscosity) with a more complex relationship that accounts for how blood viscosity changes with flow conditions or by considering the mass transfer equation to investigate aerosolized drug delivery mechanisms in bifurcated arteries.

Acknowledgment

The authors gratefully acknowledge the support of Wasit University and Mustansiriyah University for this research.

  References

[1] Kadhim, S.K., Ali, S.A.G. (2022). Numerical investigation of stenosis and degree of aneurysm on haemodynamics of blood vessels. Mathematical Modelling of Engineering Problems, 9(3): 811-818. https://doi.org/10.18280/mmep.090330

[2] Srinivasacharya, D., Rao, G.M. (2017). Micropolar fluid flow through a stenosed bifurcated artery. Nonlinear Analysis: Modelling and Control, 22(2): 147-159. https://doi.org/10.15388/NA.2017.2.1

[3] Pedley, T.J. (1980). The Fluid Mechanics of Large Blood Vessels (Cambridge Monographs on Mechanics). Cambridge: Cambridge University Press.

[4] Srinivasacharya, D., Rao, G.M. (2018). Pulsatile flow of couple stress fluid through a bifurcated artery. Ain Shams Engineering Journal, 9(4): 883-893. https://doi.org/10.1016/j.asej.2016.04.023

[5] Zain, N.M., Ismail, Z. (2023). Dynamic response of heat transfer in magnetohydrodynamic blood flow through a porous bifurcated artery with overlapping stenosis. Journal of Advanced Research in Fluid Mechanics and Thermal Sciences, 101(1): 215-235. https://doi.org/10.37934/arfmts.101.1.215235

[6] Jamali, M.S.A., Ismail, Z., Amin, N.S. (2021). Effect of different types of stenosis on generalized power law model of blood flow in a bifurcated artery. Journal of Advanced Research in Fluid Mechanics and Thermal Sciences, 87(3): 172-183. https://doi.org/10.37934/arfmts.87.3.172183

[7] Hamid, S.N.A., Thirunanasambantham, K., Ismail, Z., Jiann, L.Y. (2025). Numerical simulation of heat transfer in MHD hybrid blood nanofluid flow through a bifurcated stenosed artery. Defect and Diffusion Forum, 440: 39-45. https://doi.org/10.4028/p-qs5rJe

[8] Jamali, M.S.A., Ismail, Z. (2019). Simulation of heat transfer on blood flow through a stenosed bifurcated artery. Journal of Advanced Research in Fluid Mechanics and Thermal Sciences, 60(2): 310-323.

[9] Chakravarty, S., Mandal, P.K. (1997). An analysis of pulsatile flow in a model aortic bifurcation. International Journal of Engineering Science, 35(4): 409-422. https://doi.org/10.1016/S0020-7225(96)00081-X

[10] Dubey, A., Vasu, B., Bég, O.A., Gorla, R.S.R. (2021). Finite element computation of magneto-hemodynamic flow and heat transfer in a bifurcated artery with saccular aneurysm using the Carreau-Yasuda biorheological model. Microvascular Research, 138: 104221. https://doi.org/10.1016/j.mvr.2021.104221

[11] Meem, A.T., Hossain, M.Z., Molla, M.M., Saha, S.C., Saha, G. (2025). Biomagnetic field effects on pulsatile non-Newtonian blood flow with gold nanoparticles in a bifurcated artery in the presence of aneurysm. Pramana, 99(3): 75. https://doi.org/10.1007/s12043-025-02930-7

[12] Shaw, S., Gorla, R.S.R., Murthy, P.V.S.N., Ng, C.O. (2009). Pulsatile Casson fluid flow through a stenosed bifurcated artery. International Journal of Fluid Mechanics Research, 36(1): 43-63. https://doi.org/10.1615/interjfluidmechres.v36.i1.30

[13] Tan, Y.B., Mustapha, N. (2018). Gravitational influences on micropolar blood flow in a bifurcated artery with mild stenosis. Journal of Advanced and Applied Science, 5(11): 24-32. https://doi.org/10.21833/ijaas.2018.11.003

[14] Srinivasacharya, D., Rao, G.M. (2016). Mathematical model for blood flow through a bifurcated artery using couple stress fluid. Mathematical Biosciences, 278: 37-47. https://doi.org/10.1016/j.mbs.2016.05.003

[15] Fan, Y., Jiang, W., Zou, Y., Li, J., Chen, J., Deng, X. (2009). Numerical simulation of pulsatile non-Newtonian flow in the carotid artery bifurcation. Acta Mechanica Sinica, 25(2): 249-255. https://doi.org/10.1007/s10409-009-0227-9

[16] Ismail, Z., Jamali, M.S.A., Zain, N.M., Jiann, L. Y. (2022). Bioheat transfer of blood flow on healthy and unhealthy bifurcated artery: Stenosis. Journal of Advanced Research in Applied Sciences and Engineering Technology, 28(2): 56-79. https://doi.org/10.37934/araset.28.2.5679

[17] Thirunanasambantham, K., Ismail, Z., Jiann, L.Y., Shamjuddin, A. (2023). Numerical computational of blood flow and mass transport in stenosed bifurcated artery. Journal of Advanced Research in Fluid Mechanics and Thermal Sciences, 110(2): 79-94. https://doi.org/10.37934/arfmts.110.2.7994

[18] Zain, N.M., Ismail, Z. (2023). Numerical solution of magnetohydrodynamics effects on a generalised power law fluid model of blood flow through a bifurcated artery with an overlapping shaped stenosis, PLoS ONE, 18(2): e0276576. https://doi.org/10.1371/journal.pone.0276576

[19] Zain, N.M., Ismail, Z., Johnston, P. (2021). A stabilized finite element formulation of non-Newtonian fluid model of blood flow in a bifurcated channel with overlapping stenosis. Journal of Advanced Research in Fluid Mechanics and Thermal Sciences, 88(1): 126-138. https://doi.org/10.37934/arfmts.88.1.126139

[20] Welch, J.E., Harlow, F.H., Shannon, J.P., Daly, B.J. (1965). The MAC method-A computing technique for solving viscous, incompressible, transient fluid-flow problems involving free surfaces. Technical Report. https://doi.org/10.2172/4563173

[21] Bakheet, A., Alnussaiyri, E.A., Ismail, Z. Amin, N. (2016). Blood flow through an inclined stenosed artery. Applied Mathematical Sciences, 10(5-8): 235-254. https://doi.org/10.12988/ams.2016.511701

[22] Alnussairy, E.A., Bakheet, A. (2019). MHD micropolar blood flow model through a multiple stenosed artery. AIP Conference Proceedings, 2183(1): 090002. https://doi.org/10.1063/1.5136202

[23] Bakheet, A., Alnussairy, E.A. (2021). Numerical simulation of magnetohydrodynamic influences on Casson model for blood flow through an overlapping stenosed artery. Iraqi Journal of Science, 62(3): 1016-1024. https://doi.org/10.24996/ijs.2021.62.3.30

[24] Rabby, M.G., Shupti, S.P., Molla, M.M. (2014). Pulsatile non-Newtonian laminar blood flows through arterial double stenoses. Journal of Fluids, 2014(1): 757902 https://doi.org/10.1155/2014/757902

[25] Ahmed, S.A., Giddens, D.P. (1983). Velocity measurements in steady flow through axisymmetric stenoses at moderate Reynolds numbers. Journal of Biomechanics, 16(7): 505-516. https://doi.org/10.1016/0021-9290(83)90065-9