Modal analysis of blood flows in saccular aneurysms

IAs are typically categorized into two basic shapes: (1) saccular, or (2) fusiform Rayz and Cohen-Gadol (2020a). We focus on saccular aneurysms, which are more commonly found in the brain Frösen et al. (2012). In this work, the patient-specific imaging data of brain aneurysms was selected from the repository of Digital Imaging and Communications in Medicine (DICOM) files from the AneuRisk project Piccinelli et al. (2011a). Using a retrospective approach, the study cohort was formed from six aneurysms situated in the Internal Carotid Artery (ICA). An ICA aneurysm is the most common IA and its formation is linked to the bifurcation in the interconnected segment with the Middle Cerebral Artery (MCA) and the Anterior Communicating Artery (ACA). Based on the Bouthillier categorization Bouthillier, van Loveren, and Keller (1996), all six aneurysms occur in the vessel wall between the ophthalmic (supraclinoid) (C6) and communicative (terminal) (C7) segments of their ICAs (see Fig. 1). Specifically, the aneurysms of patients 1, 2, 3, and 4 were developed on the endothelial wall of the ICA artery at C6, while patients 5 and 6 had their aneurysms at C7.

3D Slicer Fedorov et al. (2012), an open-source software, was used for processing the DICOM images and reconstructing the 3D anatomical models. Segmentation is performed first and it yields the CoW and various brain artery segments, which are exported in stereolithography (STL) file format. Next, the model goes through several post-processing steps in another open-source software, Meshmixer (RRID:SCR_015736). In those steps, the CoW is isolated from the attached arteries and then cleaned up. A technique is carried out to smooth the arterial wall by re-meshing the model surfaces. Fig. 1 shows the final 3D models of all 6 patients from the segmentation and post-processing steps that only include the surfaces of the aneurysm and the ICA, the ACA, and the MCA. Finally, each surface mesh was converted into a solid mesh having around 50,000-100,000 vertices. The resulting models after reconstruction from the DICOM images are referred to as anatomical models.

For each aneurysm, three morphological attributes, i.e., artery diameter $D$, aneurysm diameter $W$, and aneurysm height $H$, were measured. Table 1 summarizes the anatomical attributes of the IAs. The morphological trio was used to calculate several aneurysmal indices, i.e, the Aspect Ratio $AR$:

and the Aneurysm Number $An$:

where the Pulsatility Index $PI$ is calculated from the velocity profile at the inlet such that:

where $U$ is the bulk velocity prescribed at the inlet.

The location of each aneurysm ostium was identified by following the 4-point technique from the Anatomy-Guided Multi-Level Exploration Neugebauer et al. (2011). Then, the central aneurysm axis (CAS) can be defined as a line perpendicular to the ostium at its center point. The angle of inclination Baharoglu et al. (2014), or aneurysm angle, is measured as the angle between the incoming flow and the CAS. In addition, the aneurysms can be further classified as terminal or lateral Sforza, Putman, and Cebral (2009) by measuring an angle between the directions of incoming and exiting flows.

A classification was carried out using anatomical features of the aneurysms and their parent’s ICAs. This is a prerequisite step to investigate the correlation between aneurysm size and its hemodynamics. The Unweighted Pair-Group Method with Arithmetic mean (UPGMA)Sokal and Michener (1958); Tajima (1990) is a generic method that was employed to analyze the morphological attributes in Table 1.

Assuming the attributes are of an orthogonal set that forms a three-dimensional space, then the location of each individual is a vector $\left[D_{j},W_{j},H_{j}\right]^{\intercal}$. A matrix of all positions has the following form:

with $k\in\mathbb{N}\>|\>k=[1,P]$ where $P$ is the total number of individuals. $P=6$ in our study, which makes the matrix $\mathbf{S}$ a 3-by-6 matrix.

Pairwise Euclidean distance between an arbitrary pair $(i,j)$ of the corresponding column vectors in the matrix $\mathbf{S}$ is calculated as:

The UPGMA method creates a simple agglomerative (bottom-up) hierarchical clustering by determining how objects in the dataset should be grouped into clusters. The distance between two clusters $r$ and $t$ is referred as a linkage:

where $n$ is the number of individuals in a cluster.

The dependability and viability of the DMD methods were investigated with consideration of the following parameters:

