Boundary Condition Options for Carotid Bifurcation Analysis, Using Doppler Velocity Measurements

The Common Carotid Artery plays a vital role in supplying the brain, and its bifurcation is susceptible to vascular diseases. It is often analyzed using computational fluid dynamics (CFD) simulations, but it is challenging to prescribe boundary conditions that approach patient-specific flow conditions. We examined six boundary condition (BC) groups to determine the most accurate flow conditions aligning with available measured data. We conducted CFD simulations on a stenotic carotid bifurcation, using patient-specific Doppler ultrasound sonography velocity measurements at the inlet and both outlets. Three BC methods used defined inlet flow rate and either constant pressure (Basic), Windkessel model, or constant flow ratio (Murray) at the outlets. Three other methods were defined with flow rates at two boundaries and constant pressure at the third one. Defining two boundary flow rates shows the closest results to physiologically valid data. However, the difficult Doppler measurements on the outlet branches can inaccurately amplify velocity amplitudes and may detect a false flow direction. Therefore, cross-sectional corrections were implemented to fit the outlet and inlet flow rates, while keeping the measured velocity histories. Our results show that the Murray and Basic methods, while easily available, exclude carotid-specific flow conditions by disregarding downstream flow resistances. We conclude that a Windkessel-method can produce the most accurate results without forcing outflow conditions. However, usually unavailable measurements are necessary for its application. Simulations with outlet-defined volume flow can also produce physiologically valid solutions but require the application of cross-sectional geometry correction.


Introduction
Cardiovascular diseases are the leading causes of death in western countries [1] and are therefore thoroughly researched.Atherosclerosis is one of the most common etiologies of such diseases.It causes the narrowing of the vessels, and it poses a risk of stroke through the occlusion (blockage) of arteries that supply the brain causing cerebral hypoperfusion or stenotic (narrowing) vessel section, which are prone for thromboembolus formation and causing cerebral embolization.Due to the increased susceptibility of stenosis at bifurcations and junctions [2], the carotid bifurcation stands out as a high-risk area.Consequently, numerous studies have been conducted to explore the consequences of vessel steno-occlusive disease [3][4][5].
Common carotid arteries (CCA) can be found on both sides of the neck, where it branches into the internal carotid artery (ICA) and the external carotid artery (ECA).The ICA has low resistance, since its function is to provide unidirectional flow to the brain during the entire heart cycle.In contrast, the ECA is a peripheral artery characterized by higher resistance without the need for constant flow in the entire lengths of the heart cycle.In extreme pathological cases of the CCA occlusion, a reverse (retrograde) flow in the ECA and a steal phenomenon of the ICA is documented [6], thus the cerebrovascular system main aim is to provide constant blood flow due to the ICA for the brain.
Today it is a standard method to investigate the carotid bifurcation with three-dimensional numerical flow simulations [7][8][9].In order to achieve this, a patient-specific geometry is necessary, which is usually acquired by segmenting computed tomography angiography (CTA) or magnetic resonance angiography (MRA) images.However, the segmentation method is inherently subjective.During the segmentation it is not obvious, which voxels correspond to the actual lumen geometry.This might result in the incorrect selection of additional volumes (over-segmentation).Such segmentation errors can lead to inaccuracies in the flow and simulation results [10,11], which can change the complex blood flow patterns associated with the carotid bifurcation [12].
Well-defined boundary conditions are also necessary to achieve physiologically correct, and patient-specific flow conditions.At the carotid bifurcation, it is not self-explanatory, what kind of boundary conditions should be used and what measurements are needed for their implementation.This can be seen in the wide variety of methods used in the literature.These applied boundary conditions were examined to provide a selection of methods for our investigation (Table 1).Based on the collected 38 articles [7][8][9], 11 different types of boundary condition setups were discerned from the 45 implemented cases.
Patient-specific measurement data are often not available, which restricts the range of potentially applicable boundary conditions.To gain the patient-specific data, non-invasive techniques are employed in most cases, such as the ultrasound (US) Doppler velocity measurement.However, this approach is also constrained by the difficulty of performing measurements at the ECA and ICA [48].
We investigate the applicability of six selected boundary condition groups with respect to their effect on the flow division between outlets and velocity-time traces.Furthermore, we attempt to offer simple geometry-improving methods, that can further utilize the available measurement data.
This paper is an extended and improved version of our conference paper [49].

