The effects of non-Newtonian blood modeling and pulsatility on hemodynamics in the food and drug administration’s benchmark nozzle model

Abstract

BACKGROUND:

Computational fluid dynamics (CFD) is an important tool for predicting cardiovascular device performance. The FDA developed a benchmark nozzle model in which experimental and CFD data were compared, however, the studies were limited by steady flows and Newtonian models.

OBJECTIVE:

Newtonian and non-Newtonian blood models will be compared under steady and pulsatile flows to evaluate their influence on hemodynamics in the FDA nozzle.

METHODS:

CFD simulations were validated against the FDA data for steady flow with a Newtonian model. Further simulations were performed using Newtonian and non-Newtonian models under both steady and pulsatile flows.

RESULTS:

CFD results were within the experimental standard deviations at nearly all locations and Reynolds numbers. The model differences were most evident at Re = 500, in the recirculation regions, and during diastole. The non-Newtonian model predicted blunter upstream velocity profiles, higher velocities in the throat, and differences in the recirculation flow patterns. The non-Newtonian model also predicted a greater pressure drop at Re = 500 with minimal differences observed at higher Reynolds numbers.

CONCLUSIONS:

An improved modeling framework and validation procedure were used to further investigate hemodynamics in geometries relevant to cardiovascular devices and found that accounting for blood’s non-Newtonian and pulsatile behavior can lead to large differences in predictions in hemodynamic parameters.

Keywords

FDA computational fluid dynamics non-Newtonian pulsatility benchmark

1. Introduction

Computational fluid dynamics (CFD) is an important tool for predicting the performance of cardiovascular devices and is increasingly being used to demonstrate safety and efficacy as part of submissions to the US Food and Drug Administration (FDA) [1]. In order to standardize the use of CFD for blood contacting device regulatory submissions, the FDA developed a benchmark nozzle model which contained features commonly seen in medical devices including regions of contraction and expansion, flow recirculation, and a range of flow regimes from laminar to turbulent [1]. Using this model, particle image velocimetry (PIV) flow visualization experiments were performed at three different laboratories [1,2] and round-robin CFD simulation results were submitted by 28 different groups [3–5]. The computational results showed significant variability, standard deviations were greater than 60% of the mean value, and only four of the CFD studies fell within the confidence interval of the experimental data [1].

Several other groups have used the experimental velocity and pressure data to validate their own CFD modeling efforts and evaluate improved transitional and turbulence models [6–15]. All of the aforementioned experimental and computational studies, however, were limited to steady flow conditions and Newtonian blood analogs and constitutive models, which are major limitations for extending the modeling efforts into cardiovascular applications. Trias et al. [16] used CFD to compare the hemodynamics at a Reynolds number of 500 in the FDA nozzle between a standard Newtonian and two non-Newtonian constitutive models. They found differences up to 10% in axial velocity, pressure, and shear stress between the models and that the Newtonian blood model also significantly underpredicted hemolysis compared to the non-Newtonian models. Stiehm et al. [17] also used CFD in a scaled-down version of the FDA nozzle model to mimic a coronary artery geometry and applied sinusoidal and physiological pulsatile boundary conditions. However, they only compared steady and time-averaged velocity results and did not discuss variations over the pulsatile cycle. In more cardiovascular applications, several groups have used CFD to compare hemodynamics between various constitutive blood models including between Newtonian and non-Newtonian models in straight tube [18], carotid bifurcation [19–21], and cerebral aneurysm [22] models. These studies found that the non-Newtonian models produced blunter velocity profiles and that they caused the greatest variations compared to the Newtonian model in low flow and recirculation regions. On the experimental and theoretical side, comparisons have been made between Newtonian and non-Newtonian fluids in straight pipe [23] and stenosis [24,25] models, with all studies reporting a delay in the transition to turbulence for whole blood and non-Newtonian models compared to Newtonian fluid analogs and constitutive models.

In order to provide a thorough understanding of the hemodynamics cardiovascular device geometries will experience in vivo, more accurate constitutive models for blood and pulsatile flow conditions are needed. Therefore, the FDA benchmark nozzle model and its database of available velocity and pressure data [1,5] will first be used to validate steady flow CFD simulations. A Newtonian constitutive model will first be used to match the PIV experimental working fluid and a non-Newtonian blood model will then be incorporated into the CFD solver to investigate the importance of modeling blood’s non-linear behavior. Further, pulsatile inlet conditions will be applied with average Reynolds numbers of 500, 2000, 3500, 5000, and 6500 to further elucidate the differences of the blood models at more physiological flow conditions.

2. Methods

2.1. FDA nozzle model

The FDA benchmark nozzle model (Fig. 1) shares many characteristics with cardiovascular devices including areas of contraction, expansion, flow recirculation, and local high shear stresses. As shown in Fig. 1, the model is in the “Sudden Expansion” orientation (flow from left to right) and acts as a conical concentrator followed by a sudden expansion. After the sudden expansion, flow separates from the wall and regions of recirculation are introduced that reattach further downstream. The nozzle model in this orientation also mimics a stenosis geometry and can be used to generate similar jet and recirculation flow patterns [3]. Experimental PIV data was previously collected at three laboratories using Newtonian blood analog fluids and with an acrylic nozzle model incorporated into a closed recirculating flow loop [1]. Measurements of the two-dimensional velocity fields were collected at several axial and radial sampling planes along the nozzle (Fig. 1) and five different flow conditions were investigated, based on flow rates in the nozzle throat region, spanning the range of laminar (Re = 500), transitional (Re = 2000), and turbulent (Re = 3500, 5000, and 6500) flow regimes. Further details of the experimental data collection can be found in Hariharan et al. [1] and Malinauskas et al. [5].