1. 1.

   The VOIs are subsets of the computational grids. They are set to contain the patient’s aneurysm sac. ParaView, an open-source software, was used to create the VOIs in all cases. The flow data in each case was extracted in a text file (CSV) format. The sensitivity of the DMD analysis to the choice of VOIs will be examined in detail in a later section.

2. 2.

   In order to assess the influence of spatial resolution on the DMD analysis, a spatial resampling technique was applied when extracting the VOIs with coarsening factors of two ($2\mathsf{X}$), four ($4\mathsf{X}$), and eight ($8\mathsf{X}$) times in all directions, as indicated in Table 2. At the lowest spatial resolutions (coarsening factor of $8\mathsf{X}$), the input data has a comparable spatial resolution to the ones from the 4D-flow MRI techniques ($\approx 1$ mm).

3. 3.

   In order to assess the influence of temporal resolution on the DMD analysis, a temporal resampling technique was performed by skipping the number of input instances. Three scenarios are considered, i.e., $M=50,~{}100,~{}\text{and}~{}200$, (or equivalently at the temporal resolution of $\Delta\tau=16.8,~{}8.4,~{}\text{and}~{}4.2$ milliseconds).

In summary, the 6 anatomical geometries, together with the spatial coarsening factors of $1\mathsf{X}$ (CFD resolution), $2\mathsf{X}$, $4\mathsf{X}$, and $8\mathsf{X}$ and the temporal resolutions of $M=50,~{}100,~{}\text{and}~{}200$ resulted in a total of 72 standard datasets for modal analysis.

In the Exact DMD method Tu et al. (2014); Askham and Kutz (2018), the data matrix $\mathbf{O}$ is broken into two sequences:

The DMD analysis assumes that the flow dynamics is based on the following formulation:

The Exact DMD algorithm aims to calculate the "dynamics" matrix $\mathbf{A}$ using the pseudoinverse ($\dagger$) of $\mathbf{X}$:

The Singular Value Decomposition (SVD) of the matrix $\mathbf{X}$ is carried out as a prerequisite to the calculation of the pseudoinverse matrix $\mathbf{X}^{\dagger}$. The task of finding the matrix $\mathbf{A}$ involves solving an optimization problem using the Frobenius norm as follows:

where $\mathbf{L}$ and $\mathbf{R}$ are left- and right-unitary matrices of SVD. The diagonal matrix $\mathbf{\Sigma}$ contains the singular values $\sigma$. Alg. 1 provides a step-by-step procedure of the Exact DMD method for computing the matrix $\mathbf{A}$ using non-truncated SVD on the matrix $\mathbf{X}$.

The Hankel variant Taira et al. (2017b); Kamb et al. (2020); Pan and Duraisamy (2020) involves applying the Exact DMD algorithm to a Hankelized matrix constructed from the original data, rather than directly to the data snapshots shown in the Alg. 1. This modification aims to address the pitfall of forward-time bias.

A time-delay embedding dimension $d$ can be defined when constructing the Hankelized matrix $\mathbf{X}_{H}$, given data matrix $\mathbf{O}$, as:

where $d=1$ was chosen in our modal analyses. The resulting Hankelized input matrix has the following form:

To visualize a DMD mode, a three-dimensional structure (iso-contour) is constructed at a threshold $\frac{|\mathbf{V}|}{|\mathbf{V}_{\text{max}}|}=0.5$ in the ParaView software. In addition to the 3D structures, the corresponding frequencies of these modes are calculated as $\omega_{k}=\frac{\log(\lambda_{k})}{2\pi}$ for an arbitrary $k^{th}$ DMD mode. Then, the frequency spectrum is scaled by temporal magnitude $|b_{k}\lambda^{M}_{k}|$, instead of the scalar magnitude $|b_{k}|$ Le (2021). Ultimately, the temporal magnitudes are normalized by the maximum value.

Optimized DMD Askham and Kutz (2018) is introduced to alleviate the limitations of the Exact DMD method when dealing with non-linear phenomena. The method performs the mode decomposition process on the entire data $\mathbf{O}$ with the objective is minimizing the difference between the actual data $\mathbf{O}$ and the data reconstructed from the Optimized DMD modes.

The algorithm starts by re-writing the matrix $\mathbf{O}$ as follows:

The optimization starts with initial guesses $\bm{\gamma}$ for the eigenvalues $\bm{\lambda}$ and $\bm{\beta}$ for the matrix $\mathbf{B}$. The algorithm continues to refine iteratively by solving an exponential fitting problem such that:

After the convergence, the Optimized DMD modes are recovered as:

where $\mathbf{B}^{\intercal}(:,i)$ is the $i^{th}$ column of $\mathbf{B}^{\intercal}$. And the Optimized DMD amplitudes are:

As the optimization is carried out over the entire dataset $\mathbf{O}$, this DMD algorithm significantly enhances the robustness and accuracy of mode decomposition and reduces the impact of noises which is the second disadvantage of the Exact DMD method. In addition, the Optimized DMD algorithm allows us to capture the inherent non-linearity and transient phenomena in flow dynamics. These attributes are essential in accurately assessing flows in IAs.

When comparing different variants, Askham and Kutz Askham and Kutz (2018) suggest measuring the quality of a DMD method by comparing the reconstructed snapshots $\left(\bm{\Phi}(\bm{\lambda})\bm{\Phi}(\bm{\lambda})^{\dagger}\mathbf{O}^{\intercal}\right)^{\intercal}$ with the original ones $\mathbf{O}$ such that:

The cases in Group 1 show two distinct features in the dynamics of the inflow jet: (i) the splitting of the incoming flow; and (ii) the existence of helical flow in ICA. The splitting of the incoming flow is evident at the ostium, which is indicated by elliptically highlighted areas in Figure 4. Patient 2 receives a larger portion of the incoming flow, which demonstrates that the splitting feature can depend on the morphology (refer to the $AR$ as discussed in sec. III.1) of the aneurysm. In both Patient 1 and Patient 2, the inflow jets enter the aneurysm sac at the posterior regions (note that aneurysm 1 is on the contralateral ICA). The inflow jet enters through the aneurysm ostium with a velocity of around $130$ cm/s and $60$ cm/s, respectively. The former velocity magnitude is much higher than the systolic peak of the bulk velocity at the inlet ($50$ cm/s). It is worth noting that the instantaneous streamlines in the cavernous (C4) and the C5 segments reveal a noticeable helical vortex in both cases.

In Group 2, the unique characteristic is the incoming flow traversing the aneurysm ostium near its central point. Depending on the incoming angle, the inflow jet can either approach the distal wall (Patient 3) or redirect toward the dome of the aneurysm (Patient 4). Because the morphologies of Patient 3 and Patient 4 are almost similar in terms of $AR$ (clustering results in Sec. III.1), the bend curvature of the C4 segment of the ICA, indicated by squared highlighted areas in Figure 4 for Patient 3, is responsible for the redirection. This demonstrates the influence of the incoming angle on the dynamics of the inflow jet.

The patients in Group 3 exhibit the highest velocity in the inflow jet, evidenced by regions in Fig. 4 having a magnitude exceeding $160$ cm/s. In the case of Patient 5, the high magnitude can be seen in the ICA segment before the incoming flow enters the aneurysm. Meanwhile, the inflow jet of Patient 6 achieves its highest magnitude after penetrating through the ostium. It is worth noting that the flow enters through the anterior region in this group. In addition, Group 3 possesses both differences in the aneurysmal and the ICA morphologies which have been observed separately in the previous groups as the contributor to the differences in the inflow jets’ entrances.

Fig. 4 shows that the inflow jets in patients 1, 2, 3, 5, and 6 are of the wall-impingement type. Meanwhile, Patient 4 represents the sole instance of the dome-impingement type within our cohort. It is noteworthy that within the population of wall-impingement occurrences, the locations of impingement are inconsistent.

It is evident that the impingement type draws a bidirectional relationship between the characteristics of the inflow jet and the development of the aneurysms. Patient 4 (dome-impingement) has a pointy dome in comparison to a flatter one in the aneurysm of Patient 3 (wall-impingement). Again, analyses have shown that the two patients are very morphologically similar. More noticeably in patients 5 and 6 (wall-impingement), the jets directly hit the part of the distal wall where the middle section of the crescent aneurysms is located (see Sec. III.1).