Clinical data, geometry
A series of CTA images were provided with their corresponding Doppler measurements by the Dept. of Vascular and Endovascular Surgery of Semmelweis University.For this exploratory study, one patient was chosen based on CTA image quality and the availability of US Doppler measurements at the CCA, ICA and ECA.Before further processing, all necessary patient data were anonymized.The study was conducted according to the guidelines of the Declaration of Helsinki and approved by the Institutional Review Board (or Ethics Committee) of Semmelweis University (84/2019 10 May 2019).Further uncertainties arising from the Doppler measurement technique are included in the findings of Hartley et al. [50].The CTA contains the left CCA and ICA and ECA's branches after the bifurcation.On the ICA a 68.5% stenosis is detected, which was calculated with the North American Symptomatic Endarterectomy Trial criteria [51].
A segmentation procedure was carried out to create the necessary 3D geometry from the CTA images, using the open-source software itk-SNAP.The semi-automatic method of the software was utilized to capture the relevant arteries: CCA, ECA, and ICA.A second researcher checked the quality of the created geometry to minimize the subjectivity of the segmentation method.The Vascular Modelling Toolkit (VMTK) was used to smooth the surface mesh and then generate the necessary extensions at the inlet and outlets to avoid interference with the boundaries.VMTK's extension method uses a fully automatic procedure to merge a cylindrical tube to the selected artery.The length of the extension was chosen to be 5 times the diameter at the inlet and seven times the diameter at the outlets, based on the findings of Kim et al. [52].By using these extension lengths each boundary was at least 7.5 diameter length far from any geometry change of interest.Additional boundary-corrected geometries were created, where cross-sectional changes were prescribed.In these cases, half of the extension was defined as a transition length and the other half as a straight tube.After the smoothing, the remeshing script of VMTK was applied to create a triangular surface mesh containing approximately 45 000 faces, which was later used as a basis for creating the simulation mesh.The Doppler US measurements of the CCA, ICA, and ECA branches were available as .jpegimages.The images contain the slice of the arterial geometry where the measurements were carried out and the velocity history during the measurement time (Fig. 1).US measurements were performed according to guidelines [53].From the velocity history image, the top enveloping curve was generated using the image recognition tools of MATLAB.The gained curve is the maximum velocity at the cross-section of the measurement.A section of time was chosen when four heart cycles were fully visible and an average heart cycle time was calculated (0.85 sec).Then each of the four velocity waveforms belonging to a heart cycle was adjusted to be equally long as the average periodic time.Finally, the waveforms were ensemble-averaged to gain a patient-specific velocity waveform for a single heart cycle (Fig. 2).The image of the arterial section (CCA) was used to measure the diameter.The obtained diameter was then used as a basis for the cross-sectional correction of the inlet.

