CFD Computation of Flow Fractional Reserve (FFR) in Coronary Artery Trees Using a Novel Physiologically Based Algorithm (PBA) Under 3D Steady and Pulsatile Flow Conditions

A novel physiologically based algorithm (PBA) for the computation of fractional flow reserve (FFR) in coronary artery trees (CATs) using computational fluid dynamics (CFD) is proposed and developed. The PBA was based on an extension of Murray’s law and additional inlet conditions prescribed iteratively and was implemented in OpenFOAM v1912 for testing and validation. 3D models of CATs were created using CT scans and computational meshes, and the results were compared to invasive coronary angiographic (ICA) data to validate the accuracy and effectiveness of the PBA. The discrepancy between the calculated and experimental FFR was within 2.33–5.26% in the steady-state and transient simulations, respectively, when convergence was reached. The PBA was a reliable and physiologically sound technique compared to a current lumped parameter model (LPM), which is based on empirical scaling correlations and requires nonlinear iterative computing for convergence. The accuracy of the PBA method was further confirmed using an FDA nozzle, which demonstrated good alignment with the CFD-validated values.


Introduction
The coronary artery is one of the essential blood arteries that ensure the heart receives steady blood flow. The heart muscle (myocardium) receives oxygenated blood via the coronary arteries; when these are clogged or obstructed, the myocardium starts to deteriorate, a condition known as ischemia. The most significant cause of death worldwide is coronary artery disease (CAD), also known as ischemic heart disease. According to a report by the World Health Organization, CAD causes nearly 50% of disease-related deaths in Kazakhstan [1]. An American experience a coronary episode roughly every 25 s; statistically, the pass-away rate is approximately one individual every minute [2].
Coronary angiography may offer information on anatomical stenosis, the gold standard for diagnosing coronary heart disease. Evaluations of fractional flow reserve (FFR) from CT in Europe and the United States suggest that it has the potential to decrease false diagnoses before patients are referred for invasive testing [3].
Coronary artery computerized tomography (CT) imaging is preferred over coronary angiography for prescreening asymptomatic individuals due to its noninvasive nature, ease of use, low cost, and excellent repeatability. Since FFR from CT is based on high-quality coronary artery CT angiography (CTA) image data, it does not need extra loading, scanning, or dosage. Although there is a grey area that necessitates additional functional measures for diagnosis, the diagnostic accuracy and specificity are greater than those of coronary The problem addressed in this work is the diagnosis of stenosis (obstruction or narrowing) in the coronary arteries of individual patients. Traditional methods for diagnosing stenosis, such as coronary angiography, have limitations and can be invasive and costly and can carry risks for patients. There is a need for more effective and noninvasive diagnostic methods for stenosis.
The PBA method proposed in this work is based on an extension of Murray's law and different inlet criteria [13] and is designed to be used alternately and iteratively to compute the outflow boundary conditions of coronary trees for individual patients. The PBA method is validated using FFR measurements in actual patient arteries and benchmarking with an FDA nozzle. It is intended to provide realistic and tailored outflow boundary conditions for the diagnosis of stenosis without requiring measurable data. The unique contribution of the PBA method to the field is its use of Murray's law to establish initial boundary conditions that are then modified through the addition of new inlet boundary conditions and iterations, resulting in numerical convergence. This approach may provide a more effective and noninvasive method for diagnosing stenosis in the coronary arteries of individual patients.
In the cardiovascular system, agreements between experiments and the Murray diameter model have been found only in small arteries [14] and arterioles [15,16], as present in our CAD case. The final conditions at the outlets are ultimately decided by a tree's geometry, the conservation laws built into CFD, the numerical iterations, and the additional inlet patient-specific parameters; we only utilize Murray's formula to estimate their initial conditions. The third law may not accurately predict the outcome since we use it as a preliminary approximation.
Unlike machine-learning-based data-driven simulation models, which depend on a significant quantity of data for verification, this PBA model is entirely based on physiology and physics. We are currently in the conceptual model-proof stage. It is undoubtedly desirable to undertake large-scale clinical testing, which we hope to accomplish in the next step for possible commercialization and future clinical uses. Because of this, we perform model validation using experimental benchmark data and outputs from other simulations.
The IREC approved the study related to patient data collection and ethical approval.

