# Quantifying inner-outer interactions in non-canonical wall-bounded flows

Mogeng Li1, Woutijn J. Baars1, Ivan Marusic2 and Nicholas Hutchins2 1Faculty
of Aerospace Engineering, Delft University of Technology, 2629 HS, The
Netherlands
2Department of Mechanical Engineering, University of Melbourne, Victoria 3010,
Australia

###### Abstract

We investigate the underlying physics behind the change in amplitude
modulation coefficient in non-canonical wall-bounded flows in the framework of
the inner-outer interaction model (IOIM) (Baars et al., Phys. Rev. Fluids 1
(5), 054406). The IOIM captures the amplitude modulation effect, and here we
focus on extending the model to non-canonical flows. An analytical
relationship between the amplitude modulation coefficient and IOIM parameters
is derived, which is shown to capture the increasing trend of the amplitude
modulation coefficient with an increasing Reynolds number in a smooth-wall
dataset. This relationship is then applied to classify and interpret the non-
canonical turbulent boundary layer results reported in previous works. We
further present the case study of a turbulent boundary layer after a rough-to-
smooth change. Both single-probe and two-probe hotwire measurements are
performed to acquire streamwise velocity time series in the recovering flow on
the downstream smooth wall. An increased coherence between the large-scale
motions and the small-scale envelope in the near-wall region is attributed to
the stronger footprints of the over-energetic large-scale motions in the outer
layer, whereas the near-wall cycle and its amplitude sensitivity to the
superposed structures are similar to that of a canonical smooth-wall flow.
These results indicate that the rough-wall structures above the internal layer
interact with the near-wall cycle in a similar manner as the increasingly
energetic structures in a high-Reynolds number smooth-wall boundary layer.

##  I Introduction

In turbulent boundary layers, large coherent structures are found in the
logarithmic region. They carry a high level of turbulent kinetic energy, and
make a significant contribution to Reynolds stress production [1]. These
structures can be further classified as large-scale motions (LSMs) and very-
large-scale motions (VLSMs). The former are associated with the vortex packets
formed by aligned hairpin vortices [2, 3], and typically have a streamwise
length of $\sim 3\delta$ ($\delta$ is the boundary layer thickness), while the
latter may be related to the merging of multiple LSMs [2], and can reach a
streamwise extent up to $20\delta$ with spanwise meandering [4, 5]. These
structures are observed to leave a footprint in the near-wall region [6, 4,
5]. As the Reynolds number of a turbulent boundary layer increases, the
separation between the outer-scaled motions and viscous-scaled near-wall cycle
becomes more distinct, and the strength of these large-scale footprints also
intensifies, resulting in the growth of the ‘inner-peak’ magnitude in the
broadband streamwise turbulence intensity [5, 7].

In addition to the direct superposition effect manifested as large-scale
footprints, it has also been found that the near-wall small scales are
modulated by the large scales in the outer layer [8, 5, 9]. In the near-wall
region, the amplitude and frequency of small-scale fluctuations show a
decrease when co-existing with large-scale low-speed regions, and vice versa
in the case of large-scale high-speed regions. These observations were
harnessed by the ‘two-scale’ framework, where a local (near-wall) fine-mesh
solution is coupled to the global coarse-mesh solution, in order to reduce
computational costs at high Reynolds numbers [10, 11]. In terms of the
modelling efforts, quasisteady quasihomogeneous theory [12, 13] provides an
axiomatic description of the scale interaction in near-wall turbulence. The
modulation effect was also quantified in recent works [14, 15, 16, 17, 18],
and a predictive model, termed the inner-outer interaction model (IOIM), which
outputs representative turbulence statistics in the near-wall region based on
an input signal in the logarithmic region, was developed by Marusic and
coworkers [19, 20] and by Agostini and Leschziner [21]. The former was later
revised by Baars  _et al._ [22] using spectral linear stochastic estimation.
The IOIM provides an opportunity to push the boundary of large-eddy
simulations of wall-bounded flows to very high Reynolds numbers at affordable
costs, thanks to its ability to provide representative real-time small-scale
signals in the viscous-scaled near-wall region based only on large-scale
outer-layer information [23, 24, 25].

Although similar amplitude modulation behaviours have been observed in
turbulent boundary layers, channel, and pipe flows [26], the existence of such
scale interactions remains unexplored under non-canonical conditions. A sound
understanding of how the inner- and outer-scale relationship is affected by
these conditions is essential for generalising its application to a wider
scope of flows. Enhanced amplitude modulation has been observed in various
non-canonical flows, including boundary layers over rough walls [27, 28, 29,
30, 31] or permeable surfaces [32, 33], and boundary layers with modified
outer structures, such as energetic large-scale motions injected into the flow
via freestream turbulence [34, 35], upstream dynamic roughness [36] and
synthetic large-scale signals generated by plasma actuators [37], to name a
few. It is not yet well understood how these seemingly different flow
conditions all lead to a common increase in the amplitude modulation
coefficient, and we aim to bridge this gap by establishing a physics-based
quantitative relationship between the amplitude modulation coefficient and
IOIM parameters.

Figure 1: Schematic of a turbulent boundary layer over a rough-to-smooth
change in surface condition. The roughness transition occurs at $x_{0}$, and
$\hat{x}=x-x_{0}$ denotes the fetch downstream of the transition. Reproduced
from Li  _et al._ [38].

Furthermore, here we look into another scenario where the introduced large
scales in the outer layer have an energy distribution across scales similar to
that of a canonical boundary layer, and only the amplitude is intensified.
This is achieved by a sudden rough-to-smooth surface transition occurring in
the streamwise direction, as depicted in Fig. 1. Upstream of the transition, a
turbulent boundary layer develops on a rough wall with equivalent sand grain
roughness height $k_{s}$. Here, $x$ is the streamwise direction, $x_{0}$ is
the streamwise location of the surface transition and $\hat{x}\equiv x-x_{0}$
is the distance downstream of the transition. At $x=x_{0}$, the surface
switches to a smooth wall, while the boundary layer continues to evolve and
gradually adjusts to the new surface. The effect of the new surface condition
is first felt in the near-wall region of the boundary layer and then gradually
propagates to the interior of the flow [39]. The layer that separates the
modified near-wall region from the unaffected oncoming flow farther away from
the wall is generally referred to as the internal boundary layer (IBL) with a
thickness denoted by $\delta_{i}$. For more details on the observation and
modelling of the flow recovery, we refer to the works by Elliott [40], Antonia
and Luxton [41], Hanson and Ganapathisubramani [42], Rouhi  _et al._ [43], and
Li  _et al._ [44, 45].

Turbulent boundary layers over a rough-to-smooth change in the streamwise
direction offer a new perspective to further understand the physics of inner-
outer interactions. When normalised by the friction velocity at the wall, LSMs
and VLSMs above $\delta_{i}$ are similar to their smooth-wall counterparts
[46, 47], but they are over-energised compared to the near-wall small scales.
The former retain a memory of the upstream rough wall friction velocity
whereas the latter scale on the much lower local smooth-wall friction
velocity. How these structurally similar but more energised large-scale
motions interact with the near-wall cycle will be investigated in this study.

The remainder of the paper is organised as follows: in Sec. II, we first
present a summary of the definition of the amplitude modulation coefficient
and the IOIM framework, as well as a quantitative relation between the two. We
then review the previous studies on the amplitude modulation in non-canonical
flows in Sec. III. In Sec. IV, we present a case study on the flow downstream
of a rough-to-smooth change in the surface condition.

##  II Physical underpinning of amplitude modulation

In this section, we first briefly summarise the definition of the amplitude
modulation coefficient $R$ and the IOIM formulation, and then provide a
quantitative description of the relation between the two, supported by results
computed from synthetic signals.

###  II.1 Amplitude modulation coefficient and IOIM formulation

The amplitude modulation coefficient $R$ is a commonly reported diagnostics,
largely due to the fact that it provides a straightforward quantification of
the degree of amplitude modulation within a single-point time series. An
example of the $R$ profile is shown in Fig. 2(a). It is defined as the
correlation between the low-pass-filtered envelope of small-scale fluctuations
and the large-scale fluctuations at the same location [14]:

$R(z^{+})=\frac{\left<E_{L}\left[u_{d}^{+}\left(z^{+},t^{+}\right)\right]u_{S}^{+}\left(z^{+},t^{+}\right)\right>}{\sqrt{\left<E_{L}^{2}\left[u_{d}^{+}\left(z^{+},t^{+}\right)\right]\right>}\sqrt{\left<u_{S}^{+2}\left(z^{+},t^{+}\right)\right>}}$
(1)

Here, the superscript $(\cdot)^{+}$ indicates inner scaling with the local
friction velocity as the velocity scale, the angle brackets
$\left<\cdot\right>$ denotes time average, $u$ represents the streamwise
velocity, $u_{S}^{+}$ is the zero-mean large-scale superposition signal, which
is usually obtained by low-pass filtering the time series with a threshold of
$\lambda_{x}^{+}=7000$ and $u_{d}^{+}\equiv u^{+}-u_{S}^{+}$ is the detrended
signal. $E_{L}[\cdot]$ denotes a low-pass-filtered envelope of the signal

$E_{L}\left[u_{d}^{+}\left(z^{+},t^{+}\right)\right]=\mathcal{L}\left[\sqrt{u_{d}^{+2}\left(z^{+},t^{+}\right)+\mathcal{H}^{2}\left[u_{d}^{+}\left(z^{+},t^{+}\right)\right]}\right],$
(2)

where $\mathcal{H}[\cdot]$ is the Hilbert transform,
$\sqrt{u_{d}^{+2}\left(z^{+},t^{+}\right)+\mathcal{H}^{2}\left[u_{d}^{+}\left(z^{+},t^{+}\right)\right]}$
is the analytic signal and $\mathcal{L}[\cdot]$ denotes a low-pass filter.

(a)(c)(b)(d)$\lambda_{x}^{+}=7000$ Figure 2: (a) Amplitude modulation
coefficient $R(z^{+})$ (Eq. 1) of a smooth-wall turbulent boundary layer with
$Re_{\tau}=1.3\times 10^{5}$ [48]. (b) Isocontour of the time-shifted
amplitude modulation coefficient $R_{\tau}(z^{+},\tau^{+})$ (Eq. 3) of the
same boundary layer profile. The solid contour lines are from 0.1 to 0.5 with
a step of 0.1, and the dashed contour lines are from $-0.3$ to $-0.1$ with a
step of 0.1. The solid and filled circles mark the $\tau^{+}$ values where
$R_{\tau}$ reaches its maximum and minimum at each wall-normal location,
respectively. (c) Premultiplied energy spectra of the measured velocity
$u^{+}$, universal small scales $u^{*}$ and superposition $u_{S}^{+}$ at
$z^{+}\approx 10$, and $u_{O}^{+}$, the velocity measured by the outer probe
at $z_{O}^{+}=469$. (d) The gain of the linear kernel $|\widetilde{H}_{L}|$
and the linear coherence spectrum $\gamma^{2}(f^{+})$ at $z^{+}\approx 10$.
The vertical dotted line in (c) and (d) marks the cut-off wavelength
$\lambda_{x}^{+}=7000$ (with the mean velocity at the outer probe as the
convective velocity) used by Mathis  _et al._ [20] to separate large and small
scales.

The definition can be further generalised to a time-shifted amplitude
modulation coefficient

$R_{\tau}(z^{+},\tau^{+})=\frac{\left<E_{L}\left[u_{d}^{+}\left(z^{+},t^{+}\right)\right]u_{S}^{+}\left(z^{+},t^{+}-\tau^{+}\right)\right>}{\sqrt{\left<E_{L}^{2}\left[u_{d}^{+}\left(z^{+},t^{+}\right)\right]\right>}\sqrt{\left<u_{S}^{+2}\left(z^{+},t^{+}\right)\right>}}.$
(3)

By definition, $R(z^{+})\equiv R_{\tau}(z^{+},0)$ (see Fig. 2b). The relative
shift $\tau_{a}^{+}$ is the lag between the superposition imprint $u_{S}^{+}$
and the low-pass filtered envelope $E_{L}\left[u_{d}^{+}\right]$ such that
$R_{\tau}$ reaches its maximum, i.e.
$R_{\tau}(z^{+},\tau_{a}^{+})=\text{max}\left[R_{\tau}(z^{+},\tau^{+})\right]$,
as marked by the solid circles in Fig. 2(b).

According to the IOIM [22], the statistical prediction of the fluctuating
velocity $u_{p}^{+}$ can be constructed by considering a superposition effect
of large-scale content with additively a universal signal $u^{*}$ that is
subject to an amplitude modulation

$u_{p}^{+}\left(z^{+},t^{+}\right)=\underbrace{u^{*}\left(z^{+},t^{+}\right)\left\\{1+\Gamma\left(z^{+}\right)u_{S}^{+}\left(z^{+},t^{+}-\tau_{a}^{+}\right)\right\\}}_{\text{amplitude
modulation}}+\underbrace{u_{S}^{+}\left(z^{+},t^{+}\right)}_{\text{superposition}}.$
(4)

Here, $\Gamma$ is the amplitude sensitivity, and the large-scale imprints
$u_{S}^{+}$ can be found from a given outer-layer large-scale signal
$u_{O}^{+}\left(t^{+}\right)$ via

$u_{S}^{+}\left(z^{+},t^{+}\right)=\mathcal{F}^{-1}\left[\widetilde{H}_{L}\left(z^{+},f^{+}\right)\mathcal{F}\left[u_{O}^{+}\left(t^{+}\right)\right]\right],$
(5)

