Verification Study of In Silico Computed Intracardiac Blood Flow With 4D Flow MRI

Objective: Major challenges for clinical applications of in silico medicine are limitations in time and computational resources. Computational approaches should therefore be tailored to specific applications with relatively low complexity and must be verified and validated against clinical gold standards. Methods: This study performed computational fluid dynamics simulations of left ventricular hemodynamics of different complexity based on shape reconstruction from steady state gradient echo magnetic resonance imaging (MRI) data. Computed flow results of a rigid wall model (RWM) and a prescribed motion fluid-structure interaction (PM-FSI) model were compared against phase-contrast MRI measurements for three healthy subjects. Results: Extracted boundary conditions from the steady state MRI sequences as well as computed metrics, such as flow rate, valve velocities, and kinetic energy show good agreement with in vivo flow measurements. Regional flow analysis reveals larger differences. Conclusion: Basic flow structures are well captured with RWM and PM-FSI. For the computation of further biomarkers like washout or flow efficiency, usage of PM-FSI is required. Regarding boundary-near flow, more accurate anatomical models are inevitable. Significance: These results delineate areas of application of both methods and lay a foundation for larger validation studies and sensitivity analysis for healthy and diseased cases, being an essential step upon clinical translations.


I. INTRODUCTION
T HE application of computational fluid dynamics (CFD)   to compute cardiac hemodynamics attracted increasing research over the last decade.The motivation can be found in the relevance of cardiovascular diseases to global health care systems.In addition, the potential of numerical cardiac hemodynamics has not yet been fully exploited.Although stateof-the-art non-invasive imaging modalities (i.e., phase-contrast magnetic resonance imaging (MRI) and echocardiography) for cardiac flow measurements are well-developed, they are limited in spatio-temporal resolution, entail additional patient scan times, and may be restricted by metallic implants.Regarding these aspects, image-based CFD can contribute to support diagnosis and treatment planning.Current research in the field addresses the left ventricle (LV), left atrium (LA), or aorta (AO), structures where cardiovascular diseases often originate or manifest.Subject of investigation are, e.g., intraventricular forces and their impact on remodeling processes [1], energetic efficiency of intraventricular flow [2], [3], [4], or cyclic washout and blood mixing processes in the context of thrombus formation [5], [6], [7].Furthermore, virtual anatomical manipulations allow prior estimation of post-operative outcomes [8], [9], [10], offering the potential to provide advice where clinical standard procedures are absent.
To date, various image-based CFD approaches exist to compute LV hemodynamics.They are based on different imaging modalities and numerical methods.Some approaches furthermore address high fidelity modeling in the context of multiphysic considerations [11], [12], [13] or detailed modeling of anatomical structures [3], [14], [15].Despite the quantity of methods, clinical applications are still mostly distant.Contributing factors are challenges in posing appropriate boundary conditions (BCs) and estimating associated uncertainties.Another aspect is the complexity of the methods, being accompanied by high data requirements and demanding computational as well as expert resources.In a recent study by Lesage et al. [16], a survey among clinicians associated the main barriers towards in silico medicine with computational resources and long computing times.In this context, cheaper, reduced CFD approaches, that e.g., make quasi-stationary assumptions [17], [18], may provide an adequate balance.A second possible explanation lies in the lack of verification and validation studies.Although verification and validation are essential steps prior to clinical translation, the availability of necessary high quality image-based flow data is limited.Echocardiography data can be used to validate valve velocities [19], [20], [21], [22], but does not provide all required flow field information for in-depth comparison.Phase-contrast MRI (PC-MRI), the gold standard for in vivo blood flow measurements, can provide 4D flow fields, allowing more detailed comparisons [2], [8], [21], [23], [24].Yet, PC-MRI acquisitions (also referred to as 4D flow MRI) are typically not part of clinical routine, making it challenging to obtain extensive high-quality data.Regarding the low availability of detailed verification and validation studies, more research in this context is required to gain trust of clinicians and identify optima in the trade-off between model complexity and accuracy.
In this study, we acquired an extensive MRI protocol of three healthy subjects.The protocol consists of steady state gradient echo (GE) MRI sequences for anatomical reconstruction and PC-MRI measurements for verification purposes.The geometric sequences are used as input for previously introduced CFD methods of different fidelity to compute intracardiac flow.The first method is based on a rigid wall modeling (RWM), where cardiac phases (systole and diastole) are computed separately, and the ventricular expansion is simplified to a velocity BC [17].The second method moreover models valve opening and closing processes and computes a coherent cardiac cycle being based on a prescribed-motion fluid-structure-interaction (PM-FSI) approach [25].For both methods, the intraventricular flow is computed and investigated regarding the extent to which the computed flow matches the 4D flow MRI measurements.From the verification, we aim to draw conclusions to the areas of application of both methods.I.e., where is the usage of a RWM justified and advantageous?Furthermore, we contextualize this work in state-of-the-art verification studies of intracardiac flow and give a glance on the challenges the field has to address in terms of verification and validation to target clinical translation.