Correction of the vessel diameter
The geometry was segmented from CTA images of the carotid bifurcation using itk-SNAP.The medical images provide the intensity values in Hounsfield units (HU) that can be used to decide whether a given voxel belongs to the lumen of an artery or to the surrounding area.The vessel volume is defined by choosing two voxel intensity thresholds and a seed volume that certainly belongs to the vessel.The software then expands the selected volume by adding the neighbouring voxels to it with an intensity that is between the two thresholds.Itk-snap also uses a mathematical algorithm to select the voxels near the threshold limits to create a more physiologically correct geometry.The thresholds are chosen manually, and since the voxel intensity has a gradient at the wall of the arteries, it is not possible to choose an exact threshold value objectively.The low image resolution may also introduce uncertainty in the morphological parameters of the artery.These uncertainties, that result from the subjectivity of user choices in the model, generally appear as over-segmentation [10,11], which results larger arterial diameters.
To correct the segmentation inaccuracies, the geometry has to be modified.Two types of modifications were implemented: one at the inlet and one at the outlets (Fig. 3. a.).Only the sections close to the boundaries were changed, since realistic geometry dimensions can be derived from the Doppler measurements in these places.Further from the measurement locations such modifications would not be valid since the ratio of the over-segmentation in differently sized and shaped lumen sections is different.This difference is even more significant near the stenosis, where calcification can introduce additional issues in discerning the lumen geometry.
The Doppler images were used as a basis for the correction at the CCA cross-section.Doppler measurements capture only a limited portion of the arterial geometry, but in higher resolution (8.5 pixel/mm) compared to the CTA (2.35 voxel/mm).Since the CCA has low tortuosity, it is relatively straight, its diameter can be determined with higher accuracy from the Doppler image than from the CTA.The cross-section of the chosen segmented vessel was corrected from 42.7 mm 2 to 28.0 mm 2 , using the circular shape of the geometry extensions.This means an approximately 0.7 mm reduction of the radius, which corresponds to an over-segmentation of 1.65 voxels.Based on the findings of Csippa et al. [11], this error is within the range expected from a non-medical segmentation operator.The correction was implemented into the CTA geometry by introducing a confusor (cross-section narrowing) into the CCA extension.
The second modification to the geometry aimed at achieving a similar correction of the ECA and ICA sections.Since these branches are not straight and there are significant changes in cross-section, using Doppler images proves unreliable for accurately determining the diameters.Nevertheless, the calculation of flow rate, based on measured velocity necessitates a cross-sectional correction, as continuity is not fulfilled.This error in the flow continuity was used to determine the correct cross-sectional areas of the outlets.An assumption was made that the magnitudes and the waveforms of the Dopplermeasured velocities were correct, therefore, the cross-sectional areas were adjusted to reduce the difference between the inlet and outlet volume flows.Since the Womersley number in the CCA was calculated to be 11.2, the velocity profile deviates from the parabolic shape, but as a crude approximation, the average velocity was calculated from the measured maximum velocity as where v max is the maximum velocity of the parabolic profile and v avg is the cross-sectional average of the velocity.
The time-and space-averaged velocities can be used in the time-averaged continuity equation: where Q ̅ is the cycle-averaged volumetric flow rate, v ̅ avg is the cycle and cross-sectional average of the velocities in the respective cross-sections, A * is the corrected cross-sectional area, which is to be defined.For the calculation of the unknown cross-sectional areas, the A ECA /A ICA ratio was denoted by K and assumed to be constant.This is a simplified way to handle the ratio between an ideal area and the constant segmentation error caused by image resolution.
The area ratio K of the outlets in the examined geometry is where A ECA and A ICA are the cross-sectional areas of the ECA and ICA at the location of the Doppler measurements, respectively before the area correction.The corrected cross-sectional areas can be calculated by substituting this constant into the continuity equation Eq. ( 2) as follows: The other cross-section A ECA * is calculated from Eq. 3.This correction does not change the velocity waveforms; it only fulfils the cycle-averaged continuity.Doppler measurements were not simultaneous and the change in the heart cycle between the two measurements make it impossible to correlate them with each other.In addition, because of the elasticity of the vessel wall, even in an ideal case continuity applies only in a cycle-averaged sense, not in every time instance.A difference between inflow and outflow still exists, but it is reduced using this method.
Three geometries were generated to analyze the effects of the boundary corrections.In the first geometry (RAW), no change was made after the segmentation and smoothing.In the second geometry (DCor), only the CCA was corrected since CCA Doppler measurements are most commonly available.In the third case, besides the CCA correction of the DCor geometry, ECA and ICA were also corrected (QCor).