Fig. 1.

FDA benchmark nozzle model in the Sudden Expansion orientation. Centerline and axial cross-sections (CS’s) are highlighted, as well as insets of the medium computational mesh showing the upstream contraction, throat, and downstream sudden expansion regions.

2.2. Blood modeling

While it is well established that blood is a non-Newtonian fluid, most CFD studies choose to model blood as Newtonian due to the predominance of higher shear rates ( >100 s⁻¹). However, depending on specific designs and operating conditions of cardiovascular devices, regions of low shear rate ( <100 s⁻¹) and non-steady flows often occur and will alter the hemodynamics, ultimately affecting predictions of device induced thrombosis and hemolysis.

Fluid properties from the FDA experimental and round-robin CFD studies were used in the Newtonian model simulations, with a constant viscosity of 3.5 cP and density of 1056 kg/m³ [3]. Non-Newtonian simulations were performed using the Carreau model [26], which captures the asymptomatic behavior of shear-thinning fluids at both high and low shear rates and assumes the local viscosity (𝜂) is a function of the shear rate (Eq. (1)): $\begin{eqnarray}\displaystyle {\eta}(\dot{{\gamma}})={\eta}_{\infty }+({\eta}_{0}+{\eta}_{\infty })(1+({\lambda}\dot{{\gamma}})^{2})^{\frac{m-1}{2}} & & \displaystyle\end{eqnarray}$ (1) where n _∞ is the infinite shear rate viscosity (3.5 cP), 𝜂₀ is the zero-shear rate viscosity (20 cP) [27], 𝜆 is the relaxation time for blood (3.313 s) [28], $\dot{{\gamma}}$ is the scalar shear rate, and m is the power law index that determines the transition between the low and high viscosities (0.3568) [28].

2.3. Computational modeling

OpenFOAM (OpenCFD Ltd.) [29], a C++ open-source computational software, was used for all simulations in this study and were run in parallel on the University of Tennessee’s Advanced Computing Facility high performance computers. For simulating the Newtonian cases, the OpenFOAM solver “pimpleFoam” was used and for the non-Newtonian cases, the OpenFOAM solver “nonNewtonianIcoFoam” was modified to incorporate the PIMPLE algorithm [30] and used. For all flow conditions besides Re = 500 (where no turbulence model was used), the two-equation k–𝜔 SST turbulence model [31] was applied to solve the turbulent eddy viscosity in the Reynolds-averaged Navier–Stokes (RANS) equations. This model has been shown to be optimal for blood pump applications [32] and to outperform both the standard k–𝜔 model and k–𝜖 model in previous studies of the FDA nozzle [3,9].

2.4. Computational grid refinement

Three high-quality unstructured meshes were generated using OpenFOAM’s “snappyHexMesh” utility (coarse ∼640,000; medium ∼1,860,000; and fine ∼5,720,000 cells) to perform a systematic grid refinement study. The interior cells were primarily hexahedral with nearly isotropic sizes of 1, 0.66, and 0.44 mm in the coarse, medium, and fine meshes, respectively, and refined by two levels towards the walls. The nozzle region downstream of the sudden expansion was further refined by a factor of four to improve predictions of the outlet jet’s large velocity gradients and two wall layers were added around the entire model. Images of the medium computational mesh are shown in Fig. 1 with insets highlighting the contraction, nozzle, and sudden expansion regions. As per the FDA guidelines [3], the Re = 5000 flow condition was used with the nozzle in the Sudden Expansion orientation to determine mesh independence.

2.5. Boundary conditions

Five experimental flow rates were applied as constant parabolic-velocity profiles at the FDA nozzle inlet. The velocities were set to achieve Reynolds numbers of 500, 2000, 3500, 5000, and 6500 in the narrower throat region of the nozzle model. To simulate pulsatile flow, sinusoidal velocity waveforms were generated (Fig. 2) with average Reynolds numbers of 500, 2000, 3500, 5000, and 6500 in the nozzle throat region to match the steady flow conditions (Eq. (2)): $\begin{eqnarray}\displaystyle u(t)=A\ast \sin \left(2{\pi}ft+\frac{3{\pi}}{2}\right)+A & & \displaystyle\end{eqnarray}$ (2) where u (t) is the axial velocity, t is time, A is the waveform magnitude that matches each case’s respective steady inlet velocity (0.0461 m/s for Re = 500, 0.184 m/s for Re = 2000, 0.322 m/s for Re = 3500, 0.461 m/s for Re = 5000, and 0.599 m/s for Re = 6500), and f is the waveform frequency (1 Hz for all simulations to mimic a 60 beats per minute cardiovascular flow). As shown in the Fig. 2 legend, the name for each waveform (Re = 500, 2000, 3500, 5000, and 6500) corresponds to the average Reynolds number over the entire pulsatile cycle. However, during peak systole, the actual Reynolds number in the throat region is twice as high (1000, 4000, 7000, 10,000, and 13,000, respectively). With these flow waveforms and fluid properties, the Womersley number (𝛼) (Eq. (3)), which is the ratio of the transient or oscillatory inertial forces to the viscous shear forces [33], can be calculated at different locations in the FDA nozzle model: $\begin{eqnarray}\displaystyle {\alpha}=L\left(\frac{{\omega}{\rho}}{{\mu}}\right)^{1/2} & & \displaystyle\end{eqnarray}$ (3) where L is a characteristic length scale (radius of the nozzle model), 𝜔 is the angular frequency (2𝜋), 𝜌 is the fluid density (1056 kg/m³), and 𝜇 is the fluid dynamic viscosity (3.5 cP). When 𝛼 is small ( <1), the frequency of pulsations is low enough to allow for the development of a parabolic velocity profile during each cycle. When 𝛼 is large ( >10), the frequency of pulsations is large, and the resultant velocity profiles will be relatively flat during each cycle. For the upstream and downstream portions of the nozzle model with a 6 mm radius, 𝛼 = 8.26. For the nozzle throat with a 2 mm radius, 𝛼 = 2.75. These Womersley numbers correlate well with those seen physiologically in medium sized arteries (1 < 𝛼 < 10 [34]) and specifically, the throat Womersley number of 2.75 is comparable to what has been observed in coronary arteries (𝛼 ≈ 3) [35]. Additional computational boundary conditions include no-slip velocities applied to all the nozzle model walls and the outlet boundary condition was a fixed pressure of 0 Pa. During post-processing of the simulation results, the CFD pressures were plotted relative to the pressure at the location of the sudden expansion (z = 0.0 m), following previous FDA studies [3].

