Computational Modeling of Aortic Stenosis With a Reduced Degree-of-Freedom Fluid-Structure Interaction Valve Model

An aortic stenosis refers to the abnormal narrowing of the aortic valve and is usually a result of the incomplete opening of valve leaflets [1]. This reduced valve area modifies the aortic jet and the hemodynamics in the ascending aorta, which, in turn, can generate abnormal heart sound known as bruits or murmurs [2]. The murmurs associated with aortic stenosis are thought to contain valuable information that can help detection and diagnosis [3], and many researchers have attempted to develop auscultation-based noninvasive diagnostic modalities for this and other related conditions [4–9].

The poststenotic flow that is responsible for the generation of these murmurs has been studied extensively through numerical simulations [10–12] and experiments [13–15]. In previous studies of hemodynamics across aortic stenosis, the stenosis is commonly modeled as a stationary constriction [16,17]. While such simplified models might provide some basic information about the flow, they do not account for the leaflet shapes as well as the dynamic movement of the valve leaflets, both of which, can have a significant effect on the flow patterns in the postvalvular region [18].

The simulation of the fluid-structure interaction (FSI) of the aortic valve (either natural or bioprosthetic valves) is a challenging task, and a variety of approaches have been used to model this phenomenon [18]. Given that the valves undergo large deformation during the cardiac cycle, methods using body-fitted meshes (e.g., arbitrary Lagrangian-Eulerian method [19]) struggle to maintain mesh quality. Therefore, the problem is usually solved by the immersed boundary (IB) methods [20]. Indeed, the IB method was originally developed by Peskin [21] to simulate the flow around heart valve and this method has become the defacto standard for modeling such flows [18,22]. More recently, the Sub-Grid Geometry Resolution (SGGR) method, which is closely related to the cut-cell method, has been used to investigate the FSI of prosthetic valves [23,24].

When it comes to the modeling of the FSI implicated in the aortic stenosis, another challenge lies in the structural dynamics modeling of the aortic valve. Microscopic anatomic study reveals the complex heterogeneous fibrous structure of natural aortic valve leaflets [25,26]. The literature is quite rich in this area and interested readers can refer to the reviews by Peskin [27] and Sotiropoulos et al. [18]. Depending on the purpose, models with different complexities are used to simulate the aortic valve dynamics. Relatively complex models are used when valve dynamics are a focus of the study. For example, Griffith [28] expanded the IB methods developed by Peskin to study the aortic valve dynamics under physiological driving and loading conditions wherein the aortic valve was modeled as immersed boundaries comprised of systems of elastic fibers; Sacks and colleagues [29,30] treated the leaflets as anisotropic, nonlinear hyperelastic materials to study the failure of the bioprosthetic aortic valves, while the effect of the fluid was simply modeled as a transvalvular pressure difference. Though complex material models can provide accurate prediction of the valve dynamics, accurate determination of model parameters is a highly challenging task. The Fung-type hyperelastic material model employed by Sacks and colleagues [29,30] contains 7 model parameters and requires substantial experimental effort to determine. Also, in addition to the fluid solver, a separate solver is required to solve the structure dynamics of the valve [31], increasing the difficulty of implementation. On the other hand, in problems where fluid dynamics is the interest of study, a simpler valve model can be employed. For instance, Domenichini and Pedrizzetti [32] developed a single degree-of-freedom mitral valve model whose motion was determined by a momentum balance, and the valve leaflets had no inertia or stiffness. Tagliabue et al. [33] modeled heart valves as mixed time-varying boundary conditions to resolve the fluid dynamics in the left ventricle efficiently. Laadhari and Székely [34] described the valves through level set functions, and their motions were determined by a reduced-order model governing the rotating angle of valve leaflets. The above simple models are, however, not suitable for the modeling of aortic stenosis, since the stiffness of valve leaflets cannot be controlled easily. Our objective here is to employ a simplified valve model that provides sufficient realism in terms of the effect of the valve leaflets on the aortic hemodynamics and allows us to model the effect of increasing aortic stenosis on the postvalvular flow and pressure. In the present study, we propose a novel reduced degree-of-freedom (rDOF) aortic valve model and apply it with the IB method to investigate aortic stenosis.

