Identification of vortex structures in a cohort of 204 intracranial aneurysms

Department of Mechanical and Aerospace Engineering, University at Buffalo, State University of New York, Buffalo, NY, USAToshiba Stroke and Vascular Research Center, University at Buffalo, State University of New York, Buffalo, NY, USA

Department of Mechanical and Aerospace Engineering, University at Buffalo, State University of New York, Buffalo, NY, USAToshiba Stroke and Vascular Research Center, University at Buffalo, State University of New York, Buffalo, NY, USA

Toshiba Stroke and Vascular Research Center, University at Buffalo, State University of New York, Buffalo, NY, USADepartment of Neurosurgery, University at Buffalo, State University of New York, Buffalo, NY, USA

Toshiba Stroke and Vascular Research Center, University at Buffalo, State University of New York, Buffalo, NY, USADepartment of Neurosurgery, University at Buffalo, State University of New York, Buffalo, NY, USADepartment of Radiology, University at Buffalo, State University of New York, Buffalo, NY, USA

Department of Mechanical and Aerospace Engineering, University at Buffalo, State University of New York, Buffalo, NY, USAToshiba Stroke and Vascular Research Center, University at Buffalo, State University of New York, Buffalo, NY, USADepartment of Neurosurgery, University at Buffalo, State University of New York, Buffalo, NY, USADepartment of Biomedical Engineering, University at Buffalo, State University of New York, Buffalo, NY, USA

Abstract

An intracranial aneurysm (IA) is a cerebrovascular pathology that can lead to death or disability if ruptured. Abnormal wall shear stress (WSS) has been associated with IA growth and rupture, but little is known about the underlying flow physics related to rupture-prone IAs. Previous studies, based on analysis of a few aneurysms or partial views of three-dimensional vortex structures, suggest that rupture is associated with complex vortical flow inside IAs. To further elucidate the relevance of vortical flow in aneurysm pathophysiology, we studied 204 patient IAs (56 ruptured and 148 unruptured). Using objective quantities to identify three-dimensional vortex structures, we investigated the characteristics associated with aneurysm rupture and if these features correlate with previously proposed WSS and morphological characteristics indicative of IA rupture. Based on the Q-criterion definition of a vortex, we quantified the degree of the aneurysmal region occupied by vortex structures using the volume vortex fraction (vVF) and the surface vortex fraction (sVF). Computational fluid dynamics simulations showed that the sVF, but not the vVF, discriminated ruptured from unruptured aneurysms. Furthermore, we found that the near-wall vortex structures co-localized with regions of inflow jet breakdown, and significantly correlated to previously proposed haemodynamic and morphologic characteristics of ruptured IAs.

1. Introduction

Intracranial aneurysms (IAs) are a cerebrovascular pathology present in 3–5% of the population [1,2]. Aneurysm rupture leads to subarachnoid haemorrhage, a devastating event that is associated with high rates of mortality and morbidity [3]. However, because the rupture rate is low [2] and the risk of surgical complications can be significant, the decision of whether or not to treat an unruptured intracranial aneurysm (IA) is difficult. Currently, IA size is the most commonly used metric for assessing rupture risk [4]; however, a high frequency of small aneurysms are known to rupture [5,6]. Thus, objective guidelines to determine the relative risk of rupture are warranted.

Researchers have investigated aneurysmal characteristics associated with IA pathophysiology and found that intra-aneurysmal haemodynamics play a significant role in IA remodelling, growth and rupture [7–13]. Statistical studies using image-based computational fluid dynamics (CFD) have found that time-averaged fluid shear stresses on the aneurysm wall, as well as certain morphological features, are significantly associated with IA rupture status [8,10,11,14]. However, analysing the time-averaged WSS quantities, or the aneurysm morphology itself, cannot elucidate the underlying flow physics. Therefore, there is a need to investigate the flow physics in complex aneurysm geometries, especially as it relates to the fate of the IAs.

Multiple studies in large patient IA databases have qualitatively observed two-dimensional velocity fields and reported that ruptured aneurysms tend to have more complex flows with multiple vortices when compared to unruptured IAs [8,15–17]. Additional studies based on large patient cohorts have found that vortex core lines in ruptured aneurysms exhibit more temporal and spatial fluctuations [18,19]. Because vortex core lines only show the skeleton of the entire vortex structure, and identification of two-dimensional planes can be subjective, the true association of vortical flow and IA rupture status remains unclear.