Fig. 2.

Pulsatile velocity (m/s) waveforms applied as inlet boundary conditions for the Re = 500, 2000, 3500, 5000, and 6500 cases, respectively. Four time points where pressure and velocity data were extracted are highlighted on the graph (Peak Systole, Mid-deceleration, Peak Diastole, and Mid-acceleration).

3. Results

3.1. CFD grid refinement

A systematic grid study was performed to verify that the numerical results in the CFD simulations were independent of the mesh they were solved on. Three meshes (coarse, medium, and fine) were created in OpenFOAM consisting of approximately 0.64, 1.86, and 5.72 million cells, respectively. Hemodynamic parameters of pressure and velocity were each evaluated for the steady Re = 5000 flow condition using a Newtonian blood model, with the nozzle in the Sudden Expansion orientation following FDA guidelines [3]. Centerline pressures (Fig. 3E) and velocities (Fig. 3F) and axial velocity profiles taken at nozzle CS’s 2, 3, 4, 6, 7, 8, 9, 10, 11, and 12 (Fig. 3A–D and G–J) are shown for each mesh. Along the centerline, the average difference in pressure between the coarse and medium meshes was 4.8% and between the medium and fine meshes was 2.5% (Fig. 3E). Similarly, the average difference in centerline axial velocity between the coarse and medium meshes was 1.2% and between the medium and fine meshes was 0.7% (Fig. 3F). Along the radial CS’s, the average differences in axial velocity between the coarse and medium meshes were: CS2 = 7.6%, CS3 = 2.7%, CS4 = 1.8%, CS6 = 2.8%, CS7 = 10.6%, CS8 = 5.1%, CS9 = 6.5%, CS10 = 10.5%, CS11 = 0.9%, CS12 = 0.9% (Fig. 3A–D and G–J). The average differences in axial velocity between the medium and fine meshes were: CS2 = 1.7%, CS3 = 1.3%, CS4 = 1.6%, CS6 = 2.6%, CS7 = 9.4%, CS8 = 3.8%, CS9 = 4.6%, CS10 = 6.6%, CS11 = 0.6%, CS12 = 0.7% (Fig. 3A–D and G–J). Based on the relatively small differences observed between the medium and fine mesh pressure and velocity results, the medium mesh was used for all further CFD simulations.

Fig. 3.

Grid refinement pressure and velocity results for the coarse, medium, and fine meshes at the steady Re = 5000 condition with Newtonian blood model. Pressures (Pa) were compared along the (E) nozzle centerline while axial velocities (m/s) were compared at (A) CS2, (G) CS3, (B) CS4, (H) CS6, (F) centerline, (C) CS7, (I) CS8, (D) CS9, and (J) CS12.

3.2. Validation of CFD with FDA experimental data

3.2.1. Pressure validation results

Experimental pressures [1,5] along the nozzle centerline were compared to those predicted by CFD at Re = 500 (Fig. 4E) and Re = 6500 (Fig. 5E) with a Newtonian blood model. Pressure validation results for Re = 2000, 3500, and 5000 are included in Supplemental Material section (Figs S1, S2, and S3, respectively). The multiple FDA experimental pressure data sets were averaged and plotted along with their standard deviations. Additionally, all pressures were plotted relative to the pressure at z = 0.0 m (the location of the sudden expansion).

Fig. 4.

Comparison of the FDA experimental pressure (black squares) and velocity (black lines) data with the CFD pressure (gray triangles) and velocity (gray lines) results at the steady Re = 500 condition with Newtonian blood model. Pressures (Pa) were compared along the (E) nozzle centerline while axial velocities (m/s) were compared at (A) CS2, (G) CS3, (B) CS4, (H) CS6, (F) centerline, (C) CS7, (I) CS8, (D) CS9, and (J) CS12.

In all experimental and CFD cases, a large pressure drop occurs in the nozzle contraction and continues to gradually decrease over the length of the throat region. The lowest pressure in each case occurred at the sudden expansion with a small level of pressure recovery observed downstream. For Re = 500 (Fig. 4E), the CFD pressures were within the experimental standard deviations along the centerline except the last three downstream locations, where they underpredicted pressures by approximately 30 Pa. Similarly, at Re = 6500 (Fig. 5E), the CFD pressures were also within the experimental standard deviations at all locations besides at the entrance to the throat region and the last two locations downstream of the sudden expansion, where the CFD overpredicted the experimental pressures.

Fig. 5.