Boundary conditions
Six boundary condition groups were tested in this study (Table 2).Five methods were selected based on their frequent implementation in the literature (constant flow division ratio, Windkessel, CCA-ECA, CCA-ICA, constant pressure), and one less commonly used method, based on the available measurement data (ECA-ICA).Some methods, for example, the structured-tree model and the 1D model approach, require a complex model setup.These were excluded from this study to focus on the models that are simple to implement.
Two commonly employed models were selected, where only the inlet BC requires measurement data, while the outlet BCs were defined using different rulesets.In both cases the inlet condition at the CCA was a time-varying volumetric flow rate.In one of the two cases, the outlet BCs were defined as constant 0 Pascal pressure at the outlets (Basic) with the option to allow backflow.The second chosen method describes a constant flow division ratio between the outlets, using either an approximated ratio from previous measurements or Murray's law.For our investigation, Murray's law was applied, which relies on the ratio between the outlet diameters and the flow rate of the outlets [54]: To apply this BC, the volume flow ratio had to be recalculated for the geometries where the outlet dimensions were modified.This is the most used method since it is simple to implement without measurements and has the additional benefit of creating more patient-specific flow conditions.
Three other boundary condition groups were defined utilizing the available outlet velocity measurements.Two boundaries were defined using the known CCA, ECA or ICA velocities, while constant 0 Pa pressure was used at the third boundary.This resulted in three potential configurations, each featuring at least one outlet serving as a boundary where the flow condition is enforced, based on patient-specific data (C+E, C+I, E+I).Multiple implementations can be found in the literature where outlet measurements are used to apply the known velocity waveforms (C+E, C+I), but the method to apply the measured data at both outlets is rarely employed (E+I).This is because measurements at both outlets are generally unavailable for the same patient.
Finally, a three-element Windkessel method (WK3) was implemented at the outlets, with volume flow definition at the inlet.This zero-dimensional model of the downstream system contains a proximal resistance, a distal resistance and a capacity.To calculate the Windkessel parameter for a patient-specific problem, velocity and pressure measurements are required.Since pressure measurements were unavailable, approximate values were adopted from the literature to describe the Windkesselmodel [14,27,28].This method is used more frequently in newer carotid bifurcation investigations, even if only approximate parameters are available, making it the second most common boundary condition method.
All the boundary condition groups were analyzed for the Raw, DCor, and QCor geometries.

Simulations
To obtain the volumetric flow rate values at the inlet and outlets, computational fluid dynamics (CFD) simulations were carried out.For the numerical calculations the transient laminar solver of ANSYS CFX 19.5 was utilized.The initial condition of the transient simulations was set as stationary, zero velocity fluid and zero Pa pressure throughout the volume.For each simulation the time duration was set to include three heart-cycles, but only the last one was analyzed to mitigate the influence of the initial condition.To enable transient flow structure generation, the 0.85 sec heart cycles were divided into 5000 times steps, based on the findings of Khan et al. [55].A high-resolution advection scheme and a second-order backward Euler transient scheme were used for the spatial and temporal discretization, respectively.
From the triangular surface meshes, unstructured tetrahedral meshes with 1 million elements were generated for each modified geometry.The mesh resolution was selected based on the independence analysis of Savabi et al. [56].The walls of the geometries were assumed to be rigid for the numerical simulations.The findings of Albadawi et al. [57] show that the volumetric flow rate division between the outlets is unaffected by the rigidity of the wall, while Jodko et al. [58] point out that in stenotic cases the elasticity of the wall reduces and rigid-wall approximations can have less effect on the results.The effects of the non-Newtonian properties of blood on the flow conditions are considered negligible, therefore, a Newtonian fluid with a density of 1055 kg/m 3 and dynamic viscosity of 0.0034 Pas was implemented.
In total, 18 setups were simulated for all the possible geometries and boundary condition group combinations in the present case study.

Analyzed values
Our analysis primarily focused on the volume flow division between the inlet and outlets, since it is a usual indicator of flow behavior in the carotid bifurcation [14,18,20,59].The examination was conducted by comparing the volume flow waveforms produced by the different boundary conditions and calculating the flow division between the outlets.For the flow division analysis, the volume flow ratio was calculated as: This ratio can be used as a clear indication if the flow division between the outlets is physiologically valid or not.Since the blood demand of the brain through the ICA is higher than the blood demand to the peripherals through the ECA, the flow rate ratio is between 0.5 and 1.0 in most bifurcation cases.It also indicates the backflow condition through the ECA branch if its value is higher than 1.0 (Fig. 4).
Additionally, velocity waveforms were also investigated near the stenosis.Points were chosen in the ICA branch, and the velocity values were plotted in time.The resulting velocity values were normalized with the timeand space-averaged inflow velocity: where v is the instantaneous velocity at a selected point, and A CCA is the inlet cross-sectional area of the unaffected (RAW) geometry.Velocity monitoring points were selected, where velocity oscillation was observed in the case of the Basic method in the RAW geometry (Fig. 3(b)).