To understand the complex three-dimensional flow motions that occur in sophisticated patient aneurysm geometries, identification of dynamically significant flow structures helps to classify and organize the entire flow domain. Vortex dynamics govern the evolution of coherent structures, or spatially distinct large-scale flow features [20,21]; therefore, the dynamics of vortex structures offer a key vehicle for understanding the physics of aneurysmal flows. A well-established definition of vortex structures is the Q-criterion, which quantifies the difference between the rotational and straining parts of the velocity gradient tensor [22]. Since it is straightforward to calculate, it has been frequently used for visualizing and quantifying vortical coherent structures. Recently, the Q-criterion has been used to examine the formation and evolution of vortex structures in a rabbit aneurysm model [23,24] and a few patient IA models, to demonstrate the feasibility and potential value of vortex visualization [25,26].

To further explore the relevance of vortical flow to aneurysm pathophysiology, we adopted the vortex definition based on the Q-criterion and investigated aneurysmal flow dynamics in a database of 204 patient aneurysms of diverse morphologies. Through the objective analysis of three-dimensional vortex structures, we extracted the vortex characteristics statistically relevant to aneurysm rupture, and investigated if these features were correlated to previously proposed haemodynamic and morphologic predictors of IA rupture (WSS, oscillatory shear index and size ratio (SR)).

2. Material and methods

2.1. Patient population

Under Institutional Review Board Approval, 204 IA cases (56 ruptured and 148 unruptured) were consecutively collected from 186 patients between 2006 and 2011 from Buffalo General Medical Center and Gates Vascular Institute. Cases were excluded if they did not have sufficient image quality for reconstruction. Three-dimensional angiographic images were reconstructed from digital subtraction angiography or, in the case of anterior communicating artery aneurysms, from computed tomography angiography [8]. Rupture status was determined at the time of imaging. Table 1 summarizes the location, age, gender and rupture status of the analysed IAs. This cohort is a part of a previously collected cohort analysed in past studies [8,27].

2.2. The Q-criterion definition of a vortex

To describe the flow physics in a large cohort of patient aneurysms, we examined the vortex dynamics using the Q-criterion. The Q-criterion identifies vortex regions where rotational flow dominates over straining flow, mathematically defined as when the Frobenius matrix norm of the antisymmetric vorticity tensor () is greater than the symmetric strain tensor ():2.1

Here the anti-symmetric part, , and the symmetric part, , are described by the gradient of the velocity field, , and its transpose, respectively:2.2and2.3

Traditionally, in order to visualize vortex structures using Q, arbitrary thresholds have been assigned for individual cases to best highlight the interesting structures shown by iso-surfaces. However, in order to compare structures at different thresholds across a large number of patient cases, an objective and uniform definition of Q was required. As shown in equation (2.1), when Q > 0, the magnitude of the rotational tensor dominates. Therefore, as with a previous study [25], we define a vortex region simply as Q > 0, where positive Q encompass the rotationally dominated flow domain.

2.3. Measures of prevalence of vortical flow in aneurysms

To examine the relationship between aneurysmal flow dynamics and the pathophysiological status of an aneurysm (ruptured or unruptured), we propose two parameters to measure the prevalence of vortical flow in each aneurysm: volume vortex fraction (vVF), which is the relative volume of the aneurysm with > 0, and surface vortex fraction (sVF), which is the relative area near the wall of the aneurysm with > 0.

To calculate vVF, the aneurysm was first isolated from the parent vessel and the aneurysmal volume was calculated. Then regions with positive Q were defined and their volume was calculated. To characterize the degree of vortex dominance in the aneurysm, the vortex volume is normalized by the total aneurysm volume. As shown in equation (2.4), we define vVF as the volume of the aneurysm with positive Q divided by the total aneurysm volume2.4

Preliminary qualitative analysis of the Q-criterion in ruptured and unruptured aneurysms suggested that there may be a relationship between vortex structures near the wall and the aneurysm rupture status [28]. Consequently, we also examined the degree of vortex flow dominance near the aneurysm surface. It has been mathematically established that Q at the wall must equal 0, as the magnitude of the rotational and strain tensors are equal [21]. Therefore, after the calculation of , the layer of elements adjacent to the wall were extracted. In contrast to the values of the nodes directly at the wall, the cell-centred variable values in the surface elements were used, and, therefore, velocity, velocity gradient and, consequently, Q had non-zero values. The total surface area and the surface layer area with positive Q in the aneurysm were calculated. We defined the degree of vortex flow near the aneurysm wall as the ratio of the aneurysm with positive Q near the wall normalized by the total aneurysm surface area. We quantified sVF as described in the following equation:2.5