Comparison of the FDA experimental pressure (black squares) and velocity (black lines) data with the CFD pressure (gray triangles) and velocity (gray lines) results at the steady Re = 6500 condition with Newtonian blood model. Pressures (Pa) were compared along the (E) nozzle centerline while axial velocities (m/s) were compared at (A) CS2, (G) CS3, (B) CS4, (H) CS6, (F) centerline, (C) CS7, (I) CS8, (D) CS9, and (J) CS12.

3.2.2. Velocity validation results

Experimental axial velocities [1,5] along the nozzle centerline (Figs 4F and 5F) and at CS’s 2, 3, 4, 6, 7, 8, 9, and 12 (Fig. 4A–D and G–J, Fig. 5A–D and G–J) were compared to those predicted by CFD at Re = 500 (Fig. 4) and Re = 6500 (Fig. 5) with a Newtonian blood model. Velocity validation results for Re = 2000, 3500, and 5000 are included in the Supplemental Material section (Figs S1, S2, and S3, respectively). The multiple FDA experimental velocity data sets were averaged and plotted along with their standard deviations. For all cases, centerline axial velocities showed a uniform flow upstream of the contraction region, accelerating velocities in the contraction and throat regions as the flow developed, and decreasing velocities with further axial position away from the sudden expansion.

For Re = 500, the CFD centerline velocities matched the experiments within their standard deviations at all axial locations (Fig. 4F). The simulations were also able to accurately capture the low negative velocities within the recirculation zones at CS’s 7, 8, 9, and 12 (Figs 4C, 4I, 4D, and 4J). At Re = 6500, the CFD centerline velocities matched the experiments within their standard deviations at most axial points but underpredicted the velocities in the nozzle throat (CS4 and CS5) just prior to the sudden expansion location (6.8 m/s compared to approximately 7.1 ± 0.2 m/s) (Fig. 5F). The simulation accurately predicted the axial velocity profiles upstream of the throat region (CS2 and CS3 in Figs 5A and 5G, respectively) and captured the overall velocity profile shapes downstream of the sudden expansion (CS7 and CS8 in Figs 5C and 5I, respectively) while overpredicting the amount of negative velocity (−0.76 m/s versus −0.46 m/s at CS7, and −0.37 m/s versus −0.11 m/s at CS8) in the near wall recirculation regions and predicting a slower transition out of the jet shear layer. As the flow begins to re-develop at CS9 and retains its blunt profile at CS12, (Figs 5D and 5J, respectively) the velocity profiles were accurately predicted by the CFD within the experimental standard deviations.

With all the presented pressure and velocity comparisons between the FDA experimental results and the CFD simulations at a range of laminar and turbulent Reynolds number conditions, the CFD modeling efforts were considered to be thoroughly validated for the Newtonian blood model and will be used to simulate non-Newtonian blood constitutive models and pulsatile flow conditions. Additional qualitative validation of the Newtonian and non-Newtonian models and solvers was performed by comparing to the experimental results of Biswas et al. [23] and numerical results of Tu and Deville [36] (Fig. S4 in the Supplemental Materials section). To do this, additional simulations were performed at Re = 2500 ad 3500 for comparisons to Biswas et al. [23] and at Re = 67 for comparisons to Tu and Deville [36]. The velocity profiles at CS5 in the nozzle throat region compared well with those of Biswas et al. [23] for both their whole blood and Newtonian fluids in a straight tube, increasing in bluntness with increasing Reynolds number (Fig. S4-A). Additionally, greater levels of turbulence were observed at lower Reynolds numbers with the Newtonian model compared to the non-Newtonian model, highlighted by larger levels of turbulent kinetic energy (TKE) (8 and 19% higher with the Newtonian model at Re = 2500 and 3000, respectively) in Fig. S4-A. Tu and Deville investigated pulsatile velocity profiles at a laminar Reynolds number of 67 between Newtonian and non-Newtonian constitutive models in 25% and 75% stenosis geometries [36]. An additional set of simulations were run at Re = 67 and the velocity magnitudes and profiles at the stenosis/expansion location and at locations two and five diameters downstream compared well to those of Tu and Deville (Fig. S4-B).

3.3. Comparison of Newtonian and non-Newtonian blood models

3.3.1. Steady pressure results

Centerline pressures were evaluated at all five Reynolds number flow conditions using the Newtonian and non-Newtonian blood models (Fig. 6). The greatest variations were observed at Re = 500 (Fig. 6A) where a difference of 44.5% was observed between the Newtonian and non-Newtonian models. At the higher Reynolds number conditions, however, smaller pressure differences were observed between the Newtonian and non-Newtonian models (0.1–0.6%). Both models also showed large pressure losses in the throat region while the non-Newtonian model showed greater pressure recovery following the sudden expansion compared to the Newtonian model. Pressure magnitudes at three key locations in the FDA nozzle model (upstream at CS2, in the nozzle throat at CS4, and downstream of the expansion at CS8 and CS10) are included in Table S.1 of the Supplemental Materials section to highlight the variability between constitutive models.

Fig. 6.

Steady centerline pressure (Pa) results at (A) Re = 500, (B) Re = 2000, (C) Re = 3500, (D) Re = 5000, and (E) Re = 6500 between the Newtonian (black squares)and non-Newtonian (gray triangles) blood models. [Note the different pressure scales on the y-axis for each Reynolds number.]

3.3.2. Steady velocity results

