The Impact of Cardiac Motion on Aortic Valve Flow Used in Computational Simulations of the Thoracic Aorta

Computational modeling has been used to quantify mechanical forces that correlate with disease in the abdominal aorta [1–4], coronary arteries [5,6], and cerebral vasculature [7,8]. Recent reports have extended these techniques to children with congenital cardiovascular disease and their treatments [9–12]. Advancements in this area are particularly important for patients diagnosed with coarctation of the aorta (CoA) [13–19]. Computational analysis techniques are continually being improved in an attempt to better recreate physiologic blood flow patterns in the thoracic aorta where indices of wall shear stress (WSS, defined as the tangential force per unit area exerted on a blood vessel wall as a result of flowing blood) have correlated with intimal thickening [20–23].

Recent studies have used CFD models to simulate blood flow in the AscAo using geometric [24] and patient-specific [19] representations of the AoV orifice. Accurate assessment of AoV blood flow is imperative for determining the likelihood of future complications in these patients stemming from local hemodynamic alterations and determining the optimal timing for surgical treatment. Unfortunately, the accuracy of measurements through the AoV can be compromised by the presence of cardiac motion (up to 13 mm), which occurs throughout the cardiac cycle [25–28]. Previous studies have attempted to account for the motion of the heart during blood flow measurements in a variety of ways. One study by Kayser et al. measured blood flow through the tricuspid valve by phase-contrast magnetic resonance imaging (PC-MRI) and attempted to correct for through-plane velocity by subtracting the phase shift found in pixels within 3–4 pixels of the vessel contour [29]. The reliability of this technique was highly dependent on the chosen contour, the signal-to-noise ratio of the stationary tissue, and relatively low signal from the vessel wall. In an attempt to improve on this technique, Kozerke et al. [28] elegantly described a method in which the aortic annulus could be tagged and tracked, and its location could be used to dynamically assign a PC-MRI slice throughout the cardiac cycle. The method did an excellent job of compensating for cardiac motion, but required specialized MRI sequences, software designed to run on the scanner console, and a feedback system to dynamically assign slice locations in a real-time fashion throughout the scan.

The current study describes a method analogous to that developed by Kozerke et al. [28], which can be implemented offline using standard images acquired during the course of a clinical imaging session. The organization of this investigation begins with a description of the method developed to quantify and account for the basal motion of the heart. This method was then applied to a case of congenital heart disease, using patient-specific CFD models of the thoracic aorta, to elucidate the severity of errors from basal cardiac motion that may influence calculated mechanical indices (e.g., WSS and TKE) known to be associated with cardiovascular disease.

A single patient (18-year old female with CoA repaired by resection with end-to-end anastomosis) was selected for this research from a Children's Hospital of Wisconsin Institutional Review Board (IRB) approved database of patients with congenital cardiovascular disease who had clinically indicated cardiac MRI studies. IRB approval allowed use of these anonymized MRI data for CFD research, with informed consent from the patient. The gadolinium-enhanced (0.2 m mol/kg gadodiamide, Omniscan^®, GE Healthcare, Waukesha, WI) MRI angiography (MRA) scan was used to create the CFD model. Imaging was performed using a 1.5 Tesla Siemens Symphony^® scanner (Siemens Healthcare, Erlangen, Germany) with electrocardiogram (ECG) triggering and a breath-held 3D fast gradient echo sequence. Slice thickness was 1.2 mm, with 56 oblique-sagittal slices per volume. A 384 × 224 acquisition matrix was used with an in-plane field of view (FOV) of 25 × 37 cm (in-plane spatial resolution of 1.0 × 1.0 mm). Other parameters included a repetition time (TR) of 4.3 ms, echo time (TE) of 1.4 ms, and a flip angle of 25 deg. Bi-orthogonal measurements using the MRA images were obtained throughout the thoracic aorta and branches.

Steady-state free precession cine cardiac imaging was performed through the left ventricular outflow tract (LVOT) to localize the aortic annulus and track its motion throughout the cardiac cycle. The LVOT image was prescribed through the AoV perpendicular to the three-chamber view (which in turn was acquired perpendicular to the short axis image through the AoV) to provide an accurate estimate of the plane of the aortic annulus. Images were acquired during a single breath-hold and retrospectively gated to the cardiac cycle. TR, TE, and flip angle were 2.8 ms, 1.4 ms, and 80 deg, respectively. Acquired temporal resolution was 49.8 ms, resulting in 25 interpolated images per cardiac cycle. Other parameters included a slice thickness of 5 mm, FOV of 31 × 35 cm, and acquisition matrix of 168 × 192 (in-plane spatial resolution of 1.8 × 1.8 mm).