The main vortex is formed due to the separation of the inflow jet after making an impingement. To discern unique attributes in their formation and rotation, the relative angle between the centerline of the main vortex and the CAS is qualitatively examined in Fig. 4. In the case of Patient 1, after colliding with the distal wall in the lower chamber, the inflow jet disperses and continues in a spiral pattern toward the upper chamber. This spiral motion creates the main vortex rotating clockwise around the CAS of the aneurysm. Within our study cohort, Patient 6 is the only instance showing a similar swirling pattern around the CAS, even though the size of this patient’s IA is two times larger (averaged to $12$ mm) than Patient 2’s average size of $6$ mm. The main vortices in the cases of patients 2, 3, 4, and 5 are orientated perpendicular to their CAS, creating a contrasted appearance.

The Hankel DMD modes are sorted by their scalar amplitudes. Fig. 5 shows the spatial structures (iso-contour) of the top 3 modes. The frequencies of these modes coincide with the frequency trio of the velocity waveform at the inlet (refer to Sec. II.2). In 5 out of 6 cases (patients 1, 3, 4, 5, and 6), the frequency $f_{1}=3.6$ Hz is consistently ranked as the most energetic mode (Mode 1 - yellow) followed by the energy of frequencies $f_{2}=2.4$ Hz and $f_{3}=1.2$ Hz as the second place (Mode 2 - dark green) and the third place (Mode 3 - purple), respectively. The case of Patient 2 displays a unique set of the dominant frequencies, i.e., $f_{1}=3.81$ Hz, $f_{2}=1.2$ Hz, and $f_{3}=2.65$ Hz, which differ slightly from those recorded in other patients.

The iso-contour of the three dominant modes effectively captures the visualization of the large-scale dynamics, the inflow jet. The spatial extents vary among the patient cases, as shown in Fig. 5. The traces of the inflow jets are immediately terminated after impingement in the cases of patients 1, 4, and 5. The inclusion of Patient 4 in this group indicates that the termination of the traces does not correlate with the impingement types specified in Section III.2.2. In contrast, the trajectory of the inflow jets of patients 3 and 6 can be captured after their impact with the distal walls, and the traces only terminate until they reach the domes. It is worth noting in most cases, the frequency of $1.2$ Hz possesses the highest spatial coverage and is often found in the terminal segment. The spatial extents in the case of Patient 2 possess uniqueness, again, from the study cohort, such that there exists a ring structure from Mode 2 ($f_{2}=1.2$ Hz).

In addition to the spatial extent of the iso-contour surfaces, the separation of the dominant modes can be qualitatively analyzed to classify the type of inflow jet. In this regard, patients 1, 2, and 3 exhibit "diffused" jets, while the jets are "concentrated" in the cases of patients 4, 5, and 6. One can see that the classification of the inflow jet does not correlate with the grouping scheme because both types of inflow jet are present in Group 2 (medium size). On the other hand, a comparison between two patients in Group 3 shows that two aneurysms in our cohort can belong to the same classification of inflow jet (concentrated) but differ in spatial extent.

The blood flow dynamics in IAs in one cardiac cycle are mostly contributed by the low-frequency ($<15$ Hz) modes in all the cases. Fig. 6 shows the pseudo-spectra from 4 cases, i.e., Patient 1 (Group 1), Patients 3 and 4 (Group 2), and Patient 5 (Group 3), at the temporal resolution of $M=200$ (represented by blue bars). The spectra increase uniformly from the beginning and peak at $3.6$ Hz. Observing the overall landscape, the top three modes are the ones captured in Sec. III.3.1. This quantity analysis in the temporal magnitude scale shows that the most energetic modes (ranked by scalar magnitude $|b_{k}|$) are also the most stable. After the peak at $3.6$ Hz, each spectrum exhibits a gradual decrease in normalized temporal magnitude. As noted from previous results (Fig. 5), the modes in the range $1-5$ Hz depict the inflow jet. Therefore, the dynamics in the upper range ($5-15$ Hz) are from the breakdown flows. The fact that the frequencies in this range are highly damped shows that these dynamics die out quickly in time.