Volume flow
Simulation results of the volume inflow and outflow rate were examined at the third heart cycle.The results were categorized into four groups, determined by the similarities in the boundary condition setups (Table 3).In the following subchapters, only some illustrative figures are shown, not every possible case will be presented.

Basic and Murray-law boundary condition groups
Basic and Murray-law simulations resulted in similar volumetric flow rate waveforms (Fig. 5), following the defined Q CCA characteristics, although the volume flow ratio between ICA and ECA differed.In the Murray-law case, this waveform is prescribed, while the Basic case produced these results without a forced flow at the boundaries.
The CCA geometry correction (DCor) resulted in the downscaling of the volume flows since only the Q CCA was decreased via the A CCA reduction (Table 3).
After the second geometry correction (QCor), the volume flow division of the Murray-law was recalculated,  producing a smaller difference between the outlet flow rates, consequently increasing the Q ICA ratio from 0.33 to 0.43.The QCor simulation of the Basic case showed a similar ratio increase from 0.17 to 0.32, as shown in Fig. 4.

C+E and C+I boundary condition groups
In the C+E and C+I BC groups one inlet and one outlet flow rate are defined.This creates simulation setups where the undefined second outlet volume flow can be calculated using the continuity equation.In these cases, the resulting volume flows change when applying the geometry corrections, but only because the boundary cross-section changes, while the measured maximum of the implemented velocity profile remains the same.
Ideally, both of these BC groups should yield identical results, since all the flow rates (calculated with the corresponding cross-section) are known from the Doppler velocity measurements.However, a clear difference is observable in Fig. 6.This difference is the result of a non-simultaneous measurement of the Doppler velocities at the different sites.In both cases, the volume flow ratio is higher than 0.5, and it stays so after both geometry corrections.This means that Q ICA > Q ECA , which could result in a significant velocity increase at the ICA stenosis.In contrast to the previous simulations, the flow rate ratio changes after the first geometry correction, especially in the C+I case from 1.58 to 2.41.The second correction (QCor) has the opposite effect on the flow division, resulting in a similar ratio at the C+E case (0.82 and 0.87) and a lower ratio at the C+I case (1.58 and 1.13) compared to the Raw geometry values (Fig. 4).
In the case of the C+I boundary condition group, a flow rate ratio larger than 1.0 was observed, which is equivalent to the backflow at the ECA.The Doppler measurements do not show such flow condition, but there can exist other stenosed cases, where backflow on ECA is a valid solution.

E+I boundary condition group
The E+I boundary condition group has a similar setup to the C+E and C+I cases, but the inlet is the 0 Pa boundary (opening BC) instead of one of the outlets.Defining both outlets results in the same flow rate ratio when the DCor or QCor geometry corrections are applied (Fig. 4).Thereby, the inlet velocity increases with the first correction because of the decreased inlet cross-section.The QCor correction decreases the outlet flow rates and therefore, the inlet flow rate to a similar average value that was defined in the other boundary conditions, based on the CCA Doppler measurements (Fig. 7).The average Q CCA is only 1.6% lower than that calculated with the measured velocity and the DCor corrected A CCA .Meanwhile, the flow rate ratio is only 0.96% higher than in C+E and 28.7% lower than in C+I in the QCor cases.

Windkessel method
The applied Windkessel parameters can be seen in Table 4.These parameters were chosen based on the range of values found in literature [14,27,28].Simulations with the Windkessel model provide a similar flow rate ratio as the measured values implemented in the E+I option (Fig. 4), with QICA being the higher volume flow at the outlets and QECA showing a slight amount of temporary backflow.In the case of the QCor results, the Windkessel model approximated the ECA volume flow with similar waveforms (Fig. 8).There is less than 0.2 ml/s volume flow difference between the systole values, and a similar difference between the diastole values.More noticeable are the differences at the ICA outlet, where the volume flow is 2.5 ml/s higher at the systole and 0.9 ml/s lower at the diastole.There is only a 4.6% difference in the average volume flows.This is because the aim of the QCor correction was to reduce the time-averaged difference between the measurements.
The volume flow division created by this method is not sensitive to the change in inflow rates, the flow rate ratio does not change after the first geometry corrections (DCor).Furthermore, the proportional correction of the outlets (QCor) also does not affect the flow division.