II. METHODS
For the comparison of the two image-based CFD approaches with in vivo flow measurements, MRI data of three healthy subjects were acquired.Fig. 1 illustrates how this data were employed in the simulation and comparison workflow.The steady state GE sequences were used to reconstruct the anatomical shapes and to obtain the BCs for the two CFD approaches, RWM and PM-FSI.PC-MRI acquisitions, on the other hand, were used for comparison with the simulated flow after checking data quality and identifying relevant views and time points.

A. Medical Image Acquisition
Medical image data of three healthy subjects, whose clinical and demographic data are listed in Table I, were acquired in a 1.5 T MRI scanner (Ingenia Ambition X, Koninklijke Philips N.V., Amsterdam, Netherlands) with electrocardiogram-gating.Four steady state GE sequences were used for shape reconstruction (cine short axis, 3D steady-state free-precession (SSFP), radial mitral valve (MV), cine aortic valve (AV)), and the flow was evaluated based on two 4D flow MRI sequences (4D flow threechamber view (3CH), 4D flow LV).The characteristics of the sequences are shown in Table II.Further technical details can be found in the supplementary material and in the DICOM tags of the image data that can be provided upon request.For the steady state GE sequences, the respiratory motion was compensated via breath-hold or gating, whereas for the PC-MRI sequences, no respiratory motion compensation technique was applied.The studies involving human participants were reviewed and approved by the Ethics Committee of Charité -Universitätsmedizin Berlin (EA2/084/23) and performed according to the principles of the Declaration of Helsinki.All subjects signed a written informed consent for participation.

B. Computational Fluid Dynamics
Based on the steady state GE sequences, the LV, AO, LA, MV, and projected AV area were segmented and resembled to receive shape reconstructions representing the computational domain in the CFD.The reconstructed shapes are furthermore used to determine the BCs for RWM and PM-FSI, respectively.
1) Shape Reconstruction: Detailed shape reconstructions of the lumen were obtained by segmenting SSFP, cine short axis, and cine AV sequences.This was performed manually by the same user with the software ZIBAmira (version 2015.28,Zuse Institute Berlin (ZIB), Berlin, Germany).The radial MV sequence was segmented using an in-house tool developed for MV segmentation [26] following an earlier published description [27].The annuli in the SSFP sequence were also segmented with a developed in-house tool [28].In the segmentations, chordae tendineae, papillary muscles, and trabeculae carneae were not included.Fig. 2 illustrates the schematic procedure for segmentation and shape reconstruction.
From the SSFP sequence, LA, AO, and aortic annulus were segmented in end-diastolic state, smoothed, and used in the final surface reconstruction.For alignment purposes, the mitral annulus was also segmented in the SSFP sequence.The radial MV sequence was used to segment the MV in the most open Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.