Governing Equations for Hemodynamic Flow
OpenFOAM v1912 is software used to solve fluid flow problems by numerically solving the governing equations for blood flow using the SIMPLE and PIMPLE algorithms. These algorithms can be used to perform steady-state and transient simulations, respectively, and are capable of handling 3D laminar, turbulent, steady, and unsteady flows and compressible and incompressible flows. When calculating blood flow in OpenFOAM, the software uses the incompressible 3D Navier-Stokes formula (Equation (1)) to model the flow of the fluid.
ρ-density of blood, kg/m 3 ; υ i -i component of velocity, m 3 /s; V i,t -velocity derivative with respect to time; p-pressure, Pa; τ ij -stress tensor (viscous portion). OpenFOAM uses the finite volume method (FVM) and SIMPLE/SIMPLEC/PISO pressure correction schemes to solve the Navier-Stokes equations for incompressible fluid flow. The domain of the equation is defined as Ω in three dimensions, and the boundary conditions are specified as Γ = Γ D ∪ Γ N , where Γ D represents Dirichlet boundary conditions and Γ N represents Neumann boundary conditions. The domain is discretized using n el linear, unstructured tetrahedral elements, − Ω e . The final form of the equation is given by Equation (3).
where w ∈ W k n and q ∈ P k h .

Murray's Law and PBA for Rapid Iterative Computation of Outlet Conditions
The present work proposes a method that aims to address the problems related to the Windkessel-type boundary conditions and to provide patient-specific and realistic outflow boundary conditions in the absence of observed data. This method is based on an extension of Murray's law and the use of different inlet conditions that are alternately and iteratively applied. Sumbekova et al. previously carried out experiments using commercial software [17]. It has been shown that vascular systems follow Murray's law, which states that mammalian vascular transport systems use minimum energy for blood maintenance and transport.
According to the Hagen-Poiseuille Law for laminar flow in a vessel, the power required to drive blood flow through a vessel is as follows.
Q is the volumetric flow rate, l is the length of a vessel, r is its radius, and µ is blood viscosity.
Furthermore, the power necessary for the maintenance of the blood in a vessel is proportional to the blood volume in the vessel.
The power required to maintain the metabolism in the blood is, thus, calculated as follows.
where λ is the metabolic rate of blood. The total power required to drive the blood and maintain it is calculated as follows.
The radius that meets the minimum power requirement is obtained through the differentiation of P with respect to r set to zero. Thus: Bioengineering 2023, 10, 309