Velocity-encoded 2D magnitude and through-plane PC-MRI was performed through the AoV, in the AscAo at the level of the branch pulmonary arteries, in the descending aorta (dAo) at the level of the diaphragm and orthogonal to the head and neck arteries. These blood flow data were used to prescribe inlet and outlet boundary conditions for the CFD analysis. PC-MRI parameters included a 30 × 22.5 cm FOV, 256 × 192 acquisition matrix, 46 ms TR, 3.8 ms TE, and 30 deg flip angle. The patient was breathing freely, and data were retrospectively gated to the cardiac cycle using the ECG signal (25 frames per cardiac cycle). The patient was then removed from the magnet and a single supine, bilateral upper and lower extremity blood pressure (BP) assessment was performed using Dinamap (GE Healthcare, Waukesha, WI).

A program was created in Matlab^® (Mathworks, Natick, MA) that allowed a user to identify a single landmark on the LV wall and ventricular septum in the LVOT cine image (Fig. 1). The physical coordinates of the markers were calculated using the fields in the Digital Imaging and Communications in Medicine (DICOM) header, such as Image Position and Image Orientation as well as voxel dimensions. A plane was created intersecting these points and perpendicular to the imaging plane (and thus perpendicular to the aortic annulus). The normal vector of the aortic annulus was calculated by taking the cross product of the line of intersection and the normal vector of the imaging plane. This resulted in a normal vector pointing up through the AoV. This procedure was repeated for each image in the series, and the coordinates of the basal plane in consecutive images were subtracted to calculate the distance traveled for each phase of the cardiac cycle. This information was used together with the time step between images (determined from the Trigger Delay field in the DICOM header) to calculate the velocity of the basal plane.

where n^B represents the normal vector of the aortic annulus, n^A represents the normal vector of the PC-MRI plane, v^B_through represents the velocity of the aortic annulus calculated above, and v^A_through represents the component of aortic root velocity in the through-plane direction of the AoV PC-MRI slice. This scaled, through-plane velocity was smoothed using a moving average filter with temporally adjacent points to remove any noise from the tracking procedure when applied to blood velocity through the AoV. This prevents contamination of the flow signal by the tracking procedure, especially in diastole. Finally, the area of the aortic annulus was calculated for every frame of the velocity images, and these areas were used with v^A_through to calculate the flow contribution (Q^A_through = v^A_through area) to the measurements obtained from the AoV PC-MRI scan. This flow component was removed from the AoV flow measurements, resulting in the corrected (i.e. motion-compensated) flow waveform used in this study.

CFD models of the aorta and its branches were created using several open source software packages. ITKSnap 2.2.0 (PISC Laboratory, University of Pennsylvania, Philadelphia, PA) [30] was used to create 3D volumes of the vascular regions of interest (Figs. 2(a)–2(c)). The “Snake ROI” tool was used to create the vessels of interest with “intensity regions” used to differentiate the active contour region. Three-dimensional bubbles were then placed throughout the aorta and branches contiguously and within the vessel. The propagation algorithm (i.e., region-growing) was set to “expanding” and “detailed,” and the propagation was allowed to proceed (number of iterations dependent on vessel) until the vessels were filled. Once the final volume was created, any regions which propagated outside the vessels were removed, and deficiencies in the vessels were filled in by hand. VMTK² [31,32] was then used to smooth the model and extract the vessel centerlines (Figs. 2(d)–2(f)). These centerlines were used in conjunction with the smoothed model to extract circumferential segments at distinct locations throughout the vessel using matlab (Fig. 2(g)). This was accomplished by inscribing spheres along the centerline at specified distances and extracting the (x, y, z) coordinates where these spheres intersected the smoothed model surface. Finally, these segments were refined using 20 spline points in each segment to better define the vessel and lofted into solid models using SimVascular³ (Fig. 2(h)). The individual vessels were loaded into SimVascular, and their intersections with the aorta were blended using fillets of defined radii (which were informed by branch vessel dimensions). The model originated at the aortic annulus as defined by the MRA dataset and extended the length of the thoracic aorta to the level of the diaphragm including the innominate artery (IA), right common carotid artery (RCCA), left common carotid artery (LCCA), and left subclavian artery (LSCA). Models were constructed by the same user and one additional user to determine intra- and interobserver variability in model construction, respectively. Biorthogonal measurements were performed at the same locations in the ascending, transverse, and descending thoracic aorta to determine any disparity between models.