TABLE II SPECIFICATIONS OF THE ACQUIRED MRI SEQUENCES
state, which was found to be at peak E-wave in all cases.The MV was subsequently smoothed.Furthermore, the LVOT movement was extracted landmark-based in a three-chamber view over the cardiac cycle.From the cine AV sequence, the projected orifice area of the AV in the most open state was segmented.The cine short axis sequence was used to segment the time-resolved LV shapes.Due to the slice thickness of the image data, post-processing was required, including nodal resampling, Poisson surface reconstruction, and smoothing in the open-source software MeshLab (version 2021.05,ISTI-CNR, Pisa/Rome, Italy).As the annulus region is simplified to be stiff in the CFD approaches, vertical alignment of the time-resolved LV segmentations was achieved by shifting along the long axis according to the LVOT movement identified in the radial MV sequence.Subsequently, the cine-based LV reconstructions as well as the segmented MV were aligned to the SSFP data via apex points of the end-diastolic states and the mitral annulus centers, respectively.Last, watertight LV surface meshes over the cardiac cycle are obtained as described in [25].
Based on these reconstructions, the volume curve is computed.It is spatially refined by interpolation, concatenated three times to enforce periodicity, and smoothed by removing high frequency Fourier modes.The final volume curves and the reconstructed LVs of all cases are displayed in Fig. 3.
2) Numerical Modeling: The CFD approaches, PM-FSI and RWM, are of different fidelity, thus requiring different input data.Details about both models are addressed in Sections II-B2a, and II-B2b.Numerically, both approaches solve the 3D incompressible Navier-Stokes equations for unsteady flow with the commercial flow solver STAR-CCM+ (version 2021.2.1, Siemens Industries Digital Software, Plano, TX, USA).Concomitant mass and momentum equations are shown in (1) and (2).
Therein, V denotes the control volume, Ω is the boundary of the control volume, n is the outwardly directed vector normal to dΩ, t is time, v is velocity, v g is the grid velocity, p is pressure, and ρ is density.Furthermore, I and T denote the unit tensor of second order and the viscous stress tensor, respectively.The grid velocity equals zero in the RWM approach, whereas the PM-FSI model contains a prescribed time-dependent motion inducing a mesh deformation in an arbitrary Lagrangian-Eulerian manner.Blood is modeled as non-Newtonian with a Carreau-Yasuda model using parameters according to [29]   1 mm with a refinement to 0.25 mm in the valve regions.Mesh independence studies with comparable resolution were performed for both approaches in previous studies [10], [25].Furthermore, an URANS turbulence model is applied, including solver settings for convergence tolerances of 1e-5 and time step sizes of 0.1 ms.All computations were performed on the Emmy system of the North-German Supercomputing Alliance The PM-FSI model computes the flow over the whole cardiac cycle by integrating the ventricular movement.Note, that only forces the boundary exerts on the fluid through its movement (but not vice versa) are considered.Unlike the LV, the LA, AO, and annuli remain stiff.At the pulmonary veins and the aortic outlet, generic pressure BCs are posed, and the valves are modeled via porous media theory.The PM-FSI approach allows computing cyclic flow patterns, blood washout, energetic characteristics, and deformation related phenomena.The approach is therefore suited for investigating complex cardiac flow behavior and the relation of biomarkers, rather than for being integrated into clinical routine processes.In-depths explanation and model discussion can be found in previous publications [7], [25].
In comparison to the previously published state of the modeling approach, the method was enhanced to incorporate a patient-specific volume curve, patient-specific LA and AO, as well as a patient-specific MV geometry in the opened state.The opened MV is also modeled via porous media theory: in diastole, it acts as non-permeable baffle interface and thereby as no-slip wall BC.When systole starts, it is set to a permeable internal interface behaving as a fluid region [30].In addition, the closed MV is represented by a porous plane that occludes the mitral annulus and switches between permeable and impermeable vice versa to the opened MV.The AV is modeled planar, taking the shape of the aortic annulus and scaling it to the patient-specific projected AV orifice area.Opening and closing of the valves is realized on the planes enclosing the annuli by implementing a pressure drop according to Darcy's law, as described in [25].The last model enhancement lies in the deviation of the deformation field from end-diastolic to end-systolic state.It is computed by mapping the end-systolic surface onto the end-diastolic one via Large Deformation Diffeomorphic Metric Mapping in Deformetrica [31].The established nodal correspondence is used to compute the deformation field.Subsequently, this deformation field is scaled in time to receive a time-resolved grid velocity following the specified volume curve (see [25]).Within this mapping process, Hausdorff distances from the mapped to the raw LV surface remain below 0.045 mm (mean) and 0.32 mm (max) in all cases.
b) Rigid Wall Model: The RWM builds on the approach as introduced in [17].It computes systole and diastole separately, bypassing necessity to model valve opening and closing.Furthermore, it simplifies the ventricular movement via a velocity outlet BC.The anatomical representation consists of the LA or AO at end-diastole (depending on the considered cardiac phase), the MV or AV in open state, and the LV at peak E-wave.At the pulmonary veins and the AO, the same generic pressure BCs as in PM-FSI are imposed.In comparison to PM-FSI, RWM is characterized by lower data requirements as well as pre-processing and solving times, which makes it more suitable for clinical use in terms of patient-specific therapy planning.It can be employed to compute early diastolic flow patterns around the MV [17], e.g., to estimate the outcome of interventional treatment for mitral regurgitation [10].
The vector field used for the velocity outlet is obtained by mapping the LV at peak E-wave onto the LV at end-diastole via Large Deformation Diffeomorphic Metric Mapping in Deformetrica [31].This vector field is then scaled in time such that outflow rates result in agreement with the patient-specific volume flow rates.Within the mapping process, Hausdorff distances from the mapped to the raw LV surface remain below 0.041 mm (mean) and 0.27 mm (max) in all cases.
In both approaches, the computed flow fields highly depend on the posed BCs, being a major challenge in the field.Therefore, a preliminary study on the sensitivity of BCs was performed exemplary for subject 1 and is provided as supplementary material.Therein, the impact of pulmonary vein inflow was neglected assuming the MV to have a regularizing effect on atrial inflow (see [8]).Subsequently, two known relevant BCs (MV orifice area and ejection fraction) were varied in ranges to be expected from segmentation and pre-processing.The response of PM-FSI and RWM to varied input parameters was alike.The velocities changed within the same order of magnitude as the input parameters, whereas relative pressures and kinetic energy responded disproportionately stronger with factors around 1.5-3.

