We retrospectively studied N = 6 patients. All of them were enrolled from a database of patients with existing computed tomography (CT) scans with full RR coverage over the entire LA and LAA. Selection of patients was random with the constraint of covering a wide range of LA ejection fractions as assessed by CT. Cardiac-gated cine CT scans were acquired at the National Institutes of Health (NIH), Bethesda, Maryland (N = 3) and at the University of California San Diego (UCSD), CA, United States (N = 3). The studies were approved by the Institutional Review Board at both centers.

Imaging was performed following standard clinical protocols at each participating center. Specifically, electrocardiogram (ECG) gating with inspiratory breath-hold was used to sample one entire heartbeat. The images were acquired with 3 different scanner models: 2 studies were obtained on a Siemens Force system (at NIH), 3 studies were obtained on the Toshiba Aquillion ONE (1 and NIH, 2 at UCSD) and one study was obtained on the 256-slice GE Healthcare Revolution CT at UCSD. For each image iodinated contrast was injected followed by saline flush. The doses and injection rates of contrast were chosen based on the patient’s weight using standard clinical protocols.

The images were reconstructed using the manufacturers’ standard algorithms (Toshiba reconstruction filter AIDR3D, FC08-cfa, GE Healthcare AsirV = 0.5 standard kernel, Siemens filter kernel br36d). This procedure yielded DICOM files containing z-stacks of images of 512 x 512 pixels with in-plane pixels dimensions of between 0.32 mm and 0.48 mm in the x-y axial plane. Each z-stack contained between 100 to 320 slices with resolutions of 0.5 mm to 1 mm in the z-direction. The multi-phase reconstructions were performed at regular time intervals across the heart cycle, with 5% RR and 10% RR being the minimum and maximum intervals.

The 4D CT scans described above were imported into itk-SNAP (Yushkevich et al., 2006) and a segmentation was created using the semi-automatic active contour segmentation method. The PVs were clipped by manually choosing planes for each time frame; similarly, the LV and aorta were clipped off and the LAA was identified. From each 3D segmentation, a surface triangulation was extracted in MATLAB, smoothed and resampled to the resolution required for the CFD using iso2mesh (Qianqian and Boas, 2009) (an example is shown in Supplementary Figure 1). For each subject, the point cloud of vertices of the largest-LA-volume time frame was used as a reference surface and registered to the point clouds from the other time frames using the Coherent Point Drift algorithm (Myronenko and Song, 2010). The location of the centroid of each triangle was established by Fourier interpolation with N_FOU = 6 modes to match the temporal resolution of our CFD solver. This 4D kinematic LA model provides the position and velocity of the centroids of the triangles discretizing the LA surface at all times.

Indices of atrial function were calculated from CT images using standard definitions based on chamber volume, as well as on the waveforms of PV and transmitral flow (Hoit, 2014; Table 1). LA sphericity was calculated as described by Bisbal et al. (2013). The morphology of the LAA was categorized as one of the four types previously proposed by Di Biase et al. (2012) (chicken wing, cactus, windsock or cauliflower) by a cardiologist with expertise in cardiac CT interpretation. The LAA segmentations were skeletonized using the MATLAB image processing toolbox at each reconstructed phase of the cardiac cycle, and the branches of the skeleton were pruned, segmented and labeled automatically (Figure 1). The number of LAA lobes and the length of each lobe were measured (Table 2).

The simulations were performed with the in-house code TUCAN (Moriche, 2017; Moriche et al., 2017; Flores et al., 2020), which solves the incompressible Navier-Stokes equations for a Newtonian fluid of constant kinematic viscosity, v = 0.04 cm²/s, using a fractional step method. The discretization used a staggered, Cartesian grid and centered, second-order, finite differences. Time integration was performed with a low-storage, three-stage, semi-implicit Runge–Kutta scheme.

