Prediction of Hemodynamic-Related Hemolysis in Carotid Stenosis and Aiding in Treatment Planning and Risk Stratification Using Computational Fluid Dynamics

Atherosclerosis affects human health in many ways, leading to disability or premature death due to ischemic heart disease, stroke, or limb ischemia. Poststenotic blood flow disruption may also play an essential role in artery wall impairment linked with hemolysis related to shear stress. The maximum shear stress in the atherosclerotic plaque area is the main parameter determining hemolysis risk. In our work, a 3D internal carotid artery model was built from CT scans performed on patients qualified for percutaneous angioplasty due to its symptomatic stenosis. The obtained stenosis geometries were used to conduct a series of computer simulations to identify critical parameters corresponding to the increase in shear stress in the arteries. Stenosis shape parameters responsible for the increase in shear stress were determined. The effect of changes in the carotid artery size, length, and degree of narrowing on the change in maximum shear stress was demonstrated. Then, a correlation for the quick initial diagnosis of atherosclerotic stenoses regarding the risk of hemolysis was developed. The developed relationship for rapid hemolysis risk assessment uses information from typical non-invasive tests for treated patients. Practical guidelines have been developed regarding which stenosis shape parameters pose a risk of hemolysis, which may be adapted in medical practice.


Introduction
Despite the continuous development of medicine, modern society still suffers from cardiovascular diseases [1][2][3].Today's lifestyle, combined with a lack of physical activity and an unhealthy diet, often leads to hypertension or atherosclerosis [4,5].Atherosclerosis is responsible for cardiovascular diseases leading to stroke or infarction.The increase in cholesterol deposits causes a decrease in the surface area perpendicular to the flow direction, which in turn causes a local increase in shear stresses.An increase in shear stress can destroy erythrocytes [6][7][8].This local hemolysis causes impairment of the endothelium and smooth muscle cells related to NO decrease and direct toxicity of hemoglobin and Fe ++ released from erythrocytes [9].Endothelial abrasion may expose deeper layers of the arterial wall rich with collagen to circulated platelets and may lead to its activation and clot formation.Hemolysis occurs in atherosclerosis from very early to advanced stages of artery stenosis [9].The basis of atherosclerosis treatment is to cease this pathological process as soon as possible.Guidelines related to endovascular interventions in cardiovascular diseases are focused mainly on the severity of artery stenosis before angioplasty (i.e., carotid internal artery or coronal arteries).One of the most common treatment methods for patients with advanced carotid or coronary atherosclerosis is percutaneous endovascular intervention with balloon angioplasty followed by stent placement.A better understanding of atherosclerosis hemolysis etiopathology may help to establish a new guideline with interventions performed at an earlier stage of the disease.Due to the potential risk of complications, it is worth having information about the impact of the narrowing on the hydrodynamics of blood flow in the vessel.Hence, the idea is to use non-invasive blood flow measurements combined with advanced imaging and computational fluid dynamics (CFD) for a non-invasive diagnosis of atherosclerosis [10].Appropriate diagnosis of atherosclerotic complications using CFD requires attention to blood rheology, which is a non-Newtonian fluid [11][12][13][14][15][16][17] due to the complex structure of blood.Due to the non-Newtonian nature of blood, a blood rheology model based on population balance [8] was used, which reconstructs the physiological behavior of blood related to erythrocyte agglomeration and deagglomeration depending on the shear rate.In addition to blood rheology, selecting an appropriate hemolysis model is important [6,7,[18][19][20][21][22][23][24][25].The blood rheology model based on the erythrocyte population balance has a built-in blood hemolysis model that allows for the use of information obtained from solving the size distribution of red blood cell agglomerates and the number of single erythrocytes [8].In addition, it is crucial to select an appropriate flow turbulence model that would be useful both in terms of laminar and transient flow observed in the zone after the stenosis.Printed channel models and micro particle image velocimetry (µPIV) [10,26] were used to validate the turbulence model used in the CFD simulations.The simulations were carried out for steady flows because the difference between the values of average wall shear stresses (WSSs) for average flow and time-averaged wall shear stresses (TAWSs) for pulsating flow is minimal [27] for stenosed arteries, as shown in one of the previous articles.The analysis of the influence of the shape of the blood vessel on the risk of hemolysis was started in the last article [10], where the focus was on the degree of stenosis of the vessel and its eccentricity.As part of this work, the analysis was extended to include variable vessel size, stenosis length, and variable lumen shape in the narrowest part of the stenosis.
In this article, we wanted to address the research gap in the field of rapid initial diagnosis of atherosclerotic complications regarding the risk of blood hemolysis in the space of cholesterol constrictions.The aim of the study was to identify key parameters of the shape of the atherosclerotic stenosis that influence the increase in shear stress to which red blood cells are exposed.The universal relationship between shear stresses and the geometry of the constriction fits into modern trends in personalized medicine, while increasing the efficiency of the health care system in the case of a growing number of patients.

