Introduction

The relationship between the first gradient moment m₁ and the signal phase is used in flow MRI to quantify mean velocities in each voxel. Additionally, the velocity-encoding-related signal loss can be used to characterize turbulence^1-3 and to quantify the Reynolds stress tensor (RST)⁴, a fluid-dynamic parameter that characterizes the turbulent fluctuations in fluid momentum that might be of clinical interest in pathologies with altered hemodynamics. For intravoxel velocity standard deviation (IVSD) mapping, i.e. the quantification of the RST’s diagonal elements, it has been reported that the choice of m₁ influences the quantified turbulence³, since low m₁ values lead only to insignificant signal loss, whereas large m₁ values result in strong signal dephasing such that magnitude images reflect the noise floor. The same problem exists for the quantification of the RST. Here, the normal stresses (diagonal elements) are typically substantially higher than the shear stresses (off-diagonal elements), which requires a high dynamic m₁ range to correctly quantify the full tensor.
In this work, we examined a well-known fluid-dynamic test case, the flow over periodic hills. This test case was chosen because the turbulence in this flow is strongly inhomogeneous and anisotropic, representing a challenging measurement task.

Theory

Assuming a Gaussian intra-voxel velocity distribution, the corresponding magnitude of an m₁-encoded MRI signal can be expressed as:
$$\left|S\left(\vec{k}_v\right)\right|=\left|S_0\right|\cdot\exp\left(-\frac{1}{2}\vec{k}_v~\Sigma^{-1}~\vec{k}_v\right)$$
where$$$~\vec{k}_v=\gamma\left[m_{1,x}~m_{1,y}~m_{1,x}\right]^{T}$$$, $$$\gamma$$$ the gyromagnetic ratio, $$$\Sigma$$$ the variance of the 3-dimensional Gaussian, and $$$S_0$$$ the signal of a velocity-compensated measurement.
Multiplying $$$\Sigma^{-1}~$$$by the fluid density $$$\rho~$$$is equal to the RST $$$\tau$$$. Typically, $$$\tau~$$$is obtained by fitting a 3-dimensional Gaussian to the data of six velocity-encoded measurements⁴ with the following encoding matrix $$$E$$$, according to the ICOSA6⁵ encoding scheme and an additional compensated measurement:
$$E=m_1^{enc}\cdot\begin{bmatrix}0&0&0\\1/sqrt\left(1+\psi^2\right)&0&\psi/sqrt\left(1+\psi^2\right)\\1/sqrt\left(1+\psi^2\right)&0&-\psi/sqrt\left(1+\psi^2\right)\\\psi/sqrt\left(1+\psi^2\right)&1/sqrt\left(1+\psi^2\right)&0\\-\psi/sqrt\left(1+\psi^2\right)&1/sqrt\left(1+\psi^2\right)&0\\0&\psi/sqrt\left(1+\psi^2\right)&1/sqrt\left(1+\psi^2\right)\\0&-\psi/sqrt\left(1+\psi^2\right)&1/sqrt\left(1+\psi^2\right)\end{bmatrix}$$
where $$$\psi=\left(1+sqrt\left(5\right)\right)/2~$$$and $$$m_1^{enc}~$$$is the applied encoding moment. These measurements correspond to six points distributed on a shell with radius $$$m_1^{enc}$$$.

Methods

An in-vitro flow test case as shown in Fig.1A was used that consisted of 13 identical consecutive hills with measurements performed between the 10th and 11th hill. Flow rate and temperature were kept constant at 255l/min and 22°C, resulting in a Reynolds number of Re≈60000. Data were acquired on a 3T system (Magnetom Trio, Siemens, Germany) with a custom phase-contrast sequence employing an ICOSA6 encoding with values of [2.5,$$$~$$$5,$$$~$$$10,$$$~$$$15,$$$~$$$20,$$$~$$$25,$$$~$$$30,$$$~$$$35,$$$~$$$40,$$$~$$$45, $$$~$$$50,$$$~$$$60] mT/m·ms². Acquisition parameters: resolution=1.0x1.0x6.0$$$\,$$$mm³, FOV=496x120$$$\,$$$mm², TE/TR=6.32/10.0$$$\,$$$ms, averages=256.
For velocity quantification, the phase data were corrected for background phases by additional measurements without flow and averages=4 but otherwise identical settings. The velocity vector field was then reconstructed from the data corresponding to the five lowest $$$m_1^{enc}~$$$values to minimize phase bias from higher orders of motion (e.g. acceleration). The RST was quantified by two different methods. First, by fitting a 3-dimensional Gaussian to the magnitude data of the 12 ICOSA6-encoded measurements, normalized by the velocity-compensated measurements. To avoid the offset due to the Rician nature of the magnitude noise, data points with less than 10% normalized magnitude were excluded from the fit. Second, the RST was calculated by fitting a 3-dimensional Gaussian individually for every $$$m_1^{enc}~$$$value.

Results

The three velocity vector field components shown in Fig.1C illustrate that the velocities are restricted to the x-y plane within the examined slice. Example magnitude images for all 12 shells along the second encoding direction are shown in Fig.2A, demonstrating the high range of turbulence present in the test case. This finding is supported quantitatively by Fig.2B, which displays the measured data of two regions marked in Fig.2A with corresponding Gaussian fits. Note that here, only the projection of the three-dimensional Gaussian along the second encoding direction is shown. Additionally, Fig.2C+D show differences between the Gaussian fit to the full dataset compared to individual Gaussian fits to two example $$$m_1^{enc}~$$$values. In the case of high turbulence (Fig.2D), the individual fit to $$$m_1^{enc}$$$=45$$$\,$$$mT/m·ms² leads to a substantial overestimation of $$$\Sigma~$$$and therefore to an underestimation of the RST (7-fold for$$$~\tau_{xy}$$$).
All six components of the RST, quantified by the Gaussian fit to the full dataset, are illustrated in Fig.3A-F. As expected, $$$\tau_{xx}~$$$and $$$\tau_{xy}~$$$are the highest normal and shear stresses with peak values of 6.4$$$\,$$$m²/s² and 3.8$$$\,$$$m²/s² respectively. In contrast, lower shear stresses of 0.7$$$\,$$$m²/s² and -1.0$$$\,$$$m²/s² are obtained for $$$\tau_{xz}~$$$and $$$\tau_{yz}$$$. Additionally, the turbulent kinetic energy (TKE), the trace of the RST (see Fig.3G), yielded peak values of 6$$$\,$$$kJ/m³.
Fig.4+5 display results from the individual Gaussian fits to five example $$$m_1^{enc}~$$$values, indicating amplified noise for low $$$m_1^{enc}~$$$values and peak stress underestimation for high $$$m_1^{enc}~$$$of up to 3-fold for $$$\tau_{xx}~$$$and the highest $$$m_1^{enc}$$$. Furthermore, this method fails to correctly quantify the smallest shear stress components $$$\tau_{xz}~$$$and $$$\tau_{yz}$$$.

Discussion

In this work, we successfully quantified the RST in a well-known fluid-dynamic test case that produces stresses over a high range. However, additional data either from simulations or other reference measurements would be needed to estimate the method’s accuracy. The results indicate, in analogy to IVSD mapping, that RST quantification is highly susceptible to the applied$$$~m_1^{enc}~$$$value. Furthermore, RST mapping inherently requires a higher dynamic range compared to IVSD mapping, since the shear stresses are typically much lower than the normal stresses.
We speculate that multiple$$$~m_1^{enc}~$$$values are in general necessary to correctly quantify all components of the RST, not only for this particular test case.

Acknowledgements

No acknowledgement found.