Each simulation was initialized from no flow conditions and run for 10 heartbeats at a low resolution (grid spacing Δx = 0.090 cm in all directions), with a constant time step Δt and a Courant number CFL < 0.3. The flow field was then interpolated into a finer grid (Δx = 0.051 cm, 256³ grid points) and run for another 10 heartbeats with the same CFL limit. This resolution is comparable to the finest of the resolutions employed in previous works (Vedula et al., 2015). In preliminary runs, an even finer grid was employed (Δx = 0.034 cm), showing no relevant differences in the results.

The blood residence time (T_R) was computed in the whole computational domain by solving the forced passive scalar equation (Rossini et al., 2016),

where T_R was set equal to zero at the start of the simulation and at the PVs. Eq. (1) was discretized like the Navier-Stokes equations, but spatial fluxes were evaluated with a 3rd order WENO scheme (Shu, 1998). This scheme ensures a robust solution of the residence time, avoiding Gibbs phenomena near steep gradients without sacrificing accuracy elsewhere. Given that T_R was initialized at zero, we ran the simulations for several heartbeats to let T_R rise until it converged to a quasi-periodic (i.e., recurrent albeit not exactly periodic) behavior at all points inside the LA (Figure 2). This convergence period is basically the washout time of the fluid particle with highest T_R, and can be particularly large inside the LAA. In our simulations, plotting the volume-averaged residence time inside the LAA suggested that ∼15 heartbeats were enough to achieve convergence (Figure 2). Hence, we ran our simulations for 20 heartbeats and used the last 4 heartbeats to compile statistics.

The LA surface was placed inside a cubic box of side 13 cm, using free-slip boundary conditions at the cube boundaries. The geometry and motion of the LA surface were prescribed using the patient-specific 4D kinematic LA models obtained from CT imaging as described in the previous two sections. To this end, we employed the immersed boundary method (IBM) (Mittal and Iaccarino, 2005). In these methods, the flow is solved in a fixed, Cartesian grid and the position and motion of surfaces are imposed by a localized volumetric force. In the present work, the IBM of Uhlmann (2005) was employed. Because the inlet planes of the PVs are not boundaries of the simulation domain, we also used the IBM to impose flow rates through these inlets. Specifically, the PV velocity was prescribed by applying a volumetric force f_PV in a buffer region upstream of each PV inlet plane. The target inlet velocity at each PV inlet plane was defined as vt=QiAi⁢n,where A_i and n are the cross-sectional area and the vector normal to each inlet’s plane, respectively, and Q_i (i = 1,…4) is the flow rate through each PV. The force was set proportional to the difference between the local velocity and the target inlet velocity, i.e., fP⁢V=1τ⁢(vt-v), where the constant τ is the characteristic time of the forcing. For all simulations presented here, we chose τ = 10Δt = 0.5ms, resulting deviations from the target velocity < 5% at the inlet. An analogous scheme was used to set the residence time to zero at the inlet of the PVs.

To specify the flow rates through the PVs, we calculated the combined flow rate through all the PVs, Q_PV from mass conservation inside the left atrium. We computed the difference between the LA volume change obtained from each patient’s segmentation and the flow rate exiting the LA through the mitral valve (Q_MV). The latter was obtained from the change in the ventricular volume during LV diastole, measured using each patient’s CT images. Since we did not have patient-specific measurements of the flow rate through each PV, Q_PV was split evenly among all the PVs, i.e., Q_i = Q_PV(t)/4. CFD analysis with evenly split pulmonary flow rates has been previously shown to produce reasonable LA hemodynamics when compared to phase contrast MRI data (Lantz et al., 2019). To further evaluate this choice of inflow boundary conditions, we performed preliminary simulations for one of the cases comparing this approach to splitting Q_PV proportionally to the area of the pulmonary veins, so that the bulk mean velocity was the same in all PVs. We did not find any significant differences between these two approaches in terms of the results presented in section 3. For the outflow boundary condition, the MV is considered as a plane that instantly opens (when Q_MV>0) and closes (imposing Q_MV=0). No specific model is considered for the geometry of the MV leaflets, since they are not expected to affect the flow structure inside the LA (Vedula, 2015).

