A representative anatomical model of the aorta and cerebral vasculature was used in this study (Figure 1). The model was designed and fabricated by United Biologics (Irvine, CA, United States) based on patient medical image data from several sources (e.g., the NIH Visible Human Project, patient-specific CT data). The in vitro model (Figure 1) is made of silicone with an elastic modulus of 3.1–3.4 N/mm², which is representative of human arteries. The entire model includes the aorta, common carotid arteries, internal and external carotid arteries, axillary arteries, middle cerebral arteries, and anterior cerebral arteries. In addition, a corresponding computational model of the in vitro model was reconstructed from high-resolution 3D micro computed tomography (μ-CT) scans (Figure 2A).

A closed-loop circulatory flow system was designed that included a centrifugal pump (Cole-Parmer, IL, United States), reservoir, flow meter and ultrasonic flow probes, and pressure transducers as shown in Figures 1A,B. The flow loop was connected to the cerebrovascular model, which was placed in the supine orientation. The inlet and outlet flow rates and pressures were monitored using ultrasonic flow probes (Transonic Systems, Inc., Millis, MA) and pressure transducers (Merit Medical, South Jordan, UT), respectively. The pressure transducers were connected to an analog data acquisition module (DAQ, National Instruments, Texas, United States) and recorded using LabVIEW software (National Instruments, Texas, United States).

Similar to Riley et al. (26), the working fluid consisted of a mixture of water (60% by weight) and glycerol (40% by weight) to obtain a density and dynamic viscosity that is representative of blood (1.09 ± 0.03 g/ml and 3.98 ± 0.14 cP, respectively) at an operating temperature of 22.2°C. Experiments were performed using a steady inlet flow rate of 5.17 ± 0.078 L/min, corresponding to a Reynolds number (Re) in the inlet tube of approximately 3890. This inlet flow rate was chosen to correspond to a representative mean physiological cardiac output of an adult. An extended tube of 900 mm in length was attached to the model inlet such that the flow entering the model inlet was fully developed. To study the effect of a stroke condition on the mean arterial pressure and flow rate in the cerebral arteries, nylon spherical clots of three different sizes (3.15 mm, 4.75 mm, and 6.38 mm in diameter) were manually inserted into the right MCA (RMCA) to completely block the vessel and the corresponding flow through outlet 4, as depicted in Figure 1C.

Three separate experiments were performed for both conditions (normal and stroke) to measure the flow rate and pressure at the inlet and various outlets in order to provide boundary conditions for the CFD simulations and for validation. Because the fluid heats up as it is continuously pumped through the flow loop, we waited for 30 min to allow for the fluid temperature to reach a steady state (22.2°C) before measuring the flow rate and pressure. Importantly, the viscosity of the fluid was tuned to account for the effect of this temperature rise, so as to obtain the desired value of 3.98 cP at the steady-state operating temperature. The regional distribution of the flow rate was tuned to match available literature data (30, 31) by adjusting clamps downstream of the pressure and flow rate measurement sites. This yielded a regional flow distribution such that 73.2% of the flow passed through the descending aorta and 26.8% of the flow was distributed to the remaining arteries stemming from the aortic arch. The flow rates to individual arteries were also tuned to match that reported in the literature (30, 31), which are summarized in Table 1.

The computational geometry was discretized in CF-Mesh+ (version 3.5.1; Creative Fields, Croatia) using an unstructured hexahedral mesh with 5 wall-normal inflation layers at the walls to resolve large velocity gradients. Three different resolutions of meshes were generated yielding coarse, medium, and fine meshes having 3 million, 6.5 million, and 9 million computational cells, respectively (Figures 2B-D). These generated meshes were used to perform a mesh refinement study to evaluate the sensitivity of the results to the mesh resolution. Using the topoSet utility in OpenFOAM, the mesh for the stroke condition was generated by removing the cells covering the section of the arteries that were blocked in the experiments by clots. After eliminating the cells at the clot insertion sites, the newly created surfaces were converted to walls.

Using the 3D reconstructed model, CFD simulations of the flow were performed by solving the Reynolds-averaged Navier-Stokes (RANS) equations using the open-source computational continuum mechanics library, OpenFOAM (version 2106). The gravitational force is not considered in the present simulations. The effects of turbulence were modeled using the two-equation eddy-viscosity k − ω shear stress transport (SST) turbulence model (32, 33). The k − ω SST model was chosen because of its reported good behavior in adverse pressure gradients and separated flows (34), and because it has been shown to be capable of handling multi-regime (laminar, transitional, turbulent) flows (35).

The CFD boundary conditions were specified to match the corresponding in vitro experiments. A constant, uniform inlet velocity was applied to the extended inlet to obtain fully developed flow entering the inlet to the ascending aorta at a flow rate of 5.17 L/min. A no-slip velocity boundary condition was applied on the walls, which were assumed to be rigid. The pressureInletOutletVelocity condition in OpenFOAM was applied for the velocity on all outlet boundaries, which uses a zero-gradient Neumann condition on boundary faces with outflow and an extrapolated Dirichlet condition on faces with reversed inflow. A zero-gradient pressure condition was applied at the inlet and fixed static pressure boundary conditions were prescribed on all outlets. The values of outlet pressure were specified by iteratively performing simulations to match the measured flow rate through each outlet in the experiments to within ±10%.