The rDOF valve model uses a prescribed mode-shape for the movement of the leaflets, but the instantaneous acceleration of the leaflets is determined by the equation of motion with the integrated hydrodynamic force on the leaflets as well as the inertia and the stiffness of the valve leaflets. By prescribing the mode-shape, the equation of motion is reduced to an ordinary differential equation. This simplification allows us to examine the hemodynamics of healthy as well as stenosed valves in an efficient manner.

The geometric model employed in this study consists of a canonical aorta model and the aortic valve, as shown in Fig. 1(a). The aorta is modeled as an inflexible straight tube that can be divided into three sections: the left-ventricular outflow tract, the aortic root, and the ascending aorta. Similar straight tube models are also used in other experimental and numerical studies [28,31,35–37]. The rigid tube assumption is considered reasonable in this case since results produced with such assumption agreed well with experimental studies in our previous investigation [16]. The diameter of the inflow duct (D_o) is set to $2.5cm$⁠, which is within the physiological range for both healthy subjects and patients with stenosis [38]. The dimensions of the aortic root are shown in Figs. 1(b) and 1(c), and they are based on the generic model proposed by Reul et al. [39].

The blood inside large blood vessels such as the aorta can be considered Newtonian and is governed by the Navier–Stokes equation [40]. No-slip boundary condition is enforced on the vessel wall. The flow inside the modeled aorta is driven by a pulsatile inflow. The inlet velocity profile, whose expression can be found in Appendix  A, is plotted in Fig. 2, and a uniform velocity is prescribed at the inlet. This pulsatile inflow generates a stroke volume of 70 ml. The systolic phase lasts $0.33s$⁠, and the duration of one cardiac cycle (T) is $0.9375s$ which is equivalent of 64BPM. The maximum Reynolds number (⁠$Remax=Vin,maxDo/ν$⁠) is around 6000, where ν is the kinematic viscosity of the blood and $Vin,max$ is the maximum velocity in Fig. 2. The Womersley number (⁠$Wo=Do/22π/(νT)$⁠) is around 17.73. This peak Reynolds number and the Womersley number are within the typical range for aortic flows [18]. To determine the absolute pressure inside the aorta, the after-load model proposed by Watanabe et al. [41] is adopted in the current study. The formula is given in Appendix  A, and the pressure profile in Fig. 3 is prescribed at the outlet of the modeled aorta in all of the simulations.

The complex geometry is treated by the sharp-interface IB method developed by Mittal et al. [42], and a cut-cell-based approach is used to improve the mass conservation of this IB method [43]. A Cartesian grid with $154×160×256$ (6.3 million) grid points is used to discretize the computational domain. The grid is distributed uniformly in the cross section, and clustered around the valve region in the streamwise direction. This results in around 128 grid points across the diameter of the ascending aorta. To ensure that sufficient resolution is provided to the flow, a grid convergence study is conducted and is summarized in Appendix  B. It demonstrates that this grid is capable of adequately resolving all the cases studied here. In order to reduce computational cost, a semi-adaptive time-step size is used in all the simulations. In each cardiac cycle, a small time-step size (⁠$1.25×10−4s$⁠) is used to resolve $t∈[0T,0.46T]$⁠. This time period includes the systole and part of diastole and is the most dynamic phase of the cardiac cycle. The rest of the diastolic phase is resolved by a larger time-step size (⁠$3.75×10−4s$⁠). Simulations are carried out on the Maryland Advanced Research Computing Center supercomputer with 288 cores, and it takes around 26 h of wall-time to compute 1 cardiac cycle for the most time-consuming case. 7 cardiac cycles are computed for each case, and data from the last 6 cycles are collected and analyzed.