C. 4D Flow MRI-Based Verification
To verify the CFD results, two PC-MRI sequences were acquired.The first one captures the whole LV cavity at a velocity-encoding (VENC) of 250 cm/s (4D flow LV).From the 4D flow LV sequence, the volume flow rates over MV and AV were measured in GTFlow (version 3.2.16,GyroTools LLC, Zurich, Switzerland) and used to identify the peak phases (peak systole, peak E-wave, peak A-wave).As a second PC-MRI sequence, a three-chamber view was acquired at VENCs of 50, 100, 150, and 250 cm/s (4D flow 3CH).Compared to the 4D flow LV sequence, it offers a higher spatial resolution (see Table II) and is performed in 5-15% of the time, facilitating an acquisition at different VENCs.Different VENCs, in turn, allow optimizing the signal-to-noise ratio at specific cardiac phases while preventing velocity aliasing.
To verify the simulated flow fields, focus is laid on the intraventricular kinematics.The post-processing is performed based on the derivable quantities from the 4D flow data.The three-chamber view plane (see Fig. 6) from the 4D flow 3CH sequence is equally placed in the CFD for quantitative comparison.To measure valve velocities on that plane, spherical probes are placed 3 mm below the MV leaflet and the AV center, respectively (see Fig. 3).In the CFD, the spherical probes enclose a radius of 2.5 mm, whereas in the 4D flow 3CH data, the probes are circular, with a radius of 2.5 mm.From the valve velocities, the velocity time integrals (VTIs) are computed via temporal integration of the respective valve velocity magnitude at the respective phase.The VTI is a clinically common measure, which can, e.g., be used to estimate the stroke volume [32].Furthermore, two line probes are positioned on the three-chamber view plane (see Fig. 3).The first expands along the LV long axis from the MV leaflet center towards the apex and the second along the LV short axis in the basal LV region from the LVOT to the posterior wall.Along these line probes, space-time characteristics of the velocity magnitude are investigated (see Fig. 5).The three-chamber view plane is furthermore used to qualitatively compare the velocity magnitude and the out-of-plane vorticity (see Fig. 6).The out-of-plane vorticity ω z is computed via spatial deviation of the velocity v according to (3).
The large-scale flow patterns inside the entire LV cavity are compared using streamlines (see Fig. 7).The computation of streamlines is performed equally for the 4D flow LV and the CFD data in ZIBAmira.For a regional analysis focusing on low flow regimes, five circular probes of 2.5 mm radius are placed in the three-chamber view planes close to the endocardium (see Fig. 8).The probes are moved in accordance with the ventricular movement in 4D flow measurements and PM-FSI.
The kinetic energy is computed for several points in time according to (4) with a density ρ of 1050 kgm −3 , the velocity magnitude v, and the volume V of the voxels or the computational grid for 4D flow MRI and CFD, respectively.
As last measure, relative pressures are considered (AV: MV , apex to base Δp ab = p apex − p base ), which are clinically of interest in the diagnosis of e.g., valve diseases.They are computed for PM-FSI and RWM and compared to literature values.The pressures are measured using spherical probes with a radius of 2.5 mm.For the transvalvular pressures, the probes are positioned 3 mm upstream and downstream of the respective valve.

A. Flow Rates and Valve Velocities
To analyze the agreement of volume change imposed as BC in the CFD to the 4D flow MRI measurements, the volume flow rates are compared.The qualitative temporal course matches in Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.all cases (Fig. 4).The differences in mean and max volume flow rates in systole and diastole range from 0-23% (Table III).In all cases, higher flow rates appear in A-wave in the 4D flow data.Furthermore, the peak phases appear partially shifted in time (2-60 ms).The valve velocities are displayed in Fig. 4. Qualitatively, matching trends are visible with peaks in systole, E-wave, and A-wave.Corresponding to the volume flow rates, the peak velocities are slightly shifted in time.For RWM, the A-wave peaks appear earlier at higher velocity.Comparing peak systolic

TABLE III VOLUME FLOW RATES OF 4D FLOW MRI AND RELATIVE DIFFERENCES (Δ) OF CFD COMPUTATIONS TO 4D FLOW MEASUREMENTS
and peak E-wave velocities of both CFD methods to 4D flow MRI measurements, no differences above 0.15 m/s appear.The relative differences in peak valve velocities range from 0-21% (Table III).Comparing the VTIs of both CFD methods to the 4D flow measurements, no differences above 10% are observed in subjects 1 and 2. However, for the VTI of subject 3, higher deviations occur, with a maximum of 40% for the diastolic VTI in RWM (Table III).