To explore the influence of the motion of the LA wall, we performed corresponding simulations keeping the position of the LA walls fixed at the zero mode of its Fourier temporal series. For the fixed wall-simulations, the A-wave was digitally removed from the transmitral flow, setting Q_MV = 0 for the duration of the atrial contraction.

The median age of the study subjects was 65 (range 50 – 92), and two were female. Subjects 1—3 were imaged in sinus rhythm and did not have LAA thrombus. Subjects 5 and 6 were imaged in AF with an LAA thrombus, which was digitally removed prior to segmenting the LA geometry for the CFD simulations. Subject 4 was imaged in atrial fibrillation without LAA thrombus, but had a previous history of transient brain ischemic attacks (TIAs). AF was determined based on ECG rhythm strips and lack of late diastolic LA/LAA ejection (“atrial kick”). Based on this information, we categorized subjects 1 – 3 as LAA-thrombus/TIA negative (LAAT/TIA-neg), and subjects 4 – 6 as LAAT/TIA-pos.

Subjects 1 – 3 had normal LA functions as inferred from CT-derived chamber volumes and PV flows (Table 1). Subjects 4 – 6 had enlarged LA with impaired global function. Subject 4 had decreased reservoir function but relatively normal booster function. Subjects 5 and 6 had impaired reservoir and booster function. The LA chamber of subject 5 had severely impaired function, operating mostly as a conduit between the PVs and the LV.

Figure 1 displays cardiac CT images of the 6 subjects including segmentations of their LAAs and a plane projection of the morphological skeletonization of the appendage. These images show that subject 5 had an obvious isolated LAA thrombus in the central cavity of the LAA, while subject 6 had the distal section of the LAA filled with thrombus (which was confirmed on a “late” single acquisition obtained 22 seconds after the dynamic acquisition). The LAA skeletonization was used to determine LAA lobe number (i.e., the number of skeleton branches) and LAA length (i.e., the length of the branch whose end is most distal to the LAA orifice). Table 2 contains detailed morphological information of the LAA for our study subjects. While the small size of the cohort precludes detailed statistical examination, it is worth noting that the LAAT/TIA-pos subjects had more complex LAA shapes (i.e., windsock vs chicken or cactus) and larger LAA volumes than the LAAT/TIA-neg subjects. We did not observe striking differences in LAA lobe count or LAA length (Table 2).

We first present 3D velocity vector maps for one of the subjects with normal atrial function (subject 3, Figure 3). The vectors are colored according to their proximity to four anatomical landmarks: the left PVs, the right PVs, the LAA and the mitral annulus. Three representative instants of the cardiac cycle are shown: LA diastole, LV early filling, and LA systole (Figure 3, top, center and bottom row, respectively). We juxtaposed moving-wall and fixed-wall simulations (Figure 3, left and right column, respectively) to illustrate how wall motion affects LA flow dynamics. The flow rate profiles in the PVs and mitral annulus, and the volume profiles of the LA and LAA are also plotted.

Starting at atrial diastole, the mitral valve is closed, and LA volume is increasing as blood enters the chamber in the moving-wall simulations (Figure 3A, blue arrows). The collision of the left and right PV jets redirects blood towards the mitral annulus (red arrows) and the LAA (green arrows). In contrast, the fixed-wall simulation shows almost no fluid motion at the same timepoint (Figure 3B). Given that the chamber volume does not vary in fixed-wall simulations and the mitral valve is closed, the only possible fluid motion is the remnant flow from ventricular diastole.