The four cases in Fig. 6 share a common pattern, i.e., that the contributions of the DMD modes in the high-frequency domain ($>15$ Hz) are small. None of the spectra sees a frequency peak above 40% (represented by a red dashed line) of the maxima at $3.6$ Hz. It is worth mentioning that the landscapes exhibit patient specificity in both appearance and temporal magnitude. The spectra of patients 1 and 3 witness a flat landscape, with few frequencies barely reaching the $0.1$ line (represented by a green dashed line). However, Patient 3’s landscape exhibits more fluctuation from $30$ Hz to $80$ Hz. Increasingly, the presence of fluctuation is not only more noticeable in the spectral landscape of Patient 4 but also appears instantly after $15$ Hz. High-frequency peaks around $20$ Hz, $65$ Hz and $110$ Hz reach up to 20% magnitude of the maxima. Despite the two patients being the closest in terms of morphology, high-frequency is more pronounced in Patient 4, especially when comparing the domain from $65$ to $115$ Hz. Concerning Patient 5, the spectrum is flat from $15-110$ Hz, a rather stark contrast to the smaller aneurysms in Group 2. There exists one peak of frequency at $117.54$ Hz accounting for a temporal magnitude equivalent to 25% of the maxima.

The contribution of high-frequency fluctuations can be great in certain aneurysms, those are the cases of patients 2 and 6. Fig. 7A depicts their spectra showing that there exist domains in the high-frequency region that are as stable as the large-scale dynamics. In the case of Patient 2, the majority of the DMD modes before $110$ Hz surpass the 10% line. One can clearly notice that the domain from $40$ to $70$ Hz stands out because it includes frequencies passing the 40% line (represented by a red dashed line). The landscape in this domain includes two peaks at $50.90$ Hz and $60.70$ Hz accounting for normalized magnitudes of 0.96 and 0.80, respectively. It is worth noting that Patient 2’s spectrum witnesses significant fluctuations around $10$ Hz, which is yet another uniqueness of this patient. In the case of Patient 6, the active domain is at higher frequencies ($>70$ Hz). This shift creates a frequency gap between the large-scale dynamics and the high-frequency fluctuations, in which the active domain of Patient 2 ($40-70$ Hz) sits in between. In contrast to the previous case, these frequencies in the spectrum of Patient 6 are damped. The active domain in this case begins at $72.85$ Hz and ends after $82.43$ Hz. It experiences more volatility and sees only one substantial peak at $76.93$ Hz, which has a temporal magnitude accounting for 80% of the maxima.

The iso-contours of the DMD mode in interested domains are reconstructions for Patients 2 and 6. In Fig. 7B, the spatial structures of these modes appear fractured, which is in juxtaposition to the smooth and continuous surfaces reconstructed using the dominant modes (refer to Fig. 5). The spatial extents of the high-frequency modes are located where the inflow jets (from the dominant modes, depicted in gray) separate after the impingements and form vortical structures in the aneurysm sac. In Patient 2, the fragments of high-frequency modes concentrate at the sac’s volume center at which the main vortex is formed. Meanwhile, in Patient 6, the high-frequency modes are evident after the inflow jet impinges on the aneurysm’s distal wall.

The frequency bandwidth scales proportionally with the temporal resolution. At $M=50$, all spectra are maxed out at $30$ Hz. In Fig. 6 and Fig. 7A, the landscapes of the spectra for $M=50$ (represented by yellow bars) bear resemblance to their corresponding sections at $M=200$ (highlighted by grey zones). By quantitative assessment, the frequencies are identical between the two resolutions. In addition, the similarity is qualitatively higher in the group of 4 patients having damped high-frequency modes. A consistent pattern emerges across the cohort, i.e., frequencies ranging from $5$ to $15$ Hz gain higher temporal magnitudes and frequencies beyond $20$ Hz show the opposite pattern. Ultimately at this reduced temporal resolution, high-frequency fluctuations are absent.

In Fig. 8A, the patients’ normalized CE curves appear to have no overlaps at $M=50$. The curves are separated and their order can be identified with ease especially when scaling the CE plot in the log-log scale. Here, the order of the patients from top to bottom is as follows: Patient 1 $\rightarrow$ Patient 3 $\rightarrow$ Patient 6 $\rightarrow$ Patient 4 $\rightarrow$ Patient 5 $\rightarrow$ Patient 2. It is worth noting that a significant gap is present between patients 3 and 6, consequently, creating a visual barrier that divides our cohort into 2 groups: i) the upper, and ii) the lower. The gap gains substantial evidence when the curves reach higher CE levels.