PC-MRI images were used to calculate time-resolved volumetric blood flow using specialized software (Segment)⁴ [33]. Phase errors from eddy currents and magnetic field inhomogeneities were corrected using a second-order spatially dependent algorithm [34]. Instantaneous flow rates were computed by integrating these corrected velocity values over the lumen cross sections [35]. AoV waveforms were then corrected for basal cardiac motion and interpolated using a Fourier method where the number of simulation steps was determined for optimal simulation convergence (i.e., Courant–Friedrichs–Lewy condition <1). Plug inlet velocity profiles were used to define the shape of bulk flow entering the model at each interpolated time point. While it may be desirable to use PC-MRI to sample the velocity profile downstream of the valve for direct input into the model, this requires appropriate through-plane and in-plane velocity encoding to adequately resolve flow features. This approach was difficult to implement within clinical constraints as it requires specialized sequences and obtains data that are more detailed than that commonly used in clinical imaging. The alternative approach applied here constructs CFD models with their inlet beginning at the aortic annulus, imposes the measured blood flow waveform as an assumed velocity profile at the model inlet, and allows the curvature and related geometry of the arch to influence resulting flow patterns [16]. This approach does not require specialized sequences, minimizes the introduction of noise at the model inflow due to inadequate velocity encoding, and allows for improved temporal resolution compared to three-component PC-MRI [36].

Flow waveforms from the IA, LCCA, RCCA, LSCA, and dAo were used, in conjunction with measured BP data, to prescribe outflow boundary conditions as described previously [16,37]. To replicate the impact of arteries downstream from model branches, three-element Windkessel models consisting of characteristic resistance (R_c), capacitance (C), and distal resistance (R_d) parameters were imposed using the coupled multidomain [37] and pulse pressure [38,39] methods. Briefly, the total arterial capacitance (TAC) is first determined from AscAo and upper extremity BP measurements assuming a characteristic to total resistance ratio of 6% [40]. The TAC was then distributed among outlets according to their mean PC-MRI flow distributions [41]. Terminal resistances (Rt = R_c + R_d) were then calculated from each outlet using mean BP and PC-MRI flow measurements. R_c to R_t ratios were then adjusted for each outlet using the pulse pressure method in order to replicate measured BP values.

Models were discretized using a commercially available, automatic mesh generation program (Fig. 2(i)) (MeshSim, Simmetrix, Clifton Park, NY). The same computational model and initial isotropic mesh parameters were used for simulations that were then performed independently for the original and motion-compensated AoV flow measurements. Newtonian and incompressible fluid assumptions were employed with fluid viscosity and density of 4 cP and 1.06 g/cm³, respectively. Simulations were performed with the commercially available linear flow solver, leslib (Altair Engineering Inc., Troy, MI), which uses a novel stabilized finite-element method to solve equations for conservation of mass (continuity) and balance of fluid momentum (Navier–Stokes). Final meshes contained > 3 × 10⁶ tetrahedral elements, but with localized refinement (≈8 × 10⁶ or more isotropic tetrahedral elements) that used an adaptive technique [42,43] to deposit more elements near the luminal surface and in regions prone to flow disruption specific to the simulation being performed (i.e., original or motion-compensated AoV flow). Convergence criteria included residual errors on the order of 10⁻³ and a minimum of six nonlinear iterations per time step. Simulations were run for ten cardiac cycles until the flow rate and BP fields stabilized, yielding a periodic solution.

TAWSS was calculated using an established method [44] and visualized using paraview (Kitware Inc., Clifton Park, NY). TAWSS was quantified in 3 mm circumferential bands throughout the thoracic aorta and branches as described previously [12,17]. Additionally, the surface geometry of the thoracic aorta was unwrapped using the techniques developed by Gundert et al. [45] to visualize TAWSS over the entire thoracic aorta in finer detail. Briefly, the surface of the vessel was unwrapped whereby each (x, y, z) node of the mesh was mapped to a (θ, l) coordinate system, where θ represents the circumferential location (0–360 deg) with the 0-deg location defined as the underside (i.e., inner curvature) of the aortic arch, and I represents the length down the thoracic aorta. Additionally, TAWSS was quantified at discrete locations longitudinally and circumferentially along the aorta to further elucidate regions of greatest disparity between corrected and uncorrected inflow waveforms.