During LV relaxation (E-wave, Figure 3C), blood still enters the LA through the PVs forming a jet collision pattern, but LV suction accelerates the blood towards the mitral annulus. Consequently, the flow rate exiting through the mitral annulus (Q_MV) exceeds the incoming flow rate through the PVs (Q_PVs) and LA volume decreases. The LAA also contracts during this phase, driving a relatively uniform outward flow in its interior. The fixed-wall simulation (Figure 3D) exhibits stronger flow during LV filling even if it was run with the same Q_MV. In the moving-wall simulation, passive LA contraction leads to a net difference between the outgoing and incoming flow rates Q_MV - Q_PV ≈ 400 ml/s. However, this difference must be zero in the fixed-wall simulation – thus, PV inflow jets are stronger. The LAA flow pattern in the fixed-wall simulation is reminiscent of a shear-driven cavity flow (Koseff and Street, 1984), consisting of a weak recirculating pattern driven by the external flow in the atrial body. This pattern is notably different to the one observed in the moving-wall simulations.

During atrial systole, chamber contraction drives flow into the ventricle (A-wave, Figure 3E) and simultaneously creates backflow into the PVs (Ar-wave). LAA contraction drives a relatively uniform outward flow inside the appendage, similar to the E-wave albeit somewhat slower. In contrast, the fixed-wall simulation displays a weak recirculatory flow corresponding to the late decay stage of the PV and mitral annulus jets generated during LV filling. In particular, the velocities in the LAA are much smaller than in the moving-wall simulation.

The results presented so far reveal that blood flow in the LA and LAA is sensitive to atrial wall motion. However, wall motion can be reduced in different ways depending on whether a patient has impaired reservoir and/or booster function (Hoit, 2014). Thus, we visualized LA and LAA flow patterns across the phases of the cardiac cycle for three patients: one with normal atrial function (subject 3, Figure 4A), one with impaired reservoir function (subject 4, Figure 4B), and one with both impaired reservoir and booster functions (subject 5, Figure 4C). Each snapshot shows blood velocity vectors inside the LAA and iso-surfaces of the second invariant of ∇⁡v→, which reveal the vortical structures of the atrial body. Similar to Figure 3, we juxtapose data from moving-wall and fixed-wall simulations.

The subject with normal atrial function exhibits two vortex-ring systems associated with the atrial filling jets discharging from the left and right PVs (Figure 4A). These vortex rings collide into each other and break down into a complex network of fine-scale vortices that fills the whole atrial body, and which is sustained throughout the cardiac cycle. The fixed-wall simulation presents two strong vortex-ring systems that break down into finer-scale vortices, but these decay almost completely by atrial diastole of the next cardiac cycle. As shown in Figure 3, the magnitude and dynamics of the velocity vectors in the LAA of this normal subject vary markedly between the moving-wall and fixed-wall simulations.

In patients with impaired atrial function, moving-wall and fixed-wall simulations yield qualitatively similar LA vortical structures, while the velocity fields inside the LAA exhibit interesting differences. In patient 4 (Figure 4B), who had severely impaired reservoir function (LA expansion index = 0.3, see Table 1), the LA wall was akinetic during atrial diastole and LV relaxation. Consequently, the moving-wall and fixed-wall simulations produced similar LAA velocity patterns during these phases. However, this patient had a relatively normal booster function (booster LA ejection fraction = 0.2) with normal LA wall kinetics during atrial systole. As a result, blood velocity in the LAA showed and outward flow pattern for the moving-wall simulations whereas it showed a weak recirculating motion for the fixed-wall simulations. In subjects with decreased reservoir and booster function, e.g., subject 5 in Figure 4C, the moving-wall and fixed-wall simulations produced qualitatively similar LAA velocity distributions across the whole cardiac cycle. Consistent with patient 5 having more passive than booster LA function, the most appreciable difference between fixed-wall and moving-wall simulations was observed for atrial diastole – during this phase, the moving-wall LAA has a weak inward-flowing pattern in addition to the swirl that is found in the fixed-wall LAA.

To visualize the effect of impaired atrial function on blood stasis inside the LA chamber, we mapped the blood residence time in two oblique intersecting plane sections of the LA body and the LAA. These maps are represented in Figure 5 for the same three study subjects of Figure 4, both for moving and fixed walls. Once converged (t ≥ 15 cycles, see Figure 2), the T_R maps show that blood is systematically more stagnant inside the LAA than in the LA body. Residence time gradients are also significantly sharper inside the LAA leading to elevated TR at the LAA apex, consistent with flow stirring being lower in the appendage than in the atrial body.