Simulations were performed in OpenFOAM using second-order accurate discretization schemes, including second-order upwind for the advective terms. Preliminary steady simulations were performed by solving the steady RANS equations of motion using the semi-implicit method for pressure-linked equations (SIMPLE) solver, simpleFoam. The steady simulations did not converge owing to regions of significant unsteady flow in the ascending aorta and the aortic arch. Fully transient flow simulations were performed using the hybrid pressure-implicit with splitting of operators (PISO)-SIMPLE solver, pimpleFoam. An initial unsteady simulation was performed for 4 s of physical time, which was determined to be the time that was needed for the flow to overcome the initial startup transient conditions and to obtain a stationary, quasi-steady flow state. Simulations were then run for another 4 s at these quasi-steady conditions while time-averaging the flow solution. The simulations were performed on 40 compute cores of a high-performance computing (HPC) system at the Institute for Computational and Data Sciences at the Pennsylvania State University. Post-processing and visualization of the simulation results were performed in ParaView (Kitware, Inc., Clifton Park, NY).

The values of volumetric flow rate (L/min) and pressure (mmHg) measured at various arterial outlets are shown in Figure 3A. As mentioned in Section 2.2, the arteries stemming from the aortic arch, which includes the axillary, external, middle, anterior, and other arteries, receive 26.8% of the total flow rate entering the model and 73.2% of the fluid flows through the descending aorta. Out of the total flow rate, 10.23% travels through both axillary arteries, 6.02% through both external carotid arteries (ECA), and 5.62% through the middle and anterior cerebral arteries. This measured flow distribution matches physiological values from the available literature (30, 31). However, it is important to note that these distributions are based on the steady flow measurements, whereas the values reported in the literature are based on total volumetric flow rate values over a pulsatile cardiac cycle.

Changes in mean arterial pressure and volumetric flow rate were also measured for the stroke condition when there is no flow from outlet 4. The measurements showed that the mean pressure at all outlets decreased for the stroke condition. As illustrated in Figures 3, 4, the flow rate through the right arterial outlets as well as outlets 5 and 6 increased to compensate for the flow that was obstructed due to the emulated stroke.

To verify that the k − ω SST turbulence model is behaving as expected and producing turbulent flow in regions of expected turbulence and laminar flow in the downstream cerebral arteries, we calculated the regional distribution of the turbulent viscosity ratio (νt/ν) that is defined as the ratio of turbulent viscosity (νt) to the molecular viscosity (ν). In regions where the k − ω SST model predicts turbulence, νt/ν is much larger than unity and in regions of laminar flow νt/ν is much less than 1. As illustrated in Figure 5, the k − ω SST model predicts turbulent flow in the inlet and in the aorta where the flow is expected to be turbulent. In the downstream cerebral arteries, we observe that the predicted flow laminarizes, as evidenced by the negligible turbulent viscosity ratio, indicating that the turbulence model is inactive in this region, yielding laminar flow as expected.

To evaluate the sensitivity of the CFD results to the mesh resolution, we performed a mesh refinement study to compare the results from the coarse, medium, and fine meshes that contain 3 million, 6.5 million, and 9 million computational cells, respectively. For the same inlet flow rate and with the same outlet pressures applied to each model, we compared the predicted volumetric flow rate through each outlet for the normal physiological condition. The percent differences in the volumetric flow rate at outlet 1, outlet 3, outlet 5, outlet 6, and outlet 9 were calculated to be approximately 1%, 2%, 0.3%, 0.4%, and 0.3%, respectively, between the medium and fine meshes. Between the coarse and medium meshes the percent differences for these same outlets were approximately 3%, 9%, 2%, 3%, and 1%, respectively. Thus, the results between the medium and fine meshes are in close agreement and, therefore, we chose the medium mesh for the final simulations of both the normal and stroke conditions.

Figure 6 illustrates pressure contours and velocity profiles in the cerebrovascular model from CFD for the normal condition. As shown, higher pressures are observed in the aorta compared to the cerebral arteries due to the fact that most of the pressure drop occurs in the smaller cerebral arteries. As the aortic root is attached with a smaller diameter inlet tube, there is some unsteadiness and recirculating flow in the downstream ascending aorta. Generally, the highest velocities are observed in the inlet tube. Flow speeds in the aortic arch are also much higher compared to the more distal downstream arteries. Additionally, it is interesting that the simulation revealed that the flow in the left common carotid artery has a higher velocity than that in the right common carotid.

Quantitative values of volumetric flow rate and pressure at all of the outlets from the CFD simulation at the normal condition are summarized in Figure 3B. The flow distribution in the axillary, external, and combined middle and anterior arteries are equal to 9.62%, 6.87%, and 5.68%, respectively. Because outlet boundary conditions were assigned to match the experimental flow rate measurement through each artery outlet to within ±10%, we validate the simulations by comparing the resultant outlet pressures from CFD with the experimental pressure measurements. As summarized in Tables 2, 3, all of the arterial outlet pressures from CFD are close to the experimental values to within ±7%, with many of the CFD outlet pressures within the uncertainty of the experimental measurement.