for each spatial location, $xj$⁠. Mean TAWSS values in the dAo (⁠$τmean−DAo$⁠) were used for normalization in regions of low WSS in the uncorrected simulations to prevent the overestimation of error at these points. This approach is similar to previous techniques [46]. Locations in the patient's CFD model where the influence of the inflow waveform was greater than established levels of interobserver variability in the CFD model building process [47] were then identified. Moreover, differences between TAWSS distributions for simulations using corrected and uncorrected AoV PC-MRI measurements as model inputs (see Fig. 3) were also scrutinized relative to differences from successive mesh densities (2.3 × 10⁶ versus 3.2 × 10⁶ adapted elements).

To investigate the influence of inflow boundary conditions on turbulence in the thoracic aorta, TKE was calculated using previously described methods [1,45]. The average velocity field was calculated by taking the ensemble average of spatially equivalent locations within the velocity field for the final five cardiac cycles after obtaining a converged solution. The fluctuating velocity field, defined as TKE, was determined by subtracting this averaged representation from the original velocity field. Similarly, the ensemble-averaged velocity field is defined as the KE. Finally, the TKE/KE ratio was calculated for each time point. The arch was isolated from its branches and divided into AscAo, transverse arch, and dAo sections. Average TKE, KE, and TKE/KE ratios at peak systole, mid-deceleration, and mid-diastole were then quantified. These data were used to determine the influence of basal component of inflow velocity waveform on TKE.

Measurements performed throughout the thoracic aorta and branches on the MRA data (Table 1) showed a mild-to-moderately dilated AscAo (max diameter: 3.2 cm × 3.2 cm and normal diameter: 2.5 cm × 2.5 cm), normal transverse arch (1.9 × 1.8 cm), and normal descending thoracic aorta (1.5 cm × 1.5 cm, normal: 1.6 cm × 1.6 cm). These measurements were similar to those in the simulated model. The maximum difference between models created for intra- and interobserver variability was approximately 2 mm (Table 1).

The patient's AoV waveform exhibited higher peak flow during systole and seemingly elevated antegrade flow during diastole relative to the AscAo waveform (Fig. 3(a)). Overall, cardiac output was ∼4% higher at the AoV than in the AscAo (AoV = 4.4 l/min versus AscAo = 4.2 l/min). Assessment of basal cardiac motion showed the peak through-plane velocity of the heart was 5 cm/s with a maximal long-axis contraction of 2.3 cm that occurred at 62% of the cardiac cycle. Compensating for the motion of the aortic annulus resulted in an increase in peak systolic blood flow and elimination of the majority of diastolic antegrade flow (Fig. 3(b)).

Through-plane velocity was compared between the PC-MRI data acquired in the AscAo and the corresponding slice location in CFD simulations. Comparisons between three time points during the acceleration phase of systole (t1), near peak systole (t2) and during the deceleration phase (t3) are shown in Fig. 4. During the acceleration phase, the velocity profile using the uncorrected (middle row) and corrected (bottom row) inlet conditions were similar to the PC-MRI velocity data (top row). However, near peak systole (middle column) and during the deceleration phase (left column), the profiles from the simulations differed. When compared to the PC-MRI data, the through-plane velocity profiles for motion-compensated simulations appeared to better match the PC-MRI velocity profile, specifically in the lateral extent of high velocity flow at peak systole, the anterior extent of forward flow (∼100 cm/sec), and areas of recirculation in the posterior and rightward portion of the aorta.

Velocity streamlines also obtained from CFD simulations were slightly elevated throughout the aortic arch in the corrected case as illustrated by the jet in Fig. 5. The inflow jets from corrected and uncorrected simulations, therefore, appeared to differentially impact the mid-AscAo vessel wall near a region of dilation and influenced patterns of WSS as discussed below. Peak systolic blood flow appeared laminar and fully attached to the vessel wall from the distal AscAo throughout the thoracic aorta regardless of inflow waveform.