2.4. Computational fluid dynamics methods

CFD simulations for all patient IA models were performed as described in detail previously [8]. In brief, three-dimensional angiographic images were segmented and surface meshes were generated using the open-source software Vascular Modeling Toolkit (www.vmtk.org) [29]. Unstructured volumetric meshes were generated ANSYS ICEM CFD (ANSYS Inc, Canonsburg, PA). Four additional refined prism layers were generated at the wall with a thickness of one-tenth the length of the maximum edge length and the aspect ratio was set to 1.2. To calculate sVF, we isolated the outermost wall elements and term it the surface layer.

Using laminar flow conditions, the incompressible transient Navier–Stokes equations were solved in STAR-CD (CD-adapco, Melville, NY). Because patient-specific waveforms are not routinely collected, a pulsatile waveform from a healthy volunteer was applied at the inlet and the velocity was scaled based on the inlet vessel diameter [8]. A pressure implicit with splitting of operator algorithm was used for temporal discretization. The flow split at the outlets was calculated based on the assumption of minimal work [30]. Rigid walls, a no-slip boundary condition and a density of 1056 kg m−3 were also assumed. Three cardiac cycles were run for each simulation to allow for sufficient stabilization of the solution, and the last cardiac cycle was output and analysed for all simulations. The solution was calculated at a time step of 1 millisecond and was extracted every 20 ms. Post-processing was done in Tecplot 360 (Tecplot, Inc. Bellevue, WA).

2.5. Sensitivity analysis

To examine the sensitivity of vVF and sVF to discretization, viscosity model assumptions and the cardiac cycle, we performed a series of sensitivity analyses. We chose four representative patient cases including two ruptured and two unruptured cases. The IA sizes, as previously defined [14], ranged from 6 to 10 mm. The cases were chosen because of their resemblance to the ruptured and unruptured group-averaged vVF and sVF.

First, to ensure that the fraction of vortex structures in the volume and near the surface were independent of the mesh density, we conducted a grid independence study. We refined the wall element thickness from 200 to 20 µm and monitored the change in vVF and sVF. Our criteria for convergence were a change in vVF and sVF by less than 0.01.

After grid convergence, we examined the effect of spatial discretization scheme order on the fraction of vortex structures in the aneurysms. We used a first- and second-order upwind-differencing scheme and quantified the change in vVF and sVF.

Then, to quantify the effect of a Newtonian viscosity assumption, we monitored the change in the fraction of vortex structures when a non-Newtonian viscosity model was used. We applied the Carreau–Yasuda model [31] as described in the following equation:2.6where the viscosity, μ, was dependent on the shear rate, . Similar to previous studies, we used parameters that were experimentally derived [32–34]. The zero shear rate viscosity µ0 was 0.16 Pa · s, the infinite shear rate viscosity μ∞ was 3.5 cP, the relaxation factor λ was 8.2 s, the power index n was 0.2128 and the α constant was 0.64. We compared vVF and sVF computed from a constant viscosity of 3.5 cP to vVF and sVF computed from the Carreau–Yasuda model.

To examine the prevalence of vortex structures throughout the cardiac cycle, we monitored vVF and sVF in four representative patient cases including two ruptured and two unruptured cases. In addition, to observe the temporal changes of the aneurysmal flow field, we visualized velocity streamlines and calculated vortex core line lengths throughout one cardiac cycle. Vortex core line lengths were based on the region where the second eigenvalues of Q were negative. The centres of these regions were considered the location of the core lines.

2.6. Statistical analysis

We tested our cohort of 204 aneurysms to determine if vVF, sVF, volume, surface area, mean Q or maximum Q were able to differentiate ruptured from unruptured aneurysms. Kolmogorov–Smirnov tests were performed to test if the variables were normally distributed. The aneurysms were then separated into ruptured and unruptured groups and a Mann–Whitney U-test (for non-normally distributed data) or a Student's t-test (for normally distributed data) were used to determine if there was any statistical difference between the groups. A p-value of less than 0.01 was considered statistically significant.

We also analysed the correlation of the relevant flow physics to haemodynamic forces that influence aneurysm pathophysiology. Previous studies have found an association between time-averaged WSS [8,11,35,36], OSI [8,37] and the morphologic parameter aneurysm SR (the ratio of aneurysm size to parent vessel) [8,14,38–40] to aneurysm rupture. Therefore, we tested the correlation of the flow physics significantly associated with rupture to WSS, OSI and SR using Pearson correlation analyses. All statistical analyses were performed in IBM SPSS Statistics v. 22.0 (IBM Corporation, Armonk, NY).