In the present study, we propose a simplified valve model with reduced degrees-of-freedom (rDOF) and couple it via a FSI methodology to examine the dominant (“first-order”) effects associated with the dynamics of the valve leaflets.

Here, α is a constant related to mass, which includes both the real mass and the added mass; κ is a stiffness constant; $vv$ and $dv$ are the velocity and displacement of the element; $ΔP$ is the pressure difference across the surface element; and $n$ is the unit surface normal vector of the element. The first term on the right-hand side of the equation represents the pressure force from the flow, while the second term models the elastic recoil forces on the leaflet. It should be noted that α and κ are model parameters and need to be determined.

where c₀ is the initial condition of c. Here, the initial conditions of c(t) and a(t) are all set to zero, which represents zero displacements and velocity. By imposing the mode of the valve, the partial differential equation (Eq. (1)) describing the motion of each element is reduced to a set of ordinary differential equations (Eqs. (5), (6)) that describes the collective motion of elements in each leaflet and can be solved easily. Consequently, in order to employ the proposed model, knowledge of the valve shape at fully open and closed positions is required. It is worth emphasizing that the solution procedure is carried out individually for each valve leaflet, and contacts between the leaflets are prohibited by prescribing the fully closed position.

Equations (5), (6) are coupled with the Navier–Stokes equations resulting in the FSI model. Here, an explicit coupling approach is adopted, in which the flow and the valve dynamics are solved sequentially and only once in each time-step. The flow is first solved with leaflets as internal boundaries, and then the hydrodynamic force is used in modeling the dynamics of the valve. The valve position is updated at the end of each time-step. The aforementioned small time-step size ensures the stability of this explicit coupling.

As stated before, α and κ are model parameters that need to be prescribed. The parameter α is related to the area density of the leaflets. Based on our initial parametric tests, even $O(10)kg/m2$ changes in α change valve opening and closing time by less than 10%. We, therefore, set $α=40kg/m2$ which provides stable coupling with the flow solver based on extensive numerical testing. The parameter κ determines the recoil force and this parameter has a more significant effect on the opening/closing behavior of the leaflets. To validate our model and to determine this parameter, we compare the current rDOF FSI model against a more detailed FSI model of the same valve developed by Tullio et al. [31]. In this other model, the valve leaflets are treated as a nonlinear anisotropic material described by a spring-network model, and they achieved a great agreement with the experimental study. The inlet velocity used in the test is plotted in Fig. 5. This velocity profile is the same as the one used by Tullio et al. [31], except that the negative flowrate near the end of the systole is excluded here to avoid artificial regurgitation. Figure 6 plots the normalized projected valve area (PVA) calculated over one cardiac cycle with $κ=10,600Pa/m$⁠. The result from Tullio et al. [31] is also plotted here for comparison. As can be seen in this figure, the rDOF valve model shows a rapid opening behavior similar to the more complex model, while the closing phase is not as smooth around $0.2s$ potentially due to the artificial stiffness introduced by the reduced-order motion. Despite these discrepancies, this figure clearly demonstrates that this simplified valve model can capture the first-order behavior of the more complex valve model.

Snapshots of the valve at different time instances are shown in Fig. 7(a). The phases at which these snapshots are taken are highlighted in Fig. 6. The valve is rapidly opening at phase I. Compared with the snapshot in Ref. [31] at a similar phase, the current valve profile is smoother and lacks small structures near the tip of the leaflets. Phase II and III correspond to the fully open position and early closing phase, respectively. The envelope of the valve motion is shown in Fig. 7(b), and we can see that the leaflet follows the same mode shape during the closing and opening phase due to the adoption of same trajectory $b$ in both phases.

The key factor in heart murmur generation is the poststenotic jet, whose formation is dictated by the properties of the valve opening (e.g., shape and area) [17]. Based on the above comparison, it is reasonable to conclude that the rDOF valve model is expected to generate a flow field that is close to the valve model which uses more complex constitutive relations and is sufficient for studying murmur generation mechanism. The stiffness parameter $κ=10,600Pa/m$ will be used to generate the baseline case with healthy valve motion. The incomplete opening of the valve, i.e., aortic stenosis, can be easily introduced by increasing the parameter κ.