Axial velocity profiles were evaluated at the lower three Reynolds number flow conditions using the Newtonian and non-Newtonian blood models at CS4 (throat region), CS7 (recirculation region), and CS12 (reattachment region) (Fig. 7). Velocity data for the remaining Reynolds number conditions and nozzle CS’s did not show as much variability and are included in the Supplemental Materials section (Figs S5 and S6). The greatest differences occurred at Re = 500 for each of the CS’s. At CS2, both models had fully developed velocity profiles. Initially the non-Newtonian model’s was blunter than the Newtonian model but that trend reversed as the flow moved into the converging nozzle in CS’s 3–6 (Fig. S5-ii, iii, iv). Differences at the higher Reynolds numbers were less than 1%.

Fig. 7.

Steady axial velocity (m/s) results at (A) Re = 500, (B) Re = 2000, and (C) Re = 3500 between the Newtonian (gray lines) and non-Newtonian (black dashed lines) blood models at (i) CS4, (ii) CS7, and (iii) CS12. [Note the different velocity scales on the x-axis for each Reynolds number.]

At CS4 in the nozzle throat region, the non-Newtonian model displayed larger peak velocities than the Newtonian model (0.76 m/s versus 0.68 m/s) for Re = 500 (Fig. 7A-i). At higher Reynolds numbers, however, the Newtonian and non-Newtonian models predicted similar velocities. The exact same trends were observed at CS7, where the non-Newtonian model displayed higher peak velocities at Re = 500 (Fig. 7A-ii) but nearly identical velocities at higher Reynolds numbers. Additionally, the non-Newtonian model showed slightly larger negative velocities in the recirculation regions at lower Reynolds numbers (Re = 500, 2000, and 3500) that were not observed at the highest Reynolds numbers. At CS12, the flow was re-developed for all cases besides Re = 500. In this case, the non-Newtonian model was closer to redeveloping its flow profile, displayed lower peak velocities (0.28 m/s compared to 0.51 m/s for the Newtonian case), and had only minor negative velocities remaining in the near wall regions (Fig. 7A-iii). At all other Reynolds numbers, both models predicted nearly identical velocity profiles at CS12 with differences less than 1%.

For direct biorheology applications, specifically for the modeling of blood damage, it is important to be able to accurately predict fluid shear stresses based off the velocity field. The peak stresses for all cases occurred at the entrance to the nozzle throat region, with additional high shear stress regions along the walls of the nozzle throat, and in the jet edges downstream of the sudden expansion (shear stress values reported in Table S.2 in the Supplemental Materials section). For Re = 500, the non-Newtonian model predicted lower stresses at the nozzle entrance but slightly higher stresses at the other Reynolds numbers and locations compared to the Newtonian model ( <2% differences at the nozzle entrance and in the nozzle throat but up to 5% difference in the downstream jet).

3.3.3. Pulsatile pressure results

Five sinusoidal velocity waveforms (Fig. 2) were applied to achieve average Reynolds numbers of 500, 2000, and 3500 in the nozzle throat region. Centerline pressures were evaluated at four distinct time points in each waveform (peak systole, mid-deceleration, peak diastole, and mid-acceleration) with the Newtonian and non-Newtonian blood models (Fig. 8). Pulsatile data for the other Reynolds number conditions and at the peak-diastole time point are included in the Supplemental Materials section (Fig. S7). The greatest differences were observed at Re = 500 for all time points, and at Re = 2000 and 3500 during mid-deceleration and mid-acceleration. For Re = 500, the non-Newtonian model had much higher upstream pressures compared to the Newtonian model at mid-acceleration (785 Pa versus 469 Pa) (Fig. 8A-iii) and peak systole (1273 Pa versus 947 Pa) (Fig. 8A-i). The opposite trend existed at mid-deceleration (137 Pa versus 161 Pa) (Fig. 8A-ii) and peak diastole (−7 Pa versus −16 Pa). Overall, the same trends existed at the higher Reynolds number conditions but with smaller differences observed between the blood models; (1) at mid-acceleration and peak systole, the non-Newtonian model had higher upstream pressures, (2) at mid-deceleration and peak diastole, the Newtonian model had the highest upstream pressures. At other time points, the predicted pressures were between 13 and 35%.

Fig. 8.

Pulsatile centerline pressure (Pa) results at (A) Re = 500, (B) Re = 2000, and (C) Re = 3500 between the Newtonian (black squares) and non-Newtonian (gray triangles) blood models during (i) Peak Systole, (ii) Mid-deceleration, and (iii) Mid-acceleration. [Note the different pressure scales on the y-axis for each Reynolds number.]

3.3.4. Pulsatile velocity results

Axial velocities were evaluated at four distinct time points in each waveform (peak systole, mid-deceleration, peak diastole, and mid-acceleration) with the Newtonian and non-Newtonian blood models (Figs 9–12). The greatest differences were observed at Re = 500 for all time points. At CS2 (Fig. 9A), the non-Newtonian model displayed a blunter velocity profile and predicted lower peak velocities at all four time points compared to the Newtonian model. Similar to the steady data, at CS4 in the nozzle throat (Fig. 9B), the non-Newtonian model predicted the higher velocities at peak systole (1.35 m/s versus 1.21 m/s), mid-deceleration (0.79 m/s versus 0.69 m/s), and mid-acceleration (0.7 m/s versus 0.65 m/s). The largest differences were observed at CS7 (Fig. 9C), where both blood models produced distinctly different velocity profiles. The non-Newtonian model had the highest velocities at peak systole and mid-deceleration, while the Newtonian model had the highest at mid-acceleration. At CS12 (Fig. 9D), the Newtonian model velocities were higher at peak systole but also had large negative velocity recirculation regions that were not observed in the other model which was more quickly redeveloping its parabolic velocity profile.