Visual inspection of the T_R maps suggests that LAA blood stasis in LAAT/TIA-neg patients (Figure 5A) was less prominent than in LAAT/TIA-pos patients (Figures 5B,C). To confirm this result, we compiled the statistical distributions of LAA T_R for all of our subjects using 200 time instants equally distributed over the cardiac cycle. We also compiled statistics of the kinetic energy K = 1/2 (u² + v² + w²), since this hemodynamic variable is more accessible via medical imaging than T_R. The data from moving-wall simulations show that LAA residence time increases as EF_LA worsens across different subjects (Figure 6A). Although less clear, we also found a trend for K to decrease as EF_LA worsens (Figure 6B). In the fixed-wall simulations all subjects have zero ejection fraction, but the data suggests that subjects with higher LA volume tend to have higher T_R, and lower K (Figures 6C,D).

Figures 6A,B demonstrates that blood inside the left atrial appendage of LAAT/TIA-pos patients had marked alterations in hemodynamic parameters when compared with LAAT/TIA-neg patients. Specifically, LAAT/TIA-pos patients had 57% higher mean residence time and 69% lower mean kinetic energy inside their LAA than LAAT/TIA-neg patients. Given that these trends seemed to hold for the fixed-wall simulations (94% higher mean residence time, 216% mean kinetic energy, Figures 6C,D), we examined the variation in LAA T_R and K between fixed-wall and moving-wall simulations (Figure 6E). This analysis showed that K decreased in the fixed-wall simulations across all study subjects; however, the variations in T_R were not systematic. In particular, fixed-wall simulations yielded equal or lower LAA T_R than moving-wall simulations for two of the three LAAT/TIA-neg subjects. While these results are not easy to reconcile with the dependence of T_R on EF_LA (Figure 6A), they are consistent with our finding that fixed-wall simulations tend to represent the LA and LAA flow patterns better in atria with impaired function than in normal atria (Figure 3). Finally, we examined whether fixed-wall simulations were less accurate in separating the LAAT/TIA-pos and LAAT/TIA-neg groups in our small study population. Scatter plots of LAA T_R and K (Figure 6F) showed some overlap between the LAAT/TIA-pos and LAAT/TIA-neg groups in the fixed-wall simulations. In contrast, moving-wall simulations provided good clustering of the groups for either T_R or K.

By selecting subjects with a wide range of LA ejection fractions (EF_LA), we characterized how LA and LAA flow patterns vary with atrial function. Consistent with previous MRI and CFD reports of flow in the LA (Kilner et al., 2000; Fyrenius et al., 2001; Vedula et al., 2015), we found that blood from the PVs flows mostly towards the mitral annulus during ventricular diastole and systole. In subjects with impaired reservoir function, reduced wall motion during ventricular systole caused slower flow and quicker dissipation of fine-scale vortices.

In the LAA, normal wall kinetics drove inward and outward flow jets. By contrast, flow inside the LAA was similar to a shear-driven cavity flow in subjects with impaired booster and reservoir function, consistent with Masci et al. (2019) data. Additionally, our subject with impaired reservoir function and normal booster function experienced hybrid LAA flow dynamics that oscillated between lid-driven cavity flow during LA diastole, and a weak outflow jet during LA systole. Taken together, our results indicate that flow patterns inside the LAA are not as complex as in the LA body or the LV, which we attribute to the strong geometrical confinement of the flow inside the appendage. These flow patterns led to elevated blood residence in the LAA compared to the LA body, in agreement with previous CFD results (Koizumi et al., 2015; Otani et al., 2016; Bosi et al., 2018; Garcia-Isla et al., 2018; Masci et al., 2019, 2020), and the consensus that thrombi form preferentially in the LAA. Patients with impaired atrial function also showed a trend for blood stasis to appear near the mitral annulus, since their weaker LA body vortical structures caused decreased stirring compared to normal subjects.