Blood Rheology and Hemolysis Model
The blood viscosity model based on population balance was presented in detail in various publications [8,10,27].The model assumes agglomeration of red blood cells in areas with low shear rates and deagglomerations in regions with high shear rates.As a result, a physiological distribution of the size of red blood cell agglomerates can be observed across the cross-section of the vessel so that the viscosity is greatest in the axis of the vessel.In addition, thixotropic effects of blood [28][29][30][31] related to the finite rate of deagglomeration and agglomeration of erythrocytes are observed.The hemolysis model originates from power law formulation and later linearization [6,7,32].However, population-balance-based rheology's (PBBR's) model refers to hemolysis to only fully deagglomerated single red blood cells, whose concentration is known thanks to the solution of population balance [8].
where m s is the concentration of single red blood cells in direct quadrature method of moments (DQMOM), L i is the sizes of agglomerates used in DQMOM, and w i is the weights used in DQMOM.

Numerical Settings
The calculations in this study were performed using ANSYS Fluent 2022 R2.Model implementation was presented in detail in the previous article [8].Similar to previous studies, the flow is mainly laminar; however, the turbulent flow is promoted by the rapid increase in velocity in the stenosis in the case of large flows [10].The GEKO model with the transition flow option was chosen to capture both laminar and turbulent regions.The simulations were conducted for blood flow rates in the 0.625-6.25 mL/s range, similar to the previous article [10].Simulations for different polyhedral mesh densities were performed to determine the final meshes for which the obtained results were independent of mesh size and quality.
The following boundary conditions were applied.At the INLET, the plug flow inlet was located at an appropriate distance from the stenosis to establish the fully developed velocity profile.Inlet conditions for the population balance were calculated, assuming that the hematocrit equaled 45%.It was also assumed that no hemolyzed red blood cells were at the inlet.In the case of WALLS, the no-slip condition and homogeneous Neumann conditions for user-defined scalars were used.At the OUTLET, the atmospheric pressure and homogeneous Neumann conditions for user-defined scalars were chosen.

CFD Results
The CFD analysis was carried out for an inlet velocity of 0.05 to 0.5 m/s similarly to the previous article.Within the framework of this article, the range of analyzed geometries of cholesteric stenoses was extended by examining the influence of changes in vessel diameter, stenosis length, and cross-sectional shape.The geometry model is shown in Figure 1 and detail dimensions for analyzed cases are summarized in Table 1.For artery size analysis, the base model taken from the previous article for stenosis with (A c /A 0 ) 0.5 = 0.4 [10] was scaled proportionally in three dimensions in the 50-125% range.However, in the case of the throat length analysis, the length of the throat was changed in the range of 4 to 16 mm.In addition, in the case of the shape analysis of the influence of changes in the shape of the throat in the cross-section of the throat, the baseline model had (A c /A 0 ) 0.5 = 0.3 [10], and the shape changes were obtained by asymmetric neck scaling in two directions to maintain a constant cross-sectional area.The prepared geometries were then used to conduct a detailed analysis of the impact of these three shape parameters on the maximum shear stresses and the hemolysis value.The obtained results were presented in graphs, and the tables with the results are available in Appendix A (Tables A1-A6).
single red blood cells, whose concentration is known thanks to the solution of population balance [8].
where  is the concentration of single red blood cells in direct quadrature method of moments (DQMOM),  is the sizes of agglomerates used in DQMOM, and  is the weights used in DQMOM.