Fig. 9.

Comparison of the pulsatile axial velocity (m/s) results for the Re = 500 condition at (A) CS2, (B) CS4, (C) CS7, and (D) CS12 between the Newtonian (gray lines) and non-Newtonian (black dashed lines) blood models during (i) Peak Systole, (ii) Mid-deceleration, (iii) Peak Diastole, and (iv) Mid-acceleration. [Note the different velocity scales on the x-axis for each cross-section.]

Fig. 10.

Comparison of the pulsatile axial velocity (m/s) results for the Re = 2000 condition at (A) CS2, (B) CS4, (C) CS7, and (D) CS12 between the Newtonian (gray lines) and non-Newtonian (black dashed lines) blood models during (i) Peak Systole, (ii) Mid-deceleration, (iii) Peak Diastole, and (iv) Mid-acceleration. [Note the different velocity scales on the x-axis for each cross-section.]

Fig. 11.

Comparison of the pulsatile axial velocity (m/s) results at peak diastole and (A) CS2, (B) CS4, (C) CS7, and (D) CS12 between the Newtonian (gray lines) and non-Newtonian (black dashed lines) blood models for the (i) Re = 3500, (ii) Re = 5000, and (iii) Re = 6500 flow conditions.

Fig. 12.

Comparison of the pulsatile axial velocity (m/s) results at Mid-acceleration and (A) CS2, (B) CS7 Recirculation Region, (C) CS7 Jet Region, and (D) CS12 between the Newtonian (gray lines) and non-Newtonian (black dashed lines) blood models for the (i) Re = 3500, (ii) Re = 5000, and (iii) Re = 6500 flow conditions. The velocity results at CS7 were separated and the x-axes scaled to highlight differences between the blood models. [Note the different velocity scales on the x-axis for each cross-section.]

At a higher Reynolds number of 2000 (Fig. 10), large differences between constitutive models were still observed at CS2, CS7, and CS12 throughout the pulsatile cycle while the differences at CS4 (Fig. 10B) were less than 1%. At CS2 (Fig. 10A), the non-Newtonian model had higher velocities at peak diastole and mid-acceleration but lower velocities at peak systole and mid-deceleration compared to the Newtonian model. At CS7 (Fig. 10C), both models were nearly identical at peak systole (4.3 m/s), but the non-Newtonian model predicted higher velocities at mid-deceleration (2.5 m/s versus 2.44 m/s) and peak diastole (0.3 m/s versus 0.17 m/s) and lower velocities at mid-acceleration (2.22 m/s versus 2.31 m/s). Smaller differences of 2% were observed between the models during peak systole at CS12 (Fig. 10D).

In contrast to the previous pulsatile results, the three highest Reynolds number conditions (Re = 3500, 5000, and 6500) produced velocity differences at mid-acceleration, peak systole, and mid-deceleration at CS2 that were less than 1%. For CS12, at peak systole and mid-deceleration, differences were less than 2% as well. Therefore, only the results from peak diastole (Fig. 11) and mid-acceleration (Fig. 12) are shown to highlight the constitutive model differences. For all three Reynold’s numbers and CS’s in Fig. 11, the non-Newtonian model predicted much higher velocities at the center of the nozzle during diastole compared to the Newtonian model. Additionally, the non-Newtonian model had higher negative velocities in the recirculation regions. At Re = 3500 (Fig. 12A-i) and Re = 5000 (Fig. 12A-ii), the non-Newtonian model had higher velocities while at Re = 6500 (Fig. 12A-iii) the Newtonian model had higher velocities. At CS7 (Fig. 12B, C), for Re = 3500, 5000 and 6500, the Newtonian model had higher peak velocities and greater negative recirculation velocities. At CS12 (Fig. 12D), for Re = 3500, 5000 and 6500 (Fig. 12D), the non-Newtonian model had lower velocities and a much blunter profile.

4. Discussion

This is the first study to utilize the FDA’s benchmark nozzle model to evaluate the combined effects of pulsatile flow and constitutive blood model. The results further elucidate the hemodynamics that contraction and expansion type cardiovascular device geometries experience in vivo under more physiological conditions. Trias et al. [16] previously performed CFD simulations to compare a Newtonian and two non-Newtonian constitutive models at Re = 500 and found differences up to 10% in axial velocity, pressure, and shear stress. They also found that pressure drops were larger with the non-Newtonian models at the sudden expansion compared to the Newtonian model and that centerline velocities were much lower downstream of the sudden expansion with the Casson model and Carreau models. Similar to those results, the steady Re = 500 results in Sections 3.3.1. and 3.3.2. showed larger pressure drops with the non-Newtonian model and lower velocities downstream of the sudden expansion at CS12. Trias et al. [16] also stated that incorporation of a non-Newtonian blood model added minimal computational cost to their simulations. Those same trends were observed in this study, where the non-Newtonian simulation times were approximately 14% longer compared to the Newtonian simulations, and thus should not be used as a limitation in future studies.