where $\mathcal{F}[\cdot]$ and $\mathcal{F}^{-1}[\cdot]$ represent a Fourier
transform and inverse Fourier transform, respectively, and
$\widetilde{H}_{L}\left(f^{+}\right)$ is the linear transfer kernel
incorporating the large-scale coherence between the near-wall and outer
regions.

The gain of the linear kernel $|\widetilde{H}_{L}|$ is related to the linear
coherence spectrum $\gamma^{2}(f^{+})$ via

$|\widetilde{H}_{L}(f^{+})|=\sqrt{\gamma^{2}(f^{+})\frac{\left<\left|\mathcal{F}\left[u^{+}\right]\right|^{2}\right>}{\left<\left|\mathcal{F}\left[u_{O}^{+}\right]\right|^{2}\right>}},$
(6)

where $\gamma^{2}(f^{+})$ is given by Bendat and Piersol [49]

$\gamma^{2}(f^{+})=\frac{\left|\mathcal{F}\left[u_{O}^{+}\right]\overline{\mathcal{F}\left[u^{+}\right]}\right|^{2}}{\left<\left|\mathcal{F}\left[u_{O}^{+}\right]\right|^{2}\right>\left<\left|\mathcal{F}\left[u^{+}\right]\right|^{2}\right>},$
(7)

with $\overline{\left(\cdot\right)}$ denoting the complex conjugate.

As shown in Figs. 2(c) and (d), only the large-scale energy of the outer
velocity signal $u_{O}^{+}$ is retained in the superposition signal
$u_{S}^{+}$, which contributes to the large-scale end of the energy spectrum
of $u^{+}$. The small-scale end of the spectrum, on the other hand, is mainly
from the universal small scale, $u^{*}$. The linear transfer kernel
$|\widetilde{H}_{L}|$ enables a smooth roll-off of the coherence from large to
small scales, and the scale separation here is around the commonly used cut-
off threshold of $\lambda^{+}_{x}=7000$.

###  II.2 Quantitative relationship between $R$ and IOIM parameters

In this subsection, we derive a quantitative relationship of the amplitude
modulation coefficient by expressing $R$ in the framework of IOIM. Following
Duvvuri and McKeon [36], we work with a modified expression of the amplitude
modulation coefficient

$R_{2}(z^{+})=\frac{\left<E_{2L}\left[u_{d}^{+}\left(z^{+},t^{+}\right)\right]u_{S}^{+}\left(z^{+},t^{+}\right)\right>}{\sqrt{\left<E_{2L}^{2}\left[u_{d}^{+}\left(z^{+},t^{+}\right)\right]\right>}\sqrt{\left<u_{S}^{+2}\left(z^{+},t^{+}\right)\right>}},$
(8)

where

$E_{2L}\left[v\left(t^{+}\right)\right]=\mathcal{L}\left[v^{2}\left(t^{+}\right)+\mathcal{H}^{2}\left[v\left(t^{+}\right)\right]\right]$
(9)

for an arbitrary time series $v\left(t^{+}\right)$. The modified coefficient
$R_{2}$ uses the square of the analytic signal to avoid the difficulty in
dealing with the square root in Eq. (2), leading to a simpler mathematical
expression, and given that the amplitude modulation coefficient is a
normalised measure, no significant difference is expected in the values of $R$
and $R_{2}$. This can be easily verified using experimental data.

By expressing the velocity signals $u_{d}^{+}$ and $u_{S}^{+}$ in a series of
Fourier modes, and with some trigonometric manipulations, Duvvuri and McKeon
[36] showed that

$R_{2}=\frac{2\left<u_{d}^{+2}u_{S}^{+}\right>}{\sqrt{\left<E_{2L}^{2}\left[u_{d}^{+}\right]\right>}\sqrt{\left<u_{S}^{+2}\right>}}.$
(10)

We substitute the detrended signal expressed using the notations of IOIM (i.e.
$u_{d}^{+}=u^{*}\left[1+\Gamma
u_{S}^{+}\left(t^{+}-\tau_{a}^{+}\right)\right]$) into Eq. (10) and noting
that $\left<u^{*2}u_{S}^{+}\right>=0$ because $u^{*}$ and $u_{S}^{+}$ are
uncorrelated by definition, $R_{2}$ can then be expressed as

$R_{2}=\frac{4\Gamma\left<u^{*2}\right>\left<u_{S}^{+}\left(t^{+}\right)u_{S}^{+}\left(t^{+}-\tau_{a}^{+}\right)\right>+2\Gamma^{2}\left<u^{*2}\right>\left<u_{S}^{+}\left(t^{+}\right)u_{S}^{+2}\left(t^{+}-\tau_{a}^{+}\right)\right>}{\sqrt{\left<\left\\{E_{2L}\left[u^{*}\right]+2\Gamma\mathcal{L}\left[u^{*2}u_{S}^{+}+\mathcal{H}\left[u^{*}\right]\mathcal{H}\left[u^{*}u_{S}^{+}\right]\right]+\Gamma^{2}E_{2L}\left[u^{*}u_{S}^{+}\right]\right\\}^{2}\right>}\sqrt{\left<u_{S}^{+2}\right>}}$
(11)

Comparing the orders of $u^{*}$, $u_{S}^{+}$ and $\Gamma$ in the terms in the
denominator of Eq. (11) and noting that both $\mathcal{L}[\cdot]$ and
$\mathcal{H}[\cdot]$ are linear operators give rise to

$\displaystyle E_{2L}\left[u^{*}\right]$
$\displaystyle\sim\left<u^{*2}\right>,$ (12a)
$\displaystyle\Gamma\mathcal{L}\left[u^{*2}u_{S}^{+}+\mathcal{H}\left[u^{*}\right]\mathcal{H}\left[u^{*}u_{S}^{+}\right]\right]$
$\displaystyle\sim\Gamma\left<u_{S}^{+2}\right>^{1/2}\left<u^{*2}\right>,$
(12b) $\displaystyle\Gamma^{2}E_{2L}\left[u^{*}u_{S}^{+}\right]$
$\displaystyle\sim\Gamma^{2}\left<u_{S}^{+2}\right>\left<u^{*2}\right>.$ (12c)

Typically, $\Gamma$ is a relatively small number ranging from
$\mathcal{O}(0.01)$ to $\mathcal{O}(0.10)$ [22], and
$\left<u_{S}^{+2}\right><\left<u^{*2}\right>$ especially in the near-wall
region of $z^{+}\lesssim 200$ which is the current focus. Therefore, terms
containing higher orders of $\Gamma\left<u_{S}^{+2}\right>^{1/2}$ can be
neglected, and Eq. (11) is then reduced to

$R_{2}=\frac{4\Gamma\left<u^{*2}\right>\left<u_{S}^{+}\left(t^{+}\right)u_{S}^{+}\left(t^{+}-\tau_{a}^{+}\right)\right>}{\sqrt{\left<E^{2}_{2L}\left[u^{*}\right]\right>}\sqrt{\left<u_{S}^{+2}\right>}}$
(13)

We identify the following 3 parameters that can contribute to a change in
$R_{2}$:

  1. 1.

The relative shift $\tau_{a}^{+}$.

  2. 2.

The amplitude sensitivity $\Gamma$.

  3. 3.

The amplitude of large-scale imprints $\left<u_{S}^{+2}\right>$.

Interestingly, $R_{2}$ is not affected by changes in the amplitude of
universal small scales $\left<u^{*2}\right>$. This is because both the
numerator and the denominator in Eq. (11) (and Eq. (13) as well) contain the
same order of $\left<u^{*2}\right>$, which eventually cancel out. Note that
changes in the energy distribution across scales in $u_{S}^{+}$ or $u^{*}$
have more complicated consequences: for $u_{S}^{+}$, it will affect the auto-
correlation term
$\left<u_{S}^{+}\left(t^{+}\right)u_{S}^{+}\left(t^{+}-\tau_{a}^{+}\right)\right>$
with a given $\tau_{a}^{+}$, and for $u^{*}$, it will affect how much energy
remains in $E_{2L}\left[u^{*}\right]$ after a low-pass filter is applied.
Therefore, we limit the quantitative analysis to changes in the amplitude of
$u_{S}^{+}$ and $u^{*}$ fluctuations.

In canonical smooth-wall turbulent boundary layers, the IOIM parameters
($u^{*}$, $\Gamma$, $\tau_{a}^{+}$, and $\widetilde{H}_{L}$) are Reynolds
number invariant over the range of $Re_{\tau}\approx 7350-13300$ tested in the
calibrations [50], and only a change in the large-scale imprint amplitude
$\left<u_{S}^{+2}\right>$ was observed when varying the Reynolds number.
However, $\tau^{+}_{a}$, $\Gamma$ and $\left<u_{S}^{+2}\right>$ can all be
modified under non-canonical conditions, and it is important to understand how
they individually contribute to the overall amplitude modulation coefficient
$R$.

####  II.2.1 Effect of the relative shift $\tau_{a}^{+}$

We can conclude from Eq. (13) that $R_{2}$ would increase with a decreasing
$\left|\tau_{a}^{+}\right|$, as a result of reduced time lag, or improved
alignment between modulation and superposition signals. Typically, the
relative shift $\tau_{a}^{+}$ is much smaller than the time period of
$u_{S}^{+}$ below the centre of the log region, where the predictive model is
applied. Therefore, $u_{S}^{+}\left(t^{+}\right)$ and
$u_{S}^{+}\left(t^{+}-\tau_{a}^{+}\right)$ are largely in phase. The auto-
correlation term
$\left<u_{S}^{+}\left(t^{+}\right)u_{S}^{+}\left(t^{+}-\tau_{a}^{+}\right)\right>$
is positive, and increases with a decreasing $\left|\tau_{a}^{+}\right|$.

Fig. 3(a) shows $R$ (Eq. 1) computed from the velocity signals constructed
from $u^{*}$ and $u_{S}^{+}$ with a range of different $\tau_{a}^{+}$ values
following (4). The IOIM parameters and the baseline values of $\tau_{a0}^{+}$,
$\Gamma_{0}$ and $\left<u_{S0}^{+2}\right>$ are taken from the dataset at
$Re_{\tau}=13300$ [22]. The decreasing trend of $R$ with an increasing
$\left|\tau_{a}^{+}\right|$ confirms the conclusion based on Eq. (13).

####  II.2.2 Effect of the amplitude sensitivity $\Gamma$

According to Eq. (13), $R_{2}$ will increase with an increasing $\Gamma$. This
trend is corroborated by the increasing $R$ computed from the signals
constructed following Eq. (4) with an increasing $\Gamma/\Gamma_{0}$ (see Fig.
3b).

####  II.2.3 Effect of the superposition intensity $\left<u_{S}^{+2}\right>$

Based on Eq. (13), $R_{2}$ will increase with an increasing
$\left<u_{S}^{+2}\right>$. In fact, a closer examination of the Eqs. (13) and
(11) reveals that $\Gamma$ and $\left<u_{S}^{+2}\right>$ can be grouped into a
single variable $\Gamma\left<u_{S}^{+2}\right>^{1/2}$. This is confirmed by
the same trends of $R$ with $\Gamma$ and $\left<u_{S}^{+2}\right>^{1/2}$ in
Figs. 3(b) and (c), respectively. We note that the monotonic increase of
$R_{2}$ with $\Gamma\left<u_{S}^{+2}\right>^{1/2}$ breaks down when
$\Gamma\left<u_{S}^{+2}\right>^{1/2}$ becomes comparable with or larger than
$\left<u^{*2}\right>^{1/2}$ (not shown in the figure), and the full expression
of Eq. (11) can introduce non-monotonic dependence on
$\Gamma\left<u_{S}^{+2}\right>^{1/2}$. However, given that $\Gamma$ is
typically small, and $\left<u_{S}^{+2}\right>$ is much smaller than
$\left<u^{*2}\right>$ in the near-wall region at the Reynolds number range
investigated in the experimental dataset ($Re_{\tau}\lesssim 2\times 10^{4}$),
the non-monotonic is less likely to occur.

To summarise, we have shown analytically that the amplitude modulation
coefficient will increase with (i) a reducing $\left|\tau_{a}^{+}\right|$,
(ii) an increasing $\Gamma$, and (iii) an increasing
$\left<u_{S}^{+2}\right>$, provided
$\Gamma\left<u_{S}^{+2}\right>^{1/2}\ll\left<u^{*2}\right>^{1/2}$.

(a)(b)(c) Figure 3: Amplitude modulation coefficient $R$ of the velocity
signal constructed with various (a) $\tau_{a}^{+}$, (b) $\Gamma$ and (c)
$\left<u_{S}^{+2}\right>$ values. The baseline parameters are denoted by
$\tau_{a0}^{+}$, $\Gamma_{0}$ and $\left<u_{S0}^{+2}\right>$, respectively,
and they are computed from the dataset at $Re_{\tau}=13300$ in Baars  _et al._
[22]. The arrows in the panels indicate the direction of increasing the ratio
between the varied parameter and the baseline. increasing $Re_{\tau}$ Figure
4: Amplitude modulation coefficients $R$ (grey circles) and $R_{2}$ (blue
triangles) of a canonical smooth-wall turbulent boundary layer plotted against
the corresponding superposition signal intensity $\left<u_{S}^{+2}\right>$ at
$z^{+}=50$. From pale to dark blue, the shade of blue triangles indicates
$Re_{\tau}$ increasing from 2800 to 13400. The solid black line is Eq. (13)
with a zero time lag ($\tau_{a}^{+}=0$), and the dot-dashed blue line is Eq.
(13) with the auto-correlation term
$\left<u_{S}^{+}\left(t^{+}\right)u_{S}^{+}\left(t^{+}-\tau_{a}^{+}\right)\right>$
estimated from the experimental data.