3. Results

3.1. Patient cohort and classification of flow modes in aneurysms

Using image-based CFD, we simulated time-resolved three-dimensional flow fields in each of the 204 aneurysms in our patient cohort. From the simulation results, we plotted velocity streamlines to visualize the dynamics of flow patterns, and out-of-plane vorticity on a centre plane in the aneurysm as a simplified means to visualize vortices. Based on our CFD results, we observed that the aneurysms could be classified into two flow modes: (i) Continuous Jet Mode, characterized by a single large organized vortex that is formed by flow entering into the aneurysm distally and circulating along the wall, as illustrated in figure 1a, and (ii) Jet Breakdown Mode, characterized by the inflow jet entering into the aneurysm, impinging on the wall and subsequently breaking down into counter-rotating vortices, as illustrated in figure 1b. In our aneurysm cohort, the flow in 74% of all aneurysms was characterized as that of the Continuous Jet Mode and the remaining 26% as that of the Jet Breakdown Mode.

An illustration of the classification of aneurysm flow modes based on the inflow jet. Aneurysms were classified to have either (a) Continuous Jet Mode, or (b) Jet Breakdown Mode. Both flow modes can exist in sidewall or bifurcation aneurysms.

Figures 2 and 3 show examples of the Continuous Jet Mode in a sidewall aneurysm and a bifurcation aneurysm, respectively. In both cases, as shown by streamlines in figures 2 and 3a and vorticity plots in figures 2 and 3b, flow entered into the aneurysm dome, circulated along the wall and remained organized. A single vortex was seen at five time points throughout the cardiac cycle (t/T = 0, 0.04, 0.08, 0.15 and 0.50). Additionally, as shown by the lack of positive Q contours near the aneurysm surface in figures 2 and 3c, this flow mode was first qualitatively associated with a low fraction of vortex structures next to the wall (low sVF).

A representative sidewall aneurysm case with the Continuous Jet Mode at five time points throughout the cardiac cycle (t/T = 0, 0.04, 0.08, 0.15 and 0.50). (a) Velocity streamlines and (b) out-of-plane vorticity shows a single organized vortex entering distally. This flow resulted in a low fraction of (c) vortex structures near the surface (coloured by positive Q contours). Both surface vortex structures and (d) vortex structures in the three-dimensional volume (coloured by positive Q contours) do not change throughout the cardiac cycle.

A representative bifurcation aneurysm case with the Continuous Jet Mode at five time points throughout the cardiac cycle (t/T = 0, 0.04, 0.08, 0.15 and 0.50). The single organized vortex is shown by (a) velocity streamlines and (b) out-of-plane vorticity. (c) This flow results in a low fraction of vortex structures near the surface (coloured by positive Q contours). Both vortex structures near the surface and (d) vortex structures in the three-dimensional volume (coloured by positive Q contours) do not change throughout the cardiac cycle.

Figures 4 and 5 give examples of the Jet Breakdown Mode in a sidewall aneurysm and a bifurcation aneurysm, respectively. In both cases, as shown in figures 4 and 5a,b, the more complex flow is shown by a primary jet that entered into the aneurysm, impinged on the distal aneurysm wall and resulted in the development of two counter-rotating vortices. Additionally, coloured contours in figures 4 and 5c show positive Q near the wall corresponding to the site of impingement and jet breakdown, suggesting that regions near the wall where rotational flow is dominant is more prevalent in this flow mode.

A representative sidewall aneurysm case with the Jet Breakdown Mode at five time points throughout the cardiac cycle (t/T = 0, 0.04, 0.08, 0.15 and 0.50). (a) Velocity streamlines with a solid black arrow showing the main inflow jet and dashed arrows showing the counter-rotating vortices, and (b) out-of-plane vorticity. (c) Vortex structures near the surface (coloured by positive Q contours) co-localize with the impingement jet site. Vortex structures near the surface and (d) vortex structures in the three-dimensional volume (coloured by positive Q contours) do not change throughout the cardiac cycle.

A representative bifurcation aneurysm case with the Jet Breakdown Mode at five time points throughout the cardiac cycle (t/T = 0, 0.04, 0.08, 0.15 and 0.50). (a) Velocity streamlines with a solid black arrow showing the main inflow jet and dashed arrows showing the counter-rotating vortices, and (b) out-of-plane vorticity. This impingement co-localizes with (c) vortex structures near the surface (coloured by positive Q contours). Surface vortex structures and (d) vortex structures in the three-dimensional volume (coloured by positive Q contours) do not change throughout the cardiac cycle.