of 26
The vessel radius and volumetric flow rate in individual vessels have this functional connection. To find the converged patient-specific flow outlet conditions in all the outlets of a coronary tree for a given inlet pressure and volumetric flow rate that are prescribed alternatively through iterative computation, a novel iterative scheme is, thus, proposed to couple Murray's law with CFD simulation. For a particular branch i of a coronary tree with N branches, as shown in Figure 1, its outflow is as follows.
alternatively through iterative computation, a novel iterative scheme is, thus, proposed to couple Murray's law with CFD simulation. For a particular branch i of a coronary tree with N branches, as shown in Figure 1, its outflow is as follows.
According to the conservation law, we have the following.
Equation (11) can then be obtained by dividing Equation (10) by Equation (11) and rearranging it: Our new physiologically based algorithm (PBA) aims to extract personalized, patient-specific outflow boundary conditions. It represents the interactions between a coronary tree and its microcirculation network downstream under both steady and unsteady conditions, as shown in a flowchart and illustrated in Figure 2. According to the conservation law, we have the following.
Equation (11) can then be obtained by dividing Equation (10) by Equation (11) and rearranging it: Our new physiologically based algorithm (PBA) aims to extract personalized, patientspecific outflow boundary conditions. It represents the interactions between a coronary tree and its microcirculation network downstream under both steady and unsteady conditions, as shown in a flowchart and illustrated in Figure 2.
In this work, the above-mentioned methods were implemented in OpenFOAM to perform a computational study of the hemodynamics in several patient-specific geometries with the aim of validating the methods using related ICA measurements, such as inlet flow rate, pressure, and FFR ICA . There were four vascular cases, named CT209, CHN13, CHN03, and FDA nozzle, implemented using the proposed PBA technique. Table 1 lists the workflow procedure for the PBA.
The geometries and outlets of CT209, CHN13, and CHN03 are presented in Figures 3-5, respectively. In this work, the above-mentioned methods were implemented in OpenFOAM to perform a computational study of the hemodynamics in several patient-specific geometries with the aim of validating the methods using related ICA measurements, such as inlet flow rate, pressure, and FFRICA. There were four vascular cases, named CT209, CHN13, CHN03, and FDA nozzle, implemented using the proposed PBA technique. Table  1 lists the workflow procedure for the PBA.    Tables 2-4 present the experimental results of pressure and flow rate at the inlet of a vessel (aorta) for the steady-state simulations of the three artery models: CT209, CHN13, and CHN03. Similarly, Tables 5-7 display the initial flow rates at the outlets of the three artery models, which were calculated using Equation (11) based on Murray's law with the inlet and outlet boundary conditions specified in Table 1. Figure 6 presents the inlet values for the transient case. The time-averaged value of the inlet waveform for the transient scenario was compared to the steady-state case to ensure the accuracy of the comparison.
In this study, the same set of fluid properties and blood models were used for both steady-state and transient simulations. The properties included a Newtonian dynamic viscosity of 0.0035 Ns/m 2 and a density of 1056 kg/m 3 [18]. The boundary conditions for these simulations included velocity/pressure inlets and velocity/pressure outlets, which were iteratively switched as part of the PBA scheme. The blood flow was modeled as laminar fluid flow. To ensure convergence, the criteria for the transient simulations required that the residual reduction must be lower than five orders of magnitude, ranging from 1 down to 10 −5 , while the criteria for the steady-state simulations required that the residual reduction must be lower than 10 −3 . In addition, the inlet data provided for each artery indicated that two cycles were used in the transient case.
The PBA study was replicated on an FDA benchmark model to validate the PBA approach. The nozzle geometry used in this study is shown in Figure 7 and has been previously described in detail by Stewart et al. [20,21]. The original FDA nozzle had a 12 mm diameter inlet and outlet tube, with a throat-to-inlet tube ratio of 1:3. Stiehm et al. scaled down the nozzle geometry to a 3 mm tube diameter. Normalization of the results was performed to ensure comparability with the outcomes of the FDA's round-robin trial. The FDA nozzle in this study only used one cycle, as the same boundary conditions were already tested by Stiehm et al. [18]. The input parameters for simulating an FDA nozzle using the PBA approach for steady-state and transient cases can be found in Table 8 and  Tables 5-7 display the initial flow rates at the outlets of the three artery models, which were calculated using Equation (11) based on Murray's law with the inlet and outlet boundary conditions specified in Table 1. Figure 6 presents the inlet values for the transient case. The time-averaged value of the inlet waveform for the transient scenario was compared to the steady-state case to ensure the accuracy of the comparison. In this study, the same set of fluid properties and blood models were used for both steady-state and transient simulations. The properties included a Newtonian dynamic viscosity of 0.0035 Ns/m 2 and a density of 1056 kg/m 3 [18]. The boundary conditions for these simulations included velocity/pressure inlets and velocity/pressure outlets, which were iteratively switched as part of the PBA scheme. The blood flow was modeled as laminar fluid flow. To ensure convergence, the criteria for the transient simulations required that the residual reduction must be lower than five orders of magnitude, ranging from 1 down to 10 −5 , while the criteria for the steady-state simulations required that the residual reduction must be lower than 10 −3 . In addition, the inlet data provided for each artery indicated that two cycles were used in the transient case.   Tables 5-7 display the initial flow rates at the outlets of the three artery models, which were calculated using Equation (11) based on Murray's law with the inlet and outlet boundary conditions specified in Table 1. Figure 6 presents the inlet values for the transient case. The time-averaged value of the inlet waveform for the transient scenario was compared to the steady-state case to ensure the accuracy of the comparison. In this study, the same set of fluid properties and blood models were used for both steady-state and transient simulations. The properties included a Newtonian dynamic viscosity of 0.0035 Ns/m 2 and a density of 1056 kg/m 3 [18]. The boundary conditions for these simulations included velocity/pressure inlets and velocity/pressure outlets, which were iteratively switched as part of the PBA scheme. The blood flow was modeled as laminar fluid flow. To ensure convergence, the criteria for the transient simulations required that the residual reduction must be lower than five orders of magnitude, ranging from 1 down to 10 −5 , while the criteria for the steady-state simulations required that the residual reduction must be lower than 10 −3 . In addition, the inlet data provided for each artery indicated that two cycles were used in the transient case.

Parameter Value
Experimental inlet pressure P exp 90.53 mm Hg (12,070.12 Pa) Experimental inlet flow rate Q exp 9.39944 cm 3 /s Table 3. Experimentally obtained input parameters for simulation of CHN13.