###  II.3 Verification with a smooth-wall turbulent boundary layer dataset

The analysis above is verified using an experimental smooth-wall turbulent
boundary layer dataset with $Re_{\tau}$ ranging from 2800 to 13400 [48]. In
this series of canonical smooth-wall boundary layer profiles, both
$\tau_{a}^{+}$ and $\Gamma$ are expected to remain constant, and only
$\left<u_{S}^{+2}\right>$ increases with an increasing $Re_{\tau}$ as a
consequence of more energetic large-scale motions in the outer layer,
providing an ideal test ground to examine the dependence of $R$ on
$\left<u_{S}^{+2}\right>$. Fig. 4 shows the amplitude modulation coefficients
$R$ and $R_{2}$ computed directly from the velocity time series. The large-
scale imprint intensity $\left<u_{S}^{+2}\right>$ is approximated by
$\left<u^{+2}\right>-\left<u^{*2}\right>$, where $\left<u^{*2}\right>$ is the
intensity of the universal small-scale signal from the calibration data of
Baars  _et al._ [22]. The large-scale imprint intensity
$\left<u_{S}^{+2}\right>$ increases with an increasing $Re_{\tau}$ as
expected. The two coefficients $R$ and $R_{2}$ have very similar values, and
they exhibit an increasing trend with $\left<u_{S}^{+2}\right>$, which is well
captured by Eq. (13). The auto-correlation term
$\left.\left<u_{S}^{+}\left(t^{+}\right)u_{S}^{+}\left(t^{+}-\tau_{a}^{+}\right)\right>\right/\left<u_{S}^{+2}\right>$
in Eq. (13) is smaller than 1 for a finite time shift $\tau_{a}^{+}$, but the
exact value need to be determined from the time series of $u_{S}^{+}$.
Assuming a perfect auto-correlation (which is essentially Eq. (13) with a zero
time lag, $\tau_{a}^{+}=0$) results in an overestimated $R_{2}$ as shown by
the black line. The auto-correlation term at $z^{+}\approx 50$ computed from
the experimental data increases slightly from 0.56 to 0.69 from the lowest to
the highest $Re_{\tau}$ measurements, because the relative phase shift
$\left|\tau_{a}^{+}\right|U^{+}/\lambda_{S}^{+}$ reduces with increasing
$Re_{\tau}$, with $\left|\tau_{a}^{+}\right|$ remaining constant and
$\lambda_{S}^{+}$, the most energetic wavelengths of $u_{S}^{+}$, growing with
the boundary layer thickness $\delta^{+}$ ($\equiv Re_{\tau}$). After taking
this into account, the prediction of Eq. (13) (dot-dashed blue line) captures
the trend of data points very well. In summary, the increase in $R$ and
$R_{2}$ with increasing $Re_{\tau}$ for a canonical smooth-wall turbulent
boundary layer is primarily originated from the growing
$\left<u_{S}^{+2}\right>$, and the slight increase in the auto-correlation
term also makes a small contribution to the growth.

##  III Revisiting scale interactions in previous studies

In this section, we present a summary of data from the literature and the
current work, focusing on how the IOIM parameters change in various flow types
and with increasing $Re_{\tau}$, and discuss their commonalities and
differences in the scale interaction mechanism.

Based on whether the near-wall cycle is modified from that of an impermeable
smooth wall, various flow types collected in Table 1 can be further classified
as ‘top-down’ and ‘bottom-up’ categories, which indicate whether the deviation
from a canonical flow is introduced in the outer layer (‘top-down’) or near to
the wall (‘bottom-up’). We note that although Duvvuri and McKeon [36] used a
dynamic roughness element upstream to generate the synthetic signal, the wall
remains smooth at the location where the scale interaction is quantified, and
therefore, this study is also classified as a ‘top-down’ type.

Notably, a higher amplitude modulation coefficient $R$ is found in all non-
canonical flows collected in the table. A closer examination of the IOIM
parameters would reveal that this increase in $R$ is caused by different
mechanisms in each category. For rough-wall flows, although the strength of
$\left<u_{O}^{+2}\right>$ is similar to that of a smooth-wall boundary layer,
the coherence between the inner and outer layer is reduced, presumably due to
the interruption of the shedded vortices in the roughness sublayer. The higher
$R$ value is contributed by the increased $\Gamma$, which can be explained as
follows [28]. In smooth-wall flows, small-scale fluctuations scale with
viscous units. Considering that the large-scale variation of friction velocity
$u_{\tau}$ is quasi-steady from the viewpoint of small scales, the small-scale
velocity fluctuations will then scale with $u_{\tau}$. By definition,
$u_{\tau}\propto\sqrt{\tau_{w}}$, where $\tau_{w}$ is the instantaneous wall-
shear stress and it is caused by and proportional to the large-scale velocity
signature $u_{S}$ [51, 12]. Thus, the small-scale velocity fluctuations are
proportional to $\sqrt{u_{S}}$. In fully-rough flows, however, the small-scale
fluctuations are the wake vortices originated from roughness elements, the
intensity of which is proportional to $u_{S}$. Therefore, the modulation
sensitivity in rough-wall flows is expected to be stronger than that in
smooth-wall flows. Following a similar argument, the near-wall region of a
permeable substrate is dominated by upwelling/downwelling associated with the
large-scale streamwise motions ($\propto u_{S}$) [32, 33], leading to an
increase in $R$ as well.

Regarding the ‘top-down’ flow category, more energetic large scales can be
introduced in the outer layer by increasing $Re_{\tau}$ in a canonical smooth-
wall boundary layer, through FST, or an upstream dynamic roughness element.
Another example in this category, to be discussed in Sec. IV, is that of a
rough-to-smooth change in the wall condition, where the introduced large
scales are the remnant of upstream rough-wall large-scale structures. The
near-wall cycle (represented by the intensity of $u^{*}$) is little modified
in all cases. The coherence between large scales in the inner and outer
regions, which manifests in $|\widetilde{H}_{L}|$, either remains the same
(canonical smooth wall with increasing $Re_{\tau}$), or increases (FST), and
they all lead to stronger footprints $\left<u_{S}^{+2}\right>$ close to the
wall. The stronger $\left<u_{S}^{+2}\right>$ is primarily responsible for the
increased $R$ in the ‘top-down’ flow cases.