TAWSS was generally higher in the corrected simulation (Table 2), although the distribution of TAWSS appeared similar using either inlet (Fig. 6). When quantified in 3 mm bands, the differences in TAWSS between corrected and uncorrected inflow simulations ranged from 12% to 21% (e.g., proximal AscAo 41 versus 51 dyn/cm²), with corrected TAWSS typically being higher. When the aorta is unwrapped, the largest differences between corrected and uncorrected inflow waveforms are seen along the outer wall of the AscAo. For example, the outer left wall of the AscAo showed pronounced differences through the dilated region (Fig. 6, circumferential location A: 27 dyn/cm²), while the outer right showed considerable differences between inlets just before the IA (Fig. 6, circumferential location B; difference = 25 dyn/cm²). Differences persisted, but were more modest along the inner and outer right walls of the transverse arch (e.g., 14 dyn/cm²) and in the distal dAo (e.g., 15 dyn/cm²). The differences in TAWSS between simulations were consistently above those reported for interobserver variability [47] for ∼70% of thoracic aorta (Fig. 7). Figure 8 shows the cumulative percentage distribution of area (left) exposed to values of TAWSS from 0 to 200 dyn/cm² for simulations using corrected and uncorrected AoV PC-MRI measurements as model inputs. Differences in the cumulative area between 5 and 150 dyn/cm² are more pronounced for model inlet conditions than for successive mesh densities (2.3 × 10⁶ versus 3.2 × 10⁶ elements). TAWSS values from 0 to 50 dyn/cm² are particularly interesting from the perspective of the vascular response to hemodynamics and are highlighted in the histograms for values within this range (right). A more dense mesh and correcting for cardiac motion both shift TAWSS from lower potentially adverse values to higher values traditionally through to be protective.

The distribution of TKE throughout the thoracic aorta was driven by the dilated portion of the AscAo and the curvature of the aortic arch (Fig. 9). Overall, the relative intensity of TKE was elevated for the simulation using the corrected waveform (Table 3). For example, during mid-deceleration and mid-diastole, the mean TKE was ∼16% higher for the corrected inflow (uncorr = 597 g/cm·s² versus corr = 713 g/cm·s²). The corrected waveform also had higher KE throughout the arch during systole and mid-deceleration, resulting in an overall lower TKE/KE ratio in the corrected case (AscAo, peak systole: uncorr = 0.22 versus corr = 0.17 and mid-deceleration: uncorr = 0.41 versus corr = 0.38).

Investigation into the hemodynamic basis of morbidity and treatment outcomes for adult and pediatric patients with diseases of the thoracic aorta can be aided by computational modeling tools, such as CFD. In doing so, the goal is to quantify local WSS indices or other indices known to correlate with disease in a manner that may be difficult during a routine clinical imaging session. In order for these CFD results to provide prognostic data or assist clinicians in treatment decisions, accurate inflow boundary conditions must be used. At a minimum, the inherent variability introduced by the motion of the heart should be considered.

The novel method presented here provides a way of accounting for errors inherent in blood flow measurements commonly used as inflow boundary conditions in patient-specific CFD analyses. This was achieved with several postacquisition operations conducted on data from a routine clinical imaging session. Our method includes quantification of blood flow by PC-MRI and assessment of basal cardiac motion by cine MRI images obtained through the LVOT. Compensation was accomplished offline by identifying markers surrounding the aortic annulus on the ventricular septum and LV-free wall, tracking these points throughout the cardiac cycle, and removing the velocity component resulting from this motion from subsequent blood flow assessment. Therefore, the fundamental physiologic sources of this error were accounted for without imposing any additional requirements at the time of MRI acquisition. The motion-compensated AoV flow waveform better matches the anticipated physiology (i.e., cessation of antegrade flow in diastole); however, there are still differences compared to the AscAo waveform. These differences may be due to a number of factors, including the fact that these are from separate PC-MRI acquisitions. Additionally, while the imaging plane is orthogonal to the vessel in the AscAo, the flow may not be directly though-plane. In fact, studies using 4D flow imaging show helical flow patterns in the AscAo [21,22], especially in the presence of dilation and AoV disease. This may result in an underestimation of systolic flow because this swirling flow may not be orthogonal to the imaging plane.