Parameter Value
Experimental inlet pressure P exp 90.61 mm Hg (12,870.12 Pa) Experimental inlet flow rate Q exp 7.17551 cm 3 /s Table 4. Experimentally obtained input parameters for simulation of CHN03.

Parameter Value
Experimental inlet pressure P exp 76.5 mm Hg (10,201.9 Pa) Experimental inlet flow rate Q exp 6.18 cm 3 /s Table 5. Initial calculated flow rates at each outlet by Murray's law for CT209 [19].       Table 6. Initial calculated flow rates at each outlet by Murray's law for CHN13 [19].  Table 7. Initial calculated flow rates at each outlet by Murray's law for CHN03 [19].  Figure 7 and has been previously described in detail by Stewart et al. [20,21]. The original FDA nozzle had a 12 mm diameter inlet and outlet tube, with a throat-to-inlet tube ratio of 1:3. Stiehm et al. scaled down the nozzle geometry to a 3 mm tube diameter. Normalization of the results was performed to ensure comparability with the outcomes of the FDA's round-robin trial. The FDA nozzle in this study only used one cycle, as the same boundary conditions were already tested by Stiehm et al. [18]. The input parameters for simulating an FDA nozzle using the PBA approach for steady-state and transient cases can be found in Table 8 and Figure 8, respectively.     The convergence criteria were evaluated at the end of each iteration with Equation (13). The procedure presented in Table 1 was repeated continuously until convergence was achieved. FFR was computed using Formula (13).

Results and Discussion
The mesh generation for the blood vessels was performed using the open-source cfMesh utility, which is part of OpenFOAM and generates volumetric meshes of unstructured Cartesian types. cfMesh is also able to automatically generate tetrahedral meshes The convergence criteria were evaluated at the end of each iteration with Equation (13). The procedure presented in Table 1 was repeated continuously until convergence was achieved. FFR was computed using Formula (13).

Results and Discussion
The mesh generation for the blood vessels was performed using the open-source cfMesh utility, which is part of OpenFOAM and generates volumetric meshes of unstructured Cartesian types. cfMesh is also able to automatically generate tetrahedral meshes for a variety of geometries, as shown in Figure 9. A mesh sensitivity analysis was conducted [19], and the same models were used in this study. For the CT209, CHN13, and CHN03 arteries, our research group [19] calculated mesh sizes and numbers of grid cells based on the results of the mesh sensitivity analysis.
The convergence criteria were evaluated at the end of each iteration with Equation (13). The procedure presented in Table 1 was repeated continuously until convergence was achieved. FFR was computed using Formula (13).

Results and Discussion
The mesh generation for the blood vessels was performed using the open-source cfMesh utility, which is part of OpenFOAM and generates volumetric meshes of unstructured Cartesian types. cfMesh is also able to automatically generate tetrahedral meshes for a variety of geometries, as shown in Figure 9. A mesh sensitivity analysis was conducted [19], and the same models were used in this study. For the CT209, CHN13, and CHN03 arteries, our research group [19] calculated mesh sizes and numbers of grid cells based on the results of the mesh sensitivity analysis.  The minimum mesh size was half the radius of the model's smallest artery (inlet or outlet). CT209's mesh size was 0.1702 mm, whereas CHN13's was 0.2262 mm. Mesh refinement resulted in a larger number of mesh elements without expansion layers along the vessel walls, which increased simulation computational time. Decreased mesh size from 0.20 mm to 0.1702 mm resulted in 1.7 times increase in computational time.
CHN03 had five meshes (69,681; 139,833; 286,897; 406,606; and 779,482) to examine mesh convergence. Increasing the grid number from 100,000 to 1.05 million decreased distal pressure by 2.2%. The model was discretized with 0.5 million volume cells after mesh dependency testing. At 0.4 million volume cells, the CHN03 mesh converged. In Tables 9 and 10, respectively, mesh details for the CT209 and CHN13 artery models are shown.  Figure 10 shows that, when the grid number was increased from 100,000 to about 1.05 million, the change value of the distal pressure was less than 5% (about 2.2%). After mesh dependency test, the model was discretized with a total of about 0.5 million volume cells. Further grid refinement led to <1% relative error. For CHN03, the mesh convergence results are shown in Figure 10. The number of volume cells was finally set to 0.4 million. There were five sets of meshes for CHN03 (with mesh element numbers of 69,68 139,833; 286,897; 406,606; and 779,482) to study mesh convergence. Figure 10 shows tha when the grid number was increased from 100,000 to about 1.05 million, the change valu of the distal pressure was less than 5% (about 2.2%). After mesh dependency test, th model was discretized with a total of about 0.5 million volume cells. Further grid refin ment led to < 1% relative error. For CHN03, the mesh convergence results are shown Figure 10. The number of volume cells was finally set to 0.4 million.  . FFR values in each round of iteration for CT209 [19]. Figure 11. FFR values in each round of iteration for CT209 [19]. Figure 11. FFR values in each round of iteration for CT209 [19].