| Flow type | Reference | $\left<u_{O}^{+2}\right>$ | $\left<u_{S}^{+2}\right>$ | $|\widetilde{H}_{L}|$ | $\left<u^{*2}\right>$ | $\Gamma$ |  $R$ (near wall)
---|---|---|---|---|---|---|---|---
{613mm[top-down] | R-to-S change | present study | higher | higher | higher | similar | similar | higher
|  | Smooth wall,
---
increasing $Re_{\tau}$
| Mathis  _et al._ [14, 20]
---
Baars  _et al._ [22]
higher | higher | similar | similar | similar | higher
| FST | Dogan  _et al._ [34, 35] | higher | higher | higher | similar | lower | higher
|  | Single frequency
---
large-scale input
Duvvuri and McKeon [36] | higher | higher | ? | similar | ? | ?
{515mm[bottom-up] | Rough wall |  | Squire  _et al._ [28]
---
Anderson [29]
Pathikonda and Christensen [52]
Blackman  _et al._ [31]
similar | lower | lower |  | dependent on
---
roughness
morphorlogy
higher | higher
| Permeable wall |  | Efstathiou and Luhar [53]
---
Kim  _et al._ [32]
similar/lower | ? | ? |  | dependent on
---
substrate
morphorlogy
? | higher
Table 1: Summary of IOIM parameters in various flow types. Quantities that cannot be inferred from the reported results are marked by ‘?’. Changes in the parameters are relative to a smooth-wall turbulent boundary layer at comparable Reynolds numbers. In the category of smooth-wall boundary layers with increasing Reynolds numbers, the comparison is made with regard to a lower Reynolds number. Results of the current study will be presented in details in Sec. IV. Case | Colour | $\hat{x}/\delta_{0}$ | $Re_{\tau}$ | $3.9\sqrt{Re_{\tau}}$ | $U_{\infty}$ | $U_{\tau}$ | $z_{O}^{+}$ | $\delta_{i}^{+}$ | $l_{I}^{+}$ | $l_{O}^{+}$
---|---|---|---|---|---|---|---|---|---|---
|  |  |  |  | ($\mathrm{ms}^{-1}$) | ($\mathrm{ms}^{-1}$) |  |  |  |
Smooth | black | - | 7350 | 334 | 10.0 | 0.34 | 441 | - | 22 | 22
R-to-S1 | green | 2.3 | 7200 | 330 | 22.5 | 0.69 | 331 | 890 | 21 | 21
R-to-S2 | magenta | 2.3 | 9600 | 382 | 31.0 | 0.92 | 388 | 1180 | 19 | 31
Table 2: Experimental parameters of the two-point hotwire measurements, where
the outer probe remains at a fixed location $z_{O}^{+}$ and the inner probe is
traversed between $[0,z_{O}^{+}]$. $l_{I}^{+}$ and $l_{O}^{+}$ are the
viscous-scaled hotwire filament length of the inner and outer probes,
respectively.

##  IV Case study: turbulent boundary layer following a step change in surface
roughness

In this section, we present the new experimental results of a turbulent
boundary layer downstream of a rough-to-smooth change in the surface
conditions, which is another non-canonical flow configuration in the ‘top-
down’ category. We will briefly present the experimental conditions in Sec.
IV.1, and then explore the detailed modulation behaviours using both single-
probe and two-probe results in Secs. IV.2–IV.3. The results are briefly
discussed in Sec. IV.4.

###  IV.1 Experimental setup

Simultaneous two-probe hotwire anemometry measurements are performed in the
High Reynolds Number Boundary Layer Wind Tunnel (HRNBLWT) at the University of
Melbourne. An upstream portion of the 27 m working section floor is covered by
P24 sandpaper (SP40F, Awuko Abrasives) from the inlet to $x_{0}=7.2$ m, while
the remaining length is a smooth aluminium surface. The peak-to-trough
roughness height is $k_{p}\approx 1.2$ mm, and the equivalent sand grain
roughness is $k_{s}\approx 2.43$ mm. A nominal zero-pressure gradient is
achieved by adjusting the bleeding slots on the tunnel roof. More details of
the facility can be found in [44] and [54].

Two-probe hotwire measurements are performed at two freestream velocities,
22.5 m s-1 and 31.0 m s-1, and these two cases are named as R-to-S1 and
R-to-S2, respectively. The flow conditions of these cases correspond to the
single-point dataset of cases Re10ks16 and Re14ks22 in [44], and the local
friction velocity $U_{\tau}$ of the current cases is interpolated from the
skin-friction versus $\hat{x}$ trajectory in the same study, which was
measured directly at the wall using oil-film interferometry. The procedures of
the two-probe hotwire anemometry measurements are similar to that described by
Mathis  _et al._ [20] and Baars  _et al._ [22]. The outer probe is fixed at
$z_{O}^{+}=3.9\sqrt{Re_{\tau}}$, which is the geometric centre of the
logarithmic region where the large-scale motions are highly active. The inner
probe is traversed from the wall to below the outer probe with approximately
20 logarithmically spaced points in between. Both probes are conventional
single-wire hotwire probes with a Wollaston wire etched to expose the sensing
element. The length-to-diameter ratio of the exposed filament is $l/d\geq 200$
[55]. For the case R-to-S1, both wires have a diameter of $d=2.5\upmu$m,
leading to a viscous-scaled filament length of $l_{I}^{+},l_{O}^{+}\approx
21$. For the case R-to-S2, $d=1.5\upmu$m wire is selected for the inner probe
to maintain a similar spatial resolution with $l_{I}^{+}\approx 19$ with an
increased freestream (and friction) velocity, while the outer probe filament
diameter remains $d=2.5\upmu$m because only the large-scale signal at this
location is of the interest. Both probes are conventional single-wire hotwire
probes operated by an in-house Melbourne University Constant Temperature
Anemometer (MUCTA). The hotwire sampling time $T_{s}$ is more than 20000
boundary layer turn-over time ($\delta_{99}/U_{\infty}$) to achieve a good
convergence of the statistics. A two-probe smooth-wall dataset from Mathis
_et al._ [20] and Baars  _et al._ [22] is also included for comparison in this
study. Parameters of the two-point measurements are summarised in Table 2.

(a)(b) Figure 5: Profiles of mean streamwise velocity and turbulence intensity
corresponding to the flow conditions of R-to-S1. (a) is inner scaled using the
local smooth-wall $U_{\tau}$, while (b) is outer scaled, but using $U_{\tau
0}$, the friction velocity measured just upstream of the rough-to-smooth
transition. Line colours indicate the fetch, from orange to green to blue
corresponds to $\hat{x}/\delta_{0}=0.2$, 0.5, 0.9, 1.9, 4.2, 7.4, 14.8, 29.2
and 53.1. The solid black line is a smooth-wall reference with
$Re_{\tau}=1.0\times 10^{4}$ acquired in the same facility and normalised
using the corresponding smooth-wall friction velocity [48], and dashed black
line is a rough-wall reference acquired just upstream of the rough-to-smooth
change and normalised by the rough-wall friction velocity. The white circles
represent the edge of the IBL, determined from the variance profile [44]. The
profile highlighted by triangular symbols is close to the streamwise location
where the two-probe measurements in this study are performed.

###  IV.2 Single-probe results

We first present the evolution of the mean velocity and turbulence statistics
downstream of the rough-to-smooth change in Fig. 5. Generally speaking,
immediately after the rough-to-smooth change, the rough-wall turbulent
boundary layer starts adapting to the new surface conditions first at the
wall, and as the modified region enlarges, to the interior of the flow. The
skin-friction coefficient experiences an undershoot before gradually
increasing to the smooth-wall value. Fig. 5(a) shows the inner-scaled mean
velocity profiles, where the velocity scale is selected as the friction
velocity $U_{\tau}$ measured locally with oil-film interferometry over the
smooth surface. Note that in this paper, the superscript $(\cdot)^{+}$
indicates inner scaling with the local friction velocity as the velocity
scale. Fig. 5(b) shows the outer-scaled profiles, where the velocity scale is
chosen as the friction velocity $U_{\tau 0}$ measured on the rough wall just
upstream of the rough-to-smooth transition at $\hat{x}\rightarrow 0^{-}$. The
inner-scaled mean velocity profiles of Fig. 5(a) collapse with the smooth-wall
reference first close to the wall, while the outer-scaled profiles in Fig.
5(b) agree well with the rough-wall reference above the IBL (marked by the
open circles). In the turbulence intensity profiles, a strong ‘outer-peak’
manifests at the IBL, which is a result of the remaining energetic rough-wall
structures. These structures leave a strong footprint in the near-wall region,
as evidenced by the increased inner-peak magnitude in wall units shown in Fig.
5(a). The case R-to-S2 has a higher $Re_{\tau 0}$ and $k_{s0}^{+}$ compared to
the case R-to-S1, thus, stronger large-scale motions (which scale on $U_{\tau
0}$) above the internal boundary layer are expected in the former.

Figure 6: Amplitude modulation coefficient $R$ at various streamwise locations
downstream of a rough-to-smooth change. Line colours indicate the fetch, from
red to blue corresponds to $\hat{x}/\delta_{0}=0.08,$ 0.2, 0.5, 0.9, 1.9, 4.2,
7.4, 14.8, 29.2 and 53.1. The white circles represent the edge of the IBL. The
solid black line is a smooth-wall reference with $Re_{\tau}=1.0\times 10^{4}$
[48]. (a)(b) Figure 7: (a) Isocontours of $R_{\tau}$ at
$\hat{x}/\delta_{0}=0.2$ (orange) and 0.9 (green) and the smooth-wall
reference (black). The solid contour lines are from 0.1 to 0.5 with a step of
0.1, and the dashed contour lines are from $-0.3$ to $-0.1$ with a step of
0.1. The solid and empty circles mark the $\tau^{+}$ values where $R_{\tau}$
reaches its maximum and minimum at each wall-normal location, respectively.
(b) Maximum $R_{\tau}$ values, which are essentially $R_{\tau}$ at the
locations marked by solid circles in (a). (a)(b) Figure 8: (a) Time shift
required for $R_{\tau}$ to reach its maximum. (b) Argument of the complex co-
spectrum of $u_{S}^{+}$ and $E_{L}[u_{d}^{+}]$ at $z^{+}=50$ (marked in (a) by
the vertical dashed line), premultiplied by the time period $1/f^{+}$. Legends
are the same as in Fig. 6. The solid black line is a smooth-wall reference
with $Re_{\tau}=1.3\times 10^{4}$ [48].

The amplitude modulation coefficient $R$ (Eq. 1) at a range of downstream
locations from $\hat{x}/\delta_{0}=0.08$ to 53.1 is shown in Fig. 6. The cut-
off wavelength of the large-scale filter $\mathcal{L[\cdot]}$ is computed
using the local friction velocity. Only the results of R-to-S1 are presented
here for brevity, while the R-to-S2 case shows similar behaviour. At small
$\hat{x}/\delta_{0}$, a high $R$ is observed in the near-wall region. The
coefficient $R$ decreases with an increasing $\hat{x}/\delta_{0}$, and beyond
$\hat{x}/\delta_{0}>20$, it becomes very similar to that of a smooth-wall
boundary layer.

Isocontours of $R_{\tau}$ (Eq. 3) at $\hat{x}/\delta_{0}=0.2$ (orange) and 0.9
(green) are shown in Fig. 7. In addition to the difference in the magnitude,
at small fetches, the contours of positive correlations also shifts to the
positive $\tau^{+}$ direction. In other words, close to the rough-to-smooth
change, there is a smaller time lag between the envelope of the small-scale
fluctuations and large-scale motions. The optimal positive correlation
achieved at the time shift $\tau_{a}^{+}$ is shown in Fig. 7(b). Similar to
the zero-time-shift $R$ (shown in Fig. 6), the maximum $R_{\tau}$ values are
also higher at smaller fetches, implying that the high magnitudes of $R$
observed in Fig. 6 is more than the consequence of a smaller lag between the
envelope and large-scale motions.

A comparison of $\tau_{a}^{+}$ (time shift required for the optimum positive
correlation) can be found in Fig. 8(a), where $\tau_{a}^{+}$ becomes more
negative with an increasing fetch. Fig. 8(b) shows the argument of the complex
co-spectrum $\phi$ at $z^{+}=50$, premultiplied by the time period $1/f^{+}$,
following Jacobi and McKeon [56, 57] and Deshpande  _et al._ [58]. The co-
spectrum is defined as
$\phi\equiv\left<\mathcal{F}\left[u_{S}^{+}\right]\overline{\mathcal{F}\left[E_{L}[u_{d}^{+}]\right]}\right>$,
and it can be viewed as the spectral equivalent of the amplitude modulation
coefficient $R$. The premultiplied argument
$\tau_{a}^{+}(f^{+})\equiv\text{arg}(\phi)/(2\pi f^{+})$ is essentially the
time shift between $u_{S}^{+}$ and $E_{L}[u_{d}^{+}]$ in each Fourier mode,
and it is reasonable that the overall time shift $\tau_{a}^{+}$ falls in the
same range as $\tau_{a}^{+}(f^{+})$ at each corresponding measurement
location. In addition, the increase of the time lag with increasing fetch is
also apparent here, confirming the trend of $\tau_{a}^{+}$ in Fig. 8(a).
Further, the absolute time lag $\left|\tau_{a}^{+}(f^{+})\right|$ is smaller
at higher frequencies, meaning that the amplitude of small scales are more in-
phase with the higher-frequency modes of the large scales. The increased
maximum $R_{\tau}$ values at smaller fetches can partially be explained by the
smaller scatter of $\tau_{a}^{+}(f^{+})$ across a range of frequencies,
because a single time shift $\tau_{a}^{+}$ can better align all Fourier modes
in the $u_{S}^{+}$ and $E_{L}[u_{d}^{+}]$ signals in these cases.

The increase in $R$ compared to the smooth-wall reference has been previously
observed in rough-wall flows and attributed to the stronger correlation
between the amplitude of the small-scale turbulence associated with the
roughness elements and the large-scale motions [28]. Compared to the most
downstream location, $R$ in the logarithmic region is still noticeably higher
at $\hat{x}/\delta_{0}=0.9$, which is equivalent to $\hat{x}/k_{p}\approx
110$, a fetch where we might expect a large portion of the small-scale motions
directly generated from the flow interaction with the roughness elements to
diminish. However, small scales may form through the shear between the
surviving rough-wall structures, and exhibit a stronger amplitude modulation
effect with those structures from which they originate.

In summary, based on the single-point measurements of the rough-to-smooth
cases, the increase in $R$ is contributed by both reduced $|\tau_{a}^{+}|$ and
increased $\left<u_{S}^{+2}\right>$, which are two out of the three factors
identified in the analysis in Sec. II. Further two-probe measurements are
required to quantify the effect of the modulation sensitivity $\Gamma$.

(a)(b)(c)
(d)(e)(f)

Figure 9: Premultiplied spectra of (a) $u^{+}$ and (d) $u^{*}$. The contour
levels in both plots are from 0.2 to 2 with a spacing of 0.2. Line colours
black, green and magenta represent the smooth, R-to-S1 and R-to-S2 cases,
respectively. (b) and (c) are the difference in $f^{+}\phi^{+}_{uu}$ between
the R-to-S1 and R-to-S2 cases and the smooth-wall reference, and (e) and (f)
are the corresponding difference in $f^{+}\phi^{+}_{u^{*}u^{*}}$.

###  IV.3 Two-probe results

To further understand the origin of the enhanced modulation effect observed in
Sec. IV.2, IOIM calibration following the procedure detailed in Baars  _et
al._ [22] is performed on the three cases listed in Table 2. In the results
below, we will be using black for the smooth-wall reference, green and magenta
for the R-to-S1 and R-to-S2 cases, respectively.

Premultiplied energy spectra of the measured velocity fluctuation and the
universal small-scale signal are shown in Figs. 9(a) and (d), respectively.
Figs. 9(b) and (c) are the difference between the rough-to-smooth and smooth-
wall reference. A band with excess energy at $1/f^{+}\approx 2000$ is
interpreted as the large-scale footprints in the near-wall region. The
universal small scale spectra of the three cases (Fig. 9d) are very similar,
indicating a re-establishment of the near-wall cycle after the rough-to-smooth
change.

Figure 10: Amplitude sensitivity $\Gamma(z^{+})$. Line colours black, green
and magenta represent the smooth, R-to-S1 and R-to-S2 cases, respectively.

(a)(b)

Figure 11: (a) $|\widetilde{H}_{L}|$, gain of the linear kernel and (b)
$\gamma^{2}$, linear coherence spectra between the fixed outer and moving
inner probe at $z^{+}=100$. Both quantities are filtered by a 25% bandwidth
moving filter. Line colours black, green and magenta represent the smooth,
R-to-S1 and R-to-S2 cases, respectively. The premultiplied energy spectrum
$f^{+}\phi_{uu}^{+}$ is also shown in (b) by thin lines of corresponding
colours on the left vertical axis for reference.

The sensitivity of small scales to amplitude modulation is indicated by
$\Gamma$: for a higher $\Gamma$, the universal small-scale signal $u^{*}$ will
be multiplied by a higher fraction of the superposition signal $u_{S}^{+}$ to
generate the prediction. Two calibrations of a smooth-wall turbulent boundary
layer at $Re_{\tau}\approx 7350$ and 13300 result in very similar $\Gamma$
[22]. For rough-to-smooth cases, as shown in Fig. 10, $\Gamma$ from the three
calibrations reach a good overall agreement, suggesting that the over-
energised rough-wall structures in the outer layer do not seem to alter the
amplitude modulation mechanism in the near-wall region.

Fig. 11(a) shows $|\widetilde{H}_{L}|$, the gain of the linear kernel (which
relates the superposition $u_{S}^{+}$ to the outer-layer large-scale signal
$u_{O}^{+}$ via Eq. 5), at a wall-normal position of $z^{+}=100$. Similar
trends are also seen in other wall-normal positions, and are not shown here
for brevity. The magnitude of the linear transfer kernel $|\widetilde{H}_{L}|$
is found to increase in the rough-to-smooth cases compared to the smooth-wall
reference. The linear coherence spectra are shown in Fig. 11(b). In the rough-
to-smooth cases, $\gamma^{2}$ deviates from 0 at a smaller $1/f^{+}$, and
remains higher than that of the smooth-wall reference. A higher $\gamma^{2}$
indicates a stronger correlation between the velocity fluctuations obtained by
the inner and outer probes. These results suggest that for a given structure
in the outer layer, it will leave a stronger footprint (superposition) in the
near-wall region in the rough-to-smooth case as a result of the enhanced
inner-outer coherence.

###  IV.4 Discussion

Overall, the amplitude modulation mechanism appears to be little modified
after introducing a roughness heterogeneity (at least at
$\hat{x}/\delta_{0}=2.3$ as examined here). However, there is a stronger
correlation between the large-scale velocity fluctuations obtained by the
inner and outer probes, and the gain in the linear kernel is also higher in
rough-to-smooth cases. The large-scale fluctuation $\left<u_{O}^{+2}\right>$
is already stronger in the rough-to-smooth cases, and a larger fraction of it
will contribute to the near-wall superposition signal through the increased
gain $|\widetilde{H}_{L}|$. In light of the analysis in Sec. II, the increase
of $\left<u_{S}^{+2}\right>$ is primarily responsible for the higher $R$ value
downstream of a rough-to-smooth change.

It is interesting to draw a direct comparison between the FST [34, 35] and the
current rough-to-smooth cases, as both are in the ‘top-down’ category with
broadband energetic large-scale motions imposed in the outer region. The
increased strength in the outer large scales makes them less susceptible to
the interruptions from near-wall motions, leading to an increase in the
correlation between outer- and near-wall large-scale signals, which is
eventually reflected in an higher $|\widetilde{H}_{L}|$ in both. The major
difference between the two is in the coefficient $\Gamma$: the former has a
lower $\Gamma$, while the value is unchanged in the latter, similar to the
independence of $\Gamma$ on $Re_{\tau}$ values observed in canonical smooth-
wall boundary layers [22]. We speculate that such difference is rooted in the
manner in which the outer large scales are organised, as well as the energetic
wavelengths. In the FST case, the outer structures are created by an active
grid, which are inherently different from the structures organised by hairpin
vortices in a developing boundary layer. The dissimilarity in the generation
mechanism between large- and small-scale structures in the FST might be the
reason for the reduction in $\Gamma$. The highly energetic large scales
arranged in a manner similar to that of a naturally developed boundary layer
in the current configuration makes it a good mimetic of a smooth-wall
turbulent boundary layer at very high Reynolds numbers.

##  V Concluding remarks

In this work, we first show the underlying physics behind the increased
amplitude modulation coefficient as reported in many previous studies,
utilising the framework of IOIM. An analytical relationship between the
amplitude modulation coefficient $R_{2}$ and IOIM parameters is derived and
verified using a smooth-wall turbulent boundary layer dataset. This framework
is then applied to classify and interpret the reported amplitude modulation
behaviours in previous works. We then present the case study of a turbulent
boundary layer downstream of a rough-to-smooth change with both single probe
and simultaneous two-probe measurements. A stronger amplitude modulation
effect evidenced by an increased $R$ is observed. Further analysis of the two-
probe data reveals that the modulation strength $\Gamma$ is similar to that of
a canonical smooth-wall turbulent boundary layer, and it remains the same for
the two Reynolds numbers tested. The increase in $R$ is primarily attributed
to the stronger large-scale footprints $\left<u_{S}^{+2}\right>$ in the near-
wall region, which is contributed by both the over-energetic outer layer
motions $\left<u_{O}^{+2}\right>$ and stronger coherence between the inner and
outer layers. These results and analyses offer a new perspective to interpret
the abundant literature on the the scale interactions of non-canonical
turbulent boundary layers, which can be meaningful for incorporating the IOIM
in the numerical simulation of a wide range of flow conditions.

## References

  * Balakumar and Adrian [2007] B. J. Balakumar and R. J. Adrian, Large- and very-large-scale motions in channel and boundary-layer flows, Phil. Trans. Roy. Soc. Lond. A: Math., Phys. Eng. Sci. 365, 665 (2007).
  * Kim and Adrian [1999] K. C. Kim and R. J. Adrian, Very large-scale motion in the outer layer, Phys. Fluids 11, 417 (1999).
  * Smits  _et al._ [2011] A. J. Smits, B. J. McKeon, and I. Marusic, High-Reynolds number wall turbulence, Annu. Rev. Fluid Mech. 43, 353 (2011).
  * Hutchins and Marusic [2007a] N. Hutchins and I. Marusic, Evidence of very long meandering features in the logarithmic region of turbulent boundary layers, J. Fluid Mech 579, 1 (2007a).
  * Hutchins and Marusic [2007b] N. Hutchins and I. Marusic, Large-scale influences in near-wall turbulence, Phil. Trans. Royal Soc. 365, 647 (2007b).
  * Hoyas and Jiménez [2006] S. Hoyas and J. Jiménez, Scaling of the velocity fluctuations in turbulent channels up to ${R}e_{\tau}=2003$, Phys. Fluids 18, 011702 (2006).
  * Vincenti  _et al._ [2013] P. Vincenti, J. Klewicki, C. Morrill-Winter, C. M. White, and M. Wosnik, Streamwise velocity statistics in turbulent boundary layers that spatially develop to high Reynolds number, Exp. Fluids 54, 1 (2013).
  * Bandyopadhyay and Hussain [1984] P. R. Bandyopadhyay and A. K. M. F. Hussain, The coupling between scales in shear flows, Phys. Fluids 27, 2221 (1984).
  * Agostini and Leschziner [2014] L. Agostini and M. A. Leschziner, On the influence of outer large-scale structures on near-wall turbulence in channel flow, Phys. Fluids 26, 075107 (2014).
  * He [2018] L. He, Multiscale block spectral solution for unsteady flows, Int. J. Numer. Methods Fluids 86, 655 (2018).
  * Chen and He [2023] C. Chen and L. He, Two-scale solution for tripped turbulent boundary layer, J. Fluid Mech. 955, A5 (2023).
  * Zhang and Chernyshenko [2016] C. Zhang and S. I. Chernyshenko, Quasisteady quasihomogeneous description of the scale interactions in near-wall turbulence, Phys. Rev. Fluids 1, 014401 (2016).
  * Chernyshenko [2021] S. Chernyshenko, Extension of QSQH theory of scale interaction in near-wall turbulence to all velocity components, J. Fluid Mech. 916, A52 (2021).
  * Mathis  _et al._ [2009a] R. Mathis, N. Hutchins, and I. Marusic, Large-scale amplitude modulation of the small-scale structures in turbulent boundary layers, J. Fluid Mech. 628, 311 (2009a).
  * Chung and McKeon [2010] D. Chung and B. J. McKeon, Large-eddy simulation of large-scale structures in long channel flow, J. Fluid Mech. 661, 341 (2010).
  * Guala  _et al._ [2011] M. Guala, M. Metzger, and B. J. McKeon, Interactions within the turbulent boundary layer at high Reynolds number, J. Fluid Mech. 666, 573 (2011).
  * Ganapathisubramani  _et al._ [2012] B. Ganapathisubramani, N. Hutchins, J. P. Monty, D. Chung, and I. Marusic, Amplitude and frequency modulation in wall turbulence, J. Fluid Mech. 712, 61 (2012).
  * Baars  _et al._ [2015] W. J. Baars, K. M. Talluru, N. Hutchins, and I. Marusic, Wavelet analysis of wall turbulence to study large-scale modulation of small scales, Exp. Fluids 56, 1 (2015).
  * Marusic  _et al._ [2010] I. Marusic, R. Mathis, and N. Hutchins, Predictive model for wall-bounded turbulent flow, Science 329, 193 (2010).
  * Mathis  _et al._ [2011] R. Mathis, N. Hutchins, and I. Marusic, A predictive inner–outer model for streamwise turbulence statistics in wall-bounded flows, J. Fluid Mech. 681, 537 (2011).
  * Agostini and Leschziner [2016] L. Agostini and M. Leschziner, Predicting the response of small-scale near-wall turbulence to large-scale outer motions, Phys. Fluids 28, 015107 (2016).
  * Baars  _et al._ [2016] W. J. Baars, N. Hutchins, and I. Marusic, Spectral stochastic estimation of high-Reynolds-number wall-bounded turbulence for a refined inner-outer interaction model, Phys. Rev. Fluids 1, 054406 (2016).
  * Inoue  _et al._ [2012] M. Inoue, R. Mathis, I. Marusic, and D. I. Pullin, Inner-layer intensities for the flat-plate turbulent boundary layer combining a predictive wall-model with large-eddy simulations, Phys. Fluids 24, 075102 (2012).
  * Cabrit  _et al._ [2014] O. Cabrit, R. Mathis, V. Kulandaivelu, and I. Marusic, Towards a statistically accurate wall-model for large-eddy simulation, in  _Proceedings of the $19^{th}$ Australasian Fluid Mechanics Conference_ (2014).
  * Sidebottom  _et al._ [2014] W. Sidebottom, O. Cabrit, I. Marusic, C. Meneveau, A. Ooi, and D. Jones, Modelling of wall shear-stress fluctuations for large-eddy simulation, in  _Proceedings of the $19^{th}$ Australasian Fluid Mechanics Conference_ (2014).
  * Mathis  _et al._ [2009b] R. Mathis, J. P. Monty, N. Hutchins, and I. Marusic, Comparison of large-scale amplitude modulation in turbulent boundary layers, pipes, and channel flows, Phys. Fluids 21, 111703 (2009b).
  * Monty  _et al._ [2010] J. P. Monty, M. S. Chong, R. Mathis, N. Hutchins, I. Marusic, and J. J. Allen, A high Reynolds number turbulent boundary layer with regular ‘Braille-type’ roughness, in  _IUTAM Symposium on The Physics of Wall-Bounded Turbulent Flows on Rough Walls_ (Springer, 2010) pp. 69–75.
  * Squire  _et al._ [2016a] D. T. Squire, W. J. Baars, N. Hutchins, and I. Marusic, Inner–outer interactions in rough-wall turbulence, J. Turb. 17, 1159 (2016a).
  * Anderson [2016] W. Anderson, Amplitude modulation of streamwise velocity fluctuations in the roughness sublayer: evidence from large-eddy simulations, J. Fluid Mech. 789, 567 (2016).
  * Wu  _et al._ [2019] S. Wu, K. T. Christensen, and C. Pantano, Modelling smooth-and transitionally rough-wall turbulent channel flow by leveraging inner–outer interactions and principal component analysis, J. Fluid Mech. 863, 407 (2019).
  * Blackman  _et al._ [2019] K. Blackman, L. Perret, and R. Mathis, Assessment of inner–outer interactions in the urban boundary layer using a predictive model, J. Fluid Mech. 875, 44 (2019).
  * Kim  _et al._ [2020] T. Kim, G. Blois, J. L. Best, and K. T. Christensen, Experimental evidence of amplitude modulation in permeable-wall turbulence, J. Fluid Mech. 887 (2020).
  * Khorasani  _et al._ [2022] S. M. H. Khorasani, M. Luhar, and S. Bagheri, Turbulent flows over engineered anisotropic porous substrates, arXiv preprint arXiv:2210.10140 (2022).
  * Dogan  _et al._ [2016] E. Dogan, R. E. Hanson, and B. Ganapathisubramani, Interactions of large-scale free-stream turbulence with turbulent boundary layers, J. Fluid Mech. 802, 79 (2016).
  * Dogan  _et al._ [2017] E. Dogan, R. J. Hearst, and B. Ganapathisubramani, Modelling high Reynolds number wall–turbulence interactions in laboratory experiments using large-scale free-stream turbulence, Phil. Trans. Roy. Soc. Lond. A: Math., Phys. Eng. Sci. 375, 20160091 (2017).
  * Duvvuri and McKeon [2015] S. Duvvuri and B. J. McKeon, Triadic scale interactions in a turbulent boundary layer, J. Fluid Mech. 767, R4 (2015).
  * Lozier  _et al._ [2022] M. Lozier, F. O. Thomas, and S. Gordeyev, Experimental studies of boundary layer dynamics via active manipulation of large-scale structures, in  _Proceedings of the $12^{th}$ Int. Symp. Turb. Shear Flow Phenom._ (Osaka, Japan, 2022).
  * Li  _et al._ [2019] M. Li, C. M. de Silva, R. Baidya, A. Rouhi, D. Chung, I. Marusic, and N. Hutchins, Recovery of the wall-shear stress to equilibrium flow conditions after a rough-to-smooth step-change in turbulent boundary layers, J. Fluid Mech. 872, 472 (2019).
  * Garratt [1990] J. R. Garratt, The internal boundary layer – A review, Boundary-Layer Meteorol. 50, 171 (1990).
  * Elliott [1958] W. P. Elliott, The growth of the atmospheric internal boundary layer, Trans. Am. Geophys. Union 39, 1048 (1958).
  * Antonia and Luxton [1972] R. A. Antonia and R. E. Luxton, The response of a turbulent boundary layer to a step change in surface roughness. Part 2. Rough-to-smooth, J. Fluid Mech. 53, 737 (1972).
  * Hanson and Ganapathisubramani [2016] R. E. Hanson and B. Ganapathisubramani, Development of turbulent boundary layers past a step change in wall roughness, J. Fluid Mech. 795, 494 (2016).
  * Rouhi  _et al._ [2019] A. Rouhi, D. Chung, and N. Hutchins, Direct numerical simulation of open channel flow over smooth-to-rough and rough-to-smooth step changes, J. Fluid Mech. 866, 450 (2019).
  * Li  _et al._ [2021] M. Li, C. M. de Silva, D. Chung, D. I. Pullin, I. Marusic, and N. Hutchins, Experimental study of a turbulent boundary layer with a rough-to-smooth change in surface conditions at high Reynolds numbers, J. Fluid Mech. 923, A18 (2021).
  * Li  _et al._ [2022] M. Li, C. M. de Silva, D. Chung, D. I. Pullin, I. Marusic, and N. Hutchins, Modelling the downstream development of a turbulent boundary layer following a step change of roughness, J. Fluid Mech. 949, A7 (2022).
  * Townsend [1976] A. A. Townsend,  _The structure of turbulent shear flow_ , 2nd ed. (Cambridge University Press., 1976).
  * Squire  _et al._ [2016b] D. T. Squire, C. Morrill-Winter, N. Hutchins, M. P. Schultz, J. C. Klewicki, and I. Marusic, Comparison of turbulent boundary layers over smooth and rough surfaces up to high Reynolds numbers, J. Fluid Mech. 795, 210 (2016b).
  * Marusic  _et al._ [2015] I. Marusic, K. A. Chauhan, V. Kulandaivelu, and N. Hutchins, Evolution of zero-pressure-gradient boundary layers from different tripping conditions, J. Fluid Mech. 783, 379 (2015).
  * Bendat and Piersol [2011] J. S. Bendat and A. G. Piersol,  _Random data: analysis and measurement procedures_ (John Wiley & Sons, 2011).
  * Baars  _et al._ [2017] W. J. Baars, N. Hutchins, and I. Marusic, Reynolds number trend of hierarchies and scale interactions in turbulent boundary layers, Phil. Trans. Roy. Soc. Lond. A: Math., Phys. Eng. Sci. 375, 20160077 (2017).
  * Mathis  _et al._ [2013] R. Mathis, I. Marusic, S. I. Chernyshenko, and N. Hutchins, Estimating wall-shear-stress fluctuations given an outer region input, J. Fluid Mech. 715, 163 (2013).
  * Pathikonda and Christensen [2017] G. Pathikonda and K. T. Christensen, Inner–outer interactions in a turbulent boundary layer overlying complex roughness, Phys. Rev. Fluids 2, 044603 (2017).
  * Efstathiou and Luhar [2018] C. Efstathiou and M. Luhar, Mean turbulence statistics in boundary layers over high-porosity foams, J. Fluid Mech. 841, 351 (2018).
  * Kulandaivelu [2012] V. Kulandaivelu,  _Evolution of zero pressure gradient turbulent boundary layers from different initial conditions_ , Ph.D. thesis, The University of Melbourne (2012).
  * Ligrani and Bradshaw [1987] P. M. Ligrani and P. Bradshaw, Spatial resolution and measurement of turbulence in the viscous sublayer using subminiature hot-wire probes, Exp. Fluids 5, 407 (1987).
  * Jacobi and McKeon [2013] I. Jacobi and B. J. McKeon, Phase relationships between large and small scales in the turbulent boundary layer, Exp. Fluids 54, 1 (2013).
  * Jacobi and McKeon [2017] I. Jacobi and B. J. McKeon, Phase-relationships between scales in the perturbed turbulent boundary layer, J. Turb. 18, 1120 (2017).
  * Deshpande  _et al._ [2022] R. Deshpande, D. Chandran, A. J. Smits, and I. Marusic, The relationship between manipulated inter-scale phase and energy-efficient turbulent drag reduction, arXiv preprint arXiv:2211.07176 (2022).


# Anisotropic nematic fluctuations above the ferroquadrupolar transition in
TmVO4

Z. Wang  I. Vinograd  Z. Mei  Department of Physics, University of California,
Davis, California 95616, USA  P. Menegasso  Instituto de Física “Gleb
Wataghin”, UNICAMP, Campinas-SP, 13083-859, Brazil  D. Garcia  Department of
Physics, University of California, Davis, California 95616, USA  P. Massat  I.
R. Fisher  Geballe Laboratory for Advanced Materials and Department of Applied
Physics, Stanford University, CA 94305, USA  N. J. Curro  Department of
Physics, University of California, Davis, California 95616, USA

###### Abstract

TmVO4 exhibits ferroquadrupolar order below 2.15 K with a well-isolated non-
Kramers ground state doublet, and is a model system to understand Ising
nematic order. We present 51V nuclear magnetic resonance data as a function of
field orientation in a single crystal. Although the spectra are well
understood in terms of direct dipolar hyperfine couplings, the spin lattice
relaxation rate exhibits strong anisotropy that cannot be understood in terms
of magnetic fluctuations. We find that the spin lattice relaxation rate scales
with the shear elastic constant associated with the ferroquadrupole phase
transition, suggesting that quadrupole (nematic) fluctuations dominate the
spin lattice relaxation for in-plane fields.

##  I Introduction

Recently there has been growing interest in the behavior of electronic Ising
nematicity, which may play a role in the low temperature behavior of a number
of strongly correlated electron systems. The phase diagram of the iron-based
superconductors is dominated by $C_{4}$ symmetry breaking of the spin and
orbital degrees of freedom, accompanied by a tetragonal-to-orthorhombic
structural transition [1]. The physical origin of the nematicity is yet to be
established, and might even be different in different families of compounds,
with possible contributions both from orbital effects and spin fluctuations.
The correlation between large values of the nematic susceptibility, a putative
nematic quantum critical point, and optimal superconductivity in several
materials, point to a possible role for nematic fluctuations in the pairing
interaction, as well as non-Fermi-liquid behavior in the normal state [2, 3,
4, 5]. Electronic nematic correlations are also present in the high
temperature superconducting cuprates [6, 7]. Disentangling the effects of
nematicity can be complicated by the presence of other intertwined order
parameters in such systems. It is challenging to discern whether nematic or
antiferromagnetic fluctuations dominate, especially in a conducting system
[8]. Moreover, strain couples to the nematic order, hence quenched random
disorder around dopants may give rise to inhomogeneous glassy behavior [9,
10]. It is therefore important to investigate critical fluctuations in an
Ising nematic system in the absence of metallicity or inhomogeneous strain
fields.

Figure 1:  (a) 51V spectra in TmVO4 at 11.7294T for $\mathbf{H}_{0}\perp c$ at
several different temperatures. The full-width half maximum (b) and
quadrupolar splitting (c) versus temperature for several different angles,
$\theta$. The inset shows the unit cell structure. Tm is blue, V is green
(within the pyramids) and O is red.

Ferroquadrupolar ordering of non-Kramers doublets in 4f materials offers an
important avenue to investigate Ising nematicity [11]. Specifically,
ferroquadrupole order breaks all the same symmetries as Ising nematic order
(and hence is a specific realization of nematic order), with the advantage
that the underlying effective Hamiltonian describing the low-temperature
behavior is well-understood. TmVO4 is an insulator with Tm ions in the 4f12
configuration ($L=5$, $S=1$, $J=6$) that crystallizes in space group $I41/amd$
(see inset of Fig. 1). The tetragonal crystal field splits the $J=6$ multiplet
giving rise to a $\Gamma_{5}$ ground state doublet separated by a gap of $\sim
77$ K to the lowest excited state [12, 13]. The wavefunctions of the ground
state doublet are $|\psi_{1,2}\rangle=e|\pm 5\rangle+f|\pm 1\rangle+g|\mp
3\rangle$ in the $J_{z}$ basis, where $e\approx 0.92$, $f\approx-0.37$, and
$g\approx 0.12$. The degeneracy of these ground states cannot be lifted by a
magnetic field perpendicular to the $z$-direction because
$\langle\psi_{1,2}|J_{\pm}|\psi_{1,2}\rangle=0$, hence these form a non-
Kramers doublet. The doublet can, however, be linearly split by either a
magnetic field oriented along the $c$-axis, or by lattice strains with either
a $B_{1g}$ ($x^{2}-y^{2}$) or $B_{2g}$ ($xy$) symmetry. Quadrupole-quadrupole
interactions mediated by the lattice dominate the magnetic interactions, and
the material undergoes a cooperative Jahn-Teller ferroquadrupolar ordering
with a $B_{2g}$ symmetry at $T_{Q}=2.15$K, accompanied by an orthorhombic
lattice distortion of the same symmetry [14]. The low-temperature behavior of
the Tm quadrupoles can be well-described by the transverse field Ising model,
in which pseudospins ($\tilde{S}={1}/{2}$) experience in-plane
ferroquadrupolar Ising couplings and couple to a transverse magnetic field
along the $z$-axis [11]. These fields will enhance the fluctuations of the
pseudospins and can tune the system to an Ising-nematic quantum phase
transition. This material thus offers an important platform to investigate
quantum critical nematic fluctuations in an insulator.

In order to better understand the nature of these fluctuations we have
investigated 51V ($I=7/2$, $Q=52$ mb, 99.75% abundant) NMR in a single crystal
of TmVO4 in a magnetic field $\mathbf{H}_{0}=11.72$ T oriented perpendicular
to the $c$-axis. We find that the magnetic shift and spin-lattice-relaxation
rate are strongly angular dependent (rotating the field in the $a-c$ plane)
below $\sim 80$ K, reflecting the anisotropy of the ground state doublet. The
spin-lattice-relaxation rate is non-monotonic, exhibiting a large enhancement
at low temperature that may be associated with critical nematic fluctuations.
The spectrum exhibits a strong temperature-dependent magnetic broadening at
low temperatures caused by inhomogeneous magnetic demagnetization fields.

Figure 2:  magnetic shift, $K$, versus temperature (a) and versus angle (b).
The dotted lines are fits as described in the text.

##  II NMR Spectra and magnetic shift

TmVO4 crystals were grown from a Pb2V2O7 flux using 4 mole percent of Tm2O3,
following the methods described in [15, 16]. The crystals have a rod-like
morphology with the the $c$-axis along the long axis. A single crystal of
approximate dimension 1mm $\times$ 1mm $\times$ 4mm was selected and mounted
on a cryogenic goniometer NMR probe. The magnetic susceptibility is strongly
anisotropic, reflecting that of the unusual $g$-factor ($g_{c}=10.2$,
$g_{\perp}~{}0$) of the ground state doublet. Although a crystal mounted with
$c\perp\mathbf{H}_{0}$ experiences zero torque, it is an unstable equilibrium
and there is a large torque for infinitesimal deviations from 90∘. To
alleviate this issue we secured the crystal with epoxy to a mounting plate
that itself is rotated. Spin echoes were acquired at several different
frequencies, and the Fourier transforms were summed to measure the full
spectra including all nuclear spin transitions. Fig. 1(a) shows several
representative spectra of the 51V as a function of temperature. There are
seven peaks separated by the quadrupolar interaction. Because the V has axial
symmetry, the peaks frequencies are given by:

$\nu=\gamma H_{0}\left(1+K(\theta)\right)+n\nu_{q}(\theta),$ (1)

where the magnetic, $K(\theta)$, and quadrupolar, $\nu_{q}(\theta)$, shifts
vary with the angle $\theta$ between the c-axis and $\mathbf{H}_{0}$:

$\displaystyle K(\theta)$ $\displaystyle=$ $\displaystyle
K_{cc}\cos^{2}\theta+K_{aa}\sin^{2}\theta$ (2) $\displaystyle\nu_{q}(\theta)$
$\displaystyle=$ $\displaystyle\nu_{zz}\left(3\cos^{2}\theta-1\right)/2.$ (3)

Here $\gamma=11.193$ MHz/T, $n=-3,\cdots,3$, $\nu_{zz}={eQV_{zz}}/{12h}$ and
$V_{zz}$ is the electric field gradient at the V site. We measured the spectra
for several angles $86^{\circ}\leq\theta\leq 90^{\circ}$ and fit the spectra
to a sum of Lorentzians. The temperature and angular dependence of the
linewidths, EFG, and magnetic shifts are shown in Figs. 1(b, c) and 2(a).

We find that the EFG is similar to previous measurements, [12, 17] however the
spectra are broad. For this orientation, the quadrupolar splitting is
$\nu_{zz}/2$, which is comparable to the FWHM of each resonance. As a result,
the individual peaks become difficult to resolve at low temperatures. Each of
the satellite resonances has the same linewidth, implying that that the
broadening is due to magnetic field inhomogeneity within the sample. Moreover,
we find that the spectra are narrower at lower applied fields and that the
FWHM varies approximately linear with field. Since the susceptibility is
strongly anisotropic, distortions of the internal field $\mathbf{B}$ due to
demagnetization effects within the needle-like prism of the crystal can create
a large distribution of local resonance frequencies [18]. The spectra also
display a suppression of intensity for the inner satellites, particularly at
low temperature. This phenomenon arises due to fast spin-spin decoherence
rates ($T_{2}^{-1}$) with high spin nuclei [19].

Figure 3:  $K_{aa}$ and $K_{cc}$ versus temperature (a) and versus bulk
susceptibility (b). The solid lines are fits as described in the text.

The magnetic shift shown in Fig. 2 is negative and strongly angular dependent
at low temperature. We fit the angular dependence to extract the tensor
components $K_{aa}$ and $K_{cc}$, shown as dotted lines in Fig. 2(b). This
approach enables us to extract the magnetic shift for the $c$ direction
without needing to fully align the crystal in this orientation, albeit the
error bars for $K_{cc}$ are larger than for $K_{aa}$. The temperature
dependence of $K_{aa}$ and $K_{cc}$ are shown in Fig. 3(a). $K_{cc}$ is large
and positive. Fig. 3(b) shows these shift components plotted versus the bulk
susceptibility $\chi_{aa,cc}$, which was measured independently in a SQUID
magnetometer. The shift varies linearly with susceptibility as
$K_{\alpha\alpha}=K_{\alpha\alpha}^{orb}+A_{\alpha\alpha}\chi_{\alpha\alpha}$,
where $A_{\alpha\beta}$ are the components of the hyperfine coupling tensor.
We find that $K_{aa}^{orb}=-0.315\pm 0.009\%$, $K_{cc}^{orb}=-0.4\pm 0.1\%$,
$A_{aa}=-0.32\pm 0.07$ kOe/$\mu_{B}$ and $A_{cc}=1.29\pm 0.05$ kOe/$\mu_{B}$.

These values of the hyperfine couplings are consistent with a direct dipolar
coupling mechanism between the Tm moments and the V nuclear spins. The direct
dipolar coupling is given by:
$A_{\alpha\beta}^{dip}=\sum_{i}(\nabla\times\mathbf{A}_{i})_{\alpha}/\mu_{\beta}$,
where $\mathbf{A}_{i}=\mathbf{\mu}\times\mathbf{r}_{i}/r_{i}^{3}$ is the
vector potential of a dipole moment, $\mathbf{\mu}$, located at lattice site
$\mathbf{r}_{i}$ relative to a central nucleus. For the TmVO4 lattice, we
estimate $A_{aa}^{dip}=A_{bb}^{dip}=-0.336$ kOe/$\mu_{B}$ and
$A_{cc}^{dip}=0.671$ kOe/$\mu_{B}$ at the V site. The theoretical value for
the perpendicular direction is the same as the measured value within the error
limits. For the $c$ axis, the theoretical value is within a factor of two of
the measured values, and it is likely there are larger systematic measurement
errors involved in extracting this value. Thus the anisotropic magnetic shift
tensor can be fully explained via direct dipolar interactions, as expected for
an insulator.

##  III Spin Lattice Relaxation Rate

The spin-lattice-relaxation rate, $T_{1}^{-1}$, was measured by applying
inversion pulses at the central transition $(n=0)$ and measuring the echo
intensity as a function of recovery time. The magnetization recovery was fit
to the standard expression for magnetic fluctuations:
$M(t)=M_{0}\left(1-2f\phi(t/T_{1})\right)$, where $M_{0}$ is the equilibrium
magnetization, $f$ is the inversion fraction, and

$\phi(t)=\frac{1225}{1716}e^{-28t}+\frac{75}{364}e^{-15t}+\frac{3}{44}e^{-6t}+\frac{1}{84}e^{-t}.$
(4)

This expression fits the data well without the need for a stretching exponent.
Fig. 4 shows the temperature and angular dependence of $T_{1}^{-1}$. For
$\theta=90^{\circ}$, $T_{1}^{-1}$ decreases strongly below 80K as the excited
crystal field levels are thermally depopulated. In this temperature range
$T_{1}^{-1}$ become strongly angular dependent, increasing by more than a
factor of 30 as the field $\mathbf{H}_{0}$ rotates by only 4∘ away from the
perpendicular configuration. This behavior likely reflects the anisotropy of
the $g$ factor of the ground state doublet, however the anisotropy of
$T_{1}^{-1}$ is puzzling. If the relaxation is driven by magnetic fluctuations
of the Tm ground state, then $T_{1}^{-1}$ should exhibit a _maximum_ at
$\theta=90^{\circ}$ rather than a minimum because fluctuations of the non-
Kramers doublet should lie exclusively along the $c$-axis. Therefore
$T_{1}^{-1}(0^{\circ})$ should be much smaller than $T_{1}^{-1}(90^{\circ})$,
in contrast to our observations.

Figure 4:  (a) $T_{1}^{-1}$ vs temperature for multiple angles. (b) Calculated
$T_{1}^{-1}$ versus $\theta$ for magnetic fluctuations (solid lines) and for
quadrupole fluctuations (dashed lines). (c) Calculated $T_{1}^{-1}$ versus
temperature and angle using Eq. 5, where the angle is increased in $5^{\circ}$
increments between 0 and 90∘.

###  III.1 Magnetic fluctuations

On the other hand, the hyperfine couplings can give rise to a more complicated
relationship between the direction of the Tm moments and the direction of the
hyperfine fields. To properly account for these couplings we use the Moriya
expression:

$T_{1m}^{-1}=\gamma^{2}k_{B}T\lim_{\omega\rightarrow
0}\sum\limits_{\mathbf{q},\alpha,\beta}\mathcal{F}_{\alpha\beta}(\mathbf{q})\frac{\textrm{Im}\chi_{\alpha\beta}(\mathbf{q},\omega)}{\hslash\omega},$
(5)

where the form factors $\mathcal{F}_{\alpha\beta}(\mathbf{q})$ (see Appendix A
for details) depend on the local dipolar hyperfine couplings, and
$\chi_{\alpha\beta}(\mathbf{q},\omega)$ is the dynamical magnetic
susceptibility of the Tm moments. For simplicity we only include the two
nearest neighbor and four next-nearest neighbor Tm atoms in the form factors.
Because the Tm system exhibits ferroquadrupolar order at $T_{Q}$, we assume
that the structure of the dynamical susceptibility can be modeled as:

$\chi_{\alpha\alpha}(\mathbf{q},\omega)=\frac{\chi_{\alpha\alpha}(T)}{\xi^{-2}+f(\mathbf{q})-i\omega/\Gamma
q}$ (6)

where $\xi$ is a correlation length, $\Gamma$ is a characteristic fluctuation
energy, $f(\mathbf{q})=q_{x}^{2}+q_{y}^{2}+\eta q_{z}^{2}$, $\eta$ is a
dimensionless parameter that reflects the tetragonal nature, and
$\chi_{\alpha\alpha}(T)$ is the static ($\mathbf{q}=0$) susceptibility. $\xi$
and $\eta$ are unknown parameters, but we compute the temperature and angular
dependence using $\xi=2$ and $\eta=1/2$. Fig. 4(b) shows the expected angular
dependence of $T_{1m}^{-1}$ for $\chi_{cc}/\chi_{aa}=10$ (red), close to the
experimental value, and for $\chi_{cc}/\chi_{aa}=0.1$ (blue). The former
clearly exhibits a maximum of $T_{1m}^{-1}$ at $\theta=90^{\circ}$, in
contrast to our observations. The latter exhibits a shallow minimum at
$90^{\circ}$, but the susceptibility anisotropy does not agree with
experiment. Fig. 4(c) shows the temperature dependence using the measured
values of the static susceptibility. Although there is an overall decrease in
$T_{1m}^{-1}$ at lower temperatures, the detailed temperature dependence does
not match experiment, and the calculated $T_{1}^{-1}$ still exhibits a maximum
for $\theta=90^{\circ}$ at all temperatures. Despite the complex form factors
for the direct dipolar couplings, the expected magnetic fluctuations of the Tm
ground state cannot explain the observed increase in $T_{1}^{-1}$ as the field
rotates out of the plane.

###  III.2 Quadrupolar fluctuations

An alternative explanation is that the spin-lattice-relaxation rate is
dominated quadrupolar fluctuations rather than magnetic. The Tm quadrupole
moments couple to the EFG at the V site, giving rise to a second nuclear
quadrupolar relaxation channel [8]. The enhancement of $T_{1}^{-1}$ below 20K
for $\theta=90^{\circ}$ may represent the growth of critical fluctuations near
$T_{Q}$. Note that changing $\theta$ by only 0.25∘ dramatically alters
$T_{1}^{-1}$, which is close to the limit of precision of our goniometer. Thus
it is possible that the enhancement below 20K may vanish or become smaller for
better alignment. In the presence of both magnetic and quadrupolar relaxation,
the expression for $\phi(t)$ (Eq. 4) changes, and includes three independent
rates: $T_{1m}^{-1}$, $W_{Q1}$ and $W_{Q2}$, where the latter two are
associated with $\Delta m=\pm 1$ and $\Delta m=\pm 2$ quadrupolar relaxation.
We are unable, however, to independently extract these parameters with
sufficient resolution. Moreover, the line broadening observed in Fig. 1 also
means that the magnetization relaxation at the central transition may also
include contributions from nearby satellite transitions, further complicating
any attempts to extract the independent relaxation channels. Nevertheless, it
is instructive to consider the case where quadrupole fluctuations dominate and
magnetic fluctuations can be neglected.

Quadrupolar relaxation is driven by fluctuations of the spherical tensor
components of the EFG: $V_{\pm 1}=V_{zx}\pm iV_{zy}$ and $V_{\pm
2}=\frac{1}{2}(V_{xx}-V_{yy})\pm iV_{xy}$, where the $V_{\alpha\beta}$ are the
EFG tensor components relative to the direction of $\mathbf{H}_{0}$. These
give rise to nuclear spin relaxation rates:

$W_{Q1,Q2}=(eQ/\hbar)^{2}\int_{0}^{\infty}\langle
V_{+1,2}(\tau)V_{-1,2}(0)\rangle e^{-\omega_{L}\tau}d\tau$ (7)

where $\omega_{L}$ is the Larmor frequency [20]. The nematic order in this
system has $B_{2g}$ symmetry, so $V_{xx}-V_{yy}\neq 0$, where $z$ corresponds
to the $c$-direction and $x$ and $y$ are along the principal axes of the EFG
tensor, which are rotated 45∘ relative to the tetragonal $a$-axes. Above
$T_{Q}$ fluctuations of $V_{\pm 2}$ should dominate those of $V_{\pm 1}$, and
as a result we anticipate that $W_{Q1}(\theta=0)$ can be neglected. As the
field is rotated towards the plane, the EFG tensor components change, and the
relaxation rates become angular dependent (see Appendix B for details):

$\displaystyle W_{Q2}(\theta)/W_{Q2}(0)$ $\displaystyle=$
$\displaystyle\left(\cos^{4}\theta+6\cos^{2}\theta+1\right)/8$ (8)
$\displaystyle W_{Q1}(\theta)/W_{Q2}(0)$ $\displaystyle=$
$\displaystyle\sin^{2}\theta(\cos(2\theta)+3)/4.$ (9)

These quantities are shown in Fig. 4(b) as dashed lines. $W_{Q2}$ exhibits a
minimum for $\theta=90^{\circ}$, whereas $W_{Q1}$ is nearly independent of
$\theta$ at this angle. This behavior agrees qualitatively with our
observations, but the increases we observe are in fact a much stronger
function of angle than expected for quadrupolar relaxation. Rotating $\theta$
by 1-2∘ out of the plane enhances $T_{1}^{-1}$ by an order of magnitude,
whereas $W_{Q2}$ exhibits only a quadratic minimum at this angle.

Figure 5:  The shear elastic stiffness coefficient $c_{66}$ (solid line,
reproduced from [21]) and the quantity $1/(1+(aT_{1}T)^{-1})$ as as a function
of temperature. INSET: $1/(c_{66,0}/c_{66}-1)$ versus $T_{1}T$, with
temperature implicit. The solid black line is the best linear fit, giving
$a=18.5\pm 0.4$ sec-1 K-1.

The interpretation that relaxation is driven by quadrupole fluctuations is
supported by comparisons of the temperature dependence of $T_{1}^{-1}$ with
that of the shear elastic stiffness coefficient, $c_{66}$, which softens with
decreasing temperature and vanishes at $T_{Q}$ [21]. This behavior is driven
by the nematic susceptibility:
$\chi_{nem}=c_{66,0}(1-c_{66,0}/c_{66})/\lambda^{2}$, where $\lambda$ is the
coupling between the lattice and the Tm 4f orbitals, and $c_{66,0}$ is the
stiffness coefficient in the absence of the coupling [22]. If $T_{1}^{-1}$ is
also determined by the Tm orbital fluctuations, then
$(T_{1}T)^{-1}\sim\chi_{nem}$ (see Appendix C for details) [8, 23, 24]. We
thus expect $T_{1}T\sim(c_{66,0}/c_{66}-1)^{-1}$, which is demonstrated in
Fig. 5. The main panel compares the temperature dependence of $c_{66}$ with
the measured $T_{1}^{-1}$ values, and the inset shows the scaling between the
shear modulus and $T_{1}^{-1}$ with temperature as an implicit parameter. The
scaling evident in Fig. 5 suggests that the spin lattice relaxation is driven
primarily by quadrupole fluctuations, which are reflected in the softening of
$c_{66}$.

##  IV Discussion

A slightly different scaling relationship was found in the iron pnictide
superconductors via a microscopic model that assumes that the nematicity
arises in the magnetic susceptibility, which in turn affects the nuclei
through a magnetic hyperfine interaction [22]. In TmVO4, the nematicity arises
from the Tm electronic orbitals, and the coupling to the nuclei may be through
the quadrupolar interaction. Moreover, the pnictide model assumed the presence
of Landau damping by a Fermi surface of quasiparticles, which is not the case
for insulating TmVO4. The relaxation in TmVO4 must also involve a damping
term, but the origin of this term is unknown. The fact that the scaling
relationship in Fig. 4 holds suggests that this damping term is temperature
independent.

Rotating the field away from $90^{\circ}$ can enhance quadrupole fluctuations.
A rotation of $\mathbf{H}_{0}$ by $4^{\circ}$ corresponds to a field of $0.82$
T along the $c$-axis. This is greater than the critical field of $H_{c}=0.52$
T to suppress the long-range nematic order, which naturally enhances
fluctuations of both $V_{\pm 2}$ and $V_{\pm 1}$. However, these critical
fluctuations are not likely to persist to higher temperatures beyond $\sim
10T_{Q}$, thus are unlikely to be responsible for the large anisotropy
observed up to 80K. An alternative scenario is that the higher CEF levels
cannot be ignored. Indeed, even though the in-plane g-factor of the non-
Kramers doublet vanishes in zero applied field, the excited CEF levels can be
mixed into the ground state wavefunctions by an in-plane field. As a result,
there can be an induced magnetic moment in the plane, which may also
contribute to the relaxation [21].

It is likely that the spin lattice relaxation is dominated by both magnetic
and quadrupolar fluctuations, however it is difficult to disentangle these two
relaxation channels without more detailed measurements of the relaxation at
the higher satellite transitions [20]. However, as illustrated in Fig. 1, the
satellites are magnetically broadened and cannot be well resolved, especially
at lower temperatures. This broadening is due to the demagnetization field
inhomogeneity of our crystal. In principle, it is possible to improve the
spectral resolution by removing the sharp edges and corners of the sample
and/or operating at lower applied fields, in order to better discern the
individual satellite transitions.

Nuclear spin lattice relaxation rates have also been studied in both PrAlO3
and CsCuCl3, materials that exhibit structural distortions due to the
cooperative Jahn-Teller effect with non-magnetic ground states [25, 26, 27].
In contrast to our observations in TmVO4, $T_{1}^{-1}$ did not exhibit any
enhancement above the phase transition in these cases, even though the EFG
changed below. On the other hand, unlike TmVO4, the phase transitions in these
cases are first order, thus $T_{1}^{-1}$ should not reflect any critical
slowing down. The spin-lattice-relaxation rate in the disordered state was
analyzed in terms of magnetic (hyperfine) fluctuations, although the nuclei in
question (27Al, $I=5/2$ and 133Cs, $I=7/2$) are quadrupolar and should be
sensitive to fluctuations of the EFG.

In summary, we have measured the spectra and relaxation rates in TmVO4 as a
function of field direction oriented perpendicular to the c-axis. We find that
the magnetic shift tensor agrees quantitatively with direct dipolar coupling
between the V nuclear moments and the Tm 4f moments. The spin lattice
relaxation rate exhibits a steep minimum for field oriented 90∘ to the $c$
axis, which is inconsistent with purely magnetic fluctuations. We find that
$T_{1}$ scales with the lattice constant for shear strain, $c_{66}$, which
softens and vanishes at the nematic transition. It is likely that both
quadrupolar and magnetic fluctuations are present and drive spin lattice
relaxation. However, the origin of the the steep angular dependence of
$T_{1}^{-1}$ remains an open question.

##  V Acknowledgments

We acknowledge helpful discussions with R. Fernandes. Work at UC Davis was
supported by the NSF under Grants No. DMR-1807889 and PHY-1852581, as well as
the UC Davis Seed Grant program. Crystal growth performed at Stanford
University was supported by the Air Force Office of Scientific Research under
award number FA9550-20-1-0252. P. M. was partially supported by the Gordon and
Betty Moore Foundation Emergent Phenomena in Quantum Systems Initiative
through Grant GBMF9068.

##  Appendix A Magnetic Relaxation Form Factors

We assume that the dominant hyperfine fields at the V site arise from the two
nearest neighbor and four next-nearest Tm moments, whose positions are given
in Table 1. We define:

$\mathcal{A}_{\alpha\beta}(\mathbf{q})=\sum_{i}e^{i\mathbf{q}\cdot\mathbf{r}_{i}}A_{\alpha\beta}^{dip}$
(10)

where $A_{\alpha\beta}^{dip}$ is defined in the main text. For an applied
field $\mathbf{H}_{0}$ oriented at angles $\theta$ and $\phi$ relative to the
crystalline axes, the form factors are [28]:

$\mathcal{F}_{\alpha\beta}(\mathbf{q})=\sum\limits_{\epsilon,\delta}\left[R_{x\epsilon}R_{x\delta}+R_{y\epsilon}R_{y\delta}\right]\mathcal{A}_{\epsilon\alpha}(\mathbf{q})\mathcal{A}_{\delta\beta}(-\mathbf{q}),$
(11)

where the $R_{\alpha\beta}$ are elements of the 3D rotation matrix:

$\mathbb{R}=\left(\begin{array}[]{ccc}\cos\theta\cos^{2}\phi+\sin^{2}\phi&\cos\theta\cos\phi\sin\phi-\cos\phi\sin\phi&\cos\phi\sin\theta\\\
\cos\theta\cos\phi\sin\phi-\cos\phi\sin\phi&\cos^{2}\phi+\cos\theta\sin^{2}\phi&\sin\theta\sin\phi\\\
-\cos\phi\sin\theta&-\sin\theta\sin\phi&\cos\theta\end{array}\right).$ (12)

Table 1: Position vectors for six nearest neighbor Tm sites to V, in spherical coordinates. $\mathbf{r}_{i}$ |  $r$ (Å) |  $\theta$ (∘) |  $\phi$ (∘)
---|---|---|---
1 | 3.13030 | 0 | 0
2 | 3.13030 | 180 | 0
3 | 3.86654 | 66.1218 | 0
4 | 3.86654 | 113.8782 | 90
5 | 3.86654 | 66.1218 | 180
6 | 3.86654 | 113.8782 | 270

##  Appendix B Quadrupolar Relaxation Anisotropy

Equation 7 gives the expression for quadrupolar relaxation in terms of the
spherical tensor components of the EFG tensor. The quadrupolar interaction is
only on-site, so there are no form factors. However, the EFG tensor must be
rotated properly as the field direction changes. Under a rotation the tensor
operators $V_{m}(\tau)$ transform as:

$V^{\prime}_{m}(\tau)=\sum_{m^{\prime}}D_{mm^{\prime}}^{(2)}V_{m^{\prime}}(\tau)$
(13)

where

$D_{mm^{\prime}}^{(l)}(\alpha,\beta,\gamma)=e^{-im\alpha}d^{l}_{mm^{\prime}}(\beta)e^{-im^{\prime}\gamma},$
(14)

are the Wigner $D$ matrices, and the Euler angles are
$(\alpha=\phi,\beta=\theta,\gamma=0)$. The correlation functions $\langle
V_{m}(\tau)V_{-m}(0)\rangle$ are thus given by:

$\langle
V^{\prime}_{m}(\tau)V^{\prime}_{-m}(0)\rangle=\sum_{m^{\prime},m^{\prime\prime}}D_{2m^{\prime}}^{(2)}(\phi,\theta)D_{-2m^{\prime\prime}}^{(2)}(\phi,\theta)\langle
V_{m^{\prime}}(\tau)V_{m^{\prime\prime}}(0)\rangle.$ (15)

We assume that $\langle V_{m}(\tau)V_{m^{\prime}}(0)\rangle=0$ for all
$m,m^{\prime}$ except for $m=-m^{\prime}=2$ and $m=m^{\prime}=1$. Moreover, we
assume that $\langle V_{2}(\tau)V_{-2}(0)\rangle\gg\langle
V_{1}(\tau)V_{-1}(0)\rangle$, since $\langle V_{\pm 2}\rangle\neq 0$ and
$\langle V_{\pm 1}\rangle=0$ in the nematic phase. We thus expect
$W_{Q1}(\theta=0)\approx 0$, and:

$\displaystyle W_{Q2}(\theta)/W_{Q2}(0)$ $\displaystyle=$
$\displaystyle\left(\cos^{4}\theta+6\cos^{2}\theta+1\right)/8$ (16)
$\displaystyle W_{Q1}(\theta)/W_{Q2}(0)$ $\displaystyle=$
$\displaystyle\sin^{2}\theta(\cos(2\theta)+3)/4.$ (17)

as given in the main text in Eqs. 8 and Eqs. 9.

##  Appendix C Relaxation driven by nematic fluctuations

We note that Eq. 7 can be expressed in terms of the dynamical nematic
susceptibility [8]:

$W_{Q2}(0)=\left(\frac{eQ}{\hbar}\right)^{2}k_{B}T\lim_{\omega\rightarrow
0}\sum\limits_{\mathbf{q}}\frac{\mathrm{Im}\chi_{nem}(\mathbf{q},\omega)}{\hslash\omega}.$
(18)

The dynamical susceptibility can be expressed phenomenologically as:
$\chi_{nem}(\mathbf{q},\omega)={\chi_{nem}}({1-i\omega/\omega_{n}})^{-1}$,
where where $\chi_{nem}$ is the static nematic susceptibility and $\omega_{n}$
is a damping term [23, 24]. In this case
$W_{Q2}(0)=\left({eQ}\right)^{2}k_{B}T\chi_{nem}/\hbar^{2}\omega_{n}$.

## References

  * Fernandes  _et al._ [2014] R. M. Fernandes, A. V. Chubukov, and J. Schmalian, What drives nematic order in iron-based superconductors?, Nat. Phys. 10, 97 (2014).
  * Chu  _et al._ [2012] J.-H. Chu, H.-H. Kuo, J. G. Analytis, and I. R. Fisher, Divergent nematic susceptibility in an iron arsenide superconductor, Science 337, 710 (2012).
  * Kuo  _et al._ [2016] H.-H. Kuo, J.-H. Chu, J. C. Palmstrom, S. A. Kivelson, and I. R. Fisher, Ubiquitous signatures of nematic quantum criticality in optimally doped Fe-based superconductors, Science 352, 958 (2016).
  * Lederer  _et al._ [2015] S. Lederer, Y. Schattner, E. Berg, and S. A. Kivelson, Enhancement of superconductivity near a nematic quantum critical point, Phys. Rev. Lett. 114, 097001 (2015).
  * Maier and Scalapino [2014] T. A. Maier and D. J. Scalapino, Pairing interaction near a nematic quantum critical point of a three-band CuO2 model, Phys. Rev. B 90, 174510 (2014).
  * Kivelson  _et al._ [1998] S. A. Kivelson, E. Fradkin, and V. J. Emery, Electronic liquid-crystal phases of a doped Mott insulator, Nature 393, 550 (1998).
  * Vojta [2009] M. Vojta, Lattice symmetry breaking in cuprate superconductors: stripes, nematics, and superconductivity, Advances in Physics 58, 699 (2009).
  * Dioguardi  _et al._ [2016] A. P. Dioguardi, T. Kissikov, C. H. Lin, K. R. Shirer, M. M. Lawson, H.-J. Grafe, J.-H. Chu, I. R. Fisher, R. M. Fernandes, and N. J. Curro, NMR evidence for inhomogeneous nematic fluctuations in BaFe2(As1-xPx)2, Phys. Rev. Lett. 116, 107202 (2016).
  * Dioguardi  _et al._ [2015] A. P. Dioguardi, M. M. Lawson, B. T. Bush, J. Crocker, K. R. Shirer, D. M. Nisson, T. Kissikov, S. Ran, S. L. Bud’ko, P. C. Canfield, S. Yuan, P. L. Kuhns, A. P. Reyes, H.-J. Grafe, and N. J. Curro, NMR evidence for inhomogeneous glassy behavior driven by nematic fluctuations in iron arsenide superconductors, Phys. Rev. B 92, 165116 (2015).
  * Carlson  _et al._ [2006] E. W. Carlson, K. A. Dahmen, E. Fradkin, and S. A. Kivelson, Hysteresis and noise from electronic nematicity in high-temperature superconductors, Phys Rev. Lett. 96, 097003 (2006).
  * Maharaj  _et al._ [2017] A. V. Maharaj, E. W. Rosenberg, A. T. Hristov, E. Berg, R. M. Fernandes, I. R. Fisher, and S. A. Kivelson, Transverse fields to tune an Ising-nematic quantum phase transition, Proc. Natl. Acad. Sci. 114, 13430 (2017).
  * Bleaney and Wells [1980] B. Bleaney and M. R. Wells, Radiofrequency studies of tmvo4, Proceedings of the Royal Society of London. Series A, Mathematical and Physical Sciences 370, 131 (1980).
  * Bowden [1998] G. J. Bowden, A review of the low temperature properties of the rare earth vanadates, Australian Journal of Physics 51, 201 (1998).
  * Gehring and Gehring [1975] G. A. Gehring and K. A. Gehring, Co-operative jahn-teller effects, Reports on Progress in Physics 38, 1 (1975).
  * Feigenson [1968] R. Feigenson, Flux growth of type RVO4 rare-earth vanadate crystals, Journal of the American Ceramic Society 51, 538 (1968).
  * Smith and Wanklyn [1974] S. Smith and B. Wanklyn, Flux growth of rare earth vanadates and phosphates, Journal of Crystal Growth 21, 23 (1974).
  * Bleaney  _et al._ [1983] B. Bleaney, J. Gregg, A. de Oliveira, and M. Wells, Nuclear electric quadrupole interaction in LnVO4, Journal of Magnetism and Magnetic Materials 31-34, 745 (1983).
  * Lawson  _et al._ [2018] M. Lawson, B. T. Bush, T. Kissikov, Z. Brubaker, K. R. Shirer, J. R. Jeffries, S. Ran, I. Jeon, M. B. Maple, and N. J. Curro, Measurements of the NMR Knight shift tensor and nonlinear magnetization in URu2Si2, Phys. Rev. B 97, 075138 (2018).
  * Nisson  _et al._ [2013] D. M. Nisson, A. P. Dioguardi, P. Klavins, C. H. Lin, K. Shirer, A. C. Shockley, J. Crocker, and N. J. Curro, Nuclear magnetic resonance as a probe of electronic states of Bi2Se3, Phys. Rev. B 87, 195202 (2013).
  * Suter  _et al._ [1998] A. Suter, M. Mali, J. Roos, and D. Brinkmann, Mixed magnetic and quadrupolar relaxation in the presence of a dominant static Zeeman Hamiltonian, J. Phys.: Condens. Matter 10, 5977 (1998).
  * Melcher  _et al._ [1973] R. L. Melcher, E. Pytte, and B. A. Scott, Phonon instabilities in TmV${\mathrm{O}}_{4}$, Phys. Rev. Lett. 31, 307 (1973).
  * Fernandes  _et al._ [2013] R. M. Fernandes, A. E. Böhmer, C. Meingast, and J. Schmalian, Scaling between magnetic and lattice fluctuations in iron pnictide superconductors, Phys. Rev. Lett. 111, 137001 (2013).
  * Leung and Huber [1979] K. M. Leung and D. L. Huber, Low-frequency dynamics in cooperative Jahn-Teller systems, Physical Review B 19, 5483 (1979).
  * Page and Smith [1983] J. H. Page and S. R. P. Smith, Microwave ultrasonic attenuation above the Jahn-Teller phase transition in TmVO4, Journal of Physics C: Solid State Physics 16, 309 (1983).
  * Borsa  _et al._ [1978] F. Borsa, M. Corti, and A. Rigamonti, Electronic spin-dynamics at a structural phase transition by cooperative Jahn-Teller effect: An Al27 NMR study in PrAlO3, Journal of Applied Physics 49, 1383 (1978).
  * Corti  _et al._ [1981] M. Corti, A. Rigamonti, and A. Magistris, 133cs quadrupole perturbed NMR study of jahn-teller phase transitions in CsCuCl3), physica status solidi (b) 108, 29 (1981).
  * Rigamonti [1984] A. Rigamonti, NMR-NQR studies of structural phase transitions, Advances in Physics 33, 115 (1984).
  * Smerald and Shannon [2011] A. Smerald and N. Shannon, Angle-resolved NMR: Quantitative theory of 75As ${T}_{1}$ relaxation rate in BaFe2As2, Phys. Rev. B 84, 184437 (2011).