Based on qualitative observations of positive Q contours in both flow modes, we suspected that aneurysms with the Jet Breakdown Mode had higher sVF. Therefore, we compared vVF and sVF between the two flow mode groups in our aneurysm cohort. We found no statistical difference in the group-averaged vVF, which was 0.30 ± 0.12 and 0.28 ± 0.10 for aneurysms with the Continuous Jet Mode and the Jet Breakdown Mode, respectively (p = 0.089). However, we found a significant difference in sVF between the two groups, with a group-averaged sVF of 0.09 ± 0.05 and 0.12 ± 0.05 for Continuous Jet Mode and Jet Breakdown Mode aneurysms, respectively (p < 0.001).

3.2. Sensitivity analysis

To determine the sensitivity of the solution to discretization, viscosity model assumptions and the cardiac cycle, we performed sensitivity analyses on four representative cases, including two ruptured (R1 and R2) and two unruptured (UR1 and UR2) cases.

As shown in figure 6a, we first performed a grid independence study where the wall element thickness was refined from 200 to 20 µm and vVF and sVF were monitored. We found that when the maximum wall element size was refined from 25 to 20 µm, both vVF and sVF changed by less than 0.01 in all four cases (less than 1% change for both vVF and sVF). Therefore, we determined that, for a 6–10 mm IA, a wall element thickness of 25 µm was sufficient. In our cohort of 204 IAs, where IA size ranged from 2 to 25 mm, this resulted in computational domains with 300 000–1 000 000 elements. Additionally, each case had surface layer thickness between 10 and 50 µm. Therefore, the location of the cell-centred values on the surface layer, where sVF was calculated, was between 5 and 25 µm from the wall face.

The results of the sensitivity analysis performed on two ruptured (R1 and R2) and two unruptured (UR1 and UR2) patient IA cases. The metrics vVF (left) and sVF (right) were monitored. (a) The results of a grid independence study where the wall element thickness was refined. An average change of less than 0.01 was observed at a maximum wall element size of 25 µm (dashed line). (b) The results of a spatial discretization sensitivity study that quantified the effect of first- and second-order discretization schemes. Between first- and second-order discretization schemes an average change of 0.05 and less than 0.01 were observed for the vVF and the sVF, respectively. (c) The results of a non-Newtonian viscosity model. An average change of less than 0.01 was observed when the Carreau–Yasuda (C-Y) model was used. (d) vVF and sVF are shown during one cardiac cycle. Peak systole is denoted with the dashed line. Each case showed minimal changes in regions of positive Q. The average standard deviation of the vVF and the sVF was 0.02 and 0.01, respectively.

Next, as shown in figure 6b, we monitored the effect of a first- and second-order spatial discretization scheme on vVF and sVF. When a second-order discretization scheme was used, we found that vVF increased by an average of 0.05 (14% change) and sVF decreased on average by less than 0.01 (−3% change) when compared to a first-order discretization scheme. Therefore, we justified the use of a first-order discretization scheme for this study.

Then, as shown in figure 6c, we quantified the effect of a Newtonian viscosity assumption on vVF and sVF by quantifying the change when The Carreau–Yasuda model was used. We found that vVF changed by an average of 0.03 (5% change) and sVF changed by an average of less than 0.01 (0% change) when non-Newtonian viscosity was used compared to the constant viscosity model. Therefore, we used a constant viscosity of 3.5 cP in our simulations.

Last, we monitored the prevalence of vortex structures throughout the cardiac cycle. First, we examined the change in vortex core lines, vVF and sVF. Consistent with previous computational studies [18], we found that core line lengths are dynamic and change between 33 and 56% throughout the cardiac cycle (see electronic supplementary material, figure S1). Additionally, we confirmed that maximum Q occurs at peak systole (see electronic supplementary material, figure S2). However, as shown in figures 2⇑⇑–5c,d, we found that the three-dimensional vortex structures (Q > 0) both in the volume and at the wall appeared to remain unchanged despite the changing inflow. When vVF and sVF were quantified at each time point throughout the cardiac cycle, as shown in figure 6d, the average standard deviation for vVF and sVF was 0.02 and 0.01, respectively (table 2). This shows that the vortex structures (regions of positive Q) present in the aneurysms was relatively stable throughout the cardiac cycle. Therefore, we determined it was sufficient to examine vVF and sVF at a single time point and chose to examine peak systole to capture the highest maximum value of Q. As figure 6d shows, although maximum vVF and sVF do not always occur at peak systole, for consistency between all 204 cases, and to limit the computational expense of the calculation of velocity gradients at all time points, we determined that a single time point was sufficient.