B. Space-Time Characteristic
To investigate spatio-temporal patterns of the intraventricular flow, space-time characteristics along line probes are analyzed, as exemplary shown for subject 1 in Fig. 5.The respective plots of subjects 2 and 3 are provided in the supplemental material.First, the shape of the characteristics demonstrates the differences in the spatial and temporal extent of the individual methods.While for subject 1, the cardiac cycle of the 4D flow MRI data is slightly longer, the ventricular contraction becomes visible in the characteristic of the short axis line probe.For PM-FSI, no velocities are existent at septal and lateral line probe regions as the LV contracts and the respective line probe elements are positioned outside the computational domain.In RWM, the LV is not moving and represented in the peak diastolic state, whereas in 4D flow MRI, noise is measured outside the LV.Yet, the LV contraction can still be traced in the 4D flow MRI data, e.g., by high velocity regions in systole at the LVOT.This high velocity region is also visible in similar shape and velocity in PM-FSI and RWM short axis characteristics.In the long axis line probes, higher velocities are visible at peak systole in RWM for all subjects.In diastole, the short axis characteristic captures the E-wave jet in the basal LV region.For all three cases, the width, positioning, and magnitude of the jet is comparable to the 4D flow MRI data.For subject 2, the jet is more laterally oriented and narrow than in the 4D flow measurements.In early diastole, PM-FSI and RWM characteristics are very similar, both indicating a rolling up of vortices beside the jet.After peak E-wave, higher velocities develop at the wall in RWM.After the rapid inflow has slowed in mid diastole (diastasis) of subject 1, the 4D flow MRI data indicates a partial closing of the MV, with the inflowing velocities narrowing more laterally sided towards the posterior leaflet.This partial closing is also observed in the radial MV sequence.Regarding the long axis characteristic at diastole, a notably faster penetration of the E-wave jet towards the apex in the 4D flow MRI data compared to the CFD approaches is visible in all cases.

C. Large-Scale Flow Patterns
Figs. 6 and 7 visualize the large-scale flow patterns of subject 1.The respective visualizations for subjects 2 and 3 are provided in the supplemental material.Comparing PM-FSI and RWM, only minor differences in the flow field at peak E-wave are present for all subjects.The inflowing jet is of similar shape and magnitude.In RWM, a more coherent, undisturbed ring vortex formation as well as more unidirectional flow in mid and apical LV regions is observed (Figs. 6 bottom and 7 bottom).Regarding the 4D flow data, slight differences in the MV opening state compared to the MV shape incorporated in the CFD approaches become visible (Fig. 6 top).Yet, the shape and magnitude of the jet are comparable to both CFD approaches (Figs. 6 top and 7 bottom).For subject 2, the E-wave jet is more laterally oriented and penetrates deeper into the apex in the 4D flow data.The formation of ring vortices is also indicated in the 4D flow data (Fig. 7 bottom) in all subjects.In the vorticity visualizations (Fig. 6 bottom), a shear layer at the MV leaflets becomes visible.However, a distinct rolling-up of vortices downstream of the leaflets is only observed in the CFD results.In systole (Fig. 7 top) all blood rushes towards the LVOT.Whereas the outflow velocities are well comparable in the aortic root, the PM-FSI Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.and RWM velocities in the ascending AO are higher compared to the 4D flow MRI measurements.

D. Regional Flow Analysis
The regional flow analysis reveals larger deviations between PM-FSI and RWM.Fig. 8 exemplary shows boundary-near velocities for subject 1.For all subjects, the largest differences between PM-FSI and RWM are visible in basal regions (p1 and p4), whereas velocities become more similar towards the apex (for subjects 2 and 3 see supplementary material).In the comparison to 4D flow measurements, some matching trends are visible (peaking velocities in peak phases).However, larger deviations appear in all cases (e.g., p1 and p5 in subject 1).

E. Flow Energetics
In both CFD approaches, the same trends of kinetic energy are observed (beside higher values in PM-FSI in early systole and late diastole) (Fig. 9).Comparing the 4D flow MRI data to PM-FSI and RWM, the trends of local maxima at peak systole, peak E-, and A-wave and a local minimum in diastasis appear in all three methods.In subject 2, however, the local minimum in diastasis and local maximum in A-wave are not visible in the CFD results.Yet, it can also be observed, that the reconstructed volume flow rate in diastasis of subject 2 stagnates less than measured in the 4D flow data (Fig. 4(b)).Regarding the kinetic energy magnitude, the 4D flow data matches the RWM results for subject 2. In diastole of subject 1 and in the entire cycle of subject 3 the 4D flow kinetic energy is higher than in both CFD methods, while being closer to PM-FSI.

F. Relative Pressure Differences
Fig. 10 displays relative pressures over the AV, the MV, and from apex to base.Comparing PM-FSI and RWM, only small differences appear.For the subjects, temporal course and magnitude are similar, as well, illustrating the acceleration and deceleration of flow over the cardiac cycle.At the onset of systole, blood is moved towards the LVOT (positive pressure difference from apex to base).After peak systole, the flow is decelerated until the onset of diastole (negative pressure difference from apex to base).Similar behavior is found in diastole: at the onset, blood is accelerated from base to apex (negative pressure difference from apex to base), decelerates after peak E-wave, and accelerates again in A-wave.Corresponding behavior is found at the AV in systole and the MV at diastole, respectively.

IV. DISCUSSION
In this study, we verified in silico computed intracardiac blood flow of two CFD approaches by comparison to PC-MRI measurements.In the following section, we discuss the processing of the medical image data, analyze similarities and differences of the CFD results compared to PC-MRI, and contextualize this work in terms of verification and validation of in silico models.

