1. Introduction
Turbulent flows are a canonical example of chaotic phenomena. They feature complex spatiotemporal dynamics which, albeit deterministic, cannot be predicted over long times due to the growth of the uncertainty in the initial conditions (Lorenz Reference Lorenz1969, Reference Lorenz1996). This uncertainty is unavoidable and stems from the impossibility of determining the state of a turbulent flow to arbitrary precision. The mechanisms whereby the uncertainty of the initial conditions grows and propagates in time are key to understanding the emergence of complexity in turbulence (Boffetta et al. Reference Boffetta, Cencini, Falcioni and Vulpiani2002), and to quantifying the reliability and accuracy of weather forecasts (Palmer Reference Palmer2000; DelSole Reference DelSole2004), numerical simulations (Najm Reference Najm2009), turbulence models (Lucor, Meyers & Sagaut Reference Lucor, Meyers and Sagaut2007; Cheung et al. Reference Cheung, Oliver, Prudencio, Prudhomme and Moser2011; Duraisamy, Iaccarino & Xiao Reference Duraisamy, Iaccarino and Xiao2019) or flow reconstructions from partial observations (Zaki Reference Zaki2024).
The problem of characterising the evolution of uncertainty in turbulent flows has been traditionally formulated in terms of error growth, i.e. by measuring the distance between two independent flow fields that are initially close. This approach leads naturally to the definition of the Lyapunov exponents (LEs), which measure the growth rate of infinitesimal perturbations averaged over the turbulent attractor (Eckmann & Ruelle Reference Eckmann and Ruelle1985). The LEs have been calculated in many canonical turbulent flows, including forced isotropic turbulence (Grappin & Léorat Reference Grappin and Léorat1991; Mohan, Fitzsimmons & Moser Reference Mohan, Fitzsimmons and Moser2017; Boffetta & Musacchio Reference Boffetta and Musacchio2017; Berera & Ho Reference Berera and Ho2018; Hassanaly & Raman Reference Hassanaly and Raman2019), channel flow (Keefe, Moin & Kim Reference Keefe, Moin and Kim1992; Nikitin Reference Nikitin2018) or Taylor–Couette flow (Brandstäter et al. Reference Brandstäter, Swift, Swinney, Wolf, Farmer, Jen and Crutchfield1983; Brandstater & Swinney Reference Brandstater and Swinney1987; Vastano & Moser Reference Vastano and Moser1991). The inverse of the leading LE provides the time scale in which small initial errors saturate in the nonlinear growth regime, which is a first estimate on the fundamental predictability horizon of chaotic systems. On a more fundamental level, the sum of all positive LEs, the Kolmogorov–Sinai entropy, is equal to the average growth rate of the number of possible future states of the system, a measure of information production (Latora & Baranger Reference Latora and Baranger1999; Boffetta et al. Reference Boffetta, Cencini, Falcioni and Vulpiani2002; Berera & Clark Reference Berera and Clark2019). From a practical perspective, the Kolmogorov–Sinai entropy quantifies the rate at which an observer needs to acquire information about the flow (by measuring) to determine its true state. This definition connects the dynamics of uncertainty propagation with the notion of system complexity, which is intuitively related to how difficult or costly it is to predict the future (Boffetta et al. Reference Boffetta, Cencini, Falcioni and Vulpiani2002; Crutchfield Reference Crutchfield2012; Vela-Martín Reference Vela-Martín2024a ).
It is well known that, in systems with multiple time scales, high LEs describe uncertainty propagation mainly in the fastest scales of the flow. It is possible to characterise uncertainty propagation across scales by analysing the full spectrum of exponents, some of which have a support – the Lyapunov vectors – that covers the large scales (Yamada & Ohkitani Reference Yamada and Ohkitani1987). But this approach is problematic due to the large computational effort required to sample the high-dimensional phase space of turbulence, and the numerical instability of the orthonormalisation procedure used to construct the Lyapunov vectors. As a consequence, full Lyapunov spectra have only been computed in toy models (Yamada & Ohkitani Reference Yamada and Ohkitani1987) or marginally turbulent flows (Grappin & Léorat Reference Grappin and Léorat1991; Keefe et al. Reference Keefe, Moin and Kim1992; Berera & Clark Reference Berera and Clark2019). A more fundamental limitation of the LEs stems from the use of infinitesimal perturbations to describe uncertainty growth. This restricts the information that they provide to a small vicinity (in phase space) of each turbulent state, failing to capture non-local effects important for long-term uncertainty propagation (Palmer Reference Palmer1993).
A more natural characterisation of uncertainty propagation in turbulence comes from analysing error growth beyond the linear regime. This approach was pioneered by Leith & Kraichnan (Reference Leith and Kraichnan1972) in models of turbulence and extended more recently by Aurell et al. (Reference Aurell, Boffetta, Crisanti, Paladin and Vulpiani1996a , Reference Aurell, Boffetta, Crisanti, Paladin and Vulpianib , Reference Aurell, Boffetta, Crisanti, Paladin and Vulpiani1997) and Boffetta & Musacchio (Reference Boffetta and Musacchio2017). Following perturbations beyond their linear regime allows one to define finite-size Lyapunov exponents (Aurell et al. Reference Aurell, Boffetta, Crisanti, Paladin and Vulpiani1997), which capture the propagation of disturbances from small to large scales in a way that is better suited for turbulence physics than the infinitesimal LEs. In particular, the finite-size LEs show a scaling consistent with inertial dynamics, which is key to justifying the existence of a finite predictability horizon in high-Reynolds-number turbulence. Finite-size perturbations are also better connected to the phenomenology of the energy cascade than the LEs: Métais & Lesieur (Reference Métais and Lesieur1986) showed that error growth is described in scale space by a transfer of energy from correlated, error-free scales to decorrelated, error-dominated scales, opening the possibility of describing uncertainty propagation as a cascade, in agreement with the ideas of Leith & Kraichnan (Reference Leith and Kraichnan1972).
More recent works have tried to identify fundamental mechanisms of error growth by analysing the evolution of linear (Nikitin Reference Nikitin2018) and nonlinear perturbations in physical space (Encinar & Jiménez Reference Encinar and Jiménez2023; Ge, Rolland & Vassilicos Reference Ge, Rolland and Vassilicos2023). These studies find that error growth fluctuates strongly in space and time, and that it depends mainly on the local structure and intensity of the rate-of-strain tensor. In particular, Ge et al. (Reference Ge, Rolland and Vassilicos2023) showed that uncertainty amplification is mostly related to the compression of perturbations by the rate-of-strain tensor. This is in agreement with the earlier theoretical work of Ruelle (Reference Ruelle1982), who showed the connection between perturbation growth and energy dissipation, and suggests that, in the fully nonlinear regime, uncertainty propagation inherits most of the complexity of turbulence and cannot be described by the magnitude and scaling of the LEs (or their finite-size counterparts) alone. Although error growth can be attributed locally in physical space to a single mechanism mediated by the rate-of-strain tensor, the coupling of the velocity gradients and the transport due to conservative fluxes suggest that uncertainty evolution emerges to a large extent from complex non-local interactions.
This paper builds on the research introduced above and extends the study of uncertainty propagation in turbulence to a general statistical framework, guided by the idea that uncertainty is better described by probability distributions than by pairs of flow fields. Consider a turbulent flow field as a point in phase space, and its uncertainty as a probability distribution localised around it. As the flow evolves in time, the initially compact probability distribution spreads across phase space due to chaos, reaching many possible future states consistent with the same initial condition, and reducing our knowledge – increasing the uncertainty – about the true state of the flow (Latora & Baranger Reference Latora and Baranger1999; Vela-Martín & Jiménez Reference Vela-Martín and Jiménez2021). Uncertainty propagation is fully characterised by this evolution, which is, however, extremely difficult to describe due to the high-dimensionality and complex structure of the attractors underlying turbulent flows. A way to partly simplify this problem is to partition phase space into a finite set of states and study the properties of the discrete probability distribution. This is feasible in time-discrete maps, low-dimensional chaotic systems (Latora et al. Reference Latora, Baranger, Rapisarda and Tsallis2000) or very coarse-grained partitions of the turbulence attractor (e.g. Yalnız et al. Reference Yalnız, Hof and Budanur2021; Jiménez Reference Jiménez2023; Vela-Martín & Avila Reference Vela-Martín and Avila2024; Vela-Martín Reference Vela-Martín2024a ). Another possibility is to retain a continuous description of the probability distribution but truncate the analysis to the lowest-order moments. In fact, evolution equations for the first and second moments of the distribution can be derived from the Navier–Stokes equations, but these equations need to be supplemented by a closure model (Epstein Reference Epstein1969).
This closure problem can be partly circumvented by approximating the probability distribution by a cloud of points, each representing a single realisation of the Navier–Stokes equations with initial conditions drawn from the initial distribution. This method, named Monte Carlo ensemble forecasting (Leith Reference Leith1974), is widely used in numerical weather predictions that typically imply a number of ensemble members (independent simulations) of the order of a few tens (Molteni et al. Reference Molteni, Buizza, Palmer and Petroliagis1996; Leutbecher & Palmer Reference Leutbecher and Palmer2008). From a fundamental perspective, this method is limited by the high computational cost of producing ensembles with enough members to accurately capture the evolution of the probability distribution.
Another problem of the Monte Carlo ensemble method is the difficulty of sampling the invariant probability measure of turbulent flows, which is supported on a chaotic attractor with a dimension much smaller than the number of degrees of freedom. This implies that, given an initial condition, adding random perturbations in order to sample the initial probability distribution will most likely lead to flows that are outside the attractor. This is typically solved by resorting to reconstruction or data assimilation techniques, which allow one to sample statistically consistent turbulent flows conditioned to partial observations (Molteni et al. Reference Molteni, Buizza, Palmer and Petroliagis1996; Toth & Kalnay Reference Toth and Kalnay1997; Mons, Du & Zaki Reference Mons, Du and Zaki2021; Arun, Bae & McKeon Reference Arun, Bae and McKeon2023). However, these methods are computationally very expensive, limiting their application in Monte Carlo sampling methods.
In this paper, uncertainty propagation in isotropic turbulence is studied using the Monte Carlo ensemble method described in Vela-Martín (Reference Vela-Martín2024b ), which addresses the issues described above in two ways. First, the ensembles are constructed with enough flow realisations to converge the evolution of the ensemble-averaged field and its variance, locally in space and scale, with an uncertainty that is small in Kolmogorov units. This requires a number of ensemble members in the range of a few thousands, but in turns allows one to study uncertainty propagation in a way that is consistent with the evolution of probability distributions. Second, uncertainty is introduced in the initial conditions by adding a random perturbation that is small enough so that the ensembles return sufficiently fast to the turbulent attractor (Vela-Martín & Avila Reference Vela-Martín and Avila2024). This avoids the use of reconstruction techniques for ensemble initiation, and ensures that the evolution of the ensembles becomes independent of the structure of the perturbations in a short time. The magnitude of the perturbations is smaller than the Kolmogorov scale, which represents the ideal case in which all the relevant scales of the flow are accessible to observations but uncertainty emerges from unavoidable measurement errors.
Using this approach, Vela-Martín (Reference Vela-Martín2024b ) showed that the propagation of uncertainty from small to large scales in isotropic turbulence progressively filters out the small scales of the ensemble-averaged field. This process reflects the loss of predictability in the small scales, and is akin to a filtering operation with a low-pass Gaussian kernel. The scale of this effective filter, which increases with time, separates the flow in two parts: first, the scales above, which can be predicted or reconstructed up to the lead time; and second, the scales below, which have lost all predictability and cannot be reconstructed by any means. Thus, this scale represents the minimum predictable or reconstructable scale of the flow. In agreement with the evolution of finite-size perturbations (Aurell et al. Reference Aurell, Boffetta, Crisanti, Paladin and Vulpiani1996b ; Boffetta & Musacchio Reference Boffetta and Musacchio2017; Ge et al. Reference Ge, Rolland and Vassilicos2023), the minimum predictable scale grows with time following a self-similar evolution consistent with an inertial scaling.
These results will be the starting point of this paper and lay the foundations to understand uncertainty propagation in terms of a cascade, as suggested by Leith & Kraichnan (Reference Leith and Kraichnan1972) and Métais & Lesieur (Reference Métais and Lesieur1986). This is not a vague analogy, but will be well supported here by ample evidence. First, the ensemble average operator will be applied to the Navier–Stokes equations to derive equations for the ensemble-averaged flow and its variance, where the former represents the best prediction consistent with the initial uncertainty in the least-squares sense (Leith Reference Leith1974), and the latter is a natural measure of uncertainty. It will be shown that, as anticipated by Vela-Martín (Reference Vela-Martín2024b ), the evolution of the ensemble-averaged flow is described by an equation similar to that of large-eddy simulations (LES), in which a term analogous to the subgrid scale (SGS) flux moves energy from the ensemble average to the variance. This term is interpreted as an inverse flux of uncertainty from the unpredictable scales of the flow below the minimum predictable scale, to the predictable or uncertainty-free scales above. This, together with the filtering effect that uncertainty propagation has on the ensemble-averaged flow, points to a process that resembles a cascade in a mathematical and phenomenological sense. To further support this view, the uncertainty cascade will be analysed following the standard formalism of the kinetic energy cascade. The parallelisms between uncertainty propagation and the kinetic energy cascade are many, from an inertial range with a constant uncertainty flux, to scale locality and intermittency effects.
The implications of viewing uncertainty propagation as a cascade will be discussed from a fundamental and applied viewpoint, and beyond the problem of predictability. In particular, it will be argued that uncertainty fluxes are key to small-scale reconstruction problems in turbulence, in which the uncertainty dynamics is key to determining how much information propagates from large to small scales, setting a fundamental limit to the smallest scale that can be informatively reconstructed. This limit is essential to assess the performance of machine-learning super-resolution methods or data assimilation techniques, and to quantifying the uncertainty of SGS models.
This paper is organised as follows: in § 2, the methodology is laid out with details on the construction of the ensembles and the numerical method; §§ 3 and 4 describe the propagation of uncertainty in scale space, and the structure of the ensemble-averaged flow fields. Both sections contain key results that are then employed in § 5 to justify the notion of an uncertainty cascade. Next, § 6 contains an analysis on the structure of this cascade. Finally, the conclusions are offered in § 7.
2. Massive ensembles of direct numerical simulations
2.1. Direct numerical simulations of isotropic turbulence
In this work, we will consider homogeneous isotropic turbulence in a triply periodic cubic box of volume
$(2\pi )^3$
. This flow is described by the incompressible Navier–Stokes equations,
\begin{equation} \begin{aligned} \partial _t \boldsymbol{u} + \boldsymbol{u}\boldsymbol{\cdot }\boldsymbol{\nabla }\boldsymbol{u} &= -\boldsymbol{\nabla }p + \nu {\nabla} ^2\boldsymbol{u} + \boldsymbol{f}, \\ \boldsymbol{\nabla }\boldsymbol{\cdot }\boldsymbol{u} &=0, \end{aligned} \end{equation}
where
$\boldsymbol{u}$
is the velocity vector field,
$p$
is the kinematic pressure,
$\nu$
is the kinematic viscosity and
$\boldsymbol f$
is a forcing term proportional to the velocity which acts only on the large scales. In Fourier space, the forcing reads
\begin{equation} \widehat {\boldsymbol f}=\alpha (k,t)\widehat {\boldsymbol u}, \end{equation}
where
\begin{equation} \alpha (k,t)=\frac {\gamma }{\sum \limits _{k\lt 2}\widehat {\boldsymbol{u}}\widehat {\boldsymbol{u}}^*},\quad \text{ if $k \lt 2$,} \end{equation}
and
$\alpha (k,t)=0$
otherwise. The caret denotes the Fourier transform
$\widehat {\boldsymbol{u}}=F(\boldsymbol{u})$
(
$F^{-1}$
is the inverse transform) and the asterisk the complex conjugate, and the wavenumber vector and its magnitude are
$\boldsymbol{k}$
and
$k=|\boldsymbol k|$
, respectively. In (2.3), the summation is taken only over wavenumbers
$k\lt 2$
. The forcing parameter
$\alpha$
changes in time so that the instantaneous energy injection rate is constant and equal to
$\gamma$
. This implies that the average kinetic energy dissipation is
$\varepsilon =\nu \overline {|\boldsymbol{\nabla }\boldsymbol{u}|^2}=\gamma$
, where the overline denotes the unconditional ensemble average over independent flow fields in the turbulence attractor. The Navier–Stokes equations in Fourier space are integrated on a grid of
$N^3$
collocation points in physical space using a standard fully dealiased pseudospectral method for the nonlinear terms, and a third-order Runge–Kutta scheme for time marching. Further details of the code and the forcing scheme can be found in Cardesa et al. (Reference Cardesa, Vela-Martín, Dong and Jiménez2015) and Cardesa, Vela-Martín & Jiménez (Reference Cardesa, Vela-Martín and Jiménez2017).
Two different Reynolds numbers are considered in this work,
${\textit{Re}}_\lambda =U\lambda /\nu =120$
and
${\textit{Re}}_\lambda =195$
, where
$\lambda =(15\nu /{\varepsilon })^{1/2}U$
is the Taylor microscale, and
\begin{equation} U=\sqrt {\frac {1}{3}{\overline {|\boldsymbol{u}|^2}}} \end{equation}
are the root-mean-squared velocity fluctuations. The
${\textit{Re}}_\lambda =120$
and
${\textit{Re}}_\lambda =195$
simulations are run on grids of
$N=128$
and
$N=256$
points, respectively. In all simulations, the energy injection rate,
$\gamma$
, is fixed so that the numerical resolution is
$k_{max }\eta =1$
, where
$k_{max }=\sqrt {2}/3N$
is the maximum wavenumber after dealiasing and
$\eta =(\nu ^3/{\varepsilon })^{1/4}$
is the Kolmogorov length scale. This resolution is smaller than in standard direct numerical simulations (typically
$k_{max }\eta \gt 1.5$
), but this choice is justified by the focus on inertial range dynamics and the need to achieve large Reynolds numbers at an affordable computational cost, as in Encinar & Jiménez (Reference Encinar and Jiménez2023).
The time scales of the flow are the Kolmogorov time scale,
$t_\eta =(\nu /{\varepsilon })^{1/2}$
, and the large-scale eddy-turnover time
\begin{equation} T=\frac {{L}}{U}, \end{equation}
where
\begin{equation} {L}=\frac {3\pi }{4}\frac {\sum k^{-1}\overline {E}_{\boldsymbol{u}}}{\sum \overline {E}_{\boldsymbol{u}}} \end{equation}
is the integral length scale (Pope Reference Pope2001) (the summation is taken over all wavenumbers) and
\begin{equation} E_{\boldsymbol{u}}(k,t)=2\pi k^2\langle \boldsymbol{\widehat u}\boldsymbol{\widehat u}^*\rangle _{k} \end{equation}
is the instantaneous kinetic energy spectrum, where
$\langle \boldsymbol{\cdot }\rangle _k$
represents the average over spherical wavenumber shells.
2.2. Ensembles of simulations with small perturbations
The propagation of uncertainty is studied here by means of massive Monte Carlo ensemble simulations, as described in Vela-Martín (Reference Vela-Martín2024b
). First, we consider a collection of
$\mathcal{N}_b$
independent initial conditions sampled homogeneously from the turbulence attractor. These flow fields are hereafter termed base flows and marked with a superscript,
$\boldsymbol{u}^b$
, where
$b=1,\dots ,\mathcal{N}_b$
. Around each base flow,
$\mathcal{N}_p$
perturbed initial conditions are produced by adding small perturbations
\begin{equation} \boldsymbol{u}^{b,p}=\boldsymbol{u}^b + \delta \boldsymbol{\phi }^p, \end{equation}
where
$\delta$
is a scalar and each
$\boldsymbol{\phi }^p$
, for
$p=1,\dots ,\mathcal{N}_p$
, is an independent realisation of a Gaussian random vector field with a white spectrum, zero mean and unit variance. Note that the equations describing the flow in (2.1) are strictly deterministic, including the forcing, so that the only source of uncertainty comes from the perturbations in the initial conditions. The perturbation fields are automatically projected into a divergence-free field by the code so that
$\boldsymbol{\nabla }\boldsymbol{\cdot }\boldsymbol{u}^{b,p}=0$
. The number of base and perturbed flows is
$\mathcal{N}_b=8$
and
$\mathcal{N}_p=4000$
, respectively. The number of base flows is limited by computational resources and memory storage, but the analysis by (Vela-Martín Reference Vela-Martín2024b
) suggests that it provides well-converged statistics. The ensembles of perturbed flows and the base flows are advanced in time approximately
$9T$
and statistics of the ensembles are produced on the fly approximately each
$0.5T$
. The total number of simulations for each Reynolds numbers is
$32\,000$
, which is enough to converge the average and variance of the ensemble with an uncertainty that is small in Kolmogorov units (Vela-Martín Reference Vela-Martín2024b
). Table 1 shows a summary of the database analysed in this paper.
Table 1. Main parameters of the simulations, where
$N$
is the number of grid points in each direction,
$k_{max }=\sqrt {2}/3N$
is the maximum Fourier wavenumber magnitude resolved in the simulations after dealiasing with exact phase-shift and
$T_{{\textit{total}}}$
is the total number of turnover times run in all the simulations.
In each ensemble, the perturbed fields represent a sampling of a probability distribution initially concentrated around a base flow, and
$\delta$
quantifies the spread of this distribution, i.e. it is a measure of the initial uncertainty. This parameter is set to
$\delta =7\times 10^{-3}U\lt \eta /t_\eta$
in all simulations, which ensures that the evolution of the perturbed flows with respect to the base flows is initially linear, and driven by the most unstable directions of the flow. This yields a distribution of perturbations that becomes independent of its initial structure in a short time (Vela-Martín & Avila Reference Vela-Martín and Avila2024). From an empirical perspective, the small value of
$\delta$
represents the ideal case in which every scale of the flow is observable but uncertainty emerges from unavoidable errors. The perturbations could also represent the effect of limited float-point precision in direct numerical simulations, or unavoidable thermal fluctuations at sub-Kolmogorov scales (Bandak et al. Reference Bandak, Goldenfeld, Mailybaev and Eyink2022).
The magnitude of
$\delta$
is not critical to the results presented in the following provided that it is small, and decreasing it further only increases the linear evolution regime, producing only a temporal offset in the data, as described in § 3.
2.3. The ensemble-averaged flow and its variance
In the following, the ensemble-averaged flow and its variance will be used to characterise uncertainty propagation in physical and scale space. The ensemble-averaged flow is defined as
\begin{equation} \{\boldsymbol u\}^b(\boldsymbol{x},t)=\frac {1}{\mathcal{N}_p}\sum _p \boldsymbol u^{b,p}, \end{equation}
where the summation is taken over all the perturbed flow fields of the ensemble. The superscript in the operator is used to make explicit that the average is taken over perturbations of the
$b$
th base flow. This operator is different from the standard unconditional ensemble average (denoted by an overbar) in that the average is conditional to flow fields that are initially close to
${\boldsymbol{u}}^b$
.
With this definition of the ensemble average, we obtain the variance of the velocity vector
\begin{equation} \mathcal{S}(\boldsymbol{x}, t)=\{\mathcal{E}_{\boldsymbol{u}}\}^b - \mathcal{E}_{\{\boldsymbol{u}\}^b}, \end{equation}
where
\begin{equation} \mathcal{E}_{\boldsymbol{u}}(\boldsymbol{x},t)=\frac {1}{2}u^{b,p}_iu^{b,p}_i \end{equation}
is the local kinetic energy of each perturbed flow (repeated indexes imply summation), and
\begin{equation} \mathcal{E}_{\{\boldsymbol{u}\}^b}(\boldsymbol{x},t)=\frac {1}{2}\{u_i\}^b\{u_i\}^b \end{equation}
is the kinetic energy of the ensemble-averaged field. For ease of notation, the superscript ‘b’ has been omitted in the definition of
$\mathcal{S}$
, but this quantity is a property of each ensemble, similar to the ensemble-averaged field. This definition of the variance differs from the standard by a factor of one half, and does not take into account the correlation between components of the velocity vector or between the velocity vectors at different positions. This information is contained in the full covariance matrix, which has a size
$O(N^6)$
and is not used here for practical reasons. An important aspect of
$\{\boldsymbol{u}\}^b(\boldsymbol{x},t)$
and
$\mathcal{S}(\boldsymbol{x},t)$
is that they are a vector and a scalar field, respectively, which evolve in time with dynamics that emerge directly from the Navier–Stokes equations, as we will see in § 5.
The variance has an equivalent definition in scale space, the spectral variance, defined as
\begin{equation} S(k,t)=\{E_{\boldsymbol u}\}^b - E_{\{\boldsymbol u\}^b}, \end{equation}
where the first term represents the average energy spectrum of the perturbed flows in the ensemble, and the second is the spectrum of the ensemble-averaged flow,
\begin{equation} E_{\{\boldsymbol u\}^b}(k,t)=2\pi k^2\langle \{\widehat {\boldsymbol u}\}^b\{\widehat {\boldsymbol u}^*\}^b\rangle _k. \end{equation}
By Parseval’s theorem, the spectral ensemble variance is related instantaneously to the variance by
\begin{equation} \sum _{k}S(k,t)=\langle \mathcal{S}\rangle , \end{equation}
where the brackets denote the spatial average. The variance in physical and spectral space are equivalent to the error field and the error spectrum between pairs of trajectories (Métais & Lesieur Reference Métais and Lesieur1986), but with deviations measured with respect to the ensemble-averaged flow.
The significance of the ensemble average and variance is key to understand their relation to uncertainty propagation. The ensemble-averaged flow represents the best possible prediction in the least-squares sense consistent with the initial conditions, i.e. it is the field that has the minimum average error (
$L^2$
-norm of the difference between fields) with respect to all members in the ensemble (Epstein Reference Epstein1969). Similarly, the ensemble variance represents the average squared error with respect to the best possible prediction, and is thus minimum: an average squared error defined with respect to a field other than the ensemble average is larger. This includes the average error between samples in the ensemble, as in the standard pairwise analysis of errors in turbulence, which is on average double the ensemble variance.
Figure 1. Evolution of cuts of (a–e) the enstrophy of the ensemble-averaged flow field normalised by its instantaneous spatial average,
${\varOmega }_{\{\boldsymbol{u}\}^b}/\langle \varOmega _{\{\boldsymbol{u}\}^b}\rangle$
, and (f– j) the logarithm of the variance normalised by the root-mean-squared velocity,
$\log _{10}(\mathcal{S}/U^2)$
. All snapshots correspond to the same ensemble, and time goes from left to right,
$t/T=$
: (a,f) 0; (b,g) 2.6; (c,h) 3.2; (d,i) 4.8; (e,j) 5.8. In panels (a–e), the colour scale goes from
$0$
(dark blue) to
$8$
(light yellow), and in (f–j) from
$-3$
(dark blue) to
$0$
(light yellow). The size of the domain is normalised with Kolmogorov units.
3. Uncertainty propagation in scale space
This section presents a preliminary analysis on the evolution of the ensemble-averaged flow and its variance, providing a first picture of the physics of uncertainty propagation across scales. At
$t=0$
,
$\{\boldsymbol{u}\}^b=\boldsymbol{u}^b$
because the ensemble average of all the perturbations is zero by construction. As time evolves, chaos causes the divergence of the perturbed flows, which become progressively decorrelated, leading to an increment of the ensemble variance. For sufficiently long times, typically of the order of a few eddy turnover times, the ensemble members have diverged completely and are homogeneously distributed in the turbulence attractor. At this point, any statistic of the ensemble is equal to the statistic of the attractor, except for fluctuations due to finite-sample effects (DelSole Reference DelSole2004). In particular, the ensemble-averaged velocity field relaxes asymptotically towards equilibrium, in which
$\{\boldsymbol{u}\}^b\rightarrow \overline {\boldsymbol{u}}=\boldsymbol{0}$
due to spatial homogeneity.
The relaxation of an initially concentrated probability distribution to the equilibrium distribution over the attractor happens faster in the small scales. This leads to a progressive cancellation of the ensemble-averaged field and the growth of uncertainty from small to large scales (Leith Reference Leith1974). This process is well described in figure 1(a–e), which shows the evolution of cuts of the enstrophy of the ensemble-averaged flow field,
$\varOmega _{\{\boldsymbol{u}\}^b}=\{\omega _i\}^b\{\omega _i\}^b$
, where
$\boldsymbol \omega =\boldsymbol{\nabla }\times \boldsymbol u$
is the vorticity vector. The enstrophy of the ensemble-averaged field is progressively smoothed with time, with a characteristic length scale that grows from viscous to integral scales in a few turnover times. This smoothing happens in parallel with a growth of
$\mathcal S(\boldsymbol{x},t)$
, which has a small-scale structure at short times but reaches a size comparable to the integral scales, as shown in figure 1(f–j).
The progressive smoothing due to uncertainty growth is quantified by the instantaneous kinetic energy spectra of the ensemble-averaged flow,
$E_{\{\boldsymbol u\}^b}(k,t)$
, and the spectral ensemble variance,
$S(k,t)$
. The evolution of these two quantities in a single ensemble are shown in figure 2(a,b). The gradual loss of small-scale energy in
$E_{\{\boldsymbol u\}^b}$
happens at the same time as the growth of
$S(k,t)$
from small to large scales. The decay of
$E_{\{\boldsymbol u\}^b}$
may be interpreted as the loss of predictable or correlated energy, and the growth of
$S$
as the increment of unpredictable or uncorrelated energy (Métais & Lesieur Reference Métais and Lesieur1986).
Figure 2. (a) Temporal evolution of the energy spectra of the ensemble-averaged flow and of the true flow filtered at the minimum reconstructable scale . (b) Temporal evolution of the spectral variance. The markers in both figures correspond to the same times as in figure 1,
$t/T$
: circles, 2.6; squares, 3.2; diamonds, 4.8; triangles, 5.8.
Vela-Martín (Reference Vela-Martín2024b
) showed that the separation between the predictable and unpredictable scales of the flow is well described by an exponential decay of
$E_{\{\boldsymbol u\}^b}$
towards the small scales, which implies that uncertainty propagation is akin to a Gaussian filter. To show this, let us first define a filtering operation
\begin{equation} \widetilde {\boldsymbol u}=F^{-1}[\widehat {\boldsymbol{u}}\boldsymbol{\cdot }G(k,\ell )], \end{equation}
where
\begin{equation} G(k;\ell )=\exp \left (-\frac {k^2\ell ^2}{2\pi ^2} \right ) \end{equation}
is a Gaussian kernel. Now, let us take the base flow at a time
$t$
,
$\boldsymbol u^b(\boldsymbol{x},t)$
, and find the filter scale
$\ell ^b(t)$
for which the instantaneous space-averaged enstrophy of
$\{\boldsymbol u\}^b(\boldsymbol{x}, t)$
is equal to that of the filtered base flow
$\widetilde {\boldsymbol u}^b$
, i.e.
$\langle \varOmega _{\{\boldsymbol{u}\}^b}\rangle =\langle \varOmega _{\widetilde {\boldsymbol{u}}^b}\rangle$
. This condition is formulated in terms of the energy spectra as
\begin{equation} \sum _k k^2E_{\{\boldsymbol u\}^b}=\sum _k k^2E_{\widetilde {\boldsymbol u}^b}, \end{equation}
where
\begin{equation} E_{\widetilde {\boldsymbol u}^b}=E_{{\boldsymbol u}^b}G^2(k;\ell ^b)=E_{{\boldsymbol u}^b}\exp \left (-\frac {k^2{\ell ^b}^2}{\pi ^2}\right ) \end{equation}
is the spectrum of the filtered base flow at scale
$\ell ^b(t)$
. For each base flow, an iterative solver (Newton method) is used to find
$\ell ^b(t)$
as a function of time.
As shown in figure 2(a), the condition in (3.3) provides a very good fit between
$E_{\{\boldsymbol u\}^b}$
and
$E_{\widetilde {\boldsymbol u}^b}$
. This indicates that uncertainty propagation through scale acts on the ensemble-averaged flow as a Gaussian low-pass filter. The choice of a Gaussian kernel here is heuristic, and the reasons for the good fit are unclear. The robustness of the fit has been tested by calculating
$\ell ^b(t)$
to match the energy instead of the enstrophy, and the results were quantitatively very similar. The filter scale,
$\ell ^b(t)$
, increases with time and is interpreted here as the minimum scale below which informative predictions of the velocity field are not possible after time
$t$
, given the uncertainty in the initial conditions. This is so because the ensemble-averaged flow field is the best possible reconstruction of the flow in the least-squares sense, and the decay of its energy below
$\ell ^b(t)$
implies that information about these scales cannot be recovered by any means. In other words, the error made by trying to reconstruct the scales below
$\ell ^b(t)$
is larger on average than not making a reconstruction at all, i.e. than setting the energy of these scales to zero. This interpretation is limited to the space-averaged, second-order statistics of the velocity field, as
$\ell ^b(t)$
does not take into account spatial fluctuations of predictability, and, most importantly, the information contained in higher-order moments of the probability distribution described by the ensembles.
In summary, the smoothing of the ensemble-averaged flow represents the average filtering out of unpredictable flow features. In the context of forecasting,
$\ell ^b(t)$
was termed the minimum predictable scale by Vela-Martín (Reference Vela-Martín2024b
). For the sake of generality, we adopt here the term minimum reconstructable scale (MRS).
Figure 3. (a) Temporal evolution of the average MRS as a function of time for different Reynolds numbers. The dash–dotted line is proportional to
$t^{3/2}$
. The black solid line without markers is
$\overline {\ell }_r$
, as defined in (3.7), which has been linearly scaled for comparison. (b) Average uncertainty as a function of scale normalised in Kolmogorov units. Line styles as in panel (a). The dash–dotted line is proportional to
$\ell ^{2/3}$
. In panels (a) and (b), the markers correspond to the same times (scales) as in the plots of figure 1.
The temporal evolution of
$\ell ^b(t)$
averaged over the
$\mathcal{N}_b$
ensembles,
\begin{equation} \overline {\ell }(t)=\frac {1}{\mathcal{N}_b}\sum _b\ell ^b(t), \end{equation}
is shown in figure 3(a). Hereafter we use the overline indistinctly to denote the unconditional average over the attractor or over baseflows. The fluctuations of
$\ell ^b(t)$
across different base flows is much smaller that its average, with
\begin{equation} ({\overline {\ell ^2} - {\overline {\ell }}^2})^{1/2}\approx 0.05\overline {\ell }, \end{equation}
approximately constant in the inertial range, so that the evolution of the average MRS is representative of the individual ensembles. This is an indicator that the definition of the minimum predictable scale is robust to the strong fluctuations of the average energy dissipation in each base flow, as shown by Vela-Martín (Reference Vela-Martín2024b ).
The evolution of the average MRS approaches an incipient power-law regime consistent with an inertial scaling,
$\overline {\ell }\sim t^{3/2}$
. This scaling is in agreement with the self-similar growth of errors reported in models (Leith & Kraichnan Reference Leith and Kraichnan1972; Métais & Lesieur Reference Métais and Lesieur1986; Aurell et al. Reference Aurell, Boffetta, Crisanti, Paladin and Vulpiani1997) and direct numerical simulations (Boffetta & Musacchio Reference Boffetta and Musacchio2017) of isotropic turbulence. In fact, the evolution of
$\overline {\ell }$
is proportional to the evolution of an error scale defined by Leith & Kraichnan (Reference Leith and Kraichnan1972), Métais & Lesieur (Reference Métais and Lesieur1986) and Boffetta & Musacchio (Reference Boffetta and Musacchio2017) based on a prescribed error-to-energy ratio. This scale is redefined here as
$\ell _r(t)=\pi /k_r(t)$
, where
$k_r(t)$
fulfils
\begin{equation} S(k_r(t),t)=r \{E_{\boldsymbol u}\}^b(k_r(t),t), \end{equation}
and
$r$
is a constant scalar threshold. The scale
$\ell _r(t)$
is defined as the scale in which the uncertainty is a fraction
$r$
of the energy, which are measured by the spectral variance and the energy spectrum, respectively. In figure 3(a),
$\overline {\ell }_r(t)$
for
$r=0.95$
is compared with
$\overline {\ell }(t)$
to show that the former is proportional to the later in the inertial scales.
There is a good collapse between the evolution of
$\overline {\ell }(t)$
for
${\textit{Re}}_\lambda =195$
and
${\textit{Re}}_\lambda =120$
, as shown in figure 3(a). In this plot,
$\overline {\ell }(t)$
for the lower Reynolds number has been time-shifted towards the past to account for the different magnitude of the initial perturbations and of the leading Lyapunov exponents in Kolmogorov units (Crisanti et al. Reference Crisanti, Jensen, Vulpiani and Paladin1993; Boffetta & Musacchio Reference Boffetta and Musacchio2017; Mohan et al. Reference Mohan, Fitzsimmons and Moser2017; Berera & Ho Reference Berera and Ho2018). The fact that, despite this difference, the evolution of
$\overline {\ell }(t)$
is indistinguishable across the two Reynolds numbers corroborates that the evolution of the ensemble is independent of the magnitude of the initial perturbation. In fact, decreasing the magnitude of the initial perturbation only increases the initial transient in which perturbations evolve linearly, but leads to the same long-time evolution once the magnitude of the perturbations reaches the nonlinear regime (Vela-Martín & Avila Reference Vela-Martín and Avila2024). This implies that the temporal origin of
$\overline {\ell }(t)$
can be shifted arbitrarily towards the past by reducing the magnitude of the initial perturbation, without affecting its inertial evolution. To avoid this ambiguity, hereafter the MRS is used instead of the lead time as a reference to characterise the propagation of uncertainty across scales.
The inertial scaling of uncertainty propagation is further corroborated by the evolution of the total uncertainty,
\begin{equation} \langle \overline {\mathcal S}\rangle (t)=\sum _k \overline {S}(k,t), \end{equation}
which is shown in figure 3(b). There is an incipient inertial rage scaling of the form
$\langle \overline {\mathcal S}\rangle \sim {\overline {\ell }}{}^{2/3}$
. Similar to the MRS, also the evolution of the uncertainty collapses well at different Reynolds numbers when scaled in Kolmogorov units.
4. The structure of the ensemble-averaged flow fields
Up to now, we have reviewed the main result presented in Vela-Martín (Reference Vela-Martín2024b ): uncertainty propagation acts on the ensemble-averaged flow as a low-pass filter. In this section, we show that this effective filtering preserves the structure of turbulence in the ensemble-averaged flow, which is similar to that of turbulence filtered at inertial scales. As will be shown in the next section, this is key to understanding uncertainty propagation as a cascade produced by the interaction between correlated or predictable, and uncorrelated or unpredictable scales.
To describe the structure of the ensemble-averaged flow field at scale
$\ell (t)$
, we focus on the analysis of the velocity gradients, in particular the second and third invariants of the velocity gradient tensor,
$A_{\textit{ij}}=\partial _i u_j$
, defined as
\begin{equation} \begin{aligned} Q =& -\frac {1}{2}A_{\textit{ij}}A_{\textit{ji}}, \\ R =& -\frac {1}{3}A_{\textit{ij}}A_{\textit{jk}}A_{\textit{ki}}.\\ \end{aligned} \end{equation}
These invariants describe the local topology of the flow in terms of stable–unstable saddle nodes, and compressed or stretched rotating regions (Chong, Perry & Cantwell Reference Chong, Perry and Cantwell1990). They are also related to the dynamics of the velocity gradients;
$Q$
can be expressed in terms of the squared rate-of-strain tensor,
$\varSigma =2 S_{\textit{ij}}S_{\textit{ij}}$
, and the enstrophy,
$\varOmega =\omega _i\omega _i$
,
\begin{equation} Q=\frac {1}{4}(\varOmega - \varSigma ) \end{equation}
and
$R$
in terms of vortex-stretching and strain self-amplification,
\begin{equation} R=-\dfrac {1}{4}\omega _iS_{\textit{ij}}\omega _j - \dfrac {1}{3}S_{\textit{ij}}S_{\textit{jk}}S_{\textit{ki}}. \end{equation}
The joint probability distribution of
$Q$
and
$R$
represents a fundamental feature of turbulence dynamics which is preserved across scales (Lozano-Durán et al. Reference Lozano-Durán, Holzner and Jiménez2016; Danish & Meneveau Reference Danish and Meneveau2018; Johnson & Wilczek Reference Johnson and Wilczek2024), and has been connected to the prevalent mechanisms of the turbulence energy cascade (Eyink Reference Eyink2006; Carter & Coletti Reference Carter and Coletti2018; Johnson Reference Johnson2020; Carbone & Bragg Reference Carbone and Bragg2020; Johnson Reference Johnson2021) and the generation of dissipation (Tsinober Reference Tsinober2000).
In the following, the invariants of the ensemble-averaged flow field,
\begin{equation} \begin{aligned} {Q}_{\{\boldsymbol{u}\}^b}=&-\frac {1}{2}\{A_{\textit{ij}}\}^b\{A_{\textit{ji}}\}^b, \\ {R}_{\{\boldsymbol{u}\}^b}=&-\frac {1}{3}\{A_{\textit{ij}}\}^b\{A_{\textit{jk}}\}^b\{A_{\textit{ki}}\}^b,\\ \end{aligned} \end{equation}
will be compared with those of the true filtered flow,
$\widetilde {\boldsymbol{u}}^b$
, at scale
$\ell ^b(t)$
,
\begin{equation} \begin{aligned} Q_{\widetilde {\boldsymbol{u}}^b}=&-\frac {1}{2}\widetilde {A}^b_{\textit{ij}}\widetilde {A}^b_{\textit{ji}}, \\ R_{\widetilde {\boldsymbol{u}}^b}=&-\frac {1}{3}\widetilde {A}^b_{\textit{ij}}\widetilde {A}^b_{\textit{jk}}\widetilde {A}^b_{\textit{ki}},\\ \end{aligned} \end{equation}
where
$\widetilde {A}^b_{\textit{ij}}=\partial _j \widetilde {u}^b_{i}$
. In this analysis, we will also compare the squared rate-of-strain tensor and enstrophy of the filtered and ensemble-averaged field, denoted by
$\varSigma _{\widetilde {\boldsymbol{u}}^b}$
and
$\varOmega _{\widetilde {\boldsymbol{u}}^b}$
, and
$\varSigma _{\{\boldsymbol{u}\}^b}$
and
$\varOmega _{\{\boldsymbol{u}\}^b}$
, respectively.
Figure 4(a–d) shows the joint probability distribution of the invariants of the velocity gradient tensor for the filtered true flow and the ensemble-averaged flow at four different scales,
$\overline {\ell }=19\eta$
,
$53\eta$
,
$86\eta$
and
$120\eta$
. Both the filtered and the ensemble-averaged flow show a characteristic tear-drop shape with predominant vortex stretching in the
$Q\gt 0$
and
$R\lt 0$
quadrant, and a tail in
$Q\lt 0$
and
$R\gt 0$
characterised by dominant strain self-amplification. The probability density functions (p.d.f.s) of the ensemble-averaged and filtered fields collapse well, indicating that both have statistically similar features at different scales.
Figure 4. (a–d) Joint probability density function of the invariants of the velocity gradient tensor,
$Q$
and
$R$
, in the ensemble-averaged flow (dashed lines) and in the filtered true flow (solid lines) at MRSs
$\overline {\ell }=$
: (a)
$19\eta$
; (b)
$53\eta$
; (c)
$86\eta$
; (d)
$120\eta$
. The invariants are normalised by the time scale,
$t_Q=\langle \overline {Q^2}\rangle ^{-1/4}$
, where
$Q$
is calculated with the filtered true flow or the ensemble-averaged flow as it corresponds. The contours contain
$98\,\%$
,
$90\,\%$
and
$70\,\%$
of the data. The data corresponds to
${\textit{Re}}_\lambda =192$
.
Figure 5. Joint p.d.f. of the (a) squared rate-of-strain tensor and the (b) enstrophy in the ensemble-averaged flow and the filtered true flow. Quantities are normalised by
$\langle \overline {\varSigma _{\{\boldsymbol{u}\}^b}} \rangle =\langle \overline {\varSigma _{\widetilde {\boldsymbol{u}}^b}}\rangle$
at each scale. The contours contain
$90\,\%$
,
$80\,\%$
and
$60\,\%$
of the data, which correspond to the solid, dashed, and dash–dotted lines, respectively.
In addition to having similar statistics, the gradients of the ensemble-averaged and the filtered fields are also largely correlated in space. This is shown in figure 5(a,b), in which the joint p.d.f.s of
$\varSigma _{\{\boldsymbol{u}\}^b}$
and
$\varSigma _{\widetilde {\boldsymbol{u}}^b}$
, and
$\varOmega _{\{\boldsymbol{u}\}^b}$
and
$\varOmega _{\widetilde {\boldsymbol{u}}^b}$
, are plotted at different inertial scales. To quantify the degree of spatial correlation of these fields, we define the correlation coefficient
\begin{equation} \rho (\phi ,\psi )=\frac {\langle \overline {\phi '\boldsymbol{\cdot }\psi '}\rangle }{\sqrt {\langle \overline {\phi '{}^2}\rangle \langle \overline {\psi '{}^2}\rangle }} \end{equation}
where
$\psi$
and
$\phi$
are scalar fields, the prime denotes quantities without average, and the average is taken both in space and over ensembles. At the four scales shown in figure 5, the correlation coefficient of the enstrophy and the strain are
$\rho (\varOmega _{\widetilde {\boldsymbol{u}}^b},\varOmega _{\{\boldsymbol{u}\}^b})\approx 0.7$
and
$\rho (\varSigma _{\widetilde {\boldsymbol{u}}^b},\varSigma _{\{\boldsymbol{u}\}^b})\approx 0.85$
, respectively. This corroborates the high correlation between enstrophy or strain structures between the ensemble-averaged and filtered fields, in agreement with the quantitative and visual analysis by Vela-Martín (Reference Vela-Martín2024b
).
In brief, the results presented in this section indicate that the filtering caused by uncertainty propagation yields an ensemble-averaged flow with a structure that is locally similar to that of turbulence filtered at inertial scales.
5. The uncertainty cascade in physical space
5.1. The equations describing uncertainty propagation
In this section, the dynamics of uncertainty propagation are analysed by deriving evolution equations for the ensemble-averaged flow field, and its variance, in physical space. First, the ensemble average operator is applied to the Navier–Stokes equations to produce equations for each component of the ensemble-averaged velocity vector,
\begin{equation} \begin{aligned} \partial _t \{u_i\} + \{u_j\}\partial _j \{u_i\}&= -\partial _i \{p\} + \partial _j \tau ^{u}_{\textit{ij}} + \{f_i\} + \nu \partial _{kk}\{u_i\}, \\ \partial _i \{u_i\} &=0, \\ \end{aligned} \end{equation}
where the superscript ‘b’ of the ensemble-average operator is dropped for simplicity. We have taken into account that the ensemble-average operator commutes with the spatial and temporal derivatives. The ensemble-averaged advection term,
$\{u_j\partial _j u_i\}$
, has been decomposed in a self-advection term due to the ensemble-averaged flow
$\{u_j\}\partial _j \{u_i\}$
, plus the divergence of a symmetric stress tensor,
\begin{equation} \tau ^u_{\textit{ij}}=\{u_i\}\{u_j\}-\{u_i u_j\}, \end{equation}
where the ‘
$u$
’ superscript is used to mark its relation to the evolution of uncertainty. This tensor is similar to the SGS stress tensor that appears by filtering the Navier–Stokes equations (Leonard Reference Leonard1975; Eyink Reference Eyink1995; Meneveau & Katz Reference Meneveau and Katz2000; Carati, Winckelmans & Jeanmart Reference Carati, Winckelmans and Jeanmart2001; Moser, Haering & Yalla Reference Moser, Haering and Yalla2021),
\begin{equation} \tau ^e_{\textit{ij}}=\widetilde {u}_i\widetilde {u}_j-\widetilde {u_i u_j}, \end{equation}
where the ‘
$e$
’ superscript is used as a reference to the kinetic energy cascade. In the LES framework,
$\tau ^e_{\textit{ij}}$
represents the interaction between resolved and SGS (or subfilter scales). From this viewpoint, the stress tensor
$\tau ^u_{\textit{ij}}$
may be interpreted as the interaction between the reconstructable or predictable part of the flow above
$\ell ^b(t)$
, and the uncorrelated or unpredictable scales below.
To stress this interpretation, we focus now on the evolution of the variance in physical space,
$\mathcal{S}(\boldsymbol{x},t)$
, as defined in (2.10). To construct an equation for
$\mathcal{S}$
, we first obtain the equation for the kinetic energy of the ensemble-averaged flow,
$\mathcal{E}_{\{\boldsymbol{u}\}}=\tfrac {1}{2}\{u_i\}\{u_i\}$
by contracting (5.1) with
$\{u_i\}$
,
\begin{equation} \begin{aligned} \{D\}_t \mathcal{E}_{\{\boldsymbol{u}\}}&= \partial _i\varPhi ^u_i - \{S_{\textit{ij}}\}\tau ^u_{\textit{ij}} + \{u_i\}\{f_i\} - 2\nu \{S_{\textit{ij}}\}\{S_{\textit{ij}}\}, \\ \end{aligned} \end{equation}
where
$\{D\}_t=:\partial _t + \{u_j\}\partial _j$
is the material derivative following the ensemble-averaged field, and
\begin{equation} \varPhi ^u_i= - \{u_i\}\{p\} + \{u_j\}\tau ^u_{\textit{ij}} + 2\nu \{u_j\}\{S_{\textit{ij}}\} \end{equation}
contains the spatial energy fluxes due to the pressure, the uncertainty stress tensor and viscous diffusion of momentum. Following the same procedure, we derive an equation for
$\{\mathcal{E}_{\boldsymbol{u}}\}=(1/2)\{u_iu_i\}$
by contracting the Navier–Stokes equations with
$u_i$
and applying the ensemble-average operator, leading to
\begin{equation} \begin{aligned} \{D\}_t \{\mathcal{E}_{\boldsymbol{u}}\} &= \partial_i(\{u_i\}\{u_ju_j\}/2 - \{u_iu_ju_j\}/2 - \{u_ip\} + 2\nu \{u_j S_{ij}\}) + \{f_i u_i\} - 2\nu \{S_{\textit{ij}}S_{\textit{ij}}\}. \\ \end{aligned} \end{equation}
Subtracting (5.4) from (5.6), and grouping all spatial fluxes, yields
\begin{equation} \begin{aligned} \{D\}_t \mathcal{S} &= \partial_i\varPsi^u_i + \{S_{\textit{ij}}\}\tau ^u_{\textit{ij}} + \{f_i u_i\} - \{f_i\}\{u_i\} - 2\nu (\{S_{\textit{ij}} S_{\textit{ij}}\}-\{S_{\textit{ij}}\}\{S_{\textit{ij}}\}), \\ \end{aligned} \end{equation}
where
$\varPsi ^u_i = \{u_i\}\{u_ju_j\}/2 - \{u_iu_ju_j\}/2 - \{u_ip\} + 2\nu \{u_j S_{ij}\} - \varPhi ^u_i $
. The second term on the right-hand side represents interaction between the ensemble-averaged rate-of-strain tensor and the uncertainty stress tensor on the variance, and is hereafter denoted by
\begin{equation} \varPi ^u=\{S_{\textit{ij}}\}\tau ^u_{\textit{ij}}. \end{equation}
The four last terms on the right-hand side describe the effect of viscous dissipation and the forcing, and are compacted into
\begin{equation} \mathcal D=-2\nu (\{S_{\textit{ij}} S_{\textit{ij}}\}-\{S_{\textit{ij}}\}\{S_{\textit{ij}}\}), \end{equation}
and
\begin{equation} \mathcal F=\{f_i u_i\} - \{f_i\}\{u_i\}, \end{equation}
respectively. The term describing the effect of viscosity,
$\mathcal D$
, is always negative due to a Cauchy–Schwarz inequality
$\{S_{\textit{ij}}S_{\textit{ij}}\}\geqslant \{S_{\textit{ij}}\}\{S_{\textit{ij}}\}$
, which means that the dissipation always reduces the variance, i.e. it is a sink of uncertainty. The term
$\mathcal F$
can have any sign locally in space, but it is positive on average by the linear structure of the forcing term (this is easily tested by using Parseval’s theorem to recast the expression in Fourier space). The simulations are not run for long enough to let uncertainty reach the very large scales of the flow, meaning that the term
$\mathcal{F}$
is negligible and that all simulations in each ensemble have the same large-scale forcing. Expressing (5.7) using this compact notation, we obtain
\begin{equation} \begin{aligned} \{D\}_t \mathcal{S} &= \partial _i \varPsi _i + \varPi ^u + \mathcal{F} + \mathcal{D}. \\ \end{aligned} \end{equation}
A key aspect of this equation is that
$\varPi ^u$
appears also in the evolution equation of the energy of the ensemble-averaged flow with changed sign. Thus,
$\varPi ^u$
can be interpreted as a flux of energy from the ensemble-averaged flow to the variance. From the LES viewpoint,
$\mathcal{E}_{\{\boldsymbol{u}\}}$
represents the energy that is predictable or reconstructable, analogous to the resolved energy;
$\mathcal{S}$
is the unpredictable or reconstructable energy, analogous to the SGS kinetic energy, and
$\varPi ^u$
is a flux of energy from the latter to the former, analogous to the SGS energy fluxes.
To test this idea, in the following we will compare the uncertainty fluxes
$\varPi ^u$
with the kinetic energy fluxes of the cascade (Meneveau & Lund Reference Meneveau and Lund1994; Eyink Reference Eyink1995; Borue & Orszag Reference Borue and Orszag1998; Eyink & Aluie Reference Eyink and Aluie2009; Alexakis & Chibbaro Reference Alexakis and Chibbaro2020)
\begin{equation} \varPi ^e=\widetilde {S}_{\textit{ij}}\tau ^e_{\textit{ij}}, \end{equation}
where
$\tau ^e_{\textit{ij}}=\widetilde {u}_i\widetilde {u}_j - \widetilde {u_iu_j}$
and
$\widetilde {S}_{\textit{ij}}$
are calculated filtering the true base trajectory,
$\boldsymbol{u}^b$
, at the MRS. In this definition,
$\varPi ^e\gt 0$
implies an energy flux towards the SGS or subfilter scale.
5.2. Space-averaged uncertainty fluxes
The role of
$\varPi ^u$
in the evolution of uncertainty becomes evident by averaging the equations derived in the previous section. The equation for the space-averaged variance of the ensemble reads
\begin{equation} \begin{aligned} d_t {\langle {\mathcal{S}}\rangle } &= \langle \varPi ^u\rangle + \langle \mathcal{F}\rangle + \langle \mathcal{D} \rangle , \\ \end{aligned} \end{equation}
where the spatial fluxes vanish due to periodicity. This simplification also applies to domains with net zero flux at the boundaries or that decay towards infinity, or in homogeneous flows, in which the spatial average is exchanged by the standard unconditional average over base flows. Similarly, the equation for the space-averaged ensemble energy is
\begin{equation} \begin{aligned} d_t \langle {\mathcal{E}_{\{\boldsymbol{u}\}}}\rangle &= -\langle \varPi ^u\rangle + \langle \{f_i\}\{u_i\}\rangle - 2\nu \langle \{S_{\textit{ij}}\}\{S_{\textit{ij}}\}\rangle .\\ \end{aligned} \end{equation}
These equations are similar to the budget equations derived by Ge et al. (Reference Ge, Rolland and Vassilicos2023), with the particularity that uncertainty is not characterised here by the difference between a perturbed and a reference trajectory, but between trajectories and the ensemble-averaged flow field. Interestingly, the rate-of-strain tensor appears as a key element in both frameworks.
Now, let us assume that the MRS is in the inertial range, so that the large scales are free of uncertainty,
${\{f_i u_i\}}\approx {\{f_i\}\{u_i\}}$
, and the production of uncertainty due to the forcing is negligible compared with the dissipation,
$|\mathcal F|\ll |\mathcal D|$
. As shown in the previous section, the filtering effect of uncertainty growth implies that the squared strain of the ensemble-averaged field scales as
\begin{equation} \{S_{\textit{ij}}\}\{S_{\textit{ij}}\}\sim \ell ^{-4/3}, \end{equation}
while the ensemble-averaged squared strain scales in Kolmogorov units,
\begin{equation} \{S_{\textit{ij}}S_{\textit{ij}}\}\sim \eta ^{-4/3}, \end{equation}
so that the ratio of the two terms is
\begin{equation} \frac {\{S_{\textit{ij}}S_{\textit{ij}}\}}{\{S_{\textit{ij}}\}\{S_{\textit{ij}}\}}\sim \Big (\frac {\ell }{\eta }\Big )^{4/3}. \end{equation}
This implies that, for a sufficiently large
$\ell$
in the inertial range,
$\{S_{\textit{ij}}S_{\textit{ij}}\}\gg {\{S_{\textit{ij}}\}\{S_{\textit{ij}}\}}$
and
$\langle {\mathcal{D}}\rangle \approx -2\nu \langle {\{S_{\textit{ij}}S_{\textit{ij}}\}\rangle }$
, so that (5.11) simplifies to
\begin{equation} \begin{aligned} d_t {\langle {\mathcal{S}}\rangle } &= \langle \varPi ^u\rangle - 2\nu \langle \{S_{\textit{ij}}S_{\textit{ij}}\}\rangle . \\ \end{aligned} \end{equation}
Averaging this equation over base flows leads to
\begin{equation} d_t\langle \overline {\mathcal S}\rangle =\langle \overline {\varPi ^u}\rangle - \varepsilon , \end{equation}
where the relation
$\varepsilon =2\nu \langle \{\overline {S_{\textit{ij}}S_{\textit{ij}}}\}\rangle$
holds due to the equivalence between the average over base flows and the unconditional ensemble average over the attractor.
Applying the same procedure to the equation for the energy of the ensemble-averaged flow, we reach
\begin{equation} d_t\langle \overline {\mathcal E_{\{\boldsymbol{u}\}}}\rangle =-\langle \overline {\varPi ^u}\rangle + \langle \overline {\{f_i\}\{u_i\}}\rangle , \end{equation}
where we have assumed that the decay due to viscosity is negligible,
$\langle \varPi ^u\rangle \gg 2\nu \langle \{S_{\textit{ij}}\}\{S_{\textit{ij}}\}\rangle$
, because
$\varPi ^u\sim \varepsilon$
, as will be shown in the following. For
$\ell$
in the inertial range, the lack of uncertainty in the large scales implies that the energy injected by the forcing in the ensemble-averaged flow is equal to the average energy injection,
\begin{equation} \langle \overline {\{f_i\}\{u_i\}}\rangle \approx \langle \overline {\{f_iu_i\}}\rangle . \end{equation}
In addition, due to the conservative properties of the Navier–Stokes equations, the average injected energy is equal to the kinetic energy dissipation,
$\langle \overline {\{f_iu_i\}}\rangle =\varepsilon$
, leading to
\begin{equation} d_t\langle \overline {\mathcal E_{\{\boldsymbol{u}\}}}\rangle =-\langle \overline {\varPi ^u}\rangle + \varepsilon . \end{equation}
Figure 6. (a) Average interscale energy flux,
$\langle \overline {\varPi ^e} \rangle$
, and average uncertainty flux,
$\langle \overline {\varPi ^u}\rangle$
, as a function of the MRS for
${\textit{Re}}_\lambda =195$
and
${\textit{Re}}_\lambda =120$
. (b) Ratio of average uncertainty fluxes over energy fluxes as a function of the MRS. The horizontal dashed line marks
$\langle \overline {\varPi ^u}\rangle =1.44\langle \overline {\varPi ^e}\rangle$
. The markers in panels (a) and (b) correspond to the times shown in figure 1.
Equations (5.19) and (5.22) show how the energy of the ensemble-averaged flow and the variance are connected by a flux term,
$\varPi _u$
, which is the only term causing the increments of the average variance in the inertial range. These two equations are analogous to the budget equations of the classical energy cascade in statistically steady flows, which are typically derived in scale space (Domaradzki & Rogallo Reference Domaradzki and Rogallo1990; Zhou Reference Zhou1993; Domaradzki, Liu & Brachet Reference Domaradzki, Liu and Brachet1993; Mininni, Alexakis & Pouquet Reference Mininni, Alexakis and Pouquet2006; McComb et al. Reference McComb, Berera, Yoffe and Linkmann2015) or by averaging physical-space quantities (Cardesa et al. Reference Cardesa, Vela-Martín, Dong and Jiménez2015; Doan et al. Reference Doan, Swaminathan, Davidson and Tanahashi2018). The energy of the ensemble-averaged field and the variance are equivalent to the energy of the resolved scales and the SGS energy, respectively. There is, however, an important difference: the uncertainty cascade is unsteady,
$d_t\langle \overline {\mathcal{S}}\rangle \gt 0$
, as shown in the previous sections. First, this implies that the average uncertainty flux goes from the ensemble-average to the variance,
$\langle {\overline {\varPi ^u}}\rangle \gt 0$
. Here,
$\varPi ^u\gt 0$
describes the loss of energy in the ensemble-averaged field, and thus is interpreted as an inverse flux of uncertainty from the small to the large scales. Second, that this flux is larger on average than the average kinetic energy dissipation,
$\langle \overline {\varPi ^u}\rangle \gt \varepsilon$
, and, by energy conservation, than the kinetic energy fluxes,
$\langle \overline {\varPi ^u}\rangle \gt \langle \overline {\varPi ^e}\rangle$
.
In considering the ensemble-averaged flow as the optimal reconstruction or prediction,
$\varPi ^u$
can be interpreted as a sink of information, i.e. a term that depletes the energy that can be reconstructed from partial measurements of the past. The same term, viewed from the perspective of the variance, appears as an uncertainty source, increasing the fraction of energy in the unpredictable scales of the flow.
This picture is now tested quantitatively against the data. In figure 6(a), we show the average interscale energy flux and the average uncertainty flux as a function of the MRS, where the interscale energy fluxes have been calculated by filtering the base flows at scale
$\ell (t)$
with the Gaussian kernel defined in (3.2). The uncertainty fluxes are larger than the energy fluxes, as anticipated in the previous analysis, and both display an approximately constant flux in the inertial range. The collapse of the data for the two different Reynolds numbers is good, except at large scales, in which the scaling of the lower Reynolds number breaks due to limited scale separation. Figure 6(b) shows the ratio of uncertainty to kinetic energy fluxes, which seems to reach a plateau at inertial scales in which
$\langle \overline {\varPi ^u}\rangle \approx 1.4\langle \overline {\varPi ^e}\rangle$
. Due to the limited Reynolds number of the simulations, it is unclear whether this ratio is close to its asymptotic value in the large Reynolds number limit. In any case, the inertial scaling of the MRS and the average variance (see figure 3), which provide the length and velocity scales of
$\varPi ^u$
, suggest a relationship of the type
$\langle \overline {\varPi ^u}\rangle \approx C_u\langle \overline {\varPi ^e}\rangle$
where
$C_u\gt 1$
. Except for finite Reynolds number effects (Donzis & Sreenivasan Reference Donzis and Sreenivasan2010) and large-scale unsteadiness (Vassilicos Reference Vassilicos2015; Goto & Vassilicos Reference Goto and Vassilicos2016), the quantity
$C_u$
could represent a universal ratio of uncertainty to kinetic energy fluxes, similar to the Kolmogorov constant (Sreenivasan Reference Sreenivasan1995). Further research is necessary to elucidate this question.
Figure 7. Root-mean-squared magnitude of (a) the ensemble-averaged rate-of-strain tensor, and (b) the traceless part of uncertainty stress tensor
$\tau ^u_{\textit{ij}}$
, as a function of the MRS. Quantities are normalised with integrals units,
$U$
and
$L$
. The dashed lines in panels (a) and (b) are proportional to
$\ell ^{-2/3}$
and
$\ell ^{2/3}$
, respectively.
5.3. Scale locality of the uncertainty fluxes
This section focuses on the scale locality of the uncertainty cascade. It is shown that the same arguments used to claim the scale locality of the kinetic energy flux apply to the uncertainty flux. In particular, the locality conditions proposed by Eyink (Reference Eyink2005) are fulfilled by the ensemble-averaged rate-of-strain tensor and the uncertainty stress tensor, so that
\begin{equation} |\{S_{\textit{ij}}\}|\sim \ell ^{-q}, \end{equation}
and
\begin{equation} |\tau ^{u\dagger }_{\textit{ij}}|\sim \ell ^{q}, \end{equation}
respectively, where
$q\lt 1$
is a scaling exponent. Here
$\tau ^{u\dagger }_{\textit{ij}}=\tau ^{u}_{\textit{ij}} - (1/3)\tau ^{u}_{kk}$
is the traceless stress tensor, and the modulus of a tensor is defined as
\begin{equation} |\{S_{\textit{ij}}\}|=\sqrt {S_{\textit{ij}}S_{\textit{ij}}}. \end{equation}
Figure 7(a,b) shows that these conditions, known as the infrared and ultraviolet locality conditions, are fulfilled in the inertial range, with an exponent close to
$q=2/3$
, consistent with an inertial scaling. This shows that most of the strain acting to propagate uncertainty is located in scales immediately above
$\ell (t)$
, while the stresses, on which the strain acts, are in scales immediately below. This implies that, although the rate-of-strain tensor at all scales above the MRS is expected to interact with the stress tensor at all scales below, the main contribution to the fluxes comes from scales close to the MRS.
Figure 8. The p.d.f. of (a) the local uncertainty fluxes,
$\varPi ^u$
, and (b) the local energy fluxes
$\varPi ^e$
, for different
$\overline {\ell }$
as indicated in the legend. (c) Normalised variance and (d) kurtosis as a function of the MRS of (solid line) the uncertainty fluxes and (dashed line) the energy fluxes. Data from
${\textit{Re}}_\lambda =195$
.
The locality of the uncertainty cascade is key to justify its potential universal structure. Similar to the energy cascade and the universality of small-scale statistics (Nelkin Reference Nelkin1994; Sreenivasan & Antonia Reference Sreenivasan and Antonia1997; Schumacher et al. Reference Schumacher, Scheel, Krasnov, Donzis, Yakhot and Sreenivasan2014), the scale-locality of the uncertainty cascade implies that its structure does not depend on the uncertainty the initial conditions, provided that this uncertainty is initially restricted to the small scales. This is so because after a few propagation steps from small to large scales, the particular details of the initial probability distribution should be forgotten, at least for the low-order statistics of the distribution, in agreement with the statistical universality of error growth in turbulence (Aurell et al. Reference Aurell, Boffetta, Crisanti, Paladin and Vulpiani1996b , Reference Aurell, Boffetta, Crisanti, Paladin and Vulpiani1997; Boffetta et al. Reference Boffetta, Cencini, Falcioni and Vulpiani2002; Boffetta & Musacchio Reference Boffetta and Musacchio2017). This suggests that the results presented here may be viewed as representing an asymptotic, universal regime of uncertainty propagation in turbulence. This allows us to extended the framework presented here to problems in which the initial uncertainty is large compared with the Kolmogorov scale but small compared with the integral scale.
6. The structure of the uncertainty cascade
6.1. Statistics of the uncertainty fluxes
In this section, we analyse the structure of the local uncertainty fluxes in physical space,
$\varPi ^u(\boldsymbol{x},t)$
, with focus on a comparison with the kinetic energy fluxes of the energy cascade. It is important to remark that the following analysis is subject to a gauge ambiguity of the subgrid stresses, both in the energy and the uncertainty fluxes (Vela-Martín Reference Vela-Martín2022), and should be considered only in comparative terms.
First, let us focus on the probability distribution of the kinetic energy and uncertainty fluxes. In figure 8(a,b), we show the p.d.f. of the uncertainty flux and the interscale energy flux at the same scales. Both fluxes have similarly intermittent distributions with fat tails towards large positive values which increase with decreasing scale (Cerutti & Meneveau Reference Cerutti and Meneveau1998; Yasuda & Vassilicos Reference Yasuda and Vassilicos2018). The other side of the distributions shows that negative uncertainty fluxes are substantially less likely than negative energy fluxes. These inverse uncertainty fluxes are analogous to the energy backscatter of the cascade (Piomelli et al. Reference Piomelli, Cabot, Moin and Lee1991; Carati, Ghosal & Moin Reference Carati, Ghosal and Moin1995), and could be interpreted as regions where energy goes from the variance to the ensemble-averaged flow, describing a spontaneous reduction of uncertainty. These events have also been observed by Ge et al. (Reference Ge, Rolland and Vassilicos2023), although in their case they were highly more likely. However, it is difficult to conceive that the strongly chaotic nature of turbulence allows for the spontaneous reduction of uncertainty (Vela-Martín & Jiménez Reference Vela-Martín and Jiménez2021). In this sense, these events could be an artefact of the ambiguity that exists between space-local and scale-local fluxes (Vela-Martín Reference Vela-Martín2022). Uncertainty fluxes in the direction opposite to the average cascade appear to be particularly intense close to the dissipative scales, potentially indicating intermittency effects or an indirect role of viscosity in precluding the amplification of uncertainty, as shown by synchronisation experiments (Yoshida, Yamaguchi & Kaneda Reference Yoshida, Yamaguchi and Kaneda2005; Lalescu, Meneveau & Eyink Reference Lalescu, Meneveau and Eyink2013). Further research is necessary to clarify the physical meaning of these events. In any case, this analysis shows that the uncertainty cascade is somewhat more robust than the energy cascade, with a largely predominant flow of uncertainty from small to large scales.
To further compare both probability distributions, we have computed their centred second moment normalised by the mean squared
\begin{equation} \varsigma (\varPi )=\frac {\langle \overline {\varPi '^2}\rangle }{\langle \overline {\varPi }\rangle ^{2}}, \end{equation}
which is shown as a function of the MRS in figure 8(c). The normalised second moment of the kinetic energy and uncertainty fluxes have similar scaling in the inertial range, with fluctuations decreasing with increasing scale, but the latter is substantially smaller than the former. This difference is mainly due to their different average
$\langle \overline {\varPi ^u}\rangle \gt \langle \overline {\varPi ^e}\rangle$
. To show this, we have plotted in the same figure the centred second moment of the uncertainty fluxes normalised with the average kinetic energy fluxes,
\begin{equation} \varsigma ^{\star }(\varPi ^u)=\frac {\langle \overline {{\varPi ^u}'^2}\rangle }{\langle \overline {{\varPi ^e}}\rangle ^{2}}. \end{equation}
There is a very good collapse between
$\varsigma ^{\star }(\varPi ^u)$
and
$\varsigma (\varPi ^e)$
in the inertial range, indicating that the fluctuations of both quantities are quantitatively similar.
Let us focus now on the intermittency of the energy fluxes. Figure 8(d) shows the standardised fourth-order moment (kurtosis), defined as
\begin{equation} \kappa (\varPi )=\frac {\langle \overline {\varPi '^4}\rangle }{\langle \overline {\varPi '^2}\rangle ^2}, \end{equation}
as a function of the MRS. In the inertial range, the kurtosis of the uncertainty fluxes is slightly larger than that of the kinetic energy fluxes, but both have similar scaling. This, together the with the results above, indicate that the distributions of the uncertainty and kinetic energy fluxes are quite similar, and that the main difference between the two is the larger average of the former with respect to the latter. The intermittency of the uncertainty fluxes reported here is consistent with previous works, which report the intermittency of uncertainty production mechanisms (Ge et al. Reference Ge, Rolland and Vassilicos2023) and the anomalous scaling of the finite-size Lyapunov exponents (Boffetta & Musacchio Reference Boffetta and Musacchio2017).
Figure 9. Cuts of the (a,c,e) uncertainty fluxes,
$\varPi ^u$
, and the (b,d,f) interscale energy fluxes,
$\varPi ^e$
, for the same base flow at different MRS,
$\overline {\ell }=$
: (a,b)
$19\eta$
; (c,d)
$53\eta$
; (e,f)
$86\eta$
, which correspond to the same scales as in figure 1 (see markers for reference). The fluxes are normalised by their space average, and the colour bars are the same for
$\varPi ^u$
and
$\varPi ^e$
in each scale. The size of the plane is scaled in Kolmogorov units.
Figure 10. (a) Joint p.d.f. of the interscale and uncertainty fluxes at different
$\overline {\ell }$
. (b) Average interscale energy flux conditioned to the local uncertainty fluxes. The brackets
$\langle \boldsymbol{\cdot }\rangle _c$
denote the conditional average over space and ensembles, and the different markers and colours correspond to different MRSs as described in panel (a). The dotted line marks
$\langle \varPi ^e|\varPi ^u\rangle _c=0.53\varPi ^u$
.
Figure 11. (a,b) Joint p.d.f. of the (a) enstrophy of the ensemble-averaged flow and the uncertainty fluxes, and (b) the enstrophy of the filtered true flow and the interscale energy fluxes. (c,d) Joint p.d.f. of the (c) squared rate-of-strain of the ensemble-averaged flow and the uncertainty fluxes, and (d) the squared rate-of-strain of the filtered true flow and the interscale energy fluxes. In all figures, the markers are as shown in (a). All quantities are normalised by their average.
6.2. Structure of the uncertainty fluxes in physical space
We focus now on the spatial structure of the uncertainty fluxes and their correlation in space with the kinetic energy fluxes. Figure 9(a–f) shows cuts of the uncertainty fluxes and the energy fluxes in the same ensembles. By visually comparing the two fields, it is noticeable that intense energy transfer and uncertainty transfer events are located in similar areas of the flow, suggesting that both fields are partially overlapped.
To quantify the level of spatial overlapping of the two fields, we compute the joint probability distribution of
$\varPi ^u$
and
$\varPi ^e$
, which is shown in figure 10(a). The joint p.d.f. is similar across scales, and shows that
$\varPi ^u$
and
$\varPi ^e$
are largely correlated; the correlation coefficient between the two fields is
$\rho (\varPi ^u,\varPi ^e)\approx 0.6$
. This high correlation extends to large values of the fluxes, as shown by the average kinetic energy flux conditional to the uncertainty flux in figure 10(b). There is a proportionality relation of the form
$\langle \varPi ^e|\varPi ^u\rangle _c\approx 0.53\varPi ^u$
(
$\langle \boldsymbol{\cdot }\rangle _c$
is the conditional average in space and over ensembles), which holds at different scales and up to values of
$\varPi ^u$
10 times its average.
Finally, to examine the relation between uncertainty propagation and the velocity gradients, we analyse the joint p.d.f.s of the uncertainty flux with the enstrophy and the squared rate-of-strain tensor of the ensemble-averaged fields. For comparison, we also consider the joint p.d.f.s of the energy fluxes and the filtered gradients. First, the joint p.d.f.s of the fluxes with the enstrophy are shown in figure 11(a,b). The correlation of the kinetic energy and uncertainty fluxes with the enstrophy is small and decreases with increasing scale:
$\rho (\varPi {}^u,\varOmega _{\{\boldsymbol{u}\}})=0.25$
at
$\overline {\ell }=19\eta$
and
$\rho (\varPi ^u,\varOmega _{\{\boldsymbol{u}\}})=0.15$
at
$\overline {\ell }=120\eta$
. These correlations are quantitatively very similar to those of the kinetic energy with the enstrophy of the filtered fields,
$\rho (\varPi ^e,\varOmega _{\overline {\boldsymbol{u}}})$
. Figure 11(c,d) shows the joint p.d.f.s of the fluxes and the squared rate-of-strain tensor. There is a strong correlation of these two quantities, particularly at large scales (
$\overline {\ell }=120\eta$
), in which
$\rho (\varPi ^u,\varSigma _{\{\boldsymbol{u}\}})=0.91$
and
$\rho (\varPi ^e,\varSigma _{\overline {\boldsymbol{u}}})=0.82$
. Conversely to the enstrophy, the correlation with the strain increases with increasing scale. Interestingly, the uncertainty fluxes are more correlated to the strain than the energy fluxes. This suggests a particularly active role of straining motions in propagating uncertainty up in scale. This is not surprising given the direct contribution of the strain to the fluxes, as shown by the definition in (5.8), and supports previous investigations that connect the growth of disturbances with of the rate-of-strain tensor, both in the linear (Ruelle Reference Ruelle1982; Nikitin Reference Nikitin2018) and nonlinear (Ge et al. Reference Ge, Rolland and Vassilicos2023; Encinar & Jiménez Reference Encinar and Jiménez2023) regime.
In summary, the results presented in this section strongly suggest that the uncertainty cascade and the energy cascade share a similar statistical structure. Both show intermittency, with a large kurtosis that increases with decreasing scale, and both are intense in regions of intense strain. It is unclear whether other definitions of these fluxes using different gauge representations would lead to similar conclusions (Vela-Martín Reference Vela-Martín2022). In any case, it seems likely that both cascades can be described in a similar way. Moreover, this analysis suggests that the main difference between the uncertainty cascade and the energy cascade is that the average value of the fluxes is larger in the former. Other than this, the spatial fluctuations of the two cascades are similar in magnitude and are spatially correlated.
7. Conclusions
In this paper, the evolution of initially small perturbations in isotropic turbulence has been analysed using massive ensembles of direct numerical simulations, following the previous work by Vela-Martín (Reference Vela-Martín2024b ). This approach allows one to statistically converge, locally in space and scale, the average and variance of the ensemble. These are time-evolving fields that represent the optimal prediction in a least-squares sense and its uncertainty, respectively. The evolution of these two fields shows that uncertainty propagation acts on the ensemble-averaged field as a low-pass Gaussian filter with a characteristic scale that increases with time. This represents the progressive loss of predictability in the small scales due to chaos. The decay of the energy of the ensemble-averaged field is accompanied by an increment of the variance, which accounts for the growth of uncertainty from small to large scales.
By applying the ensemble average operator to the Navier–Stokes equations, it has been shown that the filtering of the ensemble-averaged flow field is described by a flux of energy from this field to the variance. This flux is formally equivalent to the SGS energy fluxes of the turbulence energy cascade, and is interpreted as an inverse uncertainty cascade due to the interaction between the predictable or correlated scales of the flow, and the unpredictable or uncorrelated scales. The comparative analysis presented here shows that, similar to the kinetic energy cascade, the uncertainty cascade has an inertial range with an uncertainty flux which is approximately constant in scale and larger than the average kinetic energy dissipation. Despite the limited scale separation of the flows analysed here, this suggests the existence of a constant inertial ratio of average uncertainty to energy fluxes in the large Reynolds number limit.
These results open avenues to describe uncertainty growth and predictability loss due to chaos as a cascade. To stress the validity of this interpretation, the locality of the uncertainty fluxes has been positively tested using the same arguments as in the kinetic energy cascade (Eyink Reference Eyink1995). Furthermore, the space-local uncertainty fluxes have been analysed and compared with the local kinetic energy fluxes at the same scales. The uncertainty fluxes spatially overlap the kinetic energy fluxes to a large extent, and display intermittency effects, with fat tails that describe high probability of large uncertainty transfer events towards the large scales. In addition, the spatial correlations of the uncertainty fluxes with the velocity gradients of the ensemble-averaged flow are very similar to the same correlations in the energy cascade, i.e. there is a strong correlation of both fluxes with the rate-of-strain tensor.
In summary, this paper has presented ample evidence that uncertainty growth can be described using the paradigm of the cascade. It seems that the rate-of-strain is a good candidate to explain much of this cascade, mainly because it appears in the equations of error propagation (Nikitin Reference Nikitin2018; Ge et al. Reference Ge, Rolland and Vassilicos2023; Encinar & Jiménez Reference Encinar and Jiménez2023) and because it has also been given a key role in the kinetic energy cascade (Tsinober Reference Tsinober2000; Carbone & Bragg Reference Carbone and Bragg2020; Johnson Reference Johnson2021). Given the similar conclusions drawn from studying perturbation and ensembles of them, it is interesting to test whether the mechanisms identified by Ge et al. (Reference Ge, Rolland and Vassilicos2023) and Ge, Rolland & Vassilicos (Reference Ge, Rolland and Vassilicos2025), particularly those related to the compression of perturbations by the rate-of-strain tensor, can be also identified using the ensemble approach. Other mechanisms, for instance those analogous to vortex stretching (Goto Reference Goto2008; Leung, Swaminathan & Davidson Reference Leung, Swaminathan and Davidson2012), cannot be discarded due to the non-local coupling of the velocity gradients (Encinar Reference Encinar2023). In this line, research on the mechanisms of uncertainty propagation could be extended to other definitions of uncertainty beyond the ensemble variance. A promising research direction is information theory (Lozano-Durán & Arranz Reference Lozano-Durán and Arranz2022; Arranz & Lozano-Durán Reference Arranz and Lozano-Durán2024), which uses full probability distributions to quantify uncertainty and information flows. An interesting question is whether it is possible to produce evolution equations, similar to those of the ensemble-averaged field and its variance, to describe uncertainty propagation in terms of entropies or mutual information.
These are future research avenues that could be exploited to construct uncertainty propagation models for turbulence. These models could encode the dynamics of perturbation growth to accurately assess the uncertainty of numerical forecasts without the need to resort to large ensembles of simulations. In the case of isotropic turbulence, the results presented here, and in Vela-Martín (Reference Vela-Martín2024b ), indicate that a simple filtering operation with a Gaussian kernel could be used as an efficient estimator of the ensemble-averaged field and its variance. The reasons why such a simple model works for isotropic turbulence are still unclear, but suggest that much of the complexity of uncertainty propagation can be reduced to a simple diffusion process (Johnson Reference Johnson2022). How to design similar estimators that can be applied to inhomogeneous or anisotropic flows is a topic for further research.
The massive ensembles used here to characterise uncertainty propagation cannot be practically implemented in many flows of interest, particularly those with very high Reynolds numbers, or with complex geometries, which are substantially more costly to simulate than homogeneous isotropic turbulence. An important question is to determine how many ensemble members are necessary to appropriately capture uncertainty propagation (Leutbecher Reference Leutbecher2019). The answer is probably more than one ensemble member but less than the 4000 simulations employed here to ensure statistical convergence. A combination of the methods developed to study finite-size perturbation growth, e.g. Boffetta & Musacchio (Reference Boffetta and Musacchio2017), Ge et al. (Reference Ge, Rolland and Vassilicos2023) and Encinar & Jiménez (Reference Encinar and Jiménez2023), together with information on the structure of uncertainty propagation in the large-ensemble limit, offer the potential to extend the study of uncertainty propagation to flows of practical interest using an affordable number of simulations.
An important application of the uncertainty cascade is the reconstruction of turbulent flows from sparse or incomplete observations (Fukami, Fukagata & Taira Reference Fukami, Fukagata and Taira2019; Bauweraerts & Meyers Reference Bauweraerts and Meyers2021; Kim et al. Reference Kim, Kim, Won and Lee2021; Arun et al. Reference Arun, Bae and McKeon2023). In this problem, unavoidable measurement uncertainty in the large scales should propagate to the small scales, leading to the progressive loss of the information provided by the reconstruction. This would set up a fundamental limit for small scale reconstruction, similar to the fundamental predictability limit of chaotic systems, which has yet to be identified. This problem could be modelled using an information source, e.g. a nudging (Di Leoni, Mazzino & Biferale Reference Di Leoni, Mazzino and Biferale2020; Zauner et al. Reference Zauner, Mons, Marquet and Leclaire2022) or an assimilation strategy (Yoshida et al. Reference Yoshida, Yamaguchi and Kaneda2005; Lalescu et al. Reference Lalescu, Meneveau and Eyink2013), to obtain a statistically steady uncertainty cascade. It is likely that the structure of this steady cascade is similar to the unsteady uncertainty cascade described here.
A similar problem appears in LES (Meneveau & Katz Reference Meneveau and Katz2000; Moser et al. Reference Moser, Haering and Yalla2021), in which the SGS are subject to small-scale uncertainty that propagates to large scales. This uncertainty is generally not explicitly considered in the design of SGS models, which relies to a large extent on filtered direct numerical simulations. Further research on the uncertainty cascade could be leveraged to construct uncertainty-robust SGS models, such as those proposed by Langford & Moser (Reference Langford and Moser1999,Reference Langford and Moser2004).
Acknowledgements
The author would like to thank Ryo Araki for valuable discussions.
Funding
The author is grateful for the use of computational resources from the ERC grant ERC-AdG-101018287. Funding for open access publication was provided by Universidad Carlos III de Madrid (Agreement CRUE-Madroño 2025).
Declaration of interests
The authors report no conflict of interest.
Open access