To further quantify the distance between the two groups, a vertical line (dotted cyan) is taken as a reference line at Mode 10 (which is equivalent to 20% of the total number of modes). CE value at the reference line is measured for each case, and then the average value can be calculated from the members of the same group. The upper group consists the case $P_{1}-1\mathsf{X}-M_{50}$ (0.91) and the case $P_{3}-1\mathsf{X}-M_{50}$ (0.88) which averages to 0.89. In the lower group, the cases $P_{2}-1\mathsf{X}-M_{50}$, $P_{4}-1\mathsf{X}-M_{50}$, $P_{5}-1\mathsf{X}-M_{50}$, and $P_{6}-1\mathsf{X}-M_{50}$ reach, respectively, 0.78, 0.81, 0.77, and 0.82 in normalized energy after 10 modes. This group averages 0.79 resulting in a normalized energy gap of 0.1 with the upper group, or 10% of the total energy.

Taking advantage of the log scale, another quantitative approach is used employing linear fitting. The technique calculates the best-fit lines using the CE values of the first 10 modes for each group. Fig. 8 shows that our approach yields two optimal lines that differ in slope. The optimal line for the upper group (in blue) has a slope of 0.27. Whereas, the slope for the lower group’s line is 0.22. This analysis presents an index showing that the upper cases have fast decay curves, and the ones from the lower cases exhibit slower decay.

In this sensitivity study, the Hankel DMD method is applied to the patient-specific dataset at the original resolution. A combination of 6 aneurysms and 3 temporal resolutions in this section yields 18 cases with the following notations $P_{1}-1\mathsf{X}-M_{50}$, $P_{2}-1\mathsf{X}-M_{50}$, $P_{3}-1\mathsf{X}-M_{50}$, $P_{4}-1\mathsf{X}-M_{50}$, $P_{5}-1\mathsf{X}-M_{50}$, $P_{6}-1\mathsf{X}-M_{50}$, $P_{1}-1\mathsf{X}-M_{100}$, $P_{2}-1\mathsf{X}-M_{100}$, $P_{3}-1\mathsf{X}-M_{100}$, $P_{4}-1\mathsf{X}-M_{100}$, $P_{5}-1\mathsf{X}-M_{100}$, $P_{6}-1\mathsf{X}-M_{100}$, $P_{1}-1\mathsf{X}-M_{200}$, $P_{2}-1\mathsf{X}-M_{200}$, $P_{3}-1\mathsf{X}-M_{200}$, $P_{4}-1\mathsf{X}-M_{200}$, $P_{5}-1\mathsf{X}-M_{200}$, $P_{6}-1\mathsf{X}-M_{200}$.

The pattern of the two groups identified in Sec. III.4.1 remains unchanged when the temporal resolution increases to $M=100$ and $M=200$. Increasingly, the gap between the groups is consistent at around 10%. In Fig. 8B, taking the same measurements and calculations at Mode 20 (20% of the total 100 modes) shows that the upper and lower averages, respectively, are 0.94 and 0.84. In Fig. 8C, the averaged values at Mode 40 (20% of the total 200 modes) are 0.98 for the upper group and 0.89 for the lower group. In addition, comparing the optimal lines in Fig. 8A, B, and C shows that the slopes decrease when more snapshots are added to the DMD input. Specifically, the upper group exhibits a rate of 0.01. The lower group decreases at a rate that is twice as fast showing a rate of 0.02.

Next, the temporal resolution is kept constant at $M=200$ while the spatial resolution varies. The analysis in this section presents two patients, and they are the representatives of their groups (refer to Fig. 8), i.e., Patient 1 is in the fast-decay group and Patient 6 is in the slow-decay group. As Patient 6 has been shown to exhibit significant high-frequency fluctuations (refer to Fig. 7B), we perform an addition resampling case for this patient’s CFD data having a coarsening factor of $16\mathsf{X}$ to achieve a low spatial resolution ($>1.0$ mm) normally obtained from 4D-flow MRI measurements. The 8 cases are denoted as $P_{1}-1\mathsf{X}-M_{200}$, $P_{1}-2\mathsf{X}-M_{200}$, $P_{1}-4\mathsf{X}-M_{200}$, $P_{1}-8\mathsf{X}-M_{200}$, $P_{6}-1\mathsf{X}-M_{200}$, $P_{6}-4\mathsf{X}-M_{200}$, $P_{6}-8\mathsf{X}-M_{200}$, $P_{6}-16\mathsf{X}-M_{200}$. In addition, the CE values in the cases with original resolution ($P_{1}-1\mathsf{X}-M_{200}$ and $P_{6}-1\mathsf{X}-M_{200}$) are chosen to be the baselines. Consequently, error bounds are calculated using the baselines.