A. Medical Image Processing
The steady state GE sequences were chosen, segmented, and processed to obtain the best possible reconstruction of the anatomy.To do so, four sequences were used, including dedicated sequences for the heart valves.The sequences are performed quickly, and the registration can be largely automated.When drawing comparison to other studies, segmenting the LV based on cine short axis images and the MV via rotations of the long axis view passing through apex and the MV center are found to be common strategies [3], [33], [34].Comparing the shape reconstructions to PC-MRI measurements yields matching results in terms of volume flow rates.The general trend as well as mean and max difference per phase (see Fig. 4, Table III) are consistent.The ventricular contraction was also found to be similar in the space-time-characteristics (see Fig. 5).Regarding the opened MV shape, it is in reasonable agreement with the inflow jet profile (Fig. 6).The partial closing of the MV in diastasis cannot be modeled with a stiff MV.In total, these findings verify our reconstruction procedure to provide viable BCs for the CFD in terms of flow rates and anatomical representation.
To be able to verify the computed flow fields of the CFD, the PC-MRI images were acquired and processed to provide a reliable data basis for comparison.Although a VENC of 250 cm/s in the volumetric LV sequence is rather high in relation to diastolic velocities, the large-scale flow structures could still be resolved.Other studies using PC-MRI to verify CFD results use lower VENCs ( [2]: 120 m/s, [21], [23], [24]: 150 cm/s).However, considering the long acquisition time of the sequences, and the aspect that aortic velocities can easily reach 200 cm/s, we consider a VENC of 250 cm/s to be appropriate for universal application in healthy subjects without velocity aliasing.For an optimal resolution of low flow domains, three-chamber views with different VENCs are faster to acquire and facilitate the measurements.In terms of post-processing of the PC-MRI acquisitions, no actions were taken, as that would have required verification itself.This is common for verification studies of CFD of intraventricular flow using PC-MRI [2], [21], [23], [24].Yet, in a recent study investigating a rigid flow phantom, significantly improved correlations between CFD and 4D flow measurements were found when post-processing the 4D flow data [35].Such post-processing may improve noisy vorticity patterns (Fig. 6 bottom) by, e.g., introducing a divergence-free condition.We consider this an option for future investigations.As in [24], the PC-MRI sequences were acquired at free breathing and no respiratory motion compensation.These settings were chosen as they are applied in clinical routine and motivated by significantly quicker acquisition times while measuring non-significant differences in flow rates compared to 2D or 4D gated acquisitions [36].To summarize, we conclude the acquired PC-MRI sequences to provide a reliable data basis for verification of the CFD results.