Figure 12.
Relative differences between calculated and experimental FFR in each round of ite for CT209 [19].
The simulations of blood flow in the CHN13 model were performed using three sizes (coarse, fine, and very fine). Initially, the number of iterations was set to ten ever, since the relative errors for some mesh sizes showed systematic increases, the lations with 0.53 and 0.22 mesh sizes were stopped at the fifth and sixth rounds of tions, respectively. The experimental FFR for the CHN13 model was 0.68. The r showed that the finest mesh size of 0.22 mm produced the closest amount to this The FFR obtained from the simulations with the finest mesh size was equal to 0.6 shown in Figure 13. The final values of FFR obtained using the 0.35 and 0.53 mesh were 0.619 and 0.531, respectively. Additionally, the results in Figure 14 show th accuracy of the FFR obtained from simulations with the 0.35 mesh size decreased w number of iterations. Consequently, the maximum error was observed in the tenth of iterations. The simulations of blood flow in the CHN13 model were performed using three mesh sizes (coarse, fine, and very fine). Initially, the number of iterations was set to ten; however, since the relative errors for some mesh sizes showed systematic increases, the simulations with 0.53 and 0.22 mesh sizes were stopped at the fifth and sixth rounds of iterations, respectively. The experimental FFR for the CHN13 model was 0.68. The results showed that the finest mesh size of 0.22 mm produced the closest amount to this value. The FFR obtained from the simulations with the finest mesh size was equal to 0.691, as shown in Figure 13. The final values of FFR obtained using the 0.35 and 0.53 mesh sizes were 0.619 and 0.531, respectively. Additionally, the results in Figure 14 show that the accuracy of the FFR obtained from simulations with the 0.35 mesh size decreased with the number of iterations. Consequently, the maximum error was observed in the tenth round of iterations.
Bioengineering 2023, 10, x FOR PEER REVIEW Figure 13. FFR values in each round of iteration for CHN13 [19].   FFR was calculated using Equation (13), and for the CT209 model, it ha mental value of 0.76 [23]. Figure 15 shows the derived FFR values from the iteration, along with the relative percentage difference between the computed imental FFR (the degree of deviation from the experimentally obtained FFR lated noninvasive FFR value of 0.73 by Zhang et al. [23] is different from the e (ICA) FFR value of approximately 0.76 for the CT209 model. After conduc sensitivity analysis, Sumbekova et al. [19] from our research group obtained value of FFR = 0.757, which was very close to the ICA FFR. The same algor mented in OpenFOAM produced results that were similar to the ICA measur a final FFR value of FFR = 0.762.
According to Zhang et al. [23], the left anterior descending (LAD) pro was the site of arterial stenosis in the CT209 model [23]. In Figure 15, the LA region is the dark blue area of the blood artery, and it was observed that the ar was in the same location as that identified by Zhang et al. [23]. Figure 15a visualized FFR findings from the tenth iteration of the steady-state PBA. The pressure and flow rate residuals at all the outlets were defined as the averages of relative percentage differences in consecutive iterations. The convergence of the CFD simulation could be monitored by examining the reduction in the values of pressure and velocity/momentum equation residuals at the outlets as well. Figure 15 reveals the pressure and flow velocity residuals at the outlets of CT209 for all the iterations, where residual values are plotted in logarithmic scale.
FFR was calculated using Equation (13), and for the CT209 model, it had an experimental value of 0.76 [23]. Figure 15 shows the derived FFR values from the most recent iteration, along with the relative percentage difference between the computed and experimental FFR (the degree of deviation from the experimentally obtained FFR). The calculated noninvasive FFR value of 0.73 by Zhang et al. [23] is different from the experimental (ICA) FFR value of approximately 0.76 for the CT209 model. After conducting a mesh sensitivity analysis, Sumbekova et al. [19] from our research group obtained a final FFR value of FFR = 0.757, which was very close to the ICA FFR. The same algorithm implemented in OpenFOAM produced results that were similar to the ICA measurement, with a final FFR value of FFR = 0.762.
According to Zhang et al. [23], the left anterior descending (LAD) proximal artery was the site of arterial stenosis in the CT209 model [23]. In Figure 15, the LAD proximal region is the dark blue area of the blood artery, and it was observed that the artery stenosis was in the same location as that identified by Zhang et al. [23]. Figure 15a includes the visualized FFR findings from the tenth iteration of the steady-state PBA.
We compared the steady-state results obtained by the PBA method in OpenFOAM with those obtained using a traditional lumped parameter method and the current approach used by our research group [19], as well as with the ICA measurement results. It was found that the PBA method demonstrated good accuracy and efficiency, similar to the traditional methods. Figures 16 and 17 show the visualized FFR distributions for the CHN13 and CHN03 models, respectively. The calculated FFRs for CHN03 and CHN13 by the OpenFOAM PBA were 0.88 and 0.683, respectively, which are in excellent agreement with the corresponding experimental values of 0.86 and 0.68. Figures 18-20 show the FFR distributions in transient simulations for the CT209, CHN13, and CHN03 models, respectively. The PBA residuals also demonstrated convergence in the velocity and pressure equations. We compared the steady-state results obtained by the PBA method in Open with those obtained using a traditional lumped parameter method and the curr proach used by our research group [19], as well as with the ICA measurement re was found that the PBA method demonstrated good accuracy and efficiency, sim the traditional methods. Figures 16 and 17 show the visualized FFR distributions CHN13 and CHN03 models, respectively. The calculated FFRs for CHN03 and CH the OpenFOAM PBA were 0.88 and 0.683, respectively, which are in excellent agr with the corresponding experimental values of 0.86 and 0.68. Figures 18-20 show the FFR distributions in transient simulations for the CHN13, and CHN03 models, respectively. The PBA residuals also demonstrated gence in the velocity and pressure equations. Preliminary steady-state simulations were conducted to achieve convergence of the model before performing OpenFOAM transient PBA simulations using the mapFields package. The mapFields utility converts one or more fields that are specific to a geometry into their equivalents, and it is universal because there is no requirement for the geometries to be similar [17].
In Figures 18b, 19b and 20b, during iterations, the PBA method switched the boundary conditions, which could lead to high-frequency oscillations. For transient simulations, it was necessary to drive the solution to convergence at each timestep, which resulted in fluctuations that appeared to have high frequencies. The solution could exhibit wave reflections resembling high-frequency oscillations due to the complex geometry. Figure 16 compares the FFR values between steady-state and transient simulations at the locations indicated by the FFR arrows in Figures 16-20. The FFR values increased in all three scenarios and tended to decrease due to alignment with the flow rate waveform. The mean FFR values are included in Table 11 as the outcome values for the transient simulations. As shown in Figure 21c, the second half of the transient instance exhibited rapid fluctuations in FFR due to changes in the PBA boundary conditions. In Figure 21a,b, the models show good agreement with the flow rate inputs, indicating that the PBA method could converge over several flow rate cycles.
Three artery models were analyzed using both methods: the standard LPM technique was used in ANSYS [23], while the suggested PBA was implemented in Simvascular [19] and OpenFOAM for steady-state and transient cases, respectively, as described in Appendices A and B. Additionally, all abbreviations in this work are provided in Appendix C.  Table 9 provides a summary of the comparison between the computed and experimental ICA FFRs [19]. This table shows that the recommended PBA produced results that were independent of the solver. Simvascular and OpenFOAM were used for the suggested PBA, while ANSYS was used for the conventional LPM approach.
The findings indicate that the relative errors or discrepancies across steady-state simulations in ANSYS CFX, Simvascular, and OpenFOAM did not exceed 3.24%. The results also show that the relative error disparities in transient simulations in OpenFOAM did not exceed 5.26%. The smallest inaccuracy was observed between the Simvascular and OpenFOAM steady-state PBA for CT209 at 0.26%. The largest inaccuracy for CT209 was in the OpenFOAM transient PBA at 5.26%. This may be due to the highest flow rate values, which could increase the overall FFR readings in the tested artery. The error in the OpenFOAM transient PBA for CHN03 was the smallest recorded at 0.02%, while the largest inaccuracy for CHN03 was in the OpenFOAM steady-state PBA at 2.33%. The error in the OpenFOAM steady-state PBA for CHN13 was the lowest at 0.44%, while the highest error was in ANSYS CFX at 3.24%. Despite the error peaks associated with various methods, the PBA approach showed a promising performance. The classic LPM and LPNM techniques, which involve the calculation of capacitance, resistance, and inductance and the creation of a fictitious downstream capillary vessel network using fractal techniques, are fundamentally different from the PBA approach. The LPM and LPNM techniques require the calculation of resistance at the outlets and the use of additional inlet measurement conditions and numerical iterations, which are not based purely on physiology. In contrast, the PBA technique is patient-specific and physiologically based. It estimates the initial conditions at the outlets using Murray's law and then uses the geometry of the arterial tree, the conservation laws incorporated into CFD, the suggested numerical iteration scheme, and the additional measured inlet patient-specific conditions to determine the final conditions at the outlets. The suggested PBA strategy, which is completely based on physics and physiology and is patient-specific, was also shown to be computationally effective. It was integrated into the standard CFD pressure correction iterations as an iterative boundary-switching system that did not require two simulation rounds. The PBA approach was further validated using the standard LPM method published by Zhang et al. [23] for real-patient arteries. The LPM method employs a reference pressure, resistances at every outlet representing the flow resistance from the downstream microvasculature, and an overall resistance for the entire CAT that is related to outlet resistances through a population-averaged empirical scaling law [23]. The LPM is an iterative procedure that calculates the resistances and reference pressure to determine the outlet pressures in each branch outlet. It has been found that the LPM does not always ensure convergence to a unique solution for every situation. In contrast, the suggested PBA also achieved convergence to precise answers for all three analyzed scenarios. The PBA approach was found to be more reliable than the LPM because it is physiologically based and patient-specific, while the LPM does not consider these characteristics. Preliminary steady-state simulations were conducted to achieve convergenc model before performing OpenFOAM transient PBA simulations using the ma package. The mapFields utility converts one or more fields that are specific to a ge into their equivalents, and it is universal because there is no requirement for the tries to be similar [17].
In Figures 18b-20b, during iterations, the PBA method switched the bounda ditions, which could lead to high-frequency oscillations. For transient simulations necessary to drive the solution to convergence at each timestep, which resulted in ations that appeared to have high frequencies. The solution could exhibit wave ref resembling high-frequency oscillations due to the complex geometry. Figure 16    Three artery models were analyzed using both methods: the standard LPM technique was used in ANSYS [23], while the suggested PBA was implemented in Simvascular [19] and OpenFOAM for steady-state and transient cases, respectively, as described in Appendices A and B. Additionally, all abbreviations in this work are provided in Appendix C. Table 9 provides a summary of the comparison between the computed and experimental ICA FFRs [19]. This table shows that the recommended PBA produced results that were independent of the solver. Simvascular and OpenFOAM were used for the sug- Additionally, the PBA approach was verified using the benchmark model of an FDA nozzle. The goal was to replicate the study by Stiehm et al. [18] using the PBA approach. The computational results for axial velocity along the centerline of the nozzle geometry were compared to the CFD data from the FDA's round-robin study [18] for validation. The inlet conditions shown in Figures 7 and 8 were the same as those in the Stiehm et al. [18] study. The FFR distribution along the profile of the FDA nozzle can be seen in Figure 22, but we used the velocity results along the centerline for validation. tive procedure that calculates the resistances and reference pressure to determine the outlet pressures in each branch outlet. It has been found that the LPM does not always ensure convergence to a unique solution for every situation. In contrast, the suggested PBA also achieved convergence to precise answers for all three analyzed scenarios. The PBA approach was found to be more reliable than the LPM because it is physiologically based and patient-specific, while the LPM does not consider these characteristics.
Additionally, the PBA approach was verified using the benchmark model of an FDA nozzle. The goal was to replicate the study by Stiehm et al. [18] using the PBA approach The computational results for axial velocity along the centerline of the nozzle geometry were compared to the CFD data from the FDA's round-robin study [18] for validation The inlet conditions shown in Figures 7 and 8 were the same as those in the Stiehm et al [18] study. The FFR distribution along the profile of the FDA nozzle can be seen in Figure  22, but we used the velocity results along the centerline for validation.  The axial velocity results for the CFD steady-state and CFDPBA datab idealized FDA nozzle matched closely, as shown in Figure 25. Additionally, t eraged axial velocity obtained from the transient CFD simulation with wav conditions agreed well with the steady-state PBA velocity data. However, the nor differences downstream of the nozzle between 5D and 10D ( Figure 7). Th of wave pressure in the PBA algorithm caused the PBA sinus and waveform Since the pressure was constant in the steady-state example, the findings w fected, but the waveform pressure caused a slight fluctuation in the centerl values. The axial velocity results for the CFD steady-state and CFDPBA databases for the idealized FDA nozzle matched closely, as shown in Figure 25. Additionally, the timeaveraged axial velocity obtained from the transient CFD simulation with waveform inlet conditions agreed well with the steady-state PBA velocity data. However, there were minor differences downstream of the nozzle between 5D and 10D (Figure 7). The inclusion of wave pressure in the PBA algorithm caused the PBA sinus and waveform deviation. Since the pressure was constant in the steady-state example, the findings were not affected, but the waveform pressure caused a slight fluctuation in the centerline velocity values.

Conclusions
In this work, a physiologically based algorithm (PBA) was proposed and implemented in OpenFOAM CFD solver software to simulate blood flow in three coronary artery tree models and an FDA nozzle model. The performance of the proposed PBA was assessed using OpenFOAM and other CFD solvers [19]. The PBA, which was based on Murray's law and patient-specific conditions at the inlet, could accurately calculate the outlet boundary conditions iteratively, a task that is difficult for both traditional LPM methods and invasive measurements. The PBA was found to work well with all the CFD solvers, such as Simvascular and OpenFOAM, and it could accurately predict FFR values in the cases tested. The results of the PBA were also compared with those of the LPM method, which was performed in ANSYS CFX [19]. The FFR values obtained from simulations using the suggested PBA were in excellent agreement with experimental ICA measurements. Unlike traditional methods, the PBA always guaranteed quick convergence in steady-state cases to an accurate solution for every situation. Because it is signif-

Conclusions
In this work, a physiologically based algorithm (PBA) was proposed and implemented in OpenFOAM CFD solver software to simulate blood flow in three coronary artery tree models and an FDA nozzle model. The performance of the proposed PBA was assessed using OpenFOAM and other CFD solvers [19]. The PBA, which was based on Murray's law and patient-specific conditions at the inlet, could accurately calculate the outlet boundary conditions iteratively, a task that is difficult for both traditional LPM methods and invasive measurements. The PBA was found to work well with all the CFD solvers, such as Simvascular and OpenFOAM, and it could accurately predict FFR values in the cases tested. The results of the PBA were also compared with those of the LPM method, which was performed in ANSYS CFX [19]. The FFR values obtained from simulations using the suggested PBA were in excellent agreement with experimental ICA measurements. Unlike traditional methods, the PBA always guaranteed quick convergence in steady-state cases to an accurate solution for every situation. Because it is significantly more accessible and effective than the conventional Windkessel circuit analogy technique, this noninvasive FFR estimate using a PBA is a potential strategy for patient-specific detection of coronary stenosis. Our study further suggests that the suggested PBA, combined with an effective CFD solver, could be used as a low-cost, effective, and precise tool for diagnosing CAD. In individual patients, the PBA method has the potential to replace the expensive and risky invasive coronary angiography (ICA) procedure. Future research on pulsatile flows in CATs with fluid-structure interaction using the PBA is suggested.
We also used a nozzle model in this work that was based on the FDA's benchmark geometry. By comparing the estimated axial velocity profile to CFD and CFDPBA data provided by the FDA, normalized fluid mechanical quantities could be used to further validate the method and its results.
Finally, two different transient CFD simulations were performed utilizing a sinus inlet condition and a physiological waveform. We concluded that, if only time-averaged results were examined, steady-state simulations were appropriate for hemodynamic studies. This approach may reduce the need for computing resources in further hemodynamic studies. Institutional Review Board Statement: Nazarbayev University Institutional Research Ethics Committee (NU-IREC) approved the submission 616/29092022 to amend (change in research personnel) and extend the research project "Novel multidisciplinary and multiscale methods for patient-specific study of biomechanical characteristics of coronary artery disease (PCAD)". Please note that the approval for this extension will commence on 17 November 2022 and will lapse on 16 November 2023.

Informed Consent Statement:
The NU IREC approved the study related to patient data collection and ethical approval. Data Availability Statement: Upon request by emailing nursultan.alzhanov@nu.edu.kz.