While previous CFD studies suggested that reduced wall motion impairs LA blood transit (Otani et al., 2016; Masci et al., 2019, 2020), there is a lack of systematic data about how functional parameters such as EF_LA affect this transit. We demonstrated that LAA T_R has a clear dependence on EF_LA. We also showed how LAA T_R is related to the risk of thrombosis, leveraging that our database includes two subjects with LAA thrombi and one with a history of TIAs. We and others have recently shown that LV blood stasis in the setting of myocardial infarction precedes LV mural thrombosis (Martinez-Legazpi et al., 2017; Garg et al., 2019) and cerebral microembolism (Delgado-Montero et al., 2019). The present study extends this idea to the LAA showing that LAA thrombosis could occur when the mean residence time inside the LAA exceeds ∼2 s (i.e., 2 cycles). While the number of subjects in our pilot study is too low to make statistically significant conclusions, this threshold agrees with Martinez-Legazpi et al. (2017) threshold for LV thrombosis derived from N = 72 patients. This agreement is particularly notable considering the different etiologies of myocardial infarction and AF, and the different methodologies used in our two studies.

In addition to residence time, we show that kinetic energy could be potentially useful to predict the risk of LAA thrombosis. This knowledge could be useful because, while T_R is a cumulative variable whose calculation requires solving a partial differential equation, calculating K is significantly less involved. In particular, K should be accessible by medical imaging modalities with moderate time resolution, such as phase-contrast MRI. Nevertheless, the observed predictive ability of K in the case of the LAA could be in large part due to the unique geometrical characteristics of the appendage and its particular location within the LA. We note that K can yield poor predictions of blood stasis in other scenarios (Hendabadi et al., 2013).

Fixed-wall simulations provided flow fields and stasis maps that differed from those obtained in moving-wall simulations, especially in subjects with normal atrial function. This discrepancy is not surprising; in addition to neglecting LA wall motion, fixed-wall simulations have oversimplified inflow-outflow boundary conditions and ignore LV – LA interactions. On the other hand, fixed-wall simulations only require static images which are more commonly available in clinical practice, and are less involved from the point of image segmentation and CFD analysis. Previous studies considered fixed walls (Bosi et al., 2018; Garcia-Isla et al., 2018) but the effect of this simplification on the simulation results had not been evaluated before. Furthermore, given the difficulty of accurately modeling wall motion and hemodynamics in fibrillating atria (Boyle et al., 2020), comparing fixed-wall and moving-wall simulations can provide useful information about how to interpret CFD results.

In our reduced study cohort, T_R and K in the fixed-wall LAA were less accurate predictors of LAA thrombus or TIAs than their moving-wall counterparts. These findings warrant further evaluation in larger patient cohorts and reflect the influence of LA and LAA geometry on blood flow. Morphological analysis of the LA chamber in our subject cohort (Table 1) revealed that the LAAT/TIA-pos subjects had the largest LA volumes and showed a trend for higher LA sphericity, consistent with clinical reports of stroke rates vs. LA/LAA geometry in AF patients (Di Biase et al., 2012; Bisbal et al., 2013; Yamamoto et al., 2014). Moreover, our LAAT/TIA-pos subjects, which had the highest LAA residence time in both our moving-wall and fixed-wall simulations, also had LAA morphologies that are ranked more pro-thrombogenic according to Di Biase et al. (2012) categorization. Although the number of LAA lobes measured by 3D ultrasound imaging and stroke risk have been shown to correlate AF patients (Yamamoto et al., 2014), we did not observe a relationship between LAA lobe number and T_R, K or presence of LAAT/TIAs (Table 2). This discrepancy could be caused by our small cohort size and differences in spatial resolution between CT and ultrasound.

