Non-invasive estimation of relative pressure in turbulent ﬂow using virtual work-energy

Vascular pressure differences are established risk markers for a number of cardiovascular diseases. Rel- ative pressures are, however, often driven by turbulence-induced ﬂow ﬂuctuations, where conventional non-invasive methods may yield inaccurate results. Recently, we proposed a novel method for non- turbulent ﬂows, ν WERP, utilizing the concept of virtual work-energy to accurately probe relative pressure through complex branching vasculature. Here, we present an extension of this approach for turbu- lent ﬂows: ν WERP-t. We present a theoretical method derivation based on ﬂow covariance, quantifying the impact of ﬂow ﬂuctuations on relative pressure. ν WERP-t is tested on a set of in-vitro stenotic ﬂow phantoms with data acquired by 4D ﬂow MRI with six-directional ﬂow encoding, as well as on a patient-speciﬁc in-silico model of an acute aortic dissection. Over all tests ν WERP-t shows improved accuracy over alternative energy-based approaches, with excellent recovery of estimated relative pressures. In particular, the use of a guaranteed divergence-free virtual ﬁeld improves accuracy in cases where turbulent ﬂows skew the apparent divergence of the acquired ﬁeld. With the original ν WERP allowing for assessment of relative pressure into previously inaccessible vasculatures, the extended ν WERP-t further en- larges the method’s clinical scope, underlining its potential as a novel tool for assessing relative pressure in-vivo


Introduction
Flow abnormalities are typical indicators of cardiovascular disease (CVD). In the presence of valvular stenosis, the development of post-stenotic turbulence is directly related to pathologi-cal changes in cardiac workload ( Schöbel et al., 1999;Dyverfeldt et al., 2013 ), and hemodynamic alterations in heart failure patients have been linked to pathological neurohormonal activation ( Schrier and Abraham, 1999 ). With disease-related flow changes even proposed to occur prior to any detectable morphological change ( Pedrizzetti et al., 2014 ), assessing hemodynamic behaviour is of direct clinical importance.
Several hemodynamic biomarkers have been proposed for the diagnosis of CVD. In particular, pressure drops or differences in relative pressure provide valuable clinical biomarkers for a range of CVDs. Transvalvular pressure drop is the recommended measure of stenotic severity ( Baumgartner et al., 2009 ), and regional changes in relative pressure have been proposed as a determinant for interventional angioplasty ( De Bruyne et al., 2006 ) or coarc-tated artery stenting ( Jenkins and Ward, 1999 ). The pressure drop over the left ventricular outflow tract is also an established risk marker in hypertrophic cardiomyopathy patients ( Bernard et al., 2011 ), and the clinical outcome for patients with pulmonary hypertension has even been linked to the degree of transpulmonary relative pressure ( Galié et al., 2015 ). Recent studies have also evaluated the production of turbulent flow in patients with aortic stenosis ( Dyverfeldt et al., 2013;Bahlmann et al., 2010 ), indicating links between cardiovascular relative pressure and disease severity.
Current clinical assessment of pressure variations is largely based on Doppler echocardiography or intravascular catheterization. While catheterization is limited by its invasive nature ( Wyman et al., 1988;Omran et al., 2003 ), in Doppler echocardiography 1D estimates of peak velocities are linked to pressure changes through the simplified Bernoulli equation ( Stamm and Martin, 1983 ). Albeit effective for certain subsets of CVD, the simplification of the assessed fluid mechanical environment limits the method's general applicability. Extended and modified Bernoullibased methods have been proposed ( Garcia et al., 20 0 0;Falahatpisheh et al., 2016 ), however discrepancies between Bernoulli estimates and invasive measures are still frequently reported ( Baumgartner et al., 1999;Garcia et al., 2003;Donati et al., 2017;Feldman and Guerrero, 2016 ).
Recent advances in medical imaging have now enabled the fullfield measurement of cardiovascular flow through 4D flow MRI . Through refined sequencing, assessment of incoherent turbulent flow structures has also been achieved ( Ha et al., 2017;Haraldsson et al., 2018 ). With this, a more comprehensive assessment of relative pressure is enabled where the Navier-Stokes equations -theoretically describing the flow of any viscous fluid -can be utilized in its complete form. Solution of a Pressure Poisson Equation (PPE) based on acquired 4D flow data has been proposed as a method for mapping relative pressure fields ( Krittian et al., 2012;Bock et al., 2011 ); however, its accuracy depends on the defined flow domain and has shown bias when applied to stenotic flows in-silico ( Donati et al., 2017;Bertoglio et al., 2018;Casas et al., 2016 ). Other methods have also been proposed, using a mixed PPE/Stokes approach ( Švihlová et al., 2016 ), or utilizing generated flow turbulence to quantify irreversible pressure changes. That turbulence can relate to increases in pressure has been extensively covered in both theoretical ( Pope, 2001;Davidson, 2015 ) and clinical literature ( Schöbel et al., 1999;Dyverfeldt et al., 2013 ), and originates from the fact that irregular velocity fluctuations will obstruct the natural passage of flow. Proposed turbulence-based methods have evaluated either turbulent kinetic energy alone , incorporated an added shearscaling approach ( Gülan et al., 2017 ), or assessed turbulent energy dissipation directly from image data ( Ha et al., 2017 ). Even though showing initial promise, the above methods have however all been limited to relatively simplified flow scenarios, and their applicability to transient turbulent flows remains somewhat unexplored.
Recently, we derived a formulation that evaluates relative pressure using either a direct work-energy form of the Navier-Stokes equations (Work-Energy Relative Pressure, WERP ( Donati et al., 2015 )) as well as an alternative form where the work-energy of an auxiliary virtual field was evaluated (virtual Work-Energy Relative Pressure, νWERP ( Marlevi et al., 2019 )). Utilizing νWERP, we showed that accurate relative pressure estimates could be achieved in arbitrary multibranched vasculatures. Furthermore, νWERP was validated in-vivo against invasive catheterization. The cohort utilized was also such that alternative approaches were inherently obstructed by utilized spatiotemporal image resolution or challenging vascular anatomy. However, while promising, the proposed νWERP method did not include analysis of turbulent energy dissipation.
In this study, we present an extension of the νWERP formulation, νWERP-t, incorporating turbulence-driven energy dissipa-tion and expanding the original method into highly turbulent flow regimes. With the original method already shown to enable arbitrary probing of relative pressure into previously inaccessible vasculatures ( Marlevi et al., 2019 ), the incorporation of turbulent energy extends the approach for more severe clinical scenarios. In the following, we validate νWERP-t in a set of in-vitro stenotic valve phantoms, and further evaluate method performance in a complex patient-specific time-dependent in-silico data set. In all cases, method performance is compared against alternative approaches, highlighting the benefits of the proposed νWERP-t as well as providing explanations to theoretical and practical differences between evaluated methods.

