Article Text
Abstract
Objective Computational fluid dynamics (CFD) simulations of intracranial aneurysm hemodynamics usually adopt the simplification of the Newtonian blood rheology model. A study was undertaken to examine whether such a model affects the predicted hemodynamics in realistic intracranial aneurysm geometries.
Methods Pulsatile CFD simulations were carried out using the Newtonian viscosity model and two nonNewtonian models (Casson and HerschelBulkley) in three typical internal carotid artery saccular aneurysms (A, sidewall, oblongshaped with a daughter sac; B, sidewall, quasispherical; C, nearspherical bifurcation). For each aneurysm model the surface distributions of shear rate, blood viscosity and wall shear stress (WSS) predicted by the three rheology models were compared.
Results All three rheology models produced similar intraaneurysmal flow patterns: aneurysm A had a slowly recirculating secondary vortex near the dome whereas aneurysms B and C contained only a large single vortex. All models predicted similar shear rate, blood viscosity and WSS in parent vessels of all aneurysms and in the sacs of B and C. However, large discrepancies in shear rate, viscosity and WSS among predictions by the various rheology models were found in the dome area of A where the flow was relatively stagnant. Here the Newtonian model predicted higher shear rate and WSS values and lower blood viscosity than the two nonNewtonian models.
Conclusions The Newtonian fluid assumption can underestimate viscosity and overestimate shear rate and WSS in regions of stasis or slowly recirculating secondary vortices, typically found at the dome in elongated or complexshaped saccular aneurysms as well as in aneurysms following endovascular treatment. Because low shear rates and low WSS in such flow conditions indicate a high propensity for thrombus formation and rupture, CFD based on the Newtonian assumption may underestimate the propensity of these events.
 Casson model
 HerschelBulkley model
 nonNewtonian
 wall shear stress
 shear rate
 viscosity
 blood rheology
 thrombosis
 aneurysm
 artery
 blood flow
 vessel wall
 brain
 tumor
 spine
 thrombectomy
 catheter
Statistics from Altmetric.com
 Casson model
 HerschelBulkley model
 nonNewtonian
 wall shear stress
 shear rate
 viscosity
 blood rheology
 thrombosis
 aneurysm
 artery
 blood flow
 vessel wall
 brain
 tumor
 spine
 thrombectomy
 catheter