The overall disparities in the CE curves in Fig. 9A and Fig. 9B among different spatial resolutions from $2\mathsf{X}$ to $8\mathsf{X}$ fall within 1% error from the baselines for both cases of patients 1 and 6. The biggest error between any low-resolution curve and the baseline is in the first DMD mode. By examining the magnified views of this mode included in Fig. 9, it becomes apparent that the lower-resolution curves exhibit an up-shifting behavior. However, this behavior is indeterminable. In the case of Patient 1, the increments in the energy of the first mode are present (see the magnified view in Fig. 9A). However, the increment is not linearly proportional to the coarsening factor. For example, the first mode of the normalized CE curve from $P_{6}-8\mathsf{X}-M_{200}$ is closer to the baseline than the first mode from $P_{6}-4\mathsf{X}-M_{200}$ (see the magnified view in Fig. 9B). Nevertheless, at an extremely coarse data of $P_{6}-16\mathsf{X}-M_{200}$, the disparity at the first mode is less than 1.5% error from the baseline. As the CE curves progress, the errors get smaller as values from cases at different temporal resolutions converge.

Finally, this section analyses the influence of various scenarios when extracting the VOI. Using a similar approach, the temporal resolution is fixed at $M=200$, and the CFD data is kept at the original resolution. The cases are notated as $P_{4}-1\mathsf{X}-M_{200}$ to reflect the current settings.

Fig. 10A illustrates three possible options for the VOI of Patient 4. This patient is chosen out of six individuals in the cohort because they represent the intricacy of aneurysm location that influences the extracting process. Specifically, the patient’s aneurysm developed in the vicinity of the bend from the C4 to the C6 segments (refer to Fig. 2A). Therefore, a section of the internal carotid artery (ICA) can interfere with the VOI box. Three distinct VOIs have been established for Patient 4 as follows: i) Domain A encompasses the complete aneurysm sac, a portion of the ICA’s C4 in the caudal volume, and the posterior boundary truncated at the ostium; ii) Domain B comprises an incomplete sac by omitting the protrusion of the proximal wall in the caudal volume, this avoids the C4’s inclusion; iii) Domain C includes the complete sac along with arteries of the ICA in the posterior direction.

The absence of aneurysmal volume at the proximal wall (Domain B) produces a minimal influence on the normalized CE curve. Still, the curve slightly shifts upward when compared to the curve obtained from Domain A (baseline) in Fig. 10B. The behavior of the CE curve in the case of Domain B can be seen as similar to the cases with low spatial resolutions. Fig. 10B shows that there is a significant displacement in the normalized CE curves obtained from Domain C compared to the baseline. In detail, adding sections of the parent artery raises the energy of the first mode by a staggering 28% from the baseline value. The inset in Fig. 10B as a magnified view of the region of modes $30-50$ shows that after mode 40, the values from the three curves have less than 1% error. Consequently, the normalized CE curve from Domain C exhibits the most significant deviation.

The Optimized DMD method outperforms the Hankel method in terms of dynamic accuracy because it has lower reconstruction error, especially during the systole period for both cases $P_{4}-1\mathsf{X}-M_{50}$ and $P_{4}-1\mathsf{X}-M_{200}$ as shown in Fig. 11. The Optimized DMD method achieves remarkably high precision during the systole phase (error $\approx 1\times 10^{-8}$). As the temporal resolution increases to $M=200$, the quality curve of Optimized DMD (in blue) attains significantly smaller errors ($\approx 1\times 10^{-6}$) compared to the curve from the Hankel method (in purple). There exists a stark contrast as the average error for Hankel DMD is $\approx 1\times 10^{-1}$, with notably high errors during systole (from $t/T=0.15$ to $t/T=0.45$). In addition, the diastole phase at $M=200$ sees the errors between the two methods become close to each other in the range $1\times 10^{-2}$. Ultimately, the Optimized DMD still outperforms Hankel DMD in terms of reconstruction quality during diastole, even at their peak error at around $t/T=0.6$, regardless of temporal resolution.