Normalized velocity
Velocities were analyzed in different probe points near the ICA stenosis.The gained velocity waveforms were normalized to negate the effect of the changing volume flow magnitude between the different simulation setups.
The first point for observation is at the bifurcation (Bifu) point (Fig. 9(a)).In this location, the velocity magnitudes exhibit the same behavior, as the ECA volumetric flow rates (Fig. 9(c)).This also means that the fluid mostly flows toward the ECA branch.The C+I result shows an exception: the diastolic velocity magnitudes are higher and sudden changes are observable when the volume flow to the ICA increases.This phenomenon can be attributed to the fact that the C+I method creates backflow from the ECA branch, showing a clear difference from the other simulations by creating oscillating velocities at the flow direction change.In this point, the effects of the geometry corrections are neutralized by the normalization, showing highly similar waveforms for the rest of the BC-s.These observations are also valid at the outer side of the ICA before the stenosis (BefSten).However, the flow is more even and the sudden incensement of the velocity magnitude only appear when the direction of the flow changes from the ECA to the ICA (Fig. 9(b)).Fig. 10 demonstrates the velocity waveforms in the points before (StenIn) and in (Sten) the stenosis.The velocity waveforms are similar in all cases.The only difference is in the magnitudes: at the stenosis, the maximum of the normalized velocity is approximately twice as high as before the stenosis, because of the smaller cross-section.Oscillatory behavior is not present at this point, and the velocity waveforms proportionally follow the ICA volume flow rate waveform.
Shown in Fig. 11., three points were selected after the stenosis where the emerging jet flow is expected (JetFlow1, 2, 3).As the distance increases from the stenosis, an increase in velocity oscillations can be observed due to the oscillation of the jet.The first geometry correction has only slight effects on the velocities; however, QCor smoothes out the velocity oscillations, indicating a stabilization of the jet.This latter phenomenon can be attributed to the fact that narrowing the outlet creates a favorable pressure gradient in the artery section, which is known to reduce instabilities.This is an artefact, produced by the geometry correction and therefore the jet oscillation is probably present in the real artery.
At the cross-sections of the post-stenotic expansion (AftSten1, AftSten2) lower amplitudes of velocity changes occur, since close to the stenosis, the formed jet has not yet developed high-amplitude oscillation.As can be observed in Fig. 12(a) the waveforms change significantly from what is expected, based on the volume flow of the ICA outlet in Fig. 12(c).This point is in a flow separation zone where the direction of flow changes frequently, creating backflow and secondary flows.Oscillations appear at higher velocities, making the jet flow more stable in the Basic and Murray methods than in the other boundary condition cases.
Similarly to the points in the jet flow, the QCor correction has a stabilizing effect in these points as well (Fig. 12(b)).