Temporal analysis of the vVF and the sVF over one cardiac cycle for four representative cases, including two ruptured (R) and two unruptured (UR). For the cardiac cycle, the average standard deviation of the vVF and the sVF was 0.02 and 0.01, respectively.

3.3. Differences in vortex presence in ruptured and unruptured aneurysms

The statistical comparison of aneurysm volume, surface area, mean Q, maximum Q, vVF and sVF between ruptured and unruptured IAs is summarized in table 3. There was no statistical difference in volume, surface area, mean Q, maximum Q and vVF between ruptured and unruptured aneurysm groups. However, we did find a significant difference between the group-averaged sVF of ruptured and unruptured IAs. Figure 7 shows coloured contours representing (a) vVF and (b) sVF for six unruptured (top) and six ruptured (bottom) aneurysm cases. While there was no difference in relative aneurysmal volume dominated by the vortex between ruptured and unruptured IAs (vVF, ruptured: 0.29 ± 0.11, unruptured: 0.30 ± 0.12, p = 0.32), the relative area near the aneurysm wall that was dominated by vortex flow was significantly higher in ruptured compared to unruptured IAs (sVF, ruptured 0.13 ± 0.06, unruptured 0.09 ± 0.04, p < 0.001). As shown in figure 8, ruptured aneurysms had nearly 1.5 times the relative vortex structures near the surface (positive Q) as compared to unruptured aneurysms.

A boxplot showing the results of the comparison of the sVF and the vVF for the 204 aneurysms by a Mann–Whitney U-test. (a) The 148 unruptured (UR) and 56 ruptured (R) cases showed no statistical difference in the vVF, with a group average of 0.29 and 0.30, respectively (p = 0.32). (b) sVF was found to be statistically different between ruptured and unruptured aneurysms, with a group average of 0.13 and 0.09, respectively (p < 0.001). n.s., not significant; *S, statistically significant (defined as p < 0.01).

Summary of group statistical comparison of ruptured and unruptured aneurysms including the group-averaged means and standard deviations (s.d.s). Statistical significance in probability tests indicated by asterisks.

To further shed light on the flow dynamics that contributes to higher sVF and its potential association with ruptured IAs, we analysed the differences in sVF between ruptured and unruptured aneurysms with the Continuous Jet Mode and the Jet Breakdown Mode. As shown in figure 9, we found that among aneurysms with the Jet Breakdown Mode (n = 52), ruptured aneurysms had a significantly higher sVF compared to unruptured aneurysms (Jet Breakdown Mode sVF, unruptured 0.11 ± 0.04, ruptured 0.14 ± 0.06, p = 0.008). Likewise, among the aneurysms with the Continuous Jet Mode (n = 152), ruptured aneurysms had a significantly higher sVF compared to unruptured aneurysms (Continuous Jet Mode sVF, unruptured 0.08 ± 0.04, ruptured 0.12 ± 0.06, p < 0.001). Not only did aneurysms with the Jet Breakdown Mode show higher sVF compared to aneurysms with the Continuous Jet Mode, but also ruptured aneurysms showed the overall highest values of sVF.

A boxplot showing the results of the comparison of the sVF between unruptured and ruptured aneurysms with the Continuous Jet Mode and the Jet Breakdown Mode. Ruptured aneurysms had the highest sVF, which was significantly higher than that of unruptured aneurysms with the same flow mode (unruptured: 0.11 ± 0.04, ruptured: 0.14 ± 0.06, p = 0.008). *S, statistically significant (defined as p < 0.01).

Furthermore, we examined the relationship of the relevant flow physics to the haemodynamic stresses at the luminal surface that influences aneurysm pathophysiology. Figure 10 shows visualization of velocity streamlines, vortex structures near the surface (Q > 0), WSS contours and OSI contours in two representative unruptured (UR3, UR4) and two representative ruptured cases (R3, R4). As qualitatively shown in cases R3 and R4, areas of positive Q co-localized with areas of high WSS and were also associated with the region of jet breakdown at the wall. Lastly, we tested the correlation of significant aneurysm-averaged rupture-discriminating parameters found in previous studies [8,27] with sVF. As shown in figure 11, a significant, albeit weak correlation, was found between sVF and low WSS (p < 0.001, r = −0.42) and higher SR (p < 0.001, r = 0.48). However, no correlation was found between sVF and OSI (p = 0.019, r = 0.16).