In this section, results from three cases (a healthy case and two stenosed cases) are presented. To ensure the clarity of the discussion, some data reduction methods are first introduced here.

Aortic stenoses are the result of the incomplete opening of the aortic valve. One of the parameters that measure the severity of the stenosis is the area stenosis ratio (AS), which is defined as the percentage of area that is blocked by the valve at maximum opening phase. It is related to the PVA through $AS=1−PVAmax/PVA0$⁠, where PVA_max is the maximum opening area of the diseased valve and PVA₀ is the maximum opening area of the healthy valve. In the current rDOF valve model, the incomplete opening of the valve can be achieved through increasing the stiffness parameter κ. This models the change of material property in valves due to calcification and/or thrombosis [45], which are major causes for aortic stenosis [1]. Depending on whether the same κ is used for all three leaflets, rotationally symmetric or asymmetric aortic stenosis can be created easily. Only rotationally symmetric stenoses are considered in the current study, i.e., all three leaflets have the same stiffness constant. Here, three values of $κ=10,600Pa/m, 15,900Pa/m$⁠, and $42,400Pa/m$ are used, and they lead to a healthy case, and two stenosed cases with AS of 31% and 47%. These three cases correspond to maximum valve opening areas of $3.71cm2, 2.56cm2$⁠, and $1.97cm2$⁠, respectively. Based on these areas, the latter two cases are clinically classified as mild stenoses [46]. Mild stenoses are of particular interest here since we are interested in auscultation-based screening of early-stage aortic stenosis, which is usually asymptomatic and hard to detect [47].

The valve motions are summarized in Fig. 8(a). It is worth emphasizing that PVA is normalized by the maximum PVA of the healthy case, i.e., PVA₀. The inlet velocity profile is also included in the figure for reference. From this figure, we can see that the valve in all three cases opens rapidly and reaches maximum opening before peak systole. Though the leaflets in the healthy case remain fully open for approximately $0.06T$⁠, the stenosed cases immediately start to close after reaching maximum PVA. Moreover, the healthy valve is not fully closed until around $0.4T$⁠, while the other two cases close at approximately $0.35T$⁠, right after the flowrate reduces to zero. The valve remains closed during diastole in all three cases. The valve profiles at the maximum opening position for the three cases are shown in Fig. 8(b). When fully open, the valve in the healthy case creates a circular orifice and has minimum blockage of the flow, while the valve in the stenosed cases creates an orifice that is more triangular-shaped.

The inclusion of the valve model and pulsatility is expected to create complex postvalvular hemodynamics. Figure 9 shows the instantaneous vortex structures for the three cases studied here through the Q-criterion [48]. Results from selected phases are plotted here to demonstrate the effect of pulsatility. When the flow is accelerating (phase I), the valves in all three cases have already reached fully open position. Due to the orifice created by the valve, a jet is formed inside the aorta. For the healthy and 31% stenosis cases, the jet is mainly laminar at this phase and a well-defined vortex ring is clearly visible in the aorta. On the other hand, due to the strong intensity of the jet created by the narrower orifice, the flow after the 47% stenosis already shows chaotic behavior and breakup. At phase II, the inflow reaches maximum flowrate. While the leaflets in the healthy case remains fully open, they are already closing in the stenosed cases. The decreased orifice area along with the increased jet intensity are potentially responsible for the breakup of the jet for the 31% case during phase II, while the vortex structures in the healthy case are still coherent and well-defined. The postvalvular jet created by the healthy valve does not break up until the flow starts to decelerate (phase III), and the disintegration of the vortex structures is not as intense as the stenosed cases. Moreover, the extent of the jet penetration into the aorta is also affected by the narrower and higher speed jet in the cases with stenosis, and consequently, the small-scale vortex structures can reach further downstream when the stenosis becomes more severe. The small vortex structures continuously dissipate during the deceleration phase of the systole and early diastole as shown in phases IV and V.