Methods
Starting with a recapitulation of the νWERP principle, we derive νWERP-t by incorporating turbulent flow fluctuations ( Section 2.1 ).
We outline the implementation of the method onto acquired flow data ( Section 2.2 ) and present some alternative energy-based relative pressure approaches ( Section 2.3 ). Lastly, we review the performed verification and validation tests including a set of in-vitro stenotic flow phantoms and a transient patient-specific in-silico case ( Section 2.4 ).

Fluid virtual work-energy for turbulent flows
The derivation of νWERP-t originates from the conservation of mass and momentum for an isothermal viscous Newtonian fluid, described by the Navier-Stokes equations as where v is the velocity field, p the pressure, ρ fluid density, and μ dynamic viscosity.
In νWERP-t, we assume that the assessed flow behaviour has a quasi-period T , meaning that the flow field exhibits semi-periodic flow behaviour (i.e. v (t) ∼ v (t + iT ) for any integer i ∈ N ) except within incoherent turbulent flow regions where more significant variability exists. Supposing that the data is collected over n cycles, we can then define a linear expected value operator E [ ·] such that for any given function f : [0 , nT ] → R m , (m = 1 , 2 , 3) . Following from the linearity of E , along with its commutativity with spatiotemporal derivative operators, applying E on Eqs. (1) and (2) gives Assigning V = E [ v ] and P = E [ p] to be the acquired mean field quantities , Eqs. (4) and (5) can be re-expressed as Following the functionality of E , it can be shown that ] is the covariance in the observed flow. Utilizing this, we reach a final form of the Navier-Stokes equations including incoherent flow regimes reading It is worthwhile to note that the above logic is identical to the Reynolds decomposition used in Reynolds-averaged Navier-Stokes (RANS) methods, being used extensively in the fluid mechanical community to model complex, turbulent flows ( Speziale, 1998 ). As shown previously ( Marlevi et al., 2019 ), derivation of the virtual work-energy equation can then be achieved by multiplying Eq. (9) by an arbitrary virtual velocity field w , and integrating over the entire fluid domain with boundary and normal n . However, when incorporating incoherent flow behaviour, an additional turbulent energy dissipation term T e appears following these stochastic fluctuations, rendering a complete energy balance reading Each entity in Eqs. (12) -(17) represent a different virtual energy component in the assessed system. Most intuitively understood is the case where w = v , in which the energy components correspond to real work-energy: K e representing kinetic energy held within the fluid, A e the rate at which kinetic energy enters, exits, or changes within , S e the power transfer into the system due to shear, V e the rate of viscous energy dissipation, H ( p ) the input hydraulic power, and T e the turbulent energy dissipation.
With the above we can isolate the relative pressure between two arbitrary boundaries by assigning certain attributes to w . Specifically, by splitting the domain boundary into an inlet, outlet and wall domain, respectively ( = i ∪ o ∪ w ), and by choosing w to be a solenoidal field ( ∇ · w = 0 ) assigned as w = 0 at w , H ( p ) can be expressed as which, under the condition of mass conservation of w (i.e. the total inflow Q through i has to be identical to the outflow through o ), simplifies into with p k being the w · n -weighted mean pressure on boundary k , and p being the relative pressure between i and o .
Further, by letting w act in the surface normal of , and knowing that the spatial gradients of v are negligibly small at the vicinity of the surface boundary, the shear input term S e becomes effectively negligible.
Lastly, again making use of the fact that w = 0 at w , and assuming that the chosen i and o are sufficiently far away from the core turbulence regions (i.e. Cov[ v, v ] Summarizing all of the above, we end up with a final isolated expression of p representing the νWERP formulation with extension into turbulent flow regimes, specifically given as with all right-hand side terms directly derivable from the acquired flow field. The final step for a νWERP-t computation is the identification of a virtual field w abiding to all of the aforementioned assumptions. This can be achieved by setting w as the solution to a Stokes problem, specifically defined as where λ is the virtual pressure field corresponding to w . At i a parabolic inflow is defined in the normal direction n , at radial position r with a perimeter radius R . At o no constraints are defined for w .

Computation from 4D flow data
Below is a condensed outline of the processing steps for the numerical realization of νWERP-t, utilized in all verification tests of this paper. An illustrative overview is provided in Fig. 1 . For further technical or numerical details, the reader is also referred to the similar processing outlined in the original νWERP formulation ( Marlevi et al., 2019 ).

Domain segmentation and labelling
To appropriately assess the fluid domain, segmentation and domain labelling was performed through a series of sequential steps: (i) Fluid domain segmentation -The entire assessed flow domain was segmented by thresholding on generated velocity magnitude images, separating fluid domain voxels from neighbouring static background. (ii) Inlet and outlet plane selection -Inlet and outlet planes between which pressure drops were to be calculated, were manually selected within the segmented flow fields. Specifically, plane positions were manually initiated before automatic adjustments were applied to assign the inlet/outlet planes as perpendicular cross-sections to the region of interest. Overview of WERP-t estimation from 4D flow data. Processing starts with the acquisition of 4D flow mean field quantities, and corresponding flow covariance. The mean field flow data is then segmented ( S ) to isolate the flow domain. Domain subsampling and labelling are then performed, before a Stokes flow virtual field is computed by means of a finite difference method scheme. In parallel, the acquired covariance data is masked ( M ) to remove non-physical negative diagonal entries. The data is then all combined, deriving relative pressure by means of virtual work-energy assessment. With such, relative pressure between the selected inlet and outlet planes is computed as a function of time.
(iii) Fluid domain labelling -To appropriately assign boundary conditions in the virtual field computation (see Section 2.2.2 ), labelling of the segmented vascular domain was performed. Specifically, voxel-wise region growing was initiated from the defined inlet plane, with all voxels inside the assessed domain labelled as: interior (entirely within the fluid domain), exterior (entirely outside the fluid domain), inlet/outlet (inside defined inlet/outlet planes), or wall (separating interior and exterior), respectively.

Virtual field computation
Using the labelled fluid domain above, the virtual field w was computed as a Stokes problem Eqs. (25) -(27) using a staggered grid Finite Difference Method (FDM). With labelling provided in Section 2.2.1 , defined boundary conditions in Eq. (27) were projected onto the related nodes in the discretized domain. To improve numerical accuracy, data subsampling was employed by a factor of 2, subdividing a single image voxel into 8 uniformly distributed subvoxels, each inheriting the value of the original parent voxel.
Once assembled, an iterative solver based on the BFBTpreconditioning method ( Elman et al., 2006 ) was then used to solve the FDM linear equation system.

Operator discretization and numerical implementation
For spatiotemporally discretized flow field data, the derived νWERP-expression in Eq. (24) can be re-expressed as with v t + 1 2 derived as the mean of v t and v t + 1 (a similar mean was derived for Cov t + 1 2 ). With the above, each separate energy component in Eq. (24) are calculated by integration over voxelized version of and i and o , here referred to as ROI , I,ROI and O,ROI , respectively. With an image voxel identified by its indices ( i, j, k ), the energy terms can then be computed as where dS = 2 i =1 x i and dV = 3 i =1 x i are the pixel area and voxel volume, respectively, both based on the voxel length x i in each spatial dimension. Additionally, G ( ·) is the discretized gradient operator, defined as a spatial central-difference operator with one-directional derivatives employed at domain boundaries.

Data noise filtering (i) Mean field noise filtering -With the verification tests of this pa-
per relying on either pre-processed in-vitro data, or simulated in-silico data, no additional mean field noise filtering was applied. However, the original νWERP-formulation utilized kernelbased k -order polynomial fitting to spatially smooth clinically acquired flow data ( Marlevi et al., 2019 ), providing a possible option for the extension of νWERP-t into clinical data.
(ii) Covariance noise filtering -In comparison to noise filtering of the mean field quantities, application of the same to covariance data (

Alternative energy-based approaches for relative pressure estimation in turbulent flows
For a comparative evaluation, νWERP-t was assessed in conjunction with a couple of alternative or previously published energy-based methods for relative pressure estimation of turbulent flows.
Firstly, the conventional, non-extended νWERP method was used. In the original νWERP formulation, no consideration of incoherent flow fluctuations are included, therefore non-invasive relative pressures are derived by excluding the turbulence-induced T e term given in Eq. (24) . In cases of significant turbulent energy dissipation, the non-extended νWERP is thus expected to diverge from the proposed νWERP-t extension. whereas no theoretical difference is to be seen against non-turbulent flow regimes. Secondly, an extended real work-energy formulation WERP-t can be computed by assigning w = v within Eqs. (12) -(17) . With such, WERP-t is identical to Eq. (24) , however with the energy components reflecting the real work-energy of the acquired flow field. Utilizing a real work-energy approach has the benefit of not having to compute a virtual field, and the non-extended WERP has been applied to assess relative pressures in-vivo ( Donati et al., 2017 ). However, when w = v , the real flow Q in Eq. (24) will now be a function of the acquired field, causing potential divergent behaviour during late diastolic phases (where Q → 0), as well as in branching vasculatures (where Q | i = Q | o , violating assumptions in the derivation of Eq. (19) ).
Thirdly, a previously published turbulence-based relative pressure method (here denoted as turbulence production, TP) was employed ( Ha et al., 2017( Ha et al., , 2019, with the method previously tested on the data set discussed later in Section 2.4 . Just as with WERPt, TP evaluates the total real work-energy within the assessed fluid domain. In the case of a purely turbulence-driven relative pressure, TP simplifies into which is identical to the T e -term of WERP-t. However, following details outlined in ( Ha et al., 2019 ), the TP method uses an additional covariance filtering such that at each voxel ( i, j, k ). In other words, all negative integral arguments in Eq. (35) , for any velocity direction and any voxel ( i, j, k ), are removed from the TP analysis. Disregarding the masking in Eq. (36) , the TP and WERP-t approach are theoretically identical, where the real work-energy of the acquired flow field is assessed to derive relative pressure. With such, the advantages and disadvantages of WERP-t also hold for TP, and the performance of the method will be potentially affected by branching or diastolic flows. The introduced covariance masking will however cause discrepancies between TP and WERP-t, where a portion of the produced turbulent energy dissipation is neglected by Eq. (36) .
Regarding alternative approaches, it is worth noting that T e as expressed in Eq. (35) has been referred to as turbulence production rather than turbulent energy dissipation in ( Ha et al., 2017( Ha et al., , 2019. Similarly, the diagonal covariance entries (Cov[ v i , v i ]) have sometimes been used to describe turbulent kinetic energy in e.g. ( Dyverfeldt et al., 2013 ). However, with incoherent flow fluctuations inherently connected to energetic dissipation, we argue that turbulent energy dissipation is the most fitting description of the observed T e .

Method evaluation
An outline of the tests performed to verify and validate νWERPt is presented below.

Validation in a steady-state turbulent flow
The νWERP-t approach, along with other benchmark methods, were first evaluated on a series of steady-state turbulent flow cases. This was done to isolate the turbulence-driven relative pressure from transient flow behaviours, otherwise shown to impact vascular pressure drops in-vivo ( Donati et al., 2015;Lamata et al., 2014 ).
For this test, a previously published data set was utilized ( Ha et al., 2019 ) consisting of turbulent flow data acquired in a set of 3D-printed stenotic aortic valve phantoms. In short, seven different stenotic valve configurations were installed in a straight acrylic pipe, respectively, with water circulated through the setup at four different steady flow rates (6.8 -25.5 L/min). As such, the Reynolds numbers ranged from 5552 to 20,0 0 0 (calculated using generated flow rates and the assessed geometry), ensuring a fully developed turbulent flow ( Pope, 2001 ).
For a given setup, turbulent flow measurements were performed by 4D flow MRI with six-directional icosahedral flow encoding (ICOSA6). Acquisitions were performed at 1.5T (Philips Achieva, Philips Medical Systems, Best, the Netherlands) with velocity encoding: 1-5 m/s, echo time: 1.5-2.2 ms, repetition time: 3.2-3.9 ms, spatial resolution: 1.5 mm 3 , and number of signal averages: 5, resulting in an average scan time of ~16 minutes/valve. From the acquisitions, mean field flow quantities (3D velocity components) and incoherent flow fluctuations (Covariance tensor, or Reynolds stresses) were derived using previously described protocols ( Ha et al., 2019( Ha et al., , 2017, including corrections of background phase errors and median noise filtering. In addition to the acquired flow data, invasive pressure measurements were obtained for all valves and flow rates using installed pressure ports inside the custom-made flow circuit (digital pressure manometers (GMH-3161-07B, Greisinger, Germany), positioned +/-five diameters away from the stenotic narrowing).
A visualization of the seven valve configurations, along with an example of the acquired mean field, covariance, and corresponding virtual field are given in Fig. 2 . For further technical details of the setup or acquisition, see ( Ha et al., 2019 ).  In-silico data from the patient-specific AAD simulation used for the transient turbulent flow analysis. Data shown for the initialization (upper row, first three cycles), and quasi-periodic phase (lower row, last seven cycles), with renderings of mean field ( V , left), covariance (Cov[ v, v ], middle), and corresponding virtual field ( w , right), respectively. Data it shown at three time points around the systolic peak (given at fractions of the complete cycle time, T ), with volume rendering opacity set as 0 → 1 for minimun → maximum magnitude. During the initialization phase, comparably high covariance is observed through the opening of the false lumen. In contrast, during quasi-periodic stabilization, the covariance magnitude decreases substantially. Note that since w is computed solely using the segmented domain and given inlet and outlet planes, it is independent of initialization or quasi-period state. Colour bar and axis are given for each field on the far right.

Validation in a transient turbulent flow
Expanding from the steady-state stenotic phantoms, evaluations of a transient flow case was also performed. Specifically, a patientspecific in-silico model of an acute type B aortic dissection (AAD) was utilized, for which 4D flow and pressure was generated using Computational Fluid Dynamics (CFD) simulations with a stabilized variational approximation of the Navier-Stokes equation for laminar and turbulent flows ( Whiting and Jansen, 2001 ). With model and simulations details provided elsewhere ( Dillon-Murphy et al., 2016 ), data from 10 consecutive simulations cycles were generated for the specific case of turbulent energy analysis. The simulated data set was then projected onto a uniform grid of 2 mm 3 , with data sampled at 32 time slices/cardiac cycle.
From the sampled simulation data, mean field and covariance quantities were generated in two separate sets: one derived from the first 3 simulated cardiac cycles, representing an initiation phase where flow field covariance was enhanced following large inter-cycle variations, and one derived from the last 7 simulated cardiac cycles, where a more quasi-periodic flow behaviour was observed with subsequent lower covariance. For both sets of mean field and covariance data, relative pressures were estimated from the left ventricular outflow tract, down to an approximate 180 • bend of the descending false lumen of the thoracic aorta. Estimations were performed using νWERP-t along with alternative methods discussed in Section 2.3 . A visualization of mean field, covariance, and virtual field for the two different phases is provided in Fig. 3 .

Statistical analysis
Linear regression was assessed to quantify the relationship between predicted and true pressure drop in the in-vitro stenotic flow phantoms. Additionally, Bland-Altman plots were generated to assess potential method bias.
For the transient AAD-case, again linear regression analysis and Bland-Altman plots were generated with each time frame considered a separate data point. Additionally, the mean error ε p was evaluated as where p e | n + 1 2 is the estimated p at discrete time step t n + 1 2 (representing the mean of t n and t n +1 as per Eq. (28) ), and p| n + 1 2 is the corresponding true data. N is the number of temporal sample points ( N = 32 ).
For all of the above, the data analysis and method implementation was performed using MATLAB R2016a (MathWorks, Natick, MA, USA).

Steady-state turbulent flow
Estimations of relative pressure through the different stenotic flow phantoms are presented for all evaluated methods in Table 1 , together with comparable reference invasive relative pressures. Linear regression and Bland-Altman plots for each estimation method are also presented in Fig. 4 . Over all acquisitions, νWERP-t showed a mean error of −1.8 ± 3.3 mmHg, with a linear regression slope of k = 1 . 00 . The non-extended νWERP showed a mean error of −6.5 ± 6.7 mmHg, with a linear regression slope of k = 0 . 78 . For the real work-energy approaches, TP showed a mean error of 1.0 ± 4.1 mmHg, with a linear regression slope of k = 0 . 89 , whereas comparably larger errors were observed for Fig. 4. Results from the in-vitro stenotic flow phantoms, presented as linear correlation (upper row) and Bland-Altman plots (lower row) for all evaluated energy-approaches, with P = measured relative pressure, and P e = estimated relative pressure. Individual data points are given by black dots, and fitted regression by the dashed black line. For the Bland-Altman plots, 95% limits are indicated by the dashed black, and mean pm standard deviation is explicitly written. All data is presented in mmHg.
For the two extended work-energy approaches ( νWERP-t and WERP-t), energy components contributing to the complete relative pressure estimate are given in Table 1 . For the virtual νWERP-t approach, A e contributed to 76% of the total relative pressure, with 21 and 3% coming from T e and V e , respectively. For the real WERP-t scenario, T e accounted for 93% of the total relative pressure, with 11 and 2% coming from V e and A e , respectively. Note that no contribution was given from the kinetic ∂ ∂ K e , following the steady-state nature of the performed experimental mean flow.

Transient turbulent flow
Estimations of relative pressure through the patient-specific AAD as a function of the different estimation approaches are given in Figs. 5 -6 for the initialization and quasi-periodic phase, respectively. Additionally, each subfigure is provided together with a detailing plot showing the variation of different virtual energy components as a function of time.
For νWERP-t, a mean error of 6.8% was given for the initialization phase, with the error changing to 6.2% during the quasiperiodic stabilization. For the non-extended νWERP, the corresponding mean errors were 47.8% and 6.9%. Regarding the real work-energy approaches WERP-t, WERP, and TP, mean errors of 89.0%, 126.0% and 362% were seen during the initialization phase, with values changing to 115.2%, 136.1% and 115.4% during the quasi-periodic stabilization.
For νWERP-t during initialization, the relative pressure at peak systole consisted of 18.5% ∂ ∂ K e , 34.3% A e , 0.3% V e and 46.9% T e . During quasi-periodic stabilization, the same components changed to 17.6% ∂ ∂ K e , 82.6% A e , 0.3% V e and 0.5% T e . For the real-work WERP-t, the relative pressure at peak systole during initialization consisted of 43.1% ∂ ∂ K e , 18.4% A e , 3.8% V e and 34.8% T e , whereas the same entities changed to 43.4% ∂ ∂ K e , 48.8% A e , 4.2% V e and 3.6% T e during quasi-periodic stabilization.
For a cumulative summation of the above estimated traces, Fig. 7 shows linear correlation and Bland-Altman plots for the insilico estimations. Over both evaluated phases, and for all discrete time points, νWERP-t showed a mean error of 0 ± 0.3 mmHg, with a linear regression slope of k = 1 . 00 . The non-extended νWERP has a mean error of −0.2 ± 0.8 mmHg, with a linear regression slope of k = 0 . 93 . For the real work-energy approaches WERP-t, WERP, and TP, mean errors of −1.5 ± 4.5 mmHg, -1.9 ± 5.4 mmHg, and, −0.1 ± 7.5 mmHg, respectively. Corresponding linear regression slopes were k = 0.36, 0.30, and 0.57, respectively.

Discussion
In this study, we have presented an extended virtual workenergy method, νWERP-t, incorporating turbulent energy dissipation to accurately assess cardiovascular relative pressure in potentially turbulent flow directly from acquired velocity data. By including stochastic flow fluctuations in the theoretical derivation, we show that accurate estimates of relative pressure can be achieved both in-silico and in-vitro . With the extended νWERP-t comparing favourably against alternative estimation methods -especially when applied on realistic, transient cardiovascular flowsmethod novelty is highlighted.

Accuracy and comparative performance of νWERP-t in steady-state turbulent flows
For the evaluated in-vitro stenotic flow phantoms, νWERP-t showed excellent accuracy, showcasing a practically 1:1 relation against reference measures ( k = 1 . 00 , R 2 = 0 . 98 ). Shown in detail in Table 1 , deviations were slightly larger in valves with lower flow magnitudes, with the prosthetic heart valve with one blocked leaflet (PVH1) having a νWERP-t mean error of 47%, corresponding to an absolute error of −3.1 mmHg. Comparably, the tricuspid aortic valve (TAV) -the valve phantom experiencing highest relative pressure -had a νWERP-t mean error of 14% or 2.2 mmHg. As shown in the Bland-Altman plots in Fig. 4 , however, only a slight underestimation of −0.9 ± 3.3 mmHg was seen over all phantoms, again highlighting the accuracy of the proposed method extension.
The importance of incorporating turbulent energy dissipation when assessing turbulent flow fields is also highlighted when comparing against results from the non-extended νWERP approach.
Here, a systematic underestimation is apparent, and both decreased linear regression ( k = 0 . 78 , R 2 = 0 . 89 ), and increased mean error ( −6.45 ± 6.7 mmHg) is reported. Hence, using the virtual work-energy approach, in average 21% of the derived relative pres- Fig. 5. Results from the in-silico simulations, presented as relative pressure ( P ) over time for the initialization phase for all evaluated energy-approaches. Ground truth relative pressures given by the solid gray line, with corresponding estimates given by the dashed black line with black dots. The lower row shows individual energy components contributing to the total relative pressure, indicated as ∂ ∂t K e in green, A e in blue, V e in yellow, and T e in red. Relative pressures are given in mmHg, with the time given in seconds. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) Fig. 6. Results from the in-silico simulations, presented as relative pressure ( P ) over time for the quasi-periodic phase for all evaluated energy-approaches. Ground truth relative pressures given by the solid gray line, with corresponding estimates given by the dashed black line with black dots. The lower row shows individual energy components contributing to the total relative pressure, indicated as ∂ ∂t K e in green, A e in blue, V e in yellow, and T e in red. Relative pressures are given in mmHg, with the time given in seconds. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) sure is accounted for by turbulent energy dissipation, and is consequently required to accurately recover true relative pressure.
Comparing against real work-energy approaches, slightly different results were obtained for WERP-t, WERP, and TP. For WERPt, again a systematic underestimation was observed, with increasing absolute errors at increasing relative pressure ( k = 0 . 55 , R 2 = 0 . 97 , mean error of −10.9 ± 12.8 mmHg). Again the valve phantoms with lower flow magnitudes experienced highest relative errors of 56 and 71% for PVH1 and the second prosthetic valve phantom, PVH2, respectively, equalling an absolute error of −3.4 and −5.7 mmHg. Conversely, the high flow magnitude TAV showed a mean relative error of 35%, equalling an absolute error of −20 mmHg. Noteworthy is that for the real work-energy approaches, the relative pressure was almost exclusively governed by turbulent energy dissipation, contributing to in average 93% of the entire relative pressure. This is also the reason to why the non-extended real work-energy WERP shows significant output deterioration, with barely any of the valve phantoms assessed accurately.
In contrast to WERP-t, the TP method did render accurate estimations over all evaluated valves, with only slight deviations from a 1:1 correlation observed ( k = 0.89, R 2 = 0.98) (the TP analysis has also been presented in previous work ( Ha et al., 2019 )). The difference between WERP-t and TP is particularly noteworthy considering that they both originate from the same theoretical description of work-energy. Instead, the reason for the observed discrepancy lies in the additional covariance masking introduced within the TP method. As outlined in Section 2.3 , a complete masking of all negative integral arguments in Eq. (35) is applied, evidently improving accuracy in the stenotic flow phantom cohort. In fact, if employing the same type of covariance masking to the WERP-t approach, the two methods do converge, and conversely, if removing the covariance masking from the TP approach, results coincide with the described underestimation of WERP-t (see Appendix A , for details of the analysis). It should however be pointed out that there exists no theoretical reason justifying the conservative masking in Eq. (36) , and in-principle negative entries can appear at least for the off-diagonal covariance terms (i.e. Fig. 7. Results from the cumulative in-silico analysis, presented as linear correlation (upper row) and Bland-Altman plots (lower row) for all evaluated energy-approaches, with P = measured relative pressure, and P e = estimated relative pressure. Individual data points are given by red (initiation) and blue (quasi-periodic) dots, and fitted regression by the black dashed line. For the Bland-Altman plots, 95% limits are indicated by the dashed black, and mean ± standard deviation is explicitly written. All data is presented in mmHg. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) for Cov [v i , v j ] where i = j ). The mechanisms behind why accuracy improves with masking in the in-vitro data set remains to be understood, however the general applicability of such masking seem limited, as highlighted in the extended transient analysis in Section 4.2 .
For other alternative estimation approaches not relying on work-energy estimations from the complete Navier-Stokes equations (such as simplified and extended Bernoulli, or shear scaling methods), the reader is referred to previously published results from the identical phantom data set ( Ha et al., 2019 ). In short however, associated method simplifications were reflected by deterioration in estimation accuracy.

Accuracy and comparative performance of νWERP-t in transient turbulent flows
The patient-specific in-silico data set was utilized as an extension of the detailed steady-state in-vitro analysis. The simulated data output did not only allow for controlled assessment of realistic flow behaviour, but did also serve as a clinically relevant example where turbulent energy dissipation act together with advective, kinetic, and viscous energy components in contributing to the total work-energy of the assessed vasculature.
As reported, νWERP-t showed high accuracy (cumulative linear regression of k = 0.98, R 2 = 1.00), with the relative pressure traces accurately captured both during phases of significant turbulent energy dissipation (initialization) as well as during phases where other energy components dominated relative pressure behaviour (quasi-periodic stabilization). Comparing to the non-extended νWERP results, the importance of including turbulent energy dissipation is again underlined, where the mean error increased from 6.8% to 47.8% when disregarding T e in the virtual work-energy evaluation (i.e. going from νWERP-t to νWERP). The deviation was particularly evident at peak systole, again highlighting the clinical relevance of νWERP-t, where peak systolic metrics are typically the ones used to represent the hemodynamic quantification in clinics.
Comparing against real work-energy approaches, a distinct reduction in accuracy was observed over all evaluated methods. For WERP-t and WERP, mean errors of above 89.0% were seen during both the initialization and quasi-periodic phase, and as evident in Figs. 5 -6 , severe output deterioration is observed during the diastolic phase, with method divergence following as the real flow Q → 0 (comments on virtual versus real work-energy behaviour is discussed in Section 4.3 ). TP had similar problems during diastole, but did also present a pronounced overestimation of peak systolic relative pressure during the high-turbulent initiation phase. As seen in the inlet plot on energy components in Figs. 5 -6 , the TP overestimation comes from a dominant T e , where the omission of negative turbulence production causes an overshoot in retrieved turbulent energy dissipation. Thus, the constricting masking of TP does seem to skew results in cases where negative turbulence components are required to accurately capture relative pressure behaviour. Instead, νWERP-t represents a viable option where contributions from kinetic, advective, viscous, and turbulent energy components together enable accurate assessment of relative pressure.

Comparative performance between virtual and real work-energy approaches
With the extended νWERP-t the attempt was to progress beyond current state-of-the-art methods and to overcome limitations or assumptions included in previous estimation approaches (such as omission of non-advective or turbulent energy components, or limitation to single-vessel geometries). As observed, νWERP-t performs favourably compared to alternative methods, and particularly, the extension into branching, transient flows is beyond what has been presented or evaluated previously.
The novelty of νWERP-t is particularly worth highlighting in conjunction to the real work-energy equivalent of WERP-t. As pointed out previously ( Donati et al., 2015;Marlevi et al., 2019 ) real work-energy approaches will be obstructed by branching flows, where the real flow Q through the defined inlet plane i is not equivalent to the real flow through the defined outlet plane o . Along the same lines, the real flow Q will vary as a function of time, and will cause divergent behaviour in phases of low flow (such as during late diastolic phases, where Q → 0). On the contrary, in the case of an introduced virtual field, mass conservation between i and o can be enforced, and the virtual Q will no longer be affected by low real flow magnitudes. This difference becomes particularly apparent in the utilized in- Table 1 Numerical data for the comparison of relative pressure estimation methods in the set of steady-state stenotic flow phantoms. The valve phantoms are labelled as TAV = tricuspid aortic valve, Circ = circular aortic valve, BAV = bicuspid aortic valve, PVH = prosthetic heart valve. Flow speeds are given as pump revolutions per minute (RPM), and all relative pressure outputs are given in mmHg. Contributions from separate energy components (advective energy, A e , viscous energy, V e , and turbulent energy dissipation, T e ) are provided for νWERP-t and WERP-t. Note that for the non-extended νWERP and WERP the energy contributions are identical, apart from the exclusion of T e . Similarly, TP is governed exclusively by T e , following associated method results in ( Ha et al., 2019 ).

Fig. 8.
Results for both in-vitro and in-silico analysis, utilizing a divergence-free V † -equivalent to V . (Left) Linear correlation and Bland-Altman plots for the stenotic flow phantoms, with individual data points given by black dots, fitted regression by the black dashed line, and 95% Bland-Altman limits indicated by the dashed black lines. P = measured relative pressure, and P e = estimated relative pressure. (Middle) Derived relative pressure traces during initialization and the quasi-periodic phase of the in-silico data set comparing ground truth (gray continuous line) to method estimate (black dashed line with black dots). (Right) Cumulative linear correlation and Bland-Altman plots from the in-silico analysis, with data given as red (initation) and blue dots (quasi-periodic). All relative pressures are given in mmHg, and the time in seconds. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.) silico model where νWERP-t shows maintained accuracy, whereas WERP-t and the other associated real work-energy approaches diverge following the complex model anatomy and varying real flow magnitudes.
However, the differences observed for the steady-state stenotic flow phantoms is less intuitive. In principle, for a steady-state single-vessel flow, there should not be any major differences between WERP-t and the virtual equivalent of νWERP-t. Regardless, WERP-t did generate systematic underestimations, whereas νWERP-t maintained accurate method output in the in-vitro cohort.
One hypothesis explaining this discrepancy is that the observed discrepancy might originate from noise in the acquired mean field V , where inter-cycle covariance variations might shift V from its expected mean field value. Thus, the acquired V would no longer represent the true mean field quantity, deteriorating the accuracy of WERP-t. To test this, we introduced a synthetic covariance shift of the virtual field w by adding a random distribution with covariance identical to V , and subsequently assessing whether fluctuation based noise in the data would cause any systematic underestimation of νWERP-t. As shown in detail in Appendix B , the introduced noise did affect accuracy; however, it did not cause any systematic bias in the generated output. Instead, results varied around the non-shifted νWERP-t equivalent, indicating that the hypothesized flow fluctuations were not the reason to the observed discrepancy.
A secondary explanation was instead that due to turbulent flow fluctuations, V might appear as non-divergence free ( ∇ · V = 0) in the acquired data set. A fundamental requirement of WERP-t is the acquired field's solenoidal or divergence free nature ( ∇ · V = 0 ).
If this would not hold, fundamental physical behaviour would be violated, and in the case of a work-energy evaluation, arbitrary additions or subtractions of energy would appear throughout the assessed flow field. Naturally, the experimental flow of the invitro stenotic flow phantoms is divergence free; however, if sampling data over cycles with high degrees of incoherent fluctua-tions, the acquired field might still appear as non-divergence free Therefore, to test this hypothesis we chose to project the acquired mean field V onto a divergence free Stokes equivalent V † , and assess WERP-t output using this guaranteed divergence free field (for details on how V † is derived, please see Appendix C ). By doing so, we do in fact observe improved accuracy of WERPt in the in-vitro stenotic flow phantoms, as shown in Fig. 8 . The improvement is such that the V † -derived WERP-t closely coincides with νWERP-t ( k = 0.95, mean error of −3.1 ± 4.0 mmHg), thus providing a reasonable explanation to the aforementioned discrepancy: due to high turbulence in the acquired field, divergent appearance of the mean field aggravates real work-energy accuracy.
In fact, the improvement is not limited to the in-vitro data, but computing a divergence free equivalent in the in-silico set also improves accuracy, as again illustrated in Fig. 8 . Specifically, the systolic peak is retrieved with higher accuracy, with deviations primarily limited to the later diastolic phase (where variations in Q still cause method deterioration in the divergence-free WERP-t case).
These findings do not only provide an explanation to the observed in-vitro discrepancies, but do also provide further support to the versatility and robustness of the proposed virtual workenergy approach. Regardless of the turbulent nature of the acquired field, divergence free behaviour is guaranteed in νWERP-t, ensuring that fundamental constraints are withheld in the evaluated virtual work-energy setup. Along the same lines, the apparent difference in energy component weighting between virtual and real work-energy approaches also have clinical implications. As observed, even in cases of dominant real turbulent energy dissipation, the use of a virtual field seem to enhance advective energy weighting compared to its real work-energy equivalent. Consequently, νWERP-t will rely on accurate mean field quantities to a higher extent than WERP-t, where focus will rather be on maintained covariance capturing. The applicability of these findings in a clinical setting remains to be explored, however the results in-dicate advantages in using a virtual work-energy approach when assessing relative pressures from acquired 4D flow data.

Data fluctuations and flow periodicity
Even though formulated with incoherent turbulent flow in mind, it should be noted that the turbulent energy dissipation T e will not only account for purely turbulent flow features, but will also include any other behaviour that give rise to variations in data. In a clinical setting, data fluctuations might arise due to patient movements, anatomical organ motion, heart rate variability, or even intrascan acquisition variations. As a result, νWERP-t could prove advantageous even when assessing non-turbulent or transitional flow domains, by simply incorporating data uncertainties into the global work-energy balance. Whether even under-resolved flow features could propagate into the estimation of T e (by influencing assessed voxel variations) remains to be clarified, but again points to the benefit of the extended νWERP-t approach in a practical setting.
With respect to the clinical acquisition, it is also worth mentioning the semi-periodic flow behaviour assumed within the formulation of E in Eq. (3) . Importantly, this periodicity would, in theory, not need to abide to the physiological periodicity of the assessed flow, but could in fact be arbitrarily chosen as long as the mean field and covariance data are also acquired in the same arbitrarily selected way. The critical point is instead that with an arbitrary periodicity, the relative pressure estimate could itself be subject to fairly large variations. With νWERP-t generating estimates of the change in mean relative pressure P , these become less relevant if the periodicity is chosen in such a way that significant fluctuations exist around P . Instead, with semi-periodicity assumed to follow the physiological cardiac cycle, the mean field estimates of νWERP-t will be more reflective of the experienced physiological pressure variations.

Models of turbulent energy dissipation
As posed, νWERP-t relies on the RANS-like description in Eqs.
(9) and (10) to represent the assessed turbulent energy dissipation. Alternative models of turbulent energy dissipation could be envisaged, however often require covariance information on subgridscale to accurately solve the posed closure problem (i.e. relating turbulent fluctuations to global flow definitions), and are hence less applicable for an approach such as νWERP-t where mean field and covariance data is provided as user input. However, experimentally acquired mean field and covariance data could instead serve as an interesting basis for how the turbulence-based closure problem could be posed, proving potentially useful when attempting to accurately model turbulent flows (such as within the field of large-eddy simulations ( Mason, 1994 )). Further evaluations would be required to clarify the utility of acquired turbulence data for such modelling, but again highlights the potential future benefits of refined, experimental turbulence mapping.

Limitations
In this study, the extended νWERP-t method has been tested on experimental in-vitro data, and transient in-silico data. Even though conceptually identical in implementation, a validation of the extended approach in an in-vivo setting thus remains to be performed. Similarly, a valuable clinical comparison of approaches ( νWERP-t, WERP-t, PPE, etc.) has yet to be performed in-vivo . However, with the original νWERP formulation successfully tested invivo , and with some alternative approaches showing limited performance in-silico ( Bertoglio et al., 2018 ), the extended νWERP-t formulation bears promise in being able to successfully assess relative pressure even in a clinic scenario.
Additionally, using the presented formulation, νWERP and its extended νWERP-t equivalent acts on acquired full-field data. That is, refined 4D flow MRI with mapping of incoherent flow fluctuations is required to compute relative pressures. To date, this would not be considered routinely available on conventionally used MRI scanners, and current implementations require fairly long scan times, limiting its clinical scope. However, the reported benefits of turbulence mapping together with the rapid introduction of novel sequences and diagnostic tools in the area of medical technology promises extended capabilities in the near future. Accelerated 4D flow MRI sequences has been recently proposed pushing aortic full-field acquisition to below two minutes ( Bollache et al., 2018 ), and similarly turbulence-mapping has been included in accelerated MRI frameworks ( Walheim et al., 2019 ), again bearing promise for including such anaylsis in future clinical practice. In addition to this, a simplification of the νWERP-t formulation for 2D acquisitions (such as using through-plane 2D PC-MRI, or by refined vector-flow ultrasound imaging ( Pedersen et al., 2014 )) could be envisaged, such as the one explored in ( Donati et al., 2017 ). For the case of extended turbulence capabilities, further detailed investigation of method performance would however be required.
The effect of spatiotemporal sampling has not been explicitly evaluated in this study, however a convergence study was performed on νWERP in a previous publication ( Marlevi et al., 2019 ).
There, robust spatiotemporal behaviour was indicated in aortic geometries even at fairly coarse sampling and high noise levels (errors of < 30% at dx: 3 mm 3 , dt: 80 ms, SNR: 10). The spatiotemporal effect on the estimation of T e remains to be evaluated, but with no temporal or spatial gradients involved in Eq. (23) , there are reasons to believe that νWERP-t might be similar to νWERP in this regard. Instead, the appropriate choice of spatiotemporal image sampling will be dependent on the physical nature of the assessed vasculature, where increased spatial sampling will be required when assessing narrow peripheral vasculatures, and increased temporal sampling will be conversely preferred when evaluating transient cardiovascular events. With respect to accurately capturing flow covariance, a previous imaging study has indicated that high spatial sampling ( < 1.5 mm 3 ) might be required to avoid underestimation of viscous losses ( Binter et al., 2016 ), and similarly a slight dependence on velocity encoding has been indicated ( Binter et al., 2013 ). However, the use of multiple velocity encodings should mitigate such effects. Regardless further studies could be required to clarify how variations in covariance assessment propagates into image-based turbulence-driven relative pressure estimates.
Lastly on the limitations side, it should be noted that the utilized in-silico data set is, despite its patient-specific nature, a simplification of the clinical equivalent of an in-vivo scan. In a clinical scenario, external patient movements, respiratory motion, intrathoracic organ displacements, heart rate variability, vascular compliance, or limitations in scan accuracy will all enhance assessed flow covariance, in comparison to the utilized in-silico data set. However, such increased flow covariance would still be handled by the included T e of the extended νWERP-t approach, and the method therefore has specific promise even with regards to a clinical, reallife in-vivo setting.

Clinical outlook
Refined hemodynamic analysis has shown the potential to improve the management of cardiovascular disease in several studies ( Pedrizzetti et al., 2014;Dyverfeldt et al., 2013;Lamata et al., 2014 ), and the assessment of regional cardiovascular relative pressure is a critical part of clinical guidelines ( Baumgartner et al., 2009;Bernard et al., 2011 ). However, to date clinical evaluations are limited to a subset of cardiovascular disorders where simplified method assumption holds, or where the in-vivo interrogation of cardiovascular flow is possible by conventional medical imaging. In fact, despite the compelling evidence that incoherent turbulent flow is associated with several cardiovascular disorders ( Dyverfeldt et al., 2013;Bahlmann et al., 2010;Stein and Sabbah, 1976 ), such behaviour has been overlooked due to the lack of accurate imaging or assessment methods.
With this in mind, the extended νWERP-t method widens the applicability of established clinical risk markers, specifically enabling the accurate assessment of relative pressure through complex, branched, turbulence-inducing vasculatures. Importantly, νWERP-t includes assessment of complete flow behaviour, where kinetic, advective, viscous, and turbulent energy terms all contribute to generate accurate estimation outputs. By so, the method has a defined versatility, proving accurate both in high-turbulent scenarios (such as in the utilized in-vitro cohort, or during the initiliazation phase of the in-silico model), as well as in cases of comparably lower turbulent flow (such as during the quasi-periodic phase of the in-silico model). Furthermore, the fact that the turbulent energy term of νWERP-t accounts for not only turbulent flow features but any behaviour resulting in data fluctuations (as noted in Section 4.4 ), there are also practical benefits of using νWERP-t in a clinical environment.
With recent analysis of hemodynamic turbulence providing new insights into pathological development ( Zajac et al., 2015;Binter et al., 2017 ), the extended νWERP-t method thus shows particular promise for improved clinical evaluation, utilizing advanced fullfield imaging techniques to provide accurate estimates of relative pressure.

Conclusion
In conclusion, we have here presented an extension of the νWERP formulation, νWERP-t, allowing for accurate assessment of relative pressure through complex multibranched vasculatures, incorporating the effects of turbulent energy dissipation. With relative pressure being a recognized clinical biomarker and with tur-bulent flow present in several cardiovascular diseases, νWERP-t has defined clinical potential, performing favourably to alternative approaches and enabling accurate assessment of relative pressure through previously inaccessible flow domains.

Declaration of Competing Interest
The authors declare that they have no competing financial interests or personal relationships that could have appeared to influence the work reported in this paper. covariance masking, respectively. For TP, the masking was removed by simply omitting Eq. (36) , whereas in WERP-t, the masking was introduced in the computation of T e in Eq. (23) (being theoretically identical to (35) .
Linear regression and Bland-Altman plots for WERP-t and TP with and without covariance masking are presented in Appendix Fig. A.1 . As observed, if utilizing the same type of covariance handling (with or without masking), TP and WERP-t gives very similar results. With masking, WERP-t and TP shows linear regression slopes of k = 0.95 and k = 0.89, whereas the same values without masking are k = 0.55 and k = 0.56, respectively. Along the same lines, mean values of 2.96 ± 4.27 mmHg and 0.95 ± 4.09 with masking, and -10.94 ± 12.77 mmHg and -11.39 ± 12.44 without masking are given for WERP-t and TP, respectively.

Appendix B. Noise and νWERP-t
To evaluate the effect of noise in the acquired mean field V , a synthetic covariance shift of the virtual field w was introduced. This was specifically performed to see whether image noise would cause any systematic bias of derived νWERP-t output similar to that observed with WERP-t in the stenotic flow phantoms.
The covariance shift of w was introduced by sampling from a random distribution with covariance identical to the one of the acquired field (i.e. Cov[ v, v ]). Additionally, to account for magnitude differences between the real and virtual field, the virtual covariance was scaled by the ratio of peak velocity of w and V in each spatial direction, respectively.
Once generated, the virtual covariance field was then added to w in a voxel-wise manner, resulting in a spatial covariance distribution similar to that of the acquired V . To account for fluctuations between samples, virtual covariance fields were samples 100 times, with νWERP-t output presented as mean ± standard deviation over all generated data.
Results for the stenotic flow phantoms using νWERP-t on synthetically noisy w are provided in Appendix Table B. and C.1 . As observed, the added image noise did affect accuracy, however did not cause any systematic bias in retrieved output. Instead, the mean estimate of νWERP-t with noise did not deviate more than 1 mmHg from the original νWERP-t results.

Appendix C. Divergence free projection of WERP-t
To evaluate the effect of potential non divergence free appearance of the acquired mean field V , a Stokes equivalent divergence free field V † was computed. Specifically, V † was constructed as V † + ∇ 2 V † + ∇λ † = V + ∇ 2 V (C.1) with no constraints on o for V † . With the above, V † maintains the major spatial properties of V , however with the added property of being strictly divergence free. To compute the above, a staggered-grid Finite Difference Method (FDM) approach was used, being identical to the one used for the computation of the virtual field w . The only difference is the modification in Eq. (C.1) , as well as the maintained inflow boundary conditions in Eq. (C.3) .
Results from WERP-t applied on the divergence free V † are given in the main manuscript in Fig. 8 , as well as in detail for the stenotic flow phantoms in Appendix Table B. and C.1 . As observed, accuracy increases significantly in the stenotic flow phantoms (linear regression of k = 0.95, R 2 = 0.98, mean error of -3.1 ± 4.0 mmHg), with output virtually coinciding with the ones from νWERP-t. For the divergence free WERP-t, T e dominated the relative pressure, accounting for 75% of the derived output. For the remaining contribution, A e and V e accounted for around 13% each.
In the in-silico data set, the divergence free WERP-t improved estimates of the peak systolic relative pressure, as visually apparent in Fig. 8 in the main manuscript. However, V † did not reduce the large deviations observed during the diastolic phase (when real Q → 0). Table B. and C.1 Numerical data for the comparison of relative pressure estimation methods in the set of steady-state stenotic flow phantoms. The valve phantoms are labelled as TAV = tricuspid aortic valve, Circ = circular aortic valve, BAV = bicuspid aortic valve, PVH = prosthetic heart valve. Flow speeds are given as pump revolutions per minute (RPM), and all relative pressure outputs are given in mmHg. Data is given for the divergence free projected WERP-t, and the νWERP-t with added image noise (sampling noise from a distribution with covariance given by V ). Contributions from separate energy components (advective energy, A e , viscous energy, V e , and turbulent energy dissipation, T e ) are provided for the divergence free WERP-t.