In other cardiovascular applications, several groups have investigated hemodynamic differences resulting from pulsatility and non-Newtonian blood modeling. Gijsen et al. [19] conducted experimental and numerical simulations of a non-Newtonian fluid in a carotid bifurcation model and reported considerable differences compared to a Newtonian fluid model. In particular, the non-Newtonian velocity profiles were flatter with lower maximum velocities. Chen and Lu [20] and Lee and Steinman [21] also reported a pronounced influence of non-Newtonian pulsatile flow in carotid bifurcation models, with the non-Newtonian flows exhibiting flattened axial velocity profiles. Additionally, they acknowledged that while the overall rheological effects were relatively small for time-averaged fields, differences in instantaneous velocity profiles were much more pronounced at the diastolic phase when relatively low shear rates were present. The results of Gijsen et al. [19], Chen and Lu [20], and Lee and Steinman [21] upstream of the carotid bifurcations are similar to those shown in Fig. 8, where the non-Newtonian model displayed blunter velocity profiles that approached the Newtonian velocities at higher Reynolds numbers. The rheological differences were also most prominent here during diastole (Figs 10-iii, 11-iii, and 12) with lower shear rates present in the flow fields.

One of the reasons for the minor differences observed between the blood models at higher Reynolds number flow conditions was the dominance of the turbulent eddy viscosity over the fluid viscosity. For the Newtonian model (constant fluid viscosity of 3.5e⁻³ Pa ⋅ s), peak turbulent viscosities were observed in the downstream jet and ranged between 1.59e⁻¹ and 4.9e⁻¹ Pa ⋅ s for the Re = 2000 to 6500 flow conditions. In the nozzle throat, the turbulent viscosities were closer in magnitude to the fluid viscosity and were only higher for the Re = 5000 and 6500 conditions (9.96e⁻³ and 1.57e⁻² Pa ⋅ s, respectively). In all other regions of the model, the turbulent viscosities were much less than the fluid viscosity and the hemodynamic differences between constitutive blood models could be realized as a result of the variable non-Newtonian model fluid viscosities.

5. Limitations

As displayed in the analysis of the FDA PIV data [1], there was still variability among the three independent laboratories’ velocity measurements. Therefore, validation of the flow fields at certain Reynolds numbers and in certain regions of the nozzle are much more uncertain and will require further investigation in the transitional and turbulent flow regimes. The use of the laminar modeling for the Re = 500 simulations and a k–𝜔 SST turbulence model for all higher Reynolds number simulations is also a limitation of this CFD work. At both the Re = 2000 steady condition and during points in the higher Reynolds number pulsatile cycles, the flow fields are clearly in a transitional regime between fully laminar and fully turbulent and thus may not be accurately modeled. Additional work should be done with advanced transitional and turbulence models to be able to capture this range of flow fields. Additionally, the pulsatile waveforms used in this study were sinusoidal with only positive forward flow and at a single beat rate (60 bpm or 1 Hz). Further work should incorporate more physiological cardiac waveforms that cover a range of beat rates to fully elucidate the effects of blood’s non-Newtonian behavior.

6. Conclusions

In order to model hemodynamics within cardiovascular devices, CFD simulations must be thoroughly validated against experimental data. They must also be able to simulate physiologically relevant flow conditions and blood properties. Therefore, the FDA benchmark nozzle and its database of experimental data were used to thoroughly validate CFD simulations under steady flow using a Newtonian blood model. Pulsatile flow conditions and non-Newtonian blood models were then incorporated into the CFD simulations and compared to the simplified Newtonian model results to investigate their relative importance. The results suggest that accounting for blood’s non-Newtonian behavior is important at certain low velocity and recirculating flow regions, mainly at the Re = 500 flow condition, at CS7 downstream of the sudden expansion, and during diastole in the pulsatile simulations. The non-Newtonian model predicted blunter upstream velocity profiles, higher velocities in the contraction and throat regions, and differences in the downstream recirculation flow patterns that changed with Reynolds number. Overall, this improved modeling framework, and an extensive validation procedure should be used in future CFD studies of geometries relevant to blood-contacting cardiovascular devices to aid in their design and predictions of their biological responses.

Footnotes

Acknowledgements

Supported by the University of Tennessee's Mechanical, Aerospace, and Biomedical Engineering Department.

Disclosure

The author declares no conflict of interest.

Supplementary material

The supplementary material is available in the electronic version of this article: .

References

Hariharan

Giarra

Reddy

Day

Manning

Deutsch

Multilaboratory particle image velocimetry analysis of the FDA benchmark nozzle model to support validation of computational fluid dynamics simulations. Journal of Biomechanical Engineering. 2011;133(4):041002.

Stewart

Day

Burgreen

Paterson

Manning

Hariharan

Preliminary results of FDA’s “Critical Path” project to validate computational fluid dynamic methods used in medical device evaluation. ASAIO J. 2009;55(2):173.

Stewart

Paterson

Burgreen

Hariharan

Giarra

Reddy

Assessment of CFD performance in simulations of an idealized medical device: results of FDA’s first computational interlaboratory study. Cardiovascular Engineering and Technology. 2012;3(2):139–60.

Stewart

Hariharan

Paterson

Burgreen

Reddy

Day

Results of FDA’s first interlaboratory computational study of a nozzle with a sudden contraction and conical diffuser. Cardiovascular Engineering and Technology. 2013;4(4):374–91.

Malinauskas

Hariharan

Day

Herbertson

Buesen

Steinseifer

FDA benchmark medical device flow models for CFD validation. ASAIO Journal. 2017;63(2):150–60.

Zmijanovic

Mendez

Moureau

Nicoud

. About the numerical robustness of biomedical benchmark cases: interlaboratory FDA’s idealized medical device. International Journal for Numerical Methods in Biomedical Engineering. 2017;33(1):e02789.

Bhushan

Walters

Burgreen

. Laminar, turbulent, and transitional simulations in benchmark cases with cardiovascular device features. Cardiovascular Engineering and Technology. 2013;4(4):408–26.

Jain