4. Discussion

Previous statistical studies have reported a link between vortical flow and IA rupture status [8,18,23–25]. However, these studies were based on either two-dimensional or partial characterizations of vortical flow [8,18] or a few IA cases [23–25]. To ascertain the relevance of vortical flow to aneurysm rupture, we analysed a large database of 204 IAs from consecutive patients who underwent three-dimensional angiographic imaging at our centre. Using patient-specific geometries, we quantified three-dimensional vortex structures based on the Q-criterion in all the aneurysms and determined the characteristics most relevant to rupture status. We developed two new concepts quantifying the degree of dominance of vortex structures in aneurysms, namely vVF and sVF.

In a previous numerical experiment, we found that increase in the SR resulted in an increased number of vortices, where SR has been found to be an independently significant predictor of rupture [41]. In another study, we found a statistically significant association between a greater number of vortices visualized on a two-dimensional plane and the rupture status of aneurysms [8]. However, because of the complexity of aneurysm geometries and the intricacies of the resulting flow field, we believed a more complete description of vortex flow in aneurysms was required. Therefore, based on three-dimensional flow fields and vortex structures, we defined the volume vortex fraction, vVF, and the surface vortex fraction, sVF, to quantify the degree of the aneurysmal volume and surface occupied by three-dimensional vortex structures. From previous findings, our initial hypothesis was that the vVF would be a parameter that could discriminate ruptured from unruptured IAs. However, in our 204 IAs we found no statistical difference in the vVF between ruptured and unruptured aneurysms. Rather, the sVF, or positive Q near the wall, was significantly higher in ruptured aneurysms than in unruptured aneurysms.

It is not entirely surprising that ruptured aneurysms have a higher sVF compared with unruptured aneurysms. Scientific investigations of the pathological response of cells have focused, almost exclusively, on flow parameters in direct contact with the cells. Aberrant WSS, or the fluid friction at the wall, sensed by the endothelial cells, in part, mediates vascular wall remodelling in IAs [7–13]. Therefore, vortex structures near the surface could be more directly relevant to IA rupture than vortex structures in the bulk of the IA.

In a recent review, Meng et al. [13] highlighted the complex interaction between haemodynamics and aneurysm pathophysiology and theorized two destructive remodelling pathways that could potentiate IA growth and rupture. The first pathway is mediated by local high WSS and a high and positive WSS gradient, a flow condition often caused by a strong inflow jet impinging on the IA wall [13]. As demonstrated by in vivo studies, high WSS can trigger localized matrix metalloproteinase production via smooth muscle cells (instead of inflammatory infiltrates), leading to internal elastic lamina damage and cell apoptosis [42,43]. This biological cascade has been shown to cause media thinning and aneurysmal bulge formation [44]. Further, in a large computational study of 210 IA cases, high maximum WSS was found to be associated with aneurysm rupture status [9]. In the current study, we first qualitatively classified the flow in IAs and found that the local regions of jet breakdown were associated with a positive Q near the wall (figures 4 and 5) and overlapped with regions of high WSS (figure 10). This flow mode, classified by complex flow and an impingement jet, termed the Jet Breakdown mode, was further associated with ruptured aneurysms and a high sVF (figure 9). Thus, we speculate that positive Q near the wall underlies this local high shear stress environment that could activate mural-cell-mediated destructive remodelling mechanisms, potentially leading to IA rupture.

The other destructive remodelling mechanism put forth by Meng et al. is an inflammatory cell-mediated pathway driven by low WSS [13]. It has been previously shown that a low and oscillatory shear stress environment increases infiltration of inflammatory cells into the vessel wall causing wall degradation, which could potentially lead to rupture [13,45,46]. In computational studies of patient-specific IA geometries, it has been found that low aneurysm-averaged WSS and high OSI predict IA rupture status [8,10,11,27,47]. In our current study, a significant correlation was found between low aneurysm-averaged WSS and sVF (figure 11). Specifically, in aneurysms with the Jet Breakdown Mode, flow impinges on the aneurysm wall, becomes complex and forms counter-rotating vortices. This slow and recirculating flow, in turn, results in a low WSS environment in the aneurysms. We speculate that, at the luminal surface, when rotational flow is dominant (positive Q), pockets of recirculating flow could cause the endothelial cells to be pro-inflammatory and increase the residence time of leucocytes near the lumen, thereby encouraging their adhesion and transmigration to the wall.