The Hankel DMD maintains a consistent quality curve shape across resolutions, while the opposite is true for the Optimized DMD method. Increasingly, the reconstruction quality of the Hankel DMD method improves as temporal resolution increases. In contrast, the increase in temporal resolution worsens the performance of Optimized DMD. The evidence is clear when comparing the value of the cases $P_{4}-1\mathsf{X}-M_{50}$ and $P_{4}-1\mathsf{X}-M_{200}$ from Hankel DMD at any timestamp in Fig. 11. Both lines remain at low error, then sharply increase during early systole. The peaks are at around $t/T=0.17$, and then the quality curves gradually decrease throughout the late systole and diastole phases. The result from the Optimized method shows a partially similar trend at the temporal resolution $M=200$, but there are differences. First, the most noticeable is the skyrocketed period when entering the systole phase. However, the peak can be seen at an earlier time $t/T=0.13$. Second, the Optimized DMD method’s precision decreases with time in the cardiac cycle, with a significant error increase (two orders of magnitude) in systole, becoming especially clear as $t/T>0.65$ (diastole). Interestingly, with $M=50$, the increasing error theme during the systole phase in the quality curve of the Optimized method is replaced by a decreasing one. The performance of the Optimized DMD across different temporal resolutions, as evidenced in Fig. 11, exhibits an unpredictable nature in maintaining consistent reconstruction quality.

The above results show that the differences in the algorithm of the two methods influence the identification of hemodynamic features. To delve deeper into the analysis between Hankel and Optimized DMD methods, the pseudo-spectra of Patient 4’s dynamic modes is depicted in Fig. 12(A). A clear distinction can be seen between Hankel and Optimized DMD methods. In Hankel DMD, the dominant frequency ($f_{2}=3.63$ Hz) agrees well with ones from the inflow waveform (refer to Fig. 2). The highest frequencies are sustained from $65.89$ Hz to $110.92$ Hz. However, the Optimized DMD does not provide such a frequency range. Instead, the spectrum of Optimized DMD concentrates to a much lower range of $20$ Hz with a peak of $15.1$ Hz. In other words, the Optimized DMD does not provide the expected frequencies of the inflow waveform. It also misses all frequencies above $16$ Hz, which indicates that the high-frequency modes of Optimized DMD have a minimal contribution to the reconstruction of large-scale hemodynamics.

The spatial extents of the two methods’ DMD modes are depicted in Fig. 12(B). The Optimized DMD provides the dominant modes at different frequencies $f_{1}=3.48$ Hz, $f_{2}=10.33$ Hz, and $f_{3}=15.10$ Hz (see Fig. 12). Regarding the dominant frequency of $15.10$ Hz, its fragments were scattered densely near the proximal wall. Hankel DMD provides a well-defined structure for the inflow jet (shadowed - $1.2$ Hz, $2.4$ Hz, and $3.6$ Hz). The higher frequencies ($15.1$ Hz, $65.89$ Hz, and $11.92$ Hz) localize around the inflow jet. The location of frequency $3.48$ Hz coincides with the location of impingement. In comparison to the peak frequencies of the Hankel DMD approach, the accuracy to identify the location of impingement was significantly higher using the Hankel dominant mode (3.63 Hz). In addition, the traces of high frequencies from $65$ Hz to $115$ Hz concentrated around the inflow jet indicated the breakdowns of the jet at impact and its interactions with the vortices.

Fig. 13 demonstrates that Hankel DMD modes exhibit greater stability compared to Optimized DMD modes throughout the cardiac cycle across all six patients. It is most evident when comparing the eigenvalue spectra in six cases $P_{1}-1\mathsf{X}-M_{200}$, $P_{2}-1\mathsf{X}-M_{200}$, $P_{3}-1\mathsf{X}-M_{200}$, $P_{4}-1\mathsf{X}-M_{200}$, $P_{5}-1\mathsf{X}-M_{200}$ and $P_{6}-1\mathsf{X}-M_{200}$. More Hankel eigenvalues localize on the unit circle’s circumference than Optimized eigenvalues. Most Optimized DMD eigenvalues within the left semicircle are contained inside the circle, indicating damping in the corresponding dynamic modes. Consequently, stable modes from the Optimized DMD method concentrate on the right side of the circumference (first and fourth quadrants), corresponding to frequencies below 16 Hz with significantly higher contributions (Fig. 12A). The disparity in eigenvalue distribution between the two methods reflects their different frequency dispersions.

In conclusion, while the Optimized DMD method shows superior reconstruction quality, especially at lower temporal resolutions, the Hankel DMD method demonstrates better frequency range coverage, more accurate representation of flow features, and greater eigenvalue stability.