The variability in flow measurements due to cardiac motion [48] resulted in differences for hemodynamic indices, such as WSS and TKE. This finding suggests that future CFD modeling should incorporate a correction for cardiac motion. Our techniques were able to replicate the findings of the Kozerke study, but with the potential added benefit of using scans obtained during a clinically ordered MRI session. The time for these post-hoc analyses to correct for cardiac motion was similar to the time required to quantify PC-MRI data (∼10 min total). The velocity and distance traveled by the base of the heart also coincide with other values found in the previous studies [25–27]. In our study, the peak velocity was slightly lower than those found by Kozerke et al. (5 cm/s versus 8 cm/s). However, the degree of cardiac motion varies from patient-to-patient [28]. In fact, in quantifying the difference in cardiac output between the AoV and AscAo for five additional patients, the difference was ∼10%. This translates to an overestimation of AoV blood flow of ∼6%, which is greater than that which produced the results for the patient featured here. The majority of this excess flow seems to occur during diastole when the basal plane is moving toward the imaging slice.

Recent advancements in postprocessing techniques have allowed the investigation and visualization of localized sites of altered hemodynamics that may be useful in predicting disease progression in the thoracic aorta. To that end, circumferential quantification of WSS indices was used here to elucidate difference with respect to inflow waveform and anatomic location. The largest disparity in both spatially averaged TAWSS (above) and local indices of WSS (e.g., see Fig. 6) was present in the anterior wall of the AscAo, which was similar to locations of deleterious blood flow velocity and indices of WSS found by Bauer et al. [49] and Hope et al. [50,51]. Similarly, dilation and aneurysm formation tends to occur at the sinotubular junction and AscAo [49,52] in patients with AoV disease (e.g., bicuspid AoV). We and others have shown that this region experiences altered WSS in patients with bicuspid AoV [19]. The collective results suggest that it is important to account for cardiac motion in studies investigating AscAo hemodynamics in patients who have a bicuspid AoV. More recently, 4D blood flow imaging using MRI has been used to investigate flow disruptions in the AscAo of these patients [21,22,53]. However, these methods require specialized pulse sequences outside of the clinical workflow, are complex to acquire and postprocess, and prolong scan time (∼15 min for a 4D navigated flow imaging data set). Additionally, these methods suffer from low spatial/temporal resolution relative to CFD simulations, which could limit the precision of the WSS results compared to those from a properly conducted CFD simulation where uncertainties in the processes employed were appropriately considered [47,54,55].

Aortic measurements in the MRA data from the patient featured here suggest a mild-to-moderately dilated AscAo, with normal dimensions elsewhere. The three-dimensional modeling method resulted in minimal differences in intra- and interobserver variability, which suggest that these methods could be more robust than previous 2D methods at replicating patient geometry. While CFD results were not obtained for models used to assess this variability, the differences in TAWSS resulting from inlet conditions exceeded the threshold for interobserver variability mentioned elsewhere [47] throughout the majority of the aorta (Fig. 7) and those due to the computational mesh (Fig. 8). This finding may be case-specific, but the current results nonetheless suggest that cardiac motion is an important factor to consider when describing patterns of WSS in the thoracic aorta, particularly if WSS results are going to be used in clinical decision-making.

Past studies have shown flow disturbances and turbulence influence BP, indices of WSS, vascular remodeling, and inflammation [56]. While mean TKE in this case study was slightly elevated for the motion-compensated waveform throughout the arch, when normalized by KE, the TKE/KE ratio was lower for the this waveform. The TKE/KE ratio for the motion-compensated case was elevated in mid-diastole most likely due to the removal of antegrade flow in diastole, suggesting that the original simulation (using only AoV flow) may overestimate turbulence in the AscAo. Conversely, it appears that arch geometry may be the dominant factor influencing flow turbulence in the dAo. Future studies including various arch types are required to better elucidate these effects. The KE normalization may overestimate the impact of turbulence when the overall flow is low (i.e., diastole). In these cases, using raw TKE may be a more reliable measure.

From a clinical perspective, correcting blood flow measurements obtained at the AoV could more reliably demonstrate the degree of aortic insufficiency. Clinicians use indices such as aortic regurgitant fraction along with LV dimensions, LV systolic function, and a patient's clinical presentation to determine the optimal time for AoV surgery [57,58]. In addition, clinicians follow aortic diameter and its rate of change, along with the underlying etiology, to determine the optimal timing for surgical interventions in the presence of concomitant AoV and AscAo disease [59,60]. This work, combined with the inclusion of patient-specific AoVs into computational simulations of the thoracic aorta, could be used to more accurately estimate aortic regurgitant fraction or the risk for continued AscAo dilation. If hemodynamic indices from these simulations could be correlated with changes in AscAo dimension, future results could be used in combination with clinical parameters to improve patient care and outcomes.