Numerical Settings
The calculations in this study were performed using ANSYS Fluent 2022 R2.Model implementation was presented in detail in the previous article [8].Similar to previous studies, the flow is mainly laminar; however, the turbulent flow is promoted by the rapid increase in velocity in the stenosis in the case of large flows [10].The GEKO model with the transition flow option was chosen to capture both laminar and turbulent regions.The simulations were conducted for blood flow rates in the 0.625-6.25 mL/s range, similar to the previous article [10].Simulations for different polyhedral mesh densities were performed to determine the final meshes for which the obtained results were independent of mesh size and quality.
The following boundary conditions were applied.At the INLET, the plug flow inlet was located at an appropriate distance from the stenosis to establish the fully developed velocity profile.Inlet conditions for the population balance were calculated, assuming that the hematocrit equaled 45%.It was also assumed that no hemolyzed red blood cells were at the inlet.In the case of WALLS, the no-slip condition and homogeneous Neumann conditions for user-defined scalars were used.At the OUTLET, the atmospheric pressure and homogeneous Neumann conditions for user-defined scalars were chosen.

CFD Results
The CFD analysis was carried out for an inlet velocity of 0.05 to 0.5 m/s similarly to the previous article.Within the framework of this article, the range of analyzed geometries of cholesteric stenoses was extended by examining the influence of changes in vessel diameter, stenosis length, and cross-sectional shape.The geometry model is shown in Figure 1 and detail dimensions for analyzed cases are summarized in Table 1.For artery size analysis, the base model taken from the previous article for stenosis with (  ⁄ ) .= 0.4 [10] was scaled proportionally in three dimensions in the 50-125% range.However, in the case of the throat length analysis, the length of the throat was changed in the range of 4 to 16 mm.In addition, in the case of the shape analysis of the influence of changes in the shape of the throat in the cross-section of the throat, the baseline model had (  ⁄ ) .= 0.3 [10], and the shape changes were obtained by asymmetric neck scaling in two directions to maintain a constant cross-sectional area.The prepared geometries were then used to conduct a detailed analysis of the impact of these three shape parameters on the maximum shear stresses and the hemolysis value.The obtained results were presented in graphs, and the tables with the results are available in Appendix A (Tables A1-A6).Based on Figure 2 and Table A1, a strong correlation can be seen between vessel size and maximum shear stresses.The smaller the blood vessel, the greater the shear stress at the same linear velocity for the cross-section before the constriction.Consequently, the smaller the vessel, the lower the limiting velocity at which hemolysis of red blood cells can occur.2 and Table A1, a strong correlation can be seen between vessel size and maximum shear stresses.The smaller the blood vessel, the greater the shear stress at the same linear velocity for the cross-section before the constriction.Consequently, the smaller the vessel, the lower the limiting velocity at which hemolysis of red blood cells can occur.Changes in maximum shear stresses are of crucial importance in assessing the risk of hemolysis; despite the increase in maximum stresses by about 50% (Figure 2 and Table A1) for the smallest artery at maximum flow, in the case of hemolysis, the difference reaches 193% (Figure 3 and Table A2) or even higher for lower velocities.This is related to both the increase in maximum stresses and the fluid volume in which the maximum shear stresses are exceeded, above which hemolysis is observed.Changes in maximum shear stresses are of crucial importance in assessing the risk of hemolysis; despite the increase in maximum stresses by about 50% (Figure 2 and Table A1) for the smallest artery at maximum flow, in the case of hemolysis, the difference reaches 193% (Figure 3 and Table A2) or even higher for lower velocities.This is related to both the increase in maximum stresses and the fluid volume in which the maximum shear stresses are exceeded, above which hemolysis is observed.
Figure 4 and Table A2 show that the differences in maximum shear stresses are even more visible in the case of changes in the throat length.For the shortest stenosis, the increase in shear stress exceeds 100% compared to the base geometry with a stenosis length of 16 mm.Therefore, the shorter the throat, the greater the velocity acceleration and the greater the shear stress.Due to the rapid changes in velocity, apart from the increased risk of hemolysis, there is also a greater risk of damage to the vessel wall and the release of deposits, which can lead to embolism, which is very dangerous for the patient.
As in the case of changes in the size of the artery, changes in the hemolysis value are also visible in the case of changes in the length of the stenosis.However, despite much more observable changes in shear stresses for changes in the narrowing length, the changes in the hemolysis value are not so large.In the case of maximum flow and minimum length, they reach 74% (Figure 5 and Table A4).However, these changes are much greater for lower velocities, which significantly reduces the critical velocity at which hemolysis begins.Figure 4 and Table A2 show that the differences in maximum shear stresses are even more visible in the case of changes in the throat length.For the shortest stenosis, the increase in shear stress exceeds 100% compared to the base geometry with a stenosis length of 16 mm.Therefore, the shorter the throat, the greater the velocity acceleration and the greater the shear stress.Due to the rapid changes in velocity, apart from the increased risk of hemolysis, there is also a greater risk of damage to the vessel wall and the release of deposits, which can lead to embolism, which is very dangerous for the patient.As in the case of changes in the size of the artery, changes in the hemolysis value are also visible in the case of changes in the length of the stenosis.However, despite much more observable changes in shear stresses for changes in the narrowing length, the changes in the hemolysis value are not so large.In the case of maximum flow and minimum length, they reach 74% (Figure 5 and Table A4).However, these changes are much greater for lower velocities, which significantly reduces the critical velocity at which hemolysis begins.A2 show that the differences in maximum shear stresses are even more visible in the case of changes in the throat length.For the shortest stenosis, the increase in shear stress exceeds 100% compared to the base geometry with a stenosis length of 16 mm.Therefore, the shorter the throat, the greater the velocity acceleration and the greater the shear stress.Due to the rapid changes in velocity, apart from the increased risk of hemolysis, there is also a greater risk of damage to the vessel wall and the release of deposits, which can lead to embolism, which is very dangerous for the patient.As in the case of changes in the size of the artery, changes in the hemolysis value are also visible in the case of changes in the length of the stenosis.However, despite much more observable changes in shear stresses for changes in the narrowing length, the changes in the hemolysis value are not so large.In the case of maximum flow and minimum length, they reach 74% (Figure 5 and Table A4).However, these changes are much greater for lower velocities, which significantly reduces the critical velocity at which hemolysis begins.In contrast to the previous two shape parameters, the change in the shape of the cross-section has a much smaller impact on the maximum shear stresses, which can be seen both in Figure 6 and in the data from Table A5.Consequently, for most stenoses, the influence of this factor can be neglected unless the stenosis is very strongly flattened.In contrast to the previous two shape parameters, the change in the shape of the cross-section has a much smaller impact on the maximum shear stresses, which can be seen both in Figure 6 and in the data from Table A5.Consequently, for most stenoses, the influence of this factor can be neglected unless the stenosis is very strongly flattened.In contrast to the previous two shape parameters, the change in the shape of the cross-section has a much smaller impact on the maximum shear stresses, which can be seen both in Figure 6 and in the data from Table A5.Consequently, for most stenoses, the influence of this factor can be neglected unless the stenosis is very strongly flattened.Despite slight differences in the maximum shear stresses observed in Figure 6, in the case of hemolysis, the differences are slightly larger, as seen in Figure 7 and Table A6.The increase in hemolysis reaches 49% for maximum flattening at maximum flow.However, in the case of lower speeds, the differences are more significant, similar to the previous issues.Despite slight differences in the maximum shear stresses observed in Figure 6, in the case of hemolysis, the differences are slightly larger, as seen in Figure 7 and Table A6.The increase in hemolysis reaches 49% for maximum flattening at maximum flow.However, in the case of lower speeds, the differences are more significant, similar to the previous issues.Based on the obtained results, it was decided to develop tables of correction factors for the equation for maximum shear stress.The resulting equation can pre-estimate the expected maximum shear stresses for a quick initial cholesteric stenosis analysis using typical medical data from non-invasive medical procedures.For this purpose, auxiliary dimensionless parameters for reading correction factors were introduced.
The first parameter is the size factor .It is defined as the ratio of the size of the analyzed geometry  to the geometric dimensions of the base geometry  with a given degree of stenosis The second parameter is the dimensionless stenosis length , which is similarly defined as the ratio of the length of the analyzed stenosis  and the base model hydraulic diameter  , taking into account the scale parameter z =   .
(3) Based on the obtained results, it was decided to develop tables of correction factors for the equation for maximum shear stress.The resulting equation can pre-estimate the expected maximum shear stresses for a quick initial cholesteric stenosis analysis using typical medical data from non-invasive medical procedures.For this purpose, auxiliary dimensionless parameters for reading correction factors were introduced.
The first parameter is the size factor k. It is defined as the ratio of the size of the analyzed geometry x to the geometric dimensions of the base geometry X with a given degree of stenosis The second parameter is the dimensionless stenosis length z, which is similarly defined as the ratio of the length of the analyzed stenosis l and the base model hydraulic diameter D h0 , taking into account the scale parameter The last parameter is the dimensionless stenosis hydraulic diameter w, which is similarly defined as the ratio of the hydraulic diameter of the stenosis of the analyzed stenosis d h and the base model D h , taking into account the scale parameter Considering the correction factors, the modified equation for the maximum shear stress is presented in Equation ( 5).However, the correction parameter related to the shape of the neck section can be omitted in most cases: ( where τ max is maximum shear stress in the analyzed case, τ 0 is maximum shear stress for baseline artery stenosis geometry [10], and β 1 , β 2 , β 3 are correction coefficients.Correction factors are listed in Tables 2-4.The presented model has been tested in the context of application in the diagnosis of carotid artery stenoses.The geometry of the stenosis (Figure 8) was obtained as a result of the routine patient's CT scans performed before endovascular treatment (angioplasty with stent placement).The hemolysis risk was checked at rest [33] and during physical exertion [34] corresponding to physical work or physical activity.Common carotid blood flow velocity data during graded exercise on a treadmill [34] were used and then recalculated to internal carotid artery (ICA) flow for physical activity conditions.
Biomedicines 2024, 12, x FOR PEER REVIEW 9 [34] corresponding to physical work or physical activity.Common carotid blood flow locity data during graded exercise on a treadmill [34] were used and then recalculat internal carotid artery (ICA) flow for physical activity conditions.The geometry of the artery was obtained after geometry processing in MeshM and SpaceClaim to remove geometry distortions resulting from the finite resolution o apparatus.The hydraulic diameter at the vessel inlet was 7.6 mm, and the distanc tween the beginning of the stenosis and the maximum stenosis was 10 mm.The meas degree of stenosis was (  ⁄ ) .= 0.29.The length of the entire geometry in the rection (Figure 9) was 50 mm.In addition, a fully developed velocity profile at the in the stenotic carotid artery was assumed.Dimensionless geometry parameters were c lated according to Equations ( 2) and (3).Coefficients β1 and β2 were interpolated u spline function Matlab software (R2022b) with respect to velocities.Later, β1 was extr lated to adjust the scale factor outside the previously considered region; however, p interpolation presented a high R 2 , higher than 0.99 for the analyzed region.In the ca β2, another spline interpolation with dimensionless stenosis length was in the table range.All results are summed up in Table 5.The geometry of the artery was obtained after geometry processing in MeshMixer and SpaceClaim to remove geometry distortions resulting from the finite resolution of the apparatus.The hydraulic diameter at the vessel inlet was 7.6 mm, and the distance between the beginning of the stenosis and the maximum stenosis was 10 mm.The measured degree of stenosis was (A c /A 0 ) 0.5 = 0.29.The length of the entire geometry in the Z direction (Figure 9) was 50 mm.In addition, a fully developed velocity profile at the inlet to the stenotic carotid artery was assumed.Dimensionless geometry parameters were calculated according to Equations ( 2) and (3).Coefficients β 1 and β 2 were interpolated using spline function Matlab software (R2022b) with respect to velocities.Later, β 1 was extrapolated to adjust the scale factor outside the previously considered region; however, power interpolation presented a high R 2 , higher than 0.99 for the analyzed region.In the case of β 2 , another spline interpolation with dimensionless stenosis length was in the table data range.All results are summed up in Table 5.
The mean velocity in the artery before the stenosis was obtained from the Fluent software (2023 R1).The maximum shear stresses for the base model τ 0 were then calculated based on Equation ( 6) introduced in the previous article [10]: where α = 24.21,β = 1.346, γ = −3.385,and u stands for velocity before stenosis.
lated according to Equations ( 2) and ( 3).Coefficients β1 and β2 were interpolated using spline function Matlab software (R2022b) with respect to velocities.Later, β1 was extrapolated to adjust the scale factor outside the previously considered region; however, power interpolation presented a high R 2 , higher than 0.99 for the analyzed region.In the case of β2, another spline interpolation with dimensionless stenosis length was in the table data range.All results are summed up in Table 5.The correction factor β 1 has a small value due to the high value of the scale factor k because the diameter of the vessel before the stenosis has a large diameter resulting from the presence of an aneurysm.On the other hand, the correction factor β 2 is greater than unity because the obtained value of the dimensionless stenosis length z is slightly smaller than for the scaled base model.And, this shape parameter was calculated by comparing the distance between the beginning of the necking and the maximum necking due to the large differences in diameters before and after the stenosis.Based on the maximum shear stresses from the correlation τ 0 for the base model and knowing the correction factors, the expected maximum stresses τ max were calculated.The results obtained from the correlation were compared with the data obtained from the CFD simulation.In the case of flow at rest, the error was several percent; however, in the case of blood flow during exercise, this error decreased to about 10% in the case of average flow and a few percent in the case of maximum flow.
The wall shear stresses are presented in Figure 9 to highlight the difference between rest and exercise conditions and also to present the areas where shear stress exceeds the shear threshold above which hemolysis occurs.It can be seen that increased blood flow is in line with increased wall shear stress.It is worth mentioning that, for mean flow, even during exercise, wall shear stresses are lower than 150 Pa [35][36][37][38], but during systole peak, there is a vast area where hemolysis may occur.The presented results are supplemented by Figure 10, where the hemolysis value can be seen.The largest fraction of hemolyzed blood cells can be observed in the constriction, where the greatest shear stresses occur.
relation were compared with the data obtained from the CFD simulation.In the case o flow at rest, the error was several percent; however, in the case of blood flow during exe cise, this error decreased to about 10% in the case of average flow and a few percent in th case of maximum flow.
The wall shear stresses are presented in Figure 9 to highlight the difference betwee rest and exercise conditions and also to present the areas where shear stress exceeds th shear threshold above which hemolysis occurs.It can be seen that increased blood flow in line with increased wall shear stress.It is worth mentioning that, for mean flow, eve during exercise, wall shear stresses are lower than 150 Pa [35][36][37][38], but during systole peak there is a vast area where hemolysis may occur.The presented results are supplemente by Figure 10, where the hemolysis value can be seen.The largest fraction of hemolyze blood cells can be observed in the constriction, where the greatest shear stresses occur.The velocity vectors presented in Figure 11 show that, in the stenosis region, velocity increases rapidly, forming a jet structure.The velocity profile in the stenosis throat is more flattened for higher blood flows.The stagnation zone is also visible and increases with higher blood flow.The velocity vectors presented in Figure 11 show that, in the stenosis region, velocity increases rapidly, forming a jet structure.The velocity profile in the stenosis throat is more flattened for higher blood flows.The stagnation zone is also visible and increases with higher blood flow.
This study analyzed three shape parameters regarding their influence on the increase in maximum shear stresses and hemolysis.A strong correlation was observed between a decrease in vessel size and an increase in maximum shear stresses.There was a clear decrease in the critical velocity at which hemolysis occurs with a reduction in the diameter of the vessel.The increase in stress was largely reflected in the increase in hemolysis observed for smaller vessels.A similar rise in maximum shear stresses was observed for shorter constrictions, where the change in local velocity was faster.Stenosis with a sharp decrease in the vessel's lumen is particularly dangerous due to the high stresses acting on the tissue, which may cause the uncontrolled release of cholesterol deposits, causing embolism.Stenoses of this shape are particularly dangerous during physical activity, resulting from physical work and sports.Patients exposed to such cholesterol changes are often still professionally active.Therefore, rapid diagnosis is very important.As far as hemolysis is concerned, the increase in shear stress observed for shorter constrictions translates into an increase in hemolysis.However, the increase is not as large as in the case of a decrease in size.This is due to a different distribution of shear stresses.The last parameter of the shape was the degree of flattening of the vessel in the narrowing, which translated into a decrease in the hydraulic diameter.In the case of shear stresses, an increase in shear stresses is observed only for strongly flattened narrowing, which is rarely observed.As in the other cases, the increase in shear stress translates into an increase in hemolysis.However, this increase is smaller than in the other analyzed cases.Based on the obtained results, Tables 2-4 present correction factors for the correlation, allowing us to estimate the maximum shear stresses presented in the previous article [10].The obtained correlation was tested on the internal carotid artery stenosis and compared with the results of the CFD simulation; the error between the equation and the full simulation did not exceed a dozen or so percent and decreased with the increase in flow; for the maximum flow during exercise, it was only a few percent.The CFD simulations for the internal carotid artery show that increased blood flow leads to higher hemolysis, hence the need for prompt treatment, especially in physically active people, as they are more exposed to hemolysis of blood in the area of constrictions.

Conclusions
A helpful correlation was established between maximum shear stresses, a significant indicator of hemolysis risk, and blood flow and stenosis geometry.It has been shown that rapid changes in the cross-section of the vessel cause a greater risk of hemolysis than the same narrowing with a mild change in diameter on a larger section of the artery.Moreover, it has been shown that with a reduction in the vessel diameter, the critical flow velocity above which the risk of hemolysis occurs also decreases.The obtained correlation was tested on the actual geometry of a carotid atherosclerosis patient.The obtained results can be used to develop new guidelines for diagnosing atherosclerosis as a tool for a quick initial assessment of the patient's condition.

Appendix A
The main equation describing the behavior of fluids during flow is the Navier-Stokes equation.Its general form for an incompressible fluid is presented below.
where ν is kinematic viscosity.For laminar flows, this equation is sufficient for solving the flow field.However, for transient and turbulent flows, additional equations are needed due to the presence of time-varying turbulence affecting the flow field.A common solution to this problem is to average the Navier-Stokes equations based on Reynolds decompositions.In the Cartesian coordinates system, the Reynolds-averaged Navier-Stokes equations are as follows: In addition to the RANS equation, equations describing turbulence are added to close the system of equations.The two most popular models using this approach are k − ϵ and k − ω and their derivatives.
Due to the complicated rheology of blood, it was necessary to use a non-Newtonian blood model.It was decided to use the PBBR model, which uses population balance to reflect the physiological processes of the agglomeration and deagglomeration of blood cells, which influence the change in the apparent fraction of red blood cell agglomerates in the cross-section of the artery.
The viscosity model used in this article can be described as follows: where µ 0 = 6.99 • 10 −4 Pa•s is the viscosity obtained by fitting viscosity experimental data.Equation (A4) presents the relative blood viscosity M(ϕ): And, Equation (A5) shows the Peclet number formula: Furthermore, the relative blood viscosity for high shear rates M ∞ (ϕ) is defined as follows: And, the relative blood viscosity for low shear rates M 0 (ϕ) is defined in Equation (A7): where .
γ is the shear rate, a is the characteristic particle size for the whole population of agglomerates and is equal to the volume-weighted mean particle size, ϕ is the effective apparent volume fraction of red blood cells, D M (φ) is the effective diffusion coefficient [8]  for molecular diffusion, C ∼ = 2 and is an equation constant, and ϕ * is the maximum apparent volume fraction equal to 0.695 [8,94,95].
Furthermore, the hemolysis model was presented in detail in previous publications [8,10,27]; however, the main equation describing hemolysis is as follows: where D EFF is the effective diffusion coefficient for molecular and turbulent diffusion and Moreover, based on various articles [35][36][37][38], the value of shear threshold τ s was set to 150 Pa.
The numerical mesh was made in Fluent Meshing.A polyhedral mesh with a boundary layer was generated.A mesh independence test was performed for the carotid artery model for mean exercise blood flow for five different meshes ranging from several hundred thousand elements to several million cells.The results of the mesh independence test are presented in Figure A1, where the pressure drop is compared depending on the mesh size.A mesh with 2.15 million cells was chosen because the difference in pressure drop compared to the densest mesh was below 1%.The final mesh had a maximum cell size set to 0.2 mm, and the boundary layer was generated with a uniform algorithm for 12 layers with a growth rate of 1.2 and a first layer height of 0.01 mm.
The tables summarize the results of CFD calculations presented in Figures 2-7.
where  is the effective diffusion coefficient for molecular and turbulent diffusion an  = 0,   <  1,   ≥  .
(A Moreover, based on various articles [35][36][37][38], the value of shear threshold  was s to 150 Pa.The numerical mesh was made in Fluent Meshing.A polyhedral mesh with a boun ary layer was generated.A mesh independence test was performed for the carotid arter model for mean exercise blood flow for five different meshes ranging from several hu dred thousand elements to several million cells.The results of the mesh independence te are presented in Figure A1, where the pressure drop is compared depending on the mes size.A mesh with 2.15 million cells was chosen because the difference in pressure dro compared to the densest mesh was below 1%.The final mesh had a maximum cell size s to 0.2 mm, and the boundary layer was generated with a uniform algorithm for 12 laye with a growth rate of 1.2 and a first layer height of 0.01 mm.
The tables summarize the results of CFD calculations presented in Figures 2-7.

Figure 1 .
Figure 1.Geometry of an artery with stenosis.Figure 1. Geometry of an artery with stenosis.

Figure 1 .
Figure 1.Geometry of an artery with stenosis.Figure 1. Geometry of an artery with stenosis.

Figure 2 .
Figure 2. Correlation between inlet hydraulic diameter [mm] and mean velocity [m/s] for maximum shear stress [Pa].

Figure 2 .
Figure 2. Correlation between inlet hydraulic diameter [mm] and mean velocity [m/s] for maximum shear stress [Pa].

Figure 4 .
Figure 4. Correlation between stenosis length [mm] and mean velocity [m/s] for maximum shear stress [Pa].

Figure 4
Figure 4 and TableA2show that the differences in maximum shear stresses are even more visible in the case of changes in the throat length.For the shortest stenosis, the increase in shear stress exceeds 100% compared to the base geometry with a stenosis length of 16 mm.Therefore, the shorter the throat, the greater the velocity acceleration and the greater the shear stress.Due to the rapid changes in velocity, apart from the increased risk of hemolysis, there is also a greater risk of damage to the vessel wall and the release of deposits, which can lead to embolism, which is very dangerous for the patient.

Figure 4 .
Figure 4. Correlation between stenosis length [mm] and mean velocity [m/s] for maximum shear stress [Pa].

Figure 6 .
Figure 6.Correlation between stenosis hydraulic diameter [mm] and mean velocity [m/s] for maximum shear stress [Pa].

Figure 6 .
Figure 6.Correlation between stenosis hydraulic diameter [mm] and mean velocity [m/s] for maximum shear stress [Pa].

Figure 8 .
Figure 8.The geometry of internal carotid artery with stenosis.

Figure 8 .
Figure 8.The geometry of internal carotid artery with stenosis.

Figure 9 .
Figure 9. Wall shear stress [Pa] results for rest and exercise.Figure 9. Wall shear stress [Pa] results for rest and exercise.

Figure 9 .Table 5 .
Figure 9. Wall shear stress [Pa] results for rest and exercise.Figure 9. Wall shear stress [Pa] results for rest and exercise.

Figure 11 .
Figure 11.3D velocity magnitude [m/s] profile results for rest and exercise.

Figure 11 .
Figure 11.3D velocity magnitude [m/s] profile results for rest and exercise.

Table 1 .
Dimensions of analyzed geometry cases.

Table 2 .
List of β 1 correction coefficients related to scale factor.

Table 3 .
List of β 2 correction coefficients related to dimensionless stenosis length.

Table 4 .
List of β 3 correction coefficients related to dimensionless stenosis hydraulic diameter.