Figure 10 plots the fluctuation kinetic energy for the three cases. Hereafter, we do not distinguish between fluctuation and turbulent kinetic energy for ease of discussion. The plots show that for the baseline (healthy) case, the aortic jet does not exhibit a breakdown in the shear layers until quite late in the cardiac cycle (phase IV). The strongest fluctuations occur around $s=1.0Do$ (s is the distance into the aorta measured from the end of the aortic root). On the other hand, for the case with 31% stenosis, a higher magnitude of TKE is observed shortly after peak systole and the TKE continues to grow in intensity and size before it starting to diminish. By end of flow ejection (phase V), the turbulence is almost fully dissipated in a manner similar to the healthy case. For the $AS=47%$ case, turbulence appears before peak systole and a region of strong turbulence occupies a large region of the proximal aorta before rapidly diminishing at the end of flow ejection.

To investigate whether the break-up of large vortex structures leads to turbulence, the spectra associated with velocity fluctuation is plotted in Fig. 11. These spectra correspond to velocity signals measured at $s=1Do$ and at the midpoint along the radius in the healthy and two stenosed cases studied here. This postvalvular location is chosen since it lies within the shear layer of the jets and all cases exhibit high levels of fluctuations here. Figure 11 clearly demonstrates that the spectra of the diseased cases are in reasonable agreement with the classic turbulent scalings [49,50], indicating that the postvalvular jets are turbulent when there is stenosis. On the other hand, the lack of a noticeable inertial subrange in the spectra for the healthy valve case suggests that the aortic jet for this case does not fully transition to turbulence. This is consistent with in vivo observations that patients with diseased aortic valves are more likely to develop turbulence in the ascending aorta [51].

The temporal and spatial distribution of the TKE in the aorta is made clear through Fig. 12, in which the cross-sectionally averaged TKE of $AS=47%$ case is plotted along with the spatial and temporal coordinates. The figure clearly shows the spatial and temporal onset of the turbulence which occurs at about $t/T=0.1$ immediately downstream of the valve. The turbulence front tracks closely with the average transvalvular velocity and reached a peak at around $t/T=0.17$ and $s/Do=1.2$⁠. From this figure, the temporal variation and spatial distribution of the TKE can be derived and are shown in Fig. 13 for all the cases. It can be seen from Fig. 13(a) that the onset of the turbulence is delayed for decreasing severity levels of the stenosis, and the strongest turbulent fluctuation occurs during the deceleration phase of the systole. Even though the peak TKE location shifts slightly during the systole, the time-averaged value (Fig. 13(b)) shows that both stenosed cases experience most intense fluctuation around $1Do$ downstream of the valve.

The most prominent difference between the current flow behavior and the one by Zhu et al. [52], where a static stenosis model with pulsatile inflow was considered, is the lack of turbulence during the acceleration phase of systole for mild stenoses in the latter study. Without a dynamic valve model, this early appearance of turbulence is only observed in the severe stenosis (⁠$AS=75%$⁠), while the flow remains laminar in $AS=50%$ during the acceleration [52]. Nevertheless, the inclusion of the valve decreases the effective opening area during the acceleration, creating a much stronger jet, and this leads to early break-up of the jet structures. As will be shown in the following section, this phenomenon has significant implications for the generation of murmurs.

Previous studies [16,17,53,54] have already demonstrated that it is the force acting on the lumen wall that generates the systolic murmurs associated with such flows. More specifically, it is the wall pressure that serves as the source for murmur generation. Hence, in this section, we focus on the analysis of the wall pressure.