B. Verification of CFD
Both CFD approaches are intended to optimize usability, clinical feasibility, and generalizability.The required data basis is limited to four MRI acquisitions, which are partly included in clinical routine and have short acquisition times.After the segmentation, the pre-processing can be performed in 4-6 h, having the potential to be further automated.Additionally, the approaches are modality-independent and can also be applied on computed tomography (CT) data [25], echocardiography data [10], or synthetic data via statistical shape modeling [7].As a consequence, the approaches are positioned at the lower end of modeling complexity.No two-way coupled FSI, solid mechanics, or electrophysiology is included.However, the results show a good agreement with the 4D flow measurements.
Regarding the MV velocity, PM-FSI and RWM results are well in line with 4D flow measurements.The VTI, another clinically used measure, is also well reproduced for subjects 1 and 2, whereas for subject 3, especially the deviations in RWM are larger.This may be caused by heart rate variability, flow rate deviations, or the spatial evaluation area of the valve velocities.In terms of large-scale flow structures, during diastole, the peak E-wave jet in the basal region is well comparable in all cases with regard to velocities for both CFD approaches and the 4D flow measurements (Figs. 5 top, 6 top, and 7 bottom).While for subject 1 and 3 the flow direction of the peak E-wave jet aligns well between CFD and 4D flow measurements, it orientates more laterally in the 4D flow data of subject 2. This may be influenced by biases in the ventricular deformation field or the representation of the MV, as investigated by Vellguth et al. [37].Towards late diastole, RWM results show higher MV velocities and high velocities in the near wall basal regions, due to a missing no-slip wall BC.The systolic ejections in PM-FSI and RWM are in line with 4D flow measurements in terms of AV velocities (Fig. 4) and large-scale flow structures (Fig. 7).However, towards the ascending AO, the velocities exceed the 4D flow measurements.This may be influenced by the AO being segmented at enddiastole, although the AO enlarges in systole [38].As in our models the heart is regarded as an isolated organ, attaching, e.g., a Windkessel model at the aortic outlet may furthermore increase the accuracy.A regional flow analysis reveals more deviations to 4D flow measurements in boundary-near areas.This can be explained by cyclic variations in ventricular flow [39] alongside measurement errors in 4D flow acquisitions in case trabeculae carneae or papillary muscles are in the measurement area.However, PM-FSI and RWM are both methodologically unable to reproduce all details of boundary-near flow.For PM-FSI, this can be overcome by increasing the level of anatomical details in the representation of the ventricular wall.For RWM, the rigid consideration of the ventricular wall interferes with detailed flow computations near the boundaries.In the investigation of kinetic energy, it was found that RWM and PM-FSI take similar temporal courses compared to 4D flow data, while showing differences in the order of magnitude.These differences are not surprising, considering noise in the MRI data and lower spatio-temporal resolution of the PC-MRI data.Thereby, vortical structures that are below the spatio-temporal resolution cannot be resolved.As a consequence, less of the turbulence spectrum is resolved.This makes direct magnitude comparisons limited a priori.Regarding relative pressure differences, no direct comparison of CFD to 4D flow MRI measurements was performed.However, temporal courses and magnitudes of the computed relative pressures are in line with literature [40].Absolute pressure values cannot be computed due to the non-patient specific pressure BCs, marking a methodological limitation of the present study.
In total, RWM is able to adequately reproduce systolic and diastolic flow without large discrepancies to PM-FSI and has the advantages of quicker computation and less pre-processing (see Table IV).Assuming linear scaling in computational expenses (which is a conservative estimate in this consideration), the pre-processing and computation of RWM is feasible in one day with a potent workstation.Consequently, RWM provides a Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.
reasonable trade-off between accuracy and feasibility in clinical application regarding specific problems.These may, e.g., be the prediction of diastolic inflow after MV treatment or stress testing [10].PM-FSI, on the other hand, can reproduce the large-scale flow patterns of 4D flow MRI data at adequate accuracy over a coherent cardiac cycle and prospectively allows higher level of details in ventricular representations for enhanced near wall flow modeling.In contrast to RWM, PM-FSI further enables to compute additional biomarkers like washout [7], [25].In comparison to 4D flow MRI both CFD approaches allow computations at higher spatio-temporal resolution and to reduce additional patient scan times.Further model enhancements, like a moving MV, may additionally improve the accuracy of the results.Yet, rising model complexity always goes in hand with more data requirements and computational as well as labor power expenses, which always comes with a cost-benefit analysis.Likewise, other simplifications, which are to name as limitations, must be critically discussed.These are the negligence of anatomical details like papillary muscles, trabeculae carneae, chordae tendineae, and the twisting LV motion.Additionally, AO, LA, and the annuli are simplified to be stiff, and the LV motion in PM-FSI is derived based on a two-state-based motion approach.For a detailed discussion of these aspects, refer to Obermeier et al. [25].

C. Verification and Validation
Verification and validation are essential steps towards clinical translation of in silico models.For diagnosis and treatment planning, valve velocities are often measured, for example to identify valve stenosis or regurgitation.From our point of view, they should be included in all verification and validation studies considering intracardiac flow.Especially in studies based on echocardiography, a valve velocity verification is usually contained [19], [20], [21], [22], as they are easily accessible via Doppler measurements.For more detailed flow comparison, PC-MRI data is required.Caballero et al. [21] and Schenkel et al. [23] performed comparisons via the velocity field in a threechamber view.In their studies, discrepancies in the velocity magnitudes in valve regions were found, which was explained by the missing 3D MV modeling [21].Le et al. [24] compared one healthy subject to 4D flow MRI measurements, in a setup without 3D MV, LA, AO and considering diastole only.Additionally, the PC-MRI flow profile was incorporated as BC, inducing a dependence of CFD results to 4D flow data.This contradicts an independent verification.To the best of our knowledge, the most detailed comparison to 4D flow MRI data was done in the study by Lantz et al. [2], where a detailed CFD model based on CT data including valves and papillary muscles was verified based on 12 patients.Whereas their computations correlate well in systole regarding flow rate and integrated kinetic energy, discrepancies increase in diastole.To note is furthermore the comparison of different flow components (i.e., direct flow, retained inflow), which is valuable when aiming to validate washout processes.To contextualize our study, we performed an in-depth verification of our CFD approaches, considering various flow measures compared to the above-mentioned studies.We analyzed basic flow measures, such as volume flow rates, valve velocities, and flow fields at peak phases, while also investigating VTI, vorticity, the three-dimensional flow fields via streamlines, regional flow, relative pressures, and flow energetics.Especially the space-time characteristics proved useful to identify spatio-temporal flow characteristics in a reduced manner.Further analysis via, e.g., principal orthogonal decomposition may further help to pin down basic flow structures and allow for quantitative comparison.
Although verification studies can estimate the validity of modeling strategies, clinical translations require validation.Therefore, a look at the ASME V&V 20 [41] provides guidance.First, as the ASME V&V 20 [41] states, there is no validation without experimental data.From our point of view, only 4D flow MRI is currently able to provide the required in vivo data.To then validate a model, two aspects need to be addressed: uncertainty and error.The first can also be referred to as robustness, as done by Larsson et al. [19].To quantify uncertainty, it is first necessary to standardize the imaging and reconstruction procedures.Subsequently, the sensitivity of flow measures on input parameters must be investigated.From our perspective, the most important input parameters are the volume flow rates, the valve shape and positioning, as well as the ventricular movement.While we performed a preliminary analysis on the influence of boundary conditions on ventricular flow, more detailed sensitivity analyses are essential before clinical translations.The second central aspect of validation is error quantification.While considering various flow measures to get a holistic view on similarities and differences is useful, recommendations for specific measures that are suited to validate a model are lacking in literature.The same goes for error margins for these measures, which a model must ensure.To identify such measures, it is useful to distinguish between clinical measures and research measures to achieve a better understanding of LV functioning.Beside these considerations from a modeling perspective, the resource demanding continuum mechanical simulations must also be feasible in clinical reality.This means the availability of powerful workstations or easy access to high performance computing.In addition, expert knowledge on sight is crucial for the execution and interpretation of computations.
To summarize, only a few studies covering aspects of verification and validation of intracardiac flow computations exist.We performed an in-depth comparative study that partially goes beyond other verification studies.However, to target actual clinical translations, more studies, broader data basis as well as systematic quantification of uncertainty and error are challenges the research field has to address.