Introduction
Intracranial aneurysms (IAs) are pathological dilations of arterial walls that affect approximately 2–5% of the entire population.1 ,2 Ruptured IAs cause subarachnoid hemorrhage and its sequelae, resulting in significant morbidity and mortality.3–5 Hemodynamic stresses such as wall shear stress (WSS, the tangential frictional stress caused by the action of flowing blood on the vessel wall endothelium) have been shown to play an important role in IA pathophysiology of initiation and rupture.6–9 To this end, computational fluid dynamics (CFD) has been widely used to obtain patientspecific flow fields in IAs to assess potential risk of rupture.8–12
CFD simulations often involve simplifying assumptions of blood properties and boundary conditions. One of the commonly adopted simplifications in large vessels is the Newtonian fluid model which prescribes a linear shear stressstrain rate relationship (constant viscosity) for blood. Although such a linear relationship represents blood behavior at high shear rates, the nonNewtonian effect becomes appreciated at low shear rates because viscosity increases with decreasing shear rate.13 The inherent nonNewtonian characteristics of blood result from the formation of rouleaux or aggregates of red blood cells under low shear conditions.13
Although traditionally the Newtonian fluid assumption for blood has been adopted in CFD simulations of large vessels (such as arteries hosting aneurysms),14 ,15 the validity of such simplification has not been tested for calculations of aneurysmal hemodynamics in which the results are interpreted to relate to IA pathophysiology. For example, previous studies have shown that regions of low WSS in the aneurysm sac are associated with IA rupture8 ,9 ,11 as well as thrombus formation.16 In such regions the blood is relatively stagnant and is therefore subjected to lower shear rates than in the parent vessel. We suspected that, in these low WSS regions, the nonNewtonian effects might not be negligible.
The objective of this study was to test the sensitivity of blood flow field, shear rate and WSS predictions using different blood rheology models in different types of patientspecific aneurysm geometries. We aimed to determine whether modeling based on the Newtonian assumption can closely represent nonNewtonian models in typical IA geometries, particularly in complex geometries that may be associated with an increased risk of thrombosis or rupture.
Methods
Three cases of internal carotid artery (ICA) aneurysms (identified as aneurysms A, B, and C (figure 1) were selected from the study by Dhar et al 17 for the present study. Aneurysm A (ruptured, in a 45yearold woman) was a sidewall aneurysm with an irregular oblong shape and a daughter sac; aneurysm B (unruptured, in a 68yearold woman) was a sidewall aneurysm with a quasispherical shape; and aneurysm C (unruptured, in a 66yearold man) was a nearspherical bifurcation aneurysm. In the present study, sidewall aneurysms were defined as lesions emanating off the side of the parent vessel with or without a very small vascular branch, and bifurcation aneurysms were defined as those located at major bifurcations in the cerebral vessel.11 ,18 Threedimensional angiography images of the patients' aneurysms were obtained with a Toshiba Infinix VFi/BP frontal Carm system (Toshiba America Medical Systems, Tustin, California, USA). Threedimensional images of the ICA were then reconstructed in surfacetriangulation format using inhouse software based on the opensource Visualization Tool Kit libraries, as previously described.17
Two widely applied nonNewtonian fluid models for blood, the Casson and HerschelBulkley (HB) models, were used for the CFD simulations of the geometric ICA aneurysm models in this study, along with the usual Newtonian fluid model.16 The mathematical details of these models are given in the online appendix. Figure 2 shows viscosity versus shear rate for the three rheology models. All three models provide similar viscosity values at the high shear rates (>100/s) that are typical for large arteries, but exhibit very different behaviors at low shear rates. Finitevolume meshes consisting of approximately 300 000 to 1 million elements were created for the flow domain of each ICA aneurysm model using ANSYS ICEM CFD (ANSYS, Canonsburg, Pennsylvania, USA). The Navier–Stokes equations were solved numerically under pulsatile flow conditions using the CFD solver StarCD (CD Adapco, Melville, New York, USA). In all simulations a mean flow rate of 4.6 ml/s at the ICA19 was used as the inlet boundary condition, and a pulsatile velocity waveform measured from transcranial Doppler ultrasound images obtained from a normal subject was scaled to the mean flow rate. Tractionfree boundary conditions were implemented at the outlets. The mass flow rate through each outlet artery was proportional to the cube of its diameter based on the principle of optimal work.20 Detailed flowgoverning equations and numerical methods of CFD simulations are shown in the online appendix. For each aneurysm geometry we performed three pulsatile flow simulations using the Newtonian model, the Casson model and the HB model, respectively. In this study we were only concerned with the shear rate, viscosity and WSS averaged over one cardiac cycle. In this paper, mention of these parameters implies timeaveraged values over one cardiac cycle.
In our analysis of the CFD hemodynamic data we plotted the luminal distributions of shear rate, viscosity and WSS derived from each simulation. For quantitative comparison of different rheology models, we calculated the average shear rate and the average blood viscosity over the aneurysm dome volume and the average WSS over the aneurysm dome surface. Here the aneurysm sac was divided into dome, body and neck regions, each occupying onethird of the aneurysm height.21 Average values for the parent vessel parameters were calculated similarly from the reconstructed parent vessel segments, including the branches.
Results
Flow characterization
CFD simulations using all three rheology models produced very similar flow patterns in each ICA aneurysm. Figure 3 shows the timeaveraged flow patterns on representative crosssectional planes predicted by CFD using the Newtonian model. Although aneurysms B and C contained a typical large single vortex with inflow and outflow pathways, aneurysm A harbored an additional secondary vortex near the dome that was slowly recirculating or nearly stagnant. Timeresolved visualization of the pulsatile flow field (not shown) revealed that the primary vortices within all three aneurysms were stable, whereas the secondary vortex in aneurysm A was slowly oscillating. These results were similar to those obtained with the nonNewtonian models (not shown).
Shear rate
Figure 4 shows the luminal distributions of shear rate calculated using the three rheological models for each aneurysm geometry. Because the nonNewtonian effect becomes most appreciated at very low shear rates (figure 2), a logarithmic scale was used to accentuate this effect for comparison between the different blood rheology models (figure 4). In each case the Casson model and HB model predicted distributions that were similar to those associated with the Newtonian model, except for the dome of aneurysm A which featured a low shear rate zone. Figure 5 provides a quantitative comparison of shear rate values that were volumeaveraged in the aneurysm dome (top onethird of the sac) and in the parent vessel, calculated by all three rheology models.
In the dome of aneurysm A the Casson model and HB model predicted average shear rates that were only 53% and 43%, respectively, of the values obtained with the Newtonian model. Elsewhere (ie, in the rest of aneurysm A and in aneurysms B and C) the discrepancies were less than 3%. In other words, the Newtonian model overestimated the shear rate in the dome of aneurysm A, a predominantly low shear region.
Viscosity
Blood viscosity distributions at the luminal wall from the two nonNewtonian models generally showed no large variations and were similar to the Newtonian viscosity (figure 6). The only exception was the dome of aneurysm A where substantially higher viscosity values were observed when using the Casson or HB model. Here the average viscosities predicted by the Casson and HB models were 174% and 274% of the Newtonian model, respectively (figure 5)—that is, the Newtonian model did not reflect the drastically increased viscosity in the stagnant region in the dome of aneurysm A.
Wall shear stress
Figure 7 shows the luminal distributions of WSS predicted by all three rheology models. For ease of comparison, the WSS is normalized by the local WSS value from the Newtonian model prediction. Because of this normalization, the Newtonian results (left column) appear to have unity WSS distributions (ie, all values are 1). In each aneurysm case the Newtonian and Casson models predicted similar WSS, except in the dome of aneurysm A where the Casson model prediction was considerably lower—as low as 55% of the Newtonian prediction (observation from figure 7) at the tip of the dome. The HB model (right column) predicted slightly lower WSS (approximately 5% less) than the Newtonian model, except in the dome of aneurysm A where the HB model predicted as low as 60% of the Newtonian model (observation from figure 7) at the tip of the dome.
Figure 5 provides quantitative comparisons of the surfaceaveraged WSS in the aneurysm dome region and in the parent vessels predicted by the three models. The average WSS predicted by the Casson and HB models was 71% and 76%, respectively, of the Newtonian WSS, showing that the Newtonian model overestimated WSS in the dome of aneurysm A.
Discussion
Previous studies have investigated the influence of nonNewtonian blood rheology on CFD modeling in idealized geometric models of aneurysms. Those models usually consist of a sphere on a cylindrical tube (straight or curved) which possesses perfect symmetry that does not exist in real aneurysms. Fisher et al 22 virtually created different types of idealized aneurysms, including bifurcation and sidewall aneurysms on straight and curved vessels. Valencia et al 23 studied the nonNewtonian effect on two virtual saccular aneurysm models with different inclination angles and one model of a nondiseased basilar artery. Such perfect symmetry is unrealistic and may cause misleading flow patterns because realistic anatomical aneurysm geometry can have complex flow dynamics that idealized geometries cannot capture. Such complex flow dynamics could be caused by vessel curvature or tortuosity, which influence the flow entering and circulating within the aneurysm sac and thus shear stress on the aneurysm wall.15 Furthermore, patients' aneurysms are often nonspherical or have multiple lobes and, because of these irregular shapes, can harbor various disturbed flow patterns including secondary and/or unsteady vortices.17 ,24–26 Blood rheology results based on idealized aneurysm models do not therefore reflect many important features of real anatomy.23 In the current study, geometric models of typical realistic ICA aneurysms representing an oblong sidewall aneurysm, a nearspherical sidewall aneurysm and a nearspherical bifurcation aneurysm were used to investigate the influence of nonNewtonian blood rheology on CFD simulations.
Our results indicate that the Newtonian fluid assumption in general is acceptable in CFD modeling of healthy vessels and aneurysms that do not harbor regions of pronounced low shear. However, the current data suggest that the nonNewtonian properties of blood should be considered when modeling hemodynamics in aneurysms with slow and recirculating flow regions (low shear flow), which are typically found in aneurysms with complex geometry like aneurysm A with an irregular oblong shape and a daughter aneurysm.
These types of aneurysms are generally suspected to be dangerous and command attention for treatment, as low WSS has been shown to be associated with aneurysm rupture8 ,9 ,11 and low shear rates may cause clot formation.16 Xiang et al 11 demonstrated that the rupture of IAs correlated with low WSS in their study of 119 aneurysms. Ishida et al 27 identified the rupture points of eight ruptured middle cerebral artery aneurysms to be located in regions of the lowest WSS. The Newtonian model cannot capture the increased viscosity in such low shear regions and consequently overestimates the shear rate and WSS. When the nonNewtonian effects are neglected in such aneurysms, the shear rate and WSS appear falsely higher in the low shear regions and thus could underrepresent the risk of clot formation and aneurysm rupture.
In clinical practice, coil embolization or flow diverter treatment for cerebral aneurysms reduces the intraaneurysmal flow, as indicated by contrast stagnancy at the dome, thereby reducing shear rate (thus promoting clotting) and WSS (thus increasing rupture risk).24 ,28 If intraaneurysm thrombosis does not develop quickly and completely to fill the whole sac but leaves a part of the aneurysmal wall exposed to low WSS for an extended period of time, such treatments could result in deterioration of the wall via inflammatory cell infiltration (caused by low WSS and immature thrombus29) and thus an increased propensity of rupture.24 However, we believe that in most situations clotting happens much faster than the wall degradation. Thus, coil or flow diverter embolization, by massively increasing stasis, usually induces a healing response through thrombosis and subsequent cicatrization. It is important to note, however, that there have been reports of aneurysm rupture after coil or flow diverter treatments.29 The typical aneurysmal flow after coil or flow diverter implantation has high stasis with a very low shear rate and WSS. Hence, the same precaution regarding the choice of blood rheology model should be taken in hemodynamic modeling of posttreatment aneurysms.
Despite our concerns regarding the application of the Newtonian assumption in CFD of IAs, we do not suggest that all CFD analyses be conducted with nonNewtonian rheology. The Newtonian model has its merits: it is simple and easy to use, implicit in all commercial CFD software and accurate enough in most situations unless low shear regions are present. Furthermore, our data show that Newtonian models produce flow patterns roughly consistent with nonNewtonian models. Hence, if the researcher only wants to know the intraaneurysmal flow patterns but not WSS or shear rate, routine CFD with Newtonian assumption might be sufficient.
NonNewtonian models have their drawbacks. Implementation of nonNewtonian models in CFD requires writing subroutines or using specialized CFD software, which might not be available to all researchers. Moreover, nonNewtonian models sometimes may underperform. For example, in the HB model, the viscosity has no lower limit when shear rate tends to infinity. Hence, at very high shear rates where blood is expected to behave like a Newtonian fluid with a constant viscosity, the HB model would produce unrealistically lower viscosity values. This behavior causes viscosity values in the parent arteries in our study to fall below the Newtonian value (figures 2 and 6) and may unduly influence the flow field. In our study, at a shear rate of 10 000/s, the corresponding viscosity in the HB model is 2.5 cP compared with the Newtonian viscosity of 3.5 cP.
As a simplified heuristic to aid decisionmaking for future aneurysm hemodynamics modeling projects, we suggest using nonNewtonian rheology modeling in the following situations where regions of high stasis or slow recirculation are expected, especially when quantitative information of shear rate, WSS or WSSbased quantities such as oscillatory shear index are desired:

(1)Aneurysms with size ratio >2. Size ratio is defined as the maximal diameter of the aneurysm divided by the parent vessel diameter.17 Tremmel et al 25 demonstrated that, once the size ratio becomes >2, the single aneurysmal flow vortex splits into multiple vortices and the low WSS area increases drastically. Dhar et al 17 found that size ratio can significantly separate ruptured from unruptured aneurysms, and this finding has been supported by other studies.11 ,30 ,31 Furthermore, Dhar et al 17 found the optimal threshold for size ratio to be 2.05.

(2) Aneurysms with daughter aneurysms or blebs.

(3) Aneurysms after coil or flow diverter implantation.
When researchers are unsure whether a particular aneurysm case commands nonNewtonian modeling, we suggest that they first run a simple and quick steadystate CFD simulation using the Newtonian model to determine whether a low shear flow region exists. If it does, a nonNewtonian model could be used for proper CFD simulation to obtain WSS and related quantities such as oscillatory shear index.
This study has several limitations. First, only three aneurysm cases were investigated. Although this can provide an initial illustration of the effect of nonNewtonian rheology on aneurysmal hemodynamics, conclusions should be validated in followup studies that include a large number of cases. Second, because we focused on comparing hemodynamics from different rheology models, we kept the inlet boundary condition for the CFD simulations the same, which may not accurately reflect the patientspecific reality. Third, although nonNewtonian rheology can improve the accuracy of blood flow simulations, the nonNewtonian models are themselves also simplifications and cannot perfectly mimic in vivo behavior. Finally, there is no consensus as to which nonNewtonian model (Casson, HB or others) is superior. Further investigations are therefore needed that are dedicated to produce a more consistent nonNewtonian model, and care should be taken to accurately determine its model parameters experimentally.
Conclusion
The Newtonian blood rheology assumption is usually acceptable for CFD simulation of cerebral aneurysm hemodynamics but could underestimate viscosity and overestimate shear rate and WSS in the slowly recirculating flow regions that are typically found at the dome in elongated or complexshaped saccular aneurysms, as well as in aneurysms following endovascular treatment. Because low shear rates and low WSS in such flow conditions indicate a high propensity for thrombus formation and aneurysm rupture, Newtonian hemodynamics may underestimate the propensity of these events.
Acknowledgments
We thank Paul H. Dressel BFA for assistance with preparation of the illustrations and Debra J Zimmer AAS CMAA for editorial assistance.
References
Supplementary materials
Supplementary Data
This web only file has been produced by the BMJ Publishing Group from an electronic file supplied by the author(s) and has not been edited for content.
Files in this Data Supplement:
 Download Supplementary Data (PDF)  Manuscript file of format pdf
Footnotes

Additional data are published online only. To view the file please visit the journal online (http://jnis.bmj.com).

Disclosures Dr. Kolega, Dr. Natarajan, Dr. Tremmel, and Mr. Xiang have no final relationships to disclose. Dr. Meng is the principal investigator of the aforementioned NIH grant. Dr. Levy receives research grant support (principal investigator: StentAssisted Recanalization in acute Ischemic Stroke, SARIS), other research support (devices), and honoraria from Boston Scientific* and research support from Codman & Shurtleff, Inc. and ev3/Covidien Vascular Therapies; has ownership interests in Intratech Medical Ltd. and Mynx/Access Closure; serves as a consultant on the board of Scientific Advisors to Codman & Shurtleff, Inc.; serves as a consultant per project and/or per hour for Codman & Shurtleff, Inc., ev3/Covidien Vascular Therapies, and TheraSyn Sensors, Inc.; and receives fees for carotid stent training from Abbott Vascular and ev3/Covidien Vascular Therapies. Dr. Levy receives no consulting salary arrangements. All consulting is per project and/or per hour. (*Boston Scientific's neurovascular business has been acquired by Stryker.)

Funding This work was partially supported by NIH grant R01NS064592 and a grant from Toshiba Medical Systems.

Competing interests None.

Provenance and peer review Not commissioned; externally peer reviewed.