. Efficacy of the FDA nozzle benchmark and the lattice Boltzmann method for the analysis of biomedical flows in transitional regime. Medical & Biological Engineering & Computing. 2020;58:1–4.

Janiga

. Large eddy simulation of the FDA benchmark nozzle for a Reynolds number of 6500. Computers in Biology and Medicine. 2014;47:113–9.

10.

Tobin

Manning

. Large-eddy simulations of flow in the FDA benchmark nozzle geometry to predict hemolysis. Cardiovascular Engineering and Technology. 2020;11:1–4.

11.

Fehn

Wall

Kronbichler

. Modern discontinuous Galerkin methods for the simulation of transitional and turbulent flows in biomedical engineering: A comprehensive LES study of the FDA benchmark nozzle model. International Journal for Numerical Methods in Biomedical Engineering. 2019;35(12):e3228.

12.

Sánchez Abad

Vinuesa

Schlatter

Andersson

Karlsson

. Simulation strategies for the food and drug administration nozzle using Nek5000. AIP Advances. 2020;10(2):025033.

13.

Manchester

. The effect of turbulence on transitional flow in the FDA’s benchmark nozzle model using large-eddy simulation. International Journal for Numerical Methods in Biomedical Engineering. 2020;36(10):e3389.

14.

Bergersen

Mortensen

Valen-Sendstad

. The FDA nozzle benchmark: “In theory there is no difference between theory and practice, but in practice there is”. International Journal for Numerical Methods in Biomedical Engineering. 2019;35(1):e3150.

15.

Delorme

Anupindi

Frankel

. Large eddy simulation of FDA’s idealized medical device. Cardiovascular Engineering and Technology. 2013;4(4):392–407.

16.

Trias

Arbona

Massó

Miñano

Bona

. FDA’s nozzle numerical simulation challenge: Non-Newtonian fluid effects and blood damage. PLoS One. 2014;9(3):e92638.

17.

Stiehm

Wüstenhagen

Siewert

Grabow

Schmitz

. Numerical simulation of pulsatile flow through a coronary nozzle model based on FDA’s benchmark geometry. Current Directions in Biomedical Engineering. 2017;3(2):775–8.

18.

Tabakova

Nikolova

Radev

, Carreau model for oscillatory blood flow in a tube. In: AIP Conference Proceedings. Vol. 1629, No. 1, American Institute of Physics; 2014. pp. 336–43.

19.

Gijsen

van de Vosse

Janssen

. The influence of the non-Newtonian properties of blood on the flow in large arteries: steady flow in a carotid bifurcation model. Journal of Biomechanics. 1999;32(6):601–8.

20.

Chen

. Numerical investigation of the non-Newtonian pulsatile blood flow in a bifurcation model with a non-planar branch. Journal of Biomechanics. 2006;39(5):818–32.

21.

Lee

Steinman

. On the relative importance of rheology for image-based CFD models of the carotid bifurcation. Journal of Biomechanical Engineering. 2007;129(2):273–8.

22.

Cebral

Castro

Appanaboyina

Putman

Millan

Frangi

. Efficient pipeline for image-based patient-specific analysis of cerebral aneurysm hemodynamics: Technique and sensitivity. IEEE Transactions on Medical Imaging. 2005;24(4):457–67.

23.

Biswas

Casey

Crowder

Steinman

Yun

Loth

. Characterization of transition to turbulence for blood in a straight pipe under steady flow conditions. Journal of Biomechanical Engineering. 2016;138(7).

24.

Khan

Valen-Sendstad

Steinman

. Direct numerical simulation of laminar-turbulent transition in a non-axisymmetric stenosis model for Newtonian vs. shear-thinning non-Newtonian rheology’s. Flow, Turbulence and Combustion. 2019;102(1):43–72.

25.

Yao

, Transition to turbulence within an eccentric stenosis geometry under steady flow conditions for a non-Newtonian and Newtonian fluid, Master Thesis, The University of Akron, Akron, 2020.

26.

Bird

Armstrong

Hassager

, Dynamics of Polymeric Liquids. Vol. 1: Fluid Mechanics1987. John Wiley and Sons Inc; New York, NY (USA); ISBN 0-471-80245-X.

27.

Good

Deutsch

Manning

. Hemodynamics in a pediatric ascending aorta using a viscoelastic pediatric blood model. Ann Biomed Eng. 2016;7(1):1019–35.

28.

Shibeshi

Collins

. The rheology of blood flow in a branched arterial system. Applied Rheology. 2005;15(6):398–405.

29.

Greenshields

, OpenFOAM User Guide. Version 3.1. OpenFOAM Foundation Ltd.

30.

Ferziger

Perić

Street

. Computational Methods for Fluid Dynamics. Berlin: Springer; 2002.

31.

Menter

. Two-equation eddy-viscosity turbulence models for engineering applications. AIAA Journal. 1994;32(8):1598–605.

32.

Song

Wood

Day

Olsen

. Studies of turbulence models in a computational fluid dynamics model of a blood pump. Artificial Organs. 2003;27(10):935–7.

33.

Womersley

. Method for the calculation of velocity, rate of flow and viscous drag in arteries when the pressure gradient is known. The Journal of Physiology. 1955;127(3):553.

34.

Caro

Pedley

Schroter

Seed

Parker

. The Mechanics of the Circulation. 2nd ed. Cambridge: Cambridge University Press; 2011 9781139013406.

35.

. Blood flow in arteries. Annual Review of Fluid Mechanics. 1997;29(1):399–434.

36.

Deville

. Pulsatile flow of non-Newtonian fluids through arterial stenoses. Journal of Biomechanics. 1996;29(7):899–908.