V. CONCLUSION
In this study, we acquired steady state GE and PC-MRI sequences that lead to accurate shape reconstructions and provide a reasonable quality for flow quantification.Therewith, PC-MRI flow measurements were compared to two previously introduced CFD approaches of different model fidelity to compute intracardiac flow.Both CFD approaches are intended to optimize usability, clinical feasibility, and generalizability.Despite concomitant modeling simplifications, both methods could be verified within certain areas of application.They yield good agreement with 4D flow measurements regarding basic flow measures like valve velocities, kinetic energy, and large-scale flow patterns.In considerations where these quantities are of interest, RWM poses a fair balance between accuracy and complexity.The approach reaches its methodological limitations when flow measures like boundary-near flow, energetic efficiency or washout are of interest.For such investigations the modeling of valve opening and closing as well as integration of the ventricular movement are required.These aspects are covered in PM-FSI.Yet, more detailed anatomical models are required in PM-FSI for accurate computations of e.g., boundary-near flow.In the verification protocol, our investigations expand on previous studies.However, in future studies, broader database and systematic investigation of uncertainties and errors are crucial when targeting clinical translations.

VI. SUPPLEMENTAL DATA
Supplementary Material provides the preliminary sensitivity analysis, the respective visualizations of subjects 2 and 3, and technical details on the MRI acquisitions.All data necessary to reproduce the simulations, including all image data as well as reconstructed geometries and boundary conditions, can be found as open data via: https://figshare.com/s/27afb22c351dd779bace.

VII. AUTHOR CONTRIBUTIONS
LO, LG, KV: conceptualization.LO, FH, CM, KV: formal analysis and methodology.TK, LG, KV: funding acquisition and supervision.LO, MW: visualization.LO, MW, KV: investigation.LO, KV: writing-original draft.All authors: writing-review and editing.All authors contributed to the article and approved the submitted version.

Fig. 2 .
Fig. 2. Schematic illustration of the segmentation and shape reconstruction based on the four steady state GE sequences of subject 3.

Fig. 3 .
Fig. 3. Reconstructed end-diastolic (gray), and end-systolic (red) shapes of subjects 1-3 ((a)-(c)).The MV is displayed in green.The black dots and white lines indicate spherical probes at the valves, and line probes through the LV, respectively.(d)-(f) shows the reconstructed volume curves and volume flow rates.

Fig. 4 .
Fig. 4. Comparison of volume flow rates (a), and valve velocities (b) for 4D flow MRI measurements and CFD results for the three healthy subjects.The valve velocities are only displayed at points in time when the respective valve is opened.

Fig. 5 .
Fig. 5. Space-time characteristics of short axis (top), and long axis (bottom) line probes for subject 1.The plots on the right side display the velocity magnitudes at the peak phases (indicated with the black dotted line) of the respective method along the respective line probe.

Fig. 6 .
Fig. 6.Three-chamber view of subject 1, showing the velocity magnitude at peak E-wave (top), as well as the out-of-plane vorticity at peak diastole (bottom) for 4D flow MRI (left), PM-FSI (mid), and RWM (right).

Fig. 9 .
Fig. 9. Kinetic energy of all subjects.Three time points around peak systole, four around peak E-wave, two in diastasis, and three time points in A-wave were evaluated for the 4D flow data.

Fig. 10 .
Fig. 10.Relative pressure drop over AV and MV at the respective cardiac phase as well as intraventricular pressure difference from apex to base of subjects 1-3 (a)-(c) for PM-FSI and RWM.

TABLE IV ADVANTAGES
AND DISADVANTAGES OF THE METHODS