Based on the results of this study, we suggest that positive Q near the IA wall is indicative of an adverse haemodynamic environment, which can include both local high and low aneurysm-averaged WSS. As shown in several cases (figure 10), flow impingement is associated with high local WSS and local positive Q. Additionally, due to the complex aneurysm flow, including slow recirculation zones, the aneurysm-averaged WSS is low. Both of these pathways could lead to the degradation of the aneurysm wall, potentially increasing the propensity of IA rupture. However, to understand the specific mechanism by which Q impacts the wall remodelling, biological experiments are required.

While we have adopted the Q-criterion for its simplicity [22], there are other definitions of vortex including the λ2-criterion [21] and the Δ-criterion [48]. Like the Q-criterion, these Galilean invariant Eulerian methods have been developed to identify vortex regions, and they also require a threshold for vortex structure visualization [49]. We previously conducted a sensitivity study for several aneurysm cases using the λ2-criterion, which is based on the three eigenvalues of extracted quantities of the gradient of velocity. Our results showed that there were negligible differences in vortex structures based on the λ2-criterion and those based on the Q-criterion [50]. This result was consistent with the observation by Chakraborty et al., who noted that methods using a threshold ‘result in remarkably similar looking vortical structures’ [51]. On the other hand, a different class of methods identify Lagrangian coherent structures. Such methods have been applied in abdominal aortic aneurysm studies [52] and can characterize topological flow features such as separation regions, vortex boundaries and flow impingement. Therefore, in future studies they could give additional insight into the behaviour of vortex dynamics in IAs.

This study has several limitations. First, the patient cohort came from a single centre, and thus there could be a selection bias. Second, we have analysed only cross-sectional data, which cannot directly address rupture propensity of unruptured IAs and the role of vortex dynamics therein. Third, patient-specific velocity waveforms at the inlet parent vessel were unavailable for our patient cohort. To perform pulsatile CFD simulations for our IA cases, we adopted the standard practice of scaling a generic waveform based on the inlet diameter. This corresponds to a constant pulsatility index [53], while previous studies have suggested that pulsatility index with geometry could influence vortex formation in the aneurysm [23,24,54]. Although our study indicates that the vVF and the sVF may be minimally influenced by the transient waveform (figure 6d), future studies should examine the impact of different pulsatility indices on these metrics. Fourth, while the vVF and the sVF quantify the dominance of vortical flow in IAs, they do not reflect the spatial complexities of vortex structures. Lastly, to make CFD tractable, we make several simplifications in computational modelling such as rigid walls. Although these common assumptions are known limitations, they have been determined to be justified for computational studies of IAs [55,56].

5. Conclusion

To elucidate the relevance of vortical flow to aneurysm pathophysiology, we developed and implemented a novel approach to quantify vortex flow structures in 204 patient IAs. Based on the Q-criterion definition of a vortex, we quantified the aneurysmal region occupied by vortex structures volumetrically and on the surface layer by the vVF and the sVF. Computational fluid dynamic simulations showed that the sVF, but not the vVF, discriminated ruptured from unruptured aneurysms. Furthermore, we found that the near-wall vortex structures co-localized with regions of inflow jet breakdown and significantly correlated to previously proposed haemodynamic and morphological characteristics of ruptured IAs.

Authors' contributions

N.V. was responsible for the study design, CFD simulations, post-processing and statistical analysis. G.T. was responsible for the CFD analysis and processing. J.P. was responsible for the CFD simulations. K.S. was responsible for the collection of data. H.M. was responsible for overseeing the project, and study design. All authors contributed to the manuscript preparation.

Competing interests

Funding

Support for this work was partially provided by NIH grant nos. (R01 NS091075 and R03 NS090193), a grant from Toshiba Medical Systems Corp and resources from the Center for Computational Research at the University at Buffalo.

Acknowledgements

The authors would like to thank Nikhil Paliwal, Hafez Asgharzadeh and Iman Borazjani for their stimulating discussion, and Vincent Tutino for his help in the preparation of the manuscript.

2012Guidelines for the management of aneurysmal subarachnoid hemorrhage: a guideline for healthcare professionals from the American Heart Association/American Stroke Association. Stroke43, 1711–1737. (doi:10.1161/STR.0b013e3182587839)

. 2009Molecular alterations associated with aneurysmal remodeling are localized in the high hemodynamic stress region of a created carotid bifurcation. Neurosurgery65, 169–177; discussion 177–168. (doi:10.1227/01.NEU.0000343541.85713.01)