The number of patients considered in this study was small (N = 6), even considering that we performed simulations with both moving and fixed walls for each subject. However, our patient cohort was diverse enough to explore how different aspects of LA function (e.g., reservoir, booster pump, and conduit) affect the flow inside this chamber with an unprecedented level of detail. Also, while small, the number of subjects considered here is higher than in most previous CFD studies (Otani et al., 2016; Bosi et al., 2018; Masci et al., 2020). More importantly, by including patients without LAA thrombus and with LAA thrombus or history of TIAs, we were able to demonstrate that CFD analysis of LA flow could predict the risk of LAA thrombogenesis.

We did not consider a model of the coagulation cascade in this study. Biochemical modeling of thrombogenesis could be incorporated into the CFD framework (Seo et al., 2016). However, the additional models would markedly increase the cost of the simulations. Also, the accuracy of these models is unclear, given their dependence on many reaction constants whose values are difficult to obtain in a patient-specific basis.

We considered the same heart rate and blood viscosity for all patients to eliminate these parameters as independent variables in our analysis, even if these parameters may vary from patient to patient. Besides, we note that heart rate in particular can vary significantly for each individual due to physical activity, stress, circadian rhythm and, of course, rhythm disorders such as AF. Analyzing how these intra-patient variations affect CFD predictions was beyond the scope of this study. In LAAT-pos patients, we segmented the thrombus and removed it before performing the CFD analyses in an effort to recapitulate the conditions that caused thrombosis. However, the presence of a thrombus could constrain the motion of the LAA walls – thus, pre-thrombus LAA wall motion might have not been fully recovered. One of our subjects had a history of TIAs of unknown origin, which we presumed to be cardioembolic from the LA/LAA, but we cannot exclude the possibility of another etiology.

Similar to previous simulations (Koizumi et al., 2015; Otani et al., 2016; Bosi et al., 2018; Garcia-Isla et al., 2018; Masci et al., 2020), our LA model did not consider the mitral valve in the outflow of the chamber. This choice is justified by whole left heart CFD studies showing that the mitral valve does not significantly affect LA flow (Vedula, 2015; Vedula et al., 2015). Mitral regurgitation, which can affect 15-25% of AF patients (Delgado and Bax, 2017; Owens et al., 2017), would require modeling the mitral valve. Our work improves over most previous computational studies of the LA flow (Bosi et al., 2018; Garcia-Isla et al., 2018; Masci et al., 2020) by considering patient-specific inflow/outflow boundary conditions, similar to Otani et al. (2016). One limitation of our simulations is that PV flow rates were evenly distributed (Lantz et al., 2019). More patient-specific modeling of inflow boundary conditions would involve measuring the flow rate through each pulmonary vein (e.g., by transesophageal echocardiography) or modeling the flow impedance of each pulmonary vein based on imaging data such as vein diameter (Shi et al., 2011).

As discussed above, fixed-wall simulations have severe limitations as a model of AF, including but not limited to their inability to capture LA motion and deformation mediated by the chamber’s passive response to incoming PV flow and LV-LA interactions. The rationale for including fixed-wall simulations in our study was to assess how these limitations affect the prediction of patient-specific hemodynamic parameters such as residence time, by comparing with moving-wall simulations. This assessment is warranted considering that several previous studies used fixed-wall simulations (Bosi et al., 2018; Garcia-Isla et al., 2018).

Finally, in patients with paroxysmal AF, imaging performed while patients are in sinus rhythm might have to be used to predict LAA thrombosis risk. Multi-scale models of atrial electro-mechanics coupled to CFD analysis (reviewed, e.g., in Boyle et al., 2020) could tackle this problem. However, these models combine patient-specific data (e.g., anatomical images) with data that is not patient-specific (e.g., fiber orientation, LA wall mechanical properties, channel gating kinetics, and fibrosis distribution). Furthermore, we note that, while there is an unequivocal correlation between AF and embolic stroke, the long-held paradigm that thrombi form in LAA during AF episodes, embolize, and travel to the brain causing stroke is now under scrutiny, after a recent clinical trial revealed absence of AF episodes in the month preceding thromboembolic events (Brambatti et al., 2014).