Discussion
Correct definition of the boundary conditions is crucial as they have significantly impact on the flow dynamics within the studied domain.Numerous BC options have been employed in the literature at the carotid bifurcation [31].Although these BCs are primarily chosen based on the available measured data, certain variations are still possible.Such variation can be the correction of the inlet or outlet geometries utilizing Doppler imaging to decrease the effect of over-segmentation.This study investigated the capabilities of 6 different boundary condition groups and assesses the impact of applying two geometry corrections.A key objective of this research is to filter out those BC options that may be physically or physiologically not realistic.
In the case of the Murray and Basic methods, the volume flow results are determined by the geometry of the modelled arterial section.These cases presented similar waveforms and the geometry corrections had similar effects on the flow ratio.The Basic method generalizes the outlet pressures as constants, while the Murray-law method applies the relation of the cross-sectional diameter and the volumetric flow rate [51].Although the simplicity of these outlet definitions improves their usability, both methods fail to account for the significant impact of downstream artery resistances at the carotid bifurcation.In our case, these two boundary condition groups underestimated the volume flow at the ICA section compared to the volume flow calculated from the Doppler velocity measurements.This underestimation of the volume flow can result in highly different flow conditions considering the velocities, even after applying the geometry corrections.The resulting lower velocities create less oscillatory flows before and after the stenosis.This means that these methods may not be suitable for finding high Oscillatory Shear Index areas on the arterial wall, which is often linked with the formation of arterial diseases [60].As a reference, it also has to be stated that because of the approximation of the rigid-wall geometry, all velocity results might be greater than in a more realistic elastic wall simulation [58].
The second set of boundary conditions was introduced with C+E, C+I, and E+I, where the volume flow was defined at two boundaries.The volume flow ratios of these simulation align with the anticipated results, because at least one of the outlet volume flows is defined on the basis of Doppler US velocity measurements.In all three cases, the volume flow ratio exceeds 0.5, indicating a successful fulfilment of the higher blood demand to the brain through the ICA [59,61].Although these three solutions exhibit comparable flow divisions, significant differences were observed in the amplitudes and waveforms of the volume flow rate cycle.While the C+E and E+I methods yield flow rate ratios in the expected range, the C+I method leads to a value higher than 1.0, which means overestimating the ICA volume flow.Higher velocity oscillations can occur in the overestimated flow, mainly when there is backflow in the ECA.When the geometry corrections of the inlets and outlets were implemented, the overestimation was decreased.Therefore, unrealistic flow conditions may occur with the use of these methods if segmentation errors are significant and are not corrected.
Some differences between the outlet-defined BC methods could be the result of the Doppler measurements not being simultaneous.A change in the heart cycle properties between the two measurements can lead to outlet volume flows that do not inherently fulfil the continuity.Additional errors may arise from the challenging nature of the ICA and ECA Doppler measurements [48].If the measurement is conducted close to a stenotic region, then the velocities are measured higher than what would be in a healthy arterial section.The velocity measurements directly impact the derived outlet flow rate, thereby such errors render the defined flow rate invalid.Therefore, these solutions are highly dependent on the accuracy of the measurements, which again depends on the operator.However, if well-executed measurements are available, a patient-specific solution can be achieved by directly forcing the downstream physics at the outlets.While such forced boundary conditions are usually avoided in general fluid simulation practice, these methods can allow flow conditions that do not develop with the more conventional Basic or Murray-law methods.For example, backflow is possible from the outlets, which can occur in highly stenotic cases [6].
The Windkessel boundary condition applies the physical properties of the downstream arteries without imposing forced outlet conditions.The WK3 model also makes it possible to induce backflow when its parameters are appropriately selected.Therefore, this method could provide the most accurate results, but it is rarely accessible for the lack of pressure and velocity data necessary for its parameter calculations [28].However, estimated parameters can be collected from the literature, producing volume flow waveforms similar to the ones calculated from velocity measurements.Wellestimated volume flow also creates expected velocity behavior, including the formation and oscillation of jet flow after the stent.Furthermore, the Windkessel method is applicable without the inlet and outlet geometry corrections since the results show no sensitivity to the scaling of the volume flow rate.This shows that the resistances implemented in the zero-dimensional model have more influence on the volume flow division than the resistance changes originating from the outlet geometry corrections.
Inlet and outlet geometry corrections were implemented to achieve a more accurate simulation geometry.Employing novel techniques, Doppler US measurements were used to mitigate errors arising from the segmentation method.The results of the corrections are dependent on the utilized boundary conditions, but in most cases, significant rescaling or flow division change is observable.Since DCor applies the measured diameter from the CCA Doppler image, and QCor applies the outlet diameters based on the measured outlet velocities, the resulting geometry is expected to be closer to the actual geometry of the patient.Simulation results conducted with more accurate geometries should also produce better approximations of patient-specific flow conditions.The jet flow after the stenosis is also a good indicator to show the effects of the geometry corrections.In the jet flow, the velocity values show higher oscillations in the RAW cases, while this oscillation is reduced after both geometry corrections.The damping of the jet oscillation by the QCor correction is artificial, and therefore the presence of jet oscillation is better represented by the RAW geometry even if the absolute velocity values are not correct.However, this is only the limitation of the not sufficiently long extensions.By increasing the extension length of the outlets this dampening effect could be removed.The problems with these corrections are similar to those encountered with complex boundary conditions.Namely, they require additional measurements to calculate the new geometry parameters, and the applicability of these measurements is dependent on their accuracy.When CCA Doppler US images are available, a more direct approach could be the use of DCor method as a quality check for the segmentation.If a significant correction is necessary based on the findings of Csippa et al. [11], it is worth considering resegmentation, or additional measurements and application of the QCor method.