It is worth noting that for some clinical conditions, the exact index, or range of values for a given index, that will result in disease progression is not yet known. Therefore, in all the clinical examples above, accounting for cardiac motion should allow for a more accurate understanding of how hemodynamic indices are associated with disease progression by eliminating cardiac motion that masks the “true” range of hemodynamic values experienced by the thoracic aorta. For example, the methods presented here could aid in the longitudinal and serial study of patients diagnosed with a bicuspid AoV by permitting more accurate CFD identification of hemodynamic characteristics that are associated with progressive AscAo dilatation and/or aneurysm formation. If these studies are able to estimate what features of AoV disease influence the rate of disease progression in the AscAo, then key indicators could be identified as risk factors for premature AscAo dilation and may serve as targets for medical intervention. Therefore, patients at greater risk for AscAo dilation, aneurysm formation, and dissection could be followed more closely with a more accurate sense of hemodynamically significant ranges of TAWSS and TKE in order to determine the optimal time for surgical intervention.

The current results should be interpreted within the constraints of several potential limitations. The approach outlined here relies on several MRI sequences, as well as several software packages that require specific operator skills and processing steps. Although the gadolinium-enhanced 3D angiography was necessary for CFD model generation, this sequence can add an invasive variable to the steps presented here if not clinically indicated. MRI cost–benefit ratio and scanning time demand should be taken into account when designing new methods. We previously quantified CFD simulation time for patient-specific carotid artery models [61]. The ITKSnap and VMTK steps, complexity of flow patterns during simulations, and modest additional post-hoc analysis time to correct for cardiac motion are believed to increase the time per simulation 5× as compared to this prior study. There is inherent variability between users throughout the modeling process from model construction, simulation, and application of the cardiac motion compensation. For example, when selecting the basal markers, they need to remain perpendicular to the aortic annulus and identify the same anatomic locations throughout the imaging set. Future studies may investigate the robustness of this method relative to clinical variables (e.g., slice selection, patient anatomy, etc.) and incremental improvements made to the correction algorithm. A tagging sequence [62] that highlights these regions can limit this variability, but the tags tend to fade in diastole when the tracking becomes most important. The images used to perform tracking in this investigation were also of relatively low-pixel resolution (∼1.75 mm) and may therefore introduce noise in relatively stationary periods in the cardiac cycle despite the smoothing algorithm that was implemented.

Other assumptions to the CFD simulation process here included rigid vessel walls to reduce computational expense, and because material properties (thickness and stiffness of the aorta) were not clinically available. Future studies will include compliant vessel walls to elucidate possible changes in distensibility and tensile stress in the AscAo and dAo of these patients. Intercostal arteries, which can account for up to 11% of aortic flow [63,64], are not included in this model and may result in some discrepancies between MRI data and simulation results. However, these arteries are usually not resolvable on MRA data due to their small size. Future work will also focus on expanding these techniques to larger patient groups to investigate any possible connection between adverse indices of WSS and thoracic aortic disease.

In summary, accounting for cardiac motion in AoV blood flow may produce more precise measurements of hemodynamic variables known to be associated with long-term morbidity for the thoracic aorta, such as WSS and TKE. This method brings these biomedical simulations of the aorta one step closer to (patho)physiologic realism. Combining this technique with prior developments like patient-specific valves [19] may ultimately enhance the acceptance and confidence in CFD models for clinical follow up of patients with thoracic aortic disease and/or for optimal surgical planning.

The authors thank Timothy Gundert and Shawn Shadden, Ph.D., for the technical assistance with software used in the analysis of results. The funding support was provided by an AHA Predoctoral Fellowship Award No. 0810093 Z (DCW), the Alvin and Marion Birnschein Foundation, and NIH Grant No. R15HL096096-01 (JFL). This research was also funded in part by the National Science Foundation Awards Nos. OCI-0923037 “MRI: Acquisition of a Parallel Computing Cluster and Storage for the Marquette University Grid (MUGrid)” and CBET-0521602 “Acquisition of a Linux Cluster to Support College-Wide Research & Teaching Activities.”