An example of the wall pressure signal at $s=1.0Do$ for the $AS=47%$ case is plotted in Fig. 14(a). The overall trend of the signal follows the after-load profile (Fig. 3) but increased high-frequency fluctuations are observed around the peak systole. The spectrogram in Fig. 14(b) confirms that the increased high-frequency component occurs at this phase. It is noted that the total pressure (Fig. 14(a)) generates the normal heart sound as well as any murmurs associated with the stenosis. Though it is difficult to decompose the total pressure precisely into individual components responsible for the generation of normal heart sound and murmurs, previous studies suggest that the fluctuation component of the pressure (⁠$p″$⁠, see Fig. 14(c)) plays a dominant role in murmur generation [16,17].

Figure 15 further compares the spectra of the total wall pressure and the turbulent pressure fluctuations. It shows that the turbulent pressure fluctuations overlap well with the total pressure from 10 Hz to 1000 Hz. This corroborates with general assumption that the contribution from fundamental frequency (heart-beat) can be safely neglected after its tenth harmonics, and the turbulent pressure fluctuations can be treated as the main source of murmur generation [55]. Moreover, typical aortic stenosis murmurs measured in clinic have the crescendo-diminuendo pattern (diamond-shaped) that span most of the systole [3]. The current dynamic valve model is better at capturing this shape (see Fig. 14(c)) than the static stenosis model of Zhu et al. [52] since it has the capability to capture the early onset of turbulence in mild stenosis.

The source location of the murmurs is defined as the location on the lumen wall with strong $p″$⁠, which is quantified by the RMS value calculated through Eq. (9). In order to investigate how the source location changes with time, the spatial-temporal distribution of the circumferentially averaged wall pressure fluctuation is plotted in Fig. 16. It is consistent with previous observation that the murmurs are expected to occur as early as in the acceleration phase in the 47% case and reach a maximum intensity during early deceleration. Moreover, the source location shifts slightly downstream during midsystole but is mostly confined between $1Do$ and $1.5Do$⁠.

The temporal and spatial distribution of the murmur source in all three cases is summarized in Fig. 17. In Fig. 17(a), the onset of peak murmur intensity is delayed as the severity decreases. In Fig. 17(b), the high magnitude murmur source for the cases with valve stenosis is located around $1Do$ downstream of the valve. Both trends are in great agreement with the TKE trend shown in Fig. 13.

Figure 17 has great implications for developing auscultation-based automated diagnostic methods. First, Fig. 17(a) confirms clinical observation that aortic stenosis causes systolic murmurs [3], and it originates from the turbulent fluctuations caused by the jet break up (see Fig. 10). Through EKG-gated murmur measurements [7], simply classify a murmur as a systolic murmur can help to rule out heart conditions that generate diastolic murmurs, such as coronary artery stenosis and aortic regurgitation [56]. Second, Fig. 17(b) demonstrates that the source location of the murmurs is not sensitive to the two severities studied here and lies around $1Do$ after the aortic valve. This insensitivity can help us to further distinguish aortic stenosis from other heart conditions that generate systolic murmurs through source localization [55]. Last but not least, for the same patient, it is clear the onset of murmurs shifts toward early systole as the severity of the stenosis increases. This provides theoretical support that longitudinal tracking (measuring over time) of heart murmurs of individual patients can help monitor disease progression and potentially trigger early intervention.

In this study, a simple reduced-degree-of-freedom valve model has been introduced to study the fluid-structure interaction in the aortic stenosis. A canonical straight aorta model is employed here, and the valve is modeled after a bioprosthetic valve. Instead of solving a full elastic membrane model for the valve leaflets, an ordinary differential equation is used to model the flow-induced dynamics of the leaflets with a prescribed mode shape. By comparing with the results from a more complex FSI valve model, this simple model is shown to be able to capture the key features of the valve motion such as the valve-open area versus time, which is crucial in poststenotic jet formation.