Summary
A single stenotic carotid bifurcation was investigated with the aim of analyzing the differences between boundary condition types.To accomplish this, six boundary condition groups were chosen.Additionally, two geometry correction methods were established, based on data derived from Doppler measurements.
Our results show that boundary conditions without patient-specific outlets are not sufficient at the carotid bifurcation, because they do not consider the resistances of downstream arteries.To resolve this problem, alternative boundary condition configurations were implemented with defined velocity waveforms at the outlets.Although forced flow physics at the outlets is unconventional in general CFD practice, these simulations demonstrated closer alignment with patient-specific conditions.The important effects of these setups are the possibility of backflow from an outlet, and the formation and oscillatory behavior of the formed jet after the stenosis.A three-element Windkessel model was also examined, which produced similar patient-specific conditions without forcing flow physics on the outlets.The drawback of the complex boundary conditions is the additional measurements required and the necessary accuracy of these measurements.
Aside from the boundary conditions, two methods were applied to correct the segmented geometry, based on the Doppler US images and measured velocities.These corrections can significantly affect the amplitudes and waveforms of the resulting flow rates and velocities.The first geometry correction can offer a simple way to test the geometry using the CCA cross-sectional diameter since the magnitude of the correction is dependent on the inaccuracies of the segmentation.The second correction further improves the volume flow results in cases where Doppler US velocity measurements are available at the outlets.However, the QCor correction introduces an artificial damping into the jet oscillations.If the focus is on the jet behavior a longer outlet expansion should be implemented or the RAW geometry should be used with a possible rescaling.
In conclusion, the Basic and Murray methods should not be used in carotid artery simulations, since, despite their simple application, the resulting flow conditions are not physiological.The unrealistic flow rate divisions create.Therefore, velocity and WSS values are also not physiological with such methods, which might result in biased conclusions concerning the effects atherosclerosis or other analyzed diseases.The Windkessel method creates physiologically valid flow conditions even with estimated parameters and without the need for geometry corrections.However, it is difficult to implement if patient-specific results are needed, because of the necessary velocity and pressure measurements.The C+E, C+I, and E+I can be considered middle ground.Additional velocity measurements are needed at the outlets for their implementation, and they require accurately segmented geometry or the correction of the geometry.However, these methods produce patient-specific results by the nature of the forced outlets.

Fig. 1 Fig. 2
Fig. 1 Measured diameter on the Doppler image (top), and Doppler velocity measurement on the CCA, with the blue enveloping curve (bottom)

Fig. 3
Fig.3 (a) Alteration of the bifurcation geometry produced by the DCor correction (bottom) and the QCor correction (top), and (b) positions of the velocity monitoring points

Fig. 4
Fig.4 The effect of geometry corrections on volume flow ratio at the different boundary condition groups, dashed lines mark the volume flow ratios calculated from the Doppler velocity measurements

Fig. 5 Fig. 6
Fig. 5 Waveforms of the outlet flow rates between the Basic and Murray boundary condition groups before the geometry corrections

Fig. 7
Fig. 7 QCCA waveforms between E+I and Basic boundary condition groups

Fig. 8 Fig. 9
Fig. 8 Outlet volume flow waveforms between Wk3 and E+I boundary condition groups after the QCor geometry correction

Fig. 10
Fig. 10 Velocity waveforms at the stenosis (Sten, StenIn), and the volume flow waveform at the ICA outlet in the case of the C+E BC group, in the QCor geometry

Table 1
Implemented boundary condition methods in the literature

Table 2
Inlet and outlet boundary condition definitions

Table 3
Main effects of geometry corrections on the volume flow

Table 4
Estimated values of the 3-element Windkessel parameters This study was supported by the National Research, Development and Innovation Office.under Grant number NKFIH-K 129277 and NKFI-PD 146299.Project no.TKP-9-8/PALY-2021 has been implemented with the support provided by the Ministry of Culture and Innovation of Hungary from the National Research, Development and Innovation Fund, financed under the TKP2021-EGA funding scheme.