This model is not limited to weakly interacting BECs; it also applies in the classical-field limit of highly occupied modes [23, 24].

It is worthwhile noting that, as in the recent experiments [14, 15], our Letter remains in a regime $ka\ll 1$, where $k$ is the highest momentum considered. The scattering cross-section is thus momentum independent and $g=4\pi {\hslash }^{2}a/m$.

Note that the hydrodynamic kinetic energy is not equal to the total kinetic energy, as the latter also contains a contribution from density modulations, $\propto |\nabla \sqrt{n(\mathbf{r})}{|}^2$, usually referred to as quantum pressure. Though this component has been discussed before (see, e.g., Refs. [27, 31]), its role in neither vortex nor wave turbulence is well understood, and we do not focus on this term here.

Note that even though the trap shape affects large-scale properties, the far-from-equilibrium state at smaller scales is largely insensitive to the details of the confining potential.

For convenience, we refer to the wavenumber $\mathbf{k}$ as momentum, even though formally the momentum is $\mathbf{p}=\hslash \mathbf{k}$.

In End Matter, we derive an analytic expression for $t_{\mathrm{D}}$ in terms of the cascade properties.

Note that [20] speculated that the steeper effective exponent could be due to the logarithmic correction. However, this Letter assumes that there is no BEC, but the presence of the BEC is crucial for understanding the effective injection scale $k_0\propto k_{\xi }$, hence [20] cannot make quantitative comparisons to the experiment [15].

A simple geometric picture is that particles hitting the box at non-normal incidence can remain trapped for some time even if their kinetic energy ${\hslash }^2{k}^2/(2m)>U_{\mathrm{D}}$.

Even though $\alpha \ne 1$ means that the energy-density flux is not $\dot{N}{U}_{\mathrm{D}}/V$, the fact that $\alpha$ has only a very weak dependence on $a$ ($\propto a^{0.08(1)}$) and is independent of $U_{\mathrm{s}}$ suggests that the determination of the energy flux via $\dot{N}{U}_{\mathrm{D}}/V$, as done in [14, 16], will only introduce a small, presumably system-specific, systematic error.

Note that the slow decrease of the plateau of ${\mathrm{\Pi }}_{ϵ}$ is due to a progressive depletion of the condensate, which results in a (slow) decrease of ${\epsilon }_{\mathrm{in}}$ [corresponding horizontal dashed lines in left panel of Fig. 2].

Note that it is only in the limit $k_{\mathrm{D}}\to \infty$ that one recovers the expected limits ${\mathrm{\Pi }}_{ϵ}\to \epsilon$ with ${\mathrm{\Pi }}_N\to 0$ for the direct energy cascade [52].

Note that ${\xi }_t$ is not directly present in the axes of Fig. 3, but is used in the calculation of $n_0=N_k{k}^3(k{\xi }_t{)}^{0.5}$.

The independence of our EOS with $k_{\mathrm{D}}$ presumably will not hold for large enough $k_{\mathrm{D}}$ (beyond what we studied), as it is reasonable to expect that as $k\to \infty$ (with $k<k_{\mathrm{D}}$), WWT results should be recovered.

For an extensive comparison of our EOS to the prediction of 4-wave theory (as well as how the respective $N_k$ can be deduced), see Sec. VII and Fig. S3 [38].

The factor of ${\xi }^{0.5}$ on the $y$ axis is a finite-$k_{\mathrm{D}}$ effect, but even in its absence our EOS would depend on $n$; within the GP model, only the WWT prediction of a power-law EOS with an exponent of $1/3$ is independent of $n$.

As shown in [38], $n_0$ and $\epsilon /n$ are roughly constant over time, resulting in the increase of $n_0/n$ and ${\tau }_t\epsilon /(n{\zeta }_t)$. The points in Fig. 3 correspond to shaking times $\approx 2t_{\mathrm{D}}$.

One might expect $\gamma (k)$ to be universal in the range $k_{\mathrm{F}}\ll k\ll k_{\mathrm{D}}$. We instead observe universal behavior for $k\xi \gtrsim 2.5$, where $n(k)$ is isotropic [19, 60] and $\gamma (k)$ is well defined. For $k\xi \lesssim 2.5$, our numerical resolution of $\gamma (k)$ is insufficient to make a definitive claim.

At lower momenta ($k\lesssim k_{\xi }$) the particles interact more strongly, and the 4-wave WWT description is not expected to work. In the presence of a condensate, the proper quasiparticles to consider for $k\xi \lesssim 1$ are the Bogoliubov phonons, for which the appropriate WWT description is a wave-kinetic equation with 3-wave interactions. In our simulations, a description in terms of an equilibrium-like condensate is meaningless since for our strong drives, the Bogoliubov approximation is invalid.

Interestingly, a recent experiment has observed that the momentum distribution of 2D shaken Bose gases also becomes isotropic at $k\approx k_{\xi }$ [60].

Experimentally, $k_{\mathrm{D}}$ is bounded, so the separation between the injection scale $k_0\propto k_{\xi }$ and $k_{\mathrm{D}}$—and hence the region where the 4-wave WWT theory is applicable—is reduced for stronger interactions.