This valve model is then used to investigate hemodynamics associated with the aortic stenosis, and a healthy case, along with two cases with area stenosis of 31% and 47% are studied. The healthy case serves as the baseline, and it shows extremely weak murmur generation. The incomplete opening is induced by increasing the stiffness coefficient of the valve model, and only rotationally symmetric stenosis is considered here. The results show that the source quantified by wall pressure fluctuations in localized around $1Do$ downstream of the valve. The intensity of the sources increases with stenosis severity and in both stenosed cases, the strongest wall pressure fluctuations occur in the deceleration phase of the systole. Compared with the static, no-valve stenosis model, the most noticeable difference here is the onset of murmur generation in mild stenosis (⁠$AS=47%$⁠) during acceleration phase of the systole. This is due to the smaller effective orifice size during the opening of the valve, and the wall pressure fluctuation from the current model is able to capture the typical diamond shaped profiles of systolic aortic stenosis murmurs [3]. This shows that, despite the idealized model employed here, this hemodynamic study provides valuable insights into the murmur generation mechanism, especially on timing and source location, which can potentially be used to guide the development of automated auscultation based diagnostic methods.

There are several limitations in the current study. First, a straight tube model is used to represent the aorta. Although this model is widely used in many studies, it fails to include the effect of the secondary flows induced by the aortic arch [17,52]. Second, the vibration of the leaflets is not captured by the current model. Its effect on the murmur generation still needs investigation. Moreover, the stenoses are assumed to be rotationally symmetric, i.e., all three leaflets have the same stiffness coefficient. However, this is usually not the case especially for bioprosthetic valves [45]. The effects of the asymmetric stenosis have been investigated in another study [44]. Last but not least, to ensure the robustness of the ghost-cell based sharp-interface IB method [42], a small gap is required between leaflets in the current model. A strict positive inflow is used here to prevent valvular leakage through the gap. But other mechanism is needed to prevent the backflow if a physiological inlet velocity profile with negative inflow is used.

The authors acknowledge support from NSF Grants IIS-1344772 and CBET-1511200, and computational resource by the Maryland Advanced Research Computing Center.

• NSF (Grant Nos. IIS-1344772 and CBET-1511200; Funder ID: 10.13039/100000001).

• Division of Chemical, Bioengineering, Environmental, and Transport Systems (Funder ID: 10.13039/100000146).

• Division of Information and Intelligent Systems (Funder ID: 10.13039/100000145).

Here, SV represents the stroke volume and is set to 70 ml in the current study. The duration of the systole is $Ts=0.33s$⁠, while the entire cardiac cycle lasts $T=0.9375s$⁠.

Here, Q is the flowrate provided in Eq. (A1). The parameters are set to: $R1=8.0Pa·s/ml, R2=160.0Pa·s/ml$ and $C=0.02ml/Pa$ in the current study. Figure 3 plots the after load generated by theses parameters.

Given the flow is driven by the same velocity profile in all three cases, the case that generates the strongest poststenotic jet is expected to be most affected by any grid resolution issues. Hence, a grid convergence study is carried out on the $AS=47%$ stenosis. Three grids with different resolutions: coarse (⁠$120×128×200$⁠), baseline (⁠$154×160×256$⁠) and fine (⁠$190×200×320$⁠) are employed, and the ratios of the average computational cell volume are $Vcoarse/Vbaseline≈2.05$ and $Vbaseline/Vfine≈1.93$⁠. The simulations are carried out for 1 cardiac cycle and the results are compared here. Figure 18 plots the instantaneous velocity magnitude (⁠$t=0.08T$⁠) across the diameter at cross section $s=0Do$⁠. It clearly shows that the velocity profile from the baseline mesh agrees well with the fine mesh, while the result from coarse mesh deviates significantly from the other two cases, especially near the valve. Moreover, the time history of the total kinetic energy over 1 cardiac cycle is shown in Fig. 19. The kinetic energy calculated on the baseline and the fine meshes overlaps well with each other over the whole cycle. However, the results from the coarse mesh is higher than the other two cases. Based on this study, the baseline mesh (⁠$154×160×256$⁠) is used to study all three cases.