Mechanism for Sequestering Magnetic Energy at Large Scales in Shear-Flow Turbulence

Turbulence is a fundamentally multiscale process in which nonlinear dynamical motions carry energy across scales johnson2020 . One of the most robust and extensively studied mechanisms for cross-scale energy transfer is straining by shear flow batchelor1950 ; batchelor1954 , long recognized for its role in cascades to small scale townsend1976 . This process is especially active in magnetohydrodynamic (MHD) turbulence maron2004 . The stretching, squeezing, and folding of magnetic field lines by shear flow is readily discernible in visualizations, and leads to, for example, the enhancement of small-scale resistive dissipation and the intensification of magnetic energy at small scales beyond the viscous cutoff in turbulence with magnetic Prandtl number greater than unity schekochihin2002 ; maron2004 . The ubiquity of turbulence and turbulent straining is at odds with the well-established observations of magnetic fields at large scales in the universe (stellar, galactic, and beyond) neronov2010 ; kulsrud1999 and part of a conundrum as to how magnetic fields exist at such large scales brandenburg2005 ; tobias2013 ; squire2015 .

In this Letter we report the first observations that magnetic energy is, to a considerable degree, counterintuitively sequestered at large scales in straining by a paradigmatic turbulent shear flow. The result applies to a Kelvin-Helmholtz- (KH-)unstable flow maintained by an external force with a uniform magnetic field that is initially flow-aligned and therefore optimally configured to promote small-scale generation.

The process that counteracts the small scale generation of magnetic energy, which is quantified in this Letter, can be traced to the nonlinear excitation of large-scale stable linear eigenmodes hatch2011 . The dispersion relation that is associated with the shear-flow instability present in this study has a stable (damped) root, which acts to return energy and momentum to the mean profile, and which is excited to a high level by the nonlinearity fraser2017 ; fraser2021 . The nonlinear excitation of stable modes and its effect on turbulence levels and transport have been studied previously, particularly in the context of fusion-relevant gyroradius-scale turbulence terry2018 ; whelan2018 ; fraser2018 ; pueschel2016 ; makwana2014 ; hatch2013 ; makwana2012 ; hatch2011 ; makwana2011 ; terry2006 . Recent studies have shown that stable modes are excited in macroscopic shear-flow driven turbulence also and transiently affect momentum transport fraser2017 ; fraser2021 . However, critical features of counteracting straining motion of the stable modes, essential for large-scale sequestering of magnetic energy, have only come to light in the study described here.

To whit, stable-mode effects have been linked to small-scale dissipation hatch2013 , suggesting that they affect turbulence only if there is an increase in entropy. However, we show here that the reversible process by which stable modes return energy to the mean flow, in opposition to extraction by the instability, is quite effective at blunting the energy cascade to small scales. Stable modes have been connected to the shear-layer contraction (i.e., build up of mean flow) fraser2021 , but only as a transient process. Here we demonstrate that for a driven flow profile, the return of fluctuation energy to the mean profile occurs continuously at a rate that nearly matches the rate at which unstable modes attempt to flatten the flow profile.

The principal result of this Letter is that, in driven shear-flow MHD turbulence, the nonlinearly excited large-scale stable eigenmodes efficiently strain the magnetic fields, in a way that counteracts the effect of unstable modes, and thus greatly weakens the magnetic energy cascade to small scales.

We study an incompressible, two-dimensional ($2$D) system with the initial fluid velocity given as $\mathbf{v}(x,z,t=0)=U_{\mathrm{ref}}(z)\hat{\mathbf{x}}=U_{0}\mathrm{tanh}(z/a)\hat{\mathbf{x}}$ and an initially uniform flow-aligned magnetic field as $\mathbf{B}(x,z,t=0)=B_{0}\hat{\mathbf{x}}$. The parameters $a,U_{0},\mathrm{and\ }B_{0}$ represent the half-width of the flow-shear, the maximum fluid velocity, and the magnetic field strength, respectively, which are used to non-dimensionalize all the variables henceforth. Consequently, time and energy (per unit mass) have units of $a/U_{0}$ and $U_{0}^{2}$. We describe the flow and the magnetic field by a stream function $\phi$ and a flux function $\psi$, where $\mathbf{v}=\hat{\mathbf{y}}\times\nabla\phi$ and $\mathbf{B}=\hat{\mathbf{y}}\times\nabla\psi$. We study their evolution using the momentum and induction equations of the standard MHD biskamp2003 ; remark_supplementary

where $\nabla^{2}\phi$ and $\nabla^{2}\psi$ are the vorticity and the current density; Poisson bracket $\{P,Q\}\equiv\partial_{x}P\cdot\partial_{z}Q-\partial_{x}Q\cdot\partial_{z}P$; the Alfv$\acute{\text{e}}$nic Mach number is $\mbox{M}_{\rm A}\propto U_{0}/B_{0}$; the fluid and magnetic Reynolds numbers are defined as $\mathrm{Re}=U_{0}a/\nu$ and $\mathrm{Rm}=U_{0}a/\eta$. In this study, we take $\mathrm{Re}=500$. Except for simulations where $\mathrm{Rm}=50$ and $\mathrm{Rm}=5~{}000$ (which will be stated explicitly), all others have $\mathrm{Rm}=500$. In Eq. (1), $f$ represents a body force that continuously regenerates the mean-flow shear. The expression $k_{x}\mathrm{=}0$ in it indicates that the force acts only on the ($x$-averaged) mean flow, where $k_{x}$ is the Fourier wavenumber. The forcing thus prevents gradual flattening of the mean flow as the instability extracts energy from it. A similar process occurs in astrophysical flows, where forces like gravitation replenish continually the shear profile, e.g., shear layers in accretion disks, stars, and planetary atmospheres. We represent this force by a Krook-like operator pueschel2014 ; marston2008 ; smith2021 , sometimes referred to as a profile relaxation term allawala2020 ,

where $D_{\mathrm{Krook}}$ represents the rate at which the mean flow is forced towards the reference flow profile. We choose $D_{\mathrm{Krook}}=2$. Detailed analyses of a complete scan of different parameters will be reported elsewhere. The force $F_{0}$ balances the viscous diffusion of the mean shear layer, $\mathrm{Re}^{-1}\nabla^{2}U_{\mathrm{ref}}(z)+F_{0}=0$, in Eq. (1), to ensure an initial equilibrium state.

We add small-amplitude perturbations to the flow and the magnetic field remark_supplementary ; fraser2021 and perform time-integration of the above set of equations using the pseudospectral code Dedalus burns2020 for long times ($t>1000$; the $\rm{e}$-folding growth time of the instability is $\approx 6$). Apart from the simulation of $\mathrm{Pm}=10$ that uses spectral resolution of $4096$ Fourier (Chebyshev) modes along the $x$-($z$-)axis for the box size of $(L_{x},L_{z})=(6\pi,8\pi)$, all other simulations are performed at $2048\times 2048$ spectral resolution for the box size of $(L_{x},L_{z})=(10\pi,20\pi)$, while convergence checks utilize resolutions as high as $8192\times 8192$, with no substantial impact on dissipation rates and spectral energy densities. (All simulations additionally use $3/2$ dealiasing rule burns2020 ).

We also perform a corresponding eigenvalue calculation of the non-dissipative equations, derived from Eqs. (1) & (1b), by linearizing about the initial flow and magnetic field profiles fraser2021 . A priori, one does not know whether such an eigenmode basis is useful to understand turbulent features. Nevertheless, since the eigenmodes of this linear operator form a complete (albeit non-orthogonal) basis remark_B ; bender2019 , we can expand the time-evolving state space of nonlinear simulations in this basis after determining amplitude $\beta_{j}(k_{x},t)$ of each eigenmode $j$ at each $k_{x}$ as a function of time, $\widetilde{\chi}_{\mathrm{turbulent}}(x,z,t)=\sum_{k_{x}:k_{x}\neq 0}e^{ik_{x}x}\left[\sum_{j}\beta_{j}(k_{x},t)\chi_{j}(k_{x},z)\right],$ where $\widetilde{\chi}_{\mathrm{turbulent}}=(\widetilde{\phi}_{\mathrm{turbulent}},\widetilde{\psi}_{\mathrm{turbulent}})$ captures the turbulent state space and $\chi_{j}$ is the $j^{\mathrm{th}}$ eigenmode fraser2021 ; tripathi2022_sherwood ; remark_energynorm ; remark_ECS_PCA . This procedure is illustrated in Fig. 1(a) where the eigenmodes of the flow fluctuations are shown. (The magnetic fluctuations have eigenmodes fraser2021 and parity ishizawa2019 ; sato2017 related to those of the flow fluctuations.)

At each wavenumber $k_{x}$, linearly unstable to perturbations (i.e., $0<k_{x}<1$ for the present shear-flow instability), there is an unstable eigenmode $\phi_{1}$ and its conjugate stable eigenmode $\phi_{2}$, with eigenvalues complex conjugate to each other chandrashekhar1961 ; fraser2017 . All the remaining eigenmodes have purely real eigenfrequencies tripathi2022_sherwood , belong to the continuum case1960 , and have narrow structures along the $z$-axis.

Tracking the time evolution of each eigenmode amplitude $\beta_{j}(k_{x},t)$ in the nonlinear simulations, the initially exponentially decaying mode has been found to be nonlinearly excited to almost the same level as the unstable mode fraser2021 . The driving mechanism for this excitation is, at early times, the nonlinear interaction between the unstable modes, see, e.g., Refs. terry2006 ; makwana2011 ; fraser2021 . In the fully nonlinear phase, however, all modes that have large energy interact dominantly.

For a given KH-unstable wavenumber, it is found, in Fig. 1(b), that a combination of two eigenmodes from the above basis reconstructs the large-scale turbulent flow at that wavenumber. Figure 1(b) shows a reconstruction of the turbulent flow features using unstable eigenmodes at each unstable wavenumber (left column) and with their conjugate stable eigenmodes added (right column). Simply adding one additional mode drastically decreases the difference between the true stream function and its reduced representation compared to what is achieved with a reduced representation based on the unstable modes alone (such as those predicted by quasilinear models) remark_field_reconstruction .

Note that the stable and unstable eigenmodes of the KH linear operator, employed above, are related to one another by a space-time-reversal symmetry operation bender2019 ; fu2020 . Thus the unstable and stable eigenmodes transfer energy in opposite directions—the former from the mean flow to the fluctuations, while the latter in the reverse direction. The energy that is available to cascade to small scales thus depends on this competition. The magnitude of the energy transfer rates between the background shear-flow $U_{\mathrm{ref}}$ and the fluctuation scale $k_{x}\neq 0$ by the unstable and stable eigenmodes can be derived from the MHD equations remark_supplementary : $Q_{1}(k_{x})=\Gamma(k_{x})|\beta_{1}(k_{x})|^{2}$ and $Q_{2}(k_{x})=\Gamma(k_{x})|\beta_{2}(k_{x})|^{2}$ where $\Gamma(k_{x})$ is the linear growth rate of the unstable mode at that wavenumber $k_{x}$.

Shown in Fig. 2 are the time traces of $Q_{1}$ and $Q_{2}$ for different $\mathrm{Pm}$ (resistivity) over a challenging two orders of magnitude. Summation over the KH-unstable wavenumbers, $0<k_{x}<1$, is used in determining $Q_{1}$ and $Q_{2}$. Within each subplot, it is observed that $Q_{1}$ and $Q_{2}$ are nearly in equipartition in the nonlinear phase.

Stable-mode-excitation process may also be impacted by magnetic fields. To analyze such, we first observe the visco-resistive dissipation with varying field strengths in Figs. 3(a) and (b). The small-scale dissipation—a proxy of small-scale cascade—enhances with stronger fields, despite the fields reducing the linear growth rate of the instability, which sometimes is argued to lead to subdued turbulent fluctuations as the energy extraction rate by the unstable mode becomes lower.

To understand such enhanced dissipation rate, it may be instructive to compute the energy transfer rates between the background flow and the fluctuations due to both the unstable and stable eigenmodes. As seen in Fig. 3(c), the (time-averaged) energy injection rate by the instability, $Q_{1}$, decreases with stronger magnetic fields, as anticipated. However, the energy return rate by the stable modes, $Q_{2}$, is impacted more by the stronger fields, allowing larger energy cascade to small scales ($\propto Q_{1}-Q_{2}$) that manifests as enhanced dissipation rate at such scales remark_instantaneousmean . Thus a stronger magnetic field, until it nearly stabilizes the system, allows more net energy extraction from the free-energy source in the nonlinear state. This questions models that assume the linear growth rate of an instability as a surrogate for the cascade rate or a proxy for nonlinear time, see, e.g., fuller2019 ; pessah2006 ; garaud2018 ; goodman1994 . Predictions based upon conventional wisdom of instability-saturation are thus challenged.

The ratio $Q_{2}/Q_{1}$ of the stable-eigenmode efficiency to the unstable-eigenmode efficiency in transferring energy between the fluctuations and the background flow is more than $60\%$ even for the strongest magnetic fields in Fig. 3(c), see the red dashed line. Such a strong Lorentz back-reaction is most potent in disrupting the vortex mak2017 and creating more small scales, leading to more reduced $Q_{2}$ than the corresponding reduction in $Q_{1}$.

We now show visualizations of highly effective role of stable modes in structuring the magnetic fields. After carrying out an initial-value simulation up to a time when the dynamics are nonlinear and turbulent, we perform two distinct simulation continuations—one unchanged with the stable modes untouched, and another with the stable modes removed. For the latter, we set only the large-scale ($0<|k_{x}|<1$) stable-eigenmode amplitudes to zero at the instant when the simulation is restarted. The differences immediately afterward are significant: Removing stable eigenmodes rapidly enhances dissipation rates, see Figs. 4(a) and (b); (Multimedia view). To understand this rise in the dissipation, we plot the squared current density $|\nabla^{2}\psi|^{2}$ at small scales ($|k_{x}|>1$) in Figs. 4(c)–(f). Observe in Figs. 4(d) and (e) that the magnetic fields get rapidly distorted by the unstable modes of the flow, creating spirals as the counter-straining motion of the stable modes is removed. Note that the (un)stable modes have the largest strain around $(x,z)\approx(10,0)$ for the shown time, as can also be seen in the figure in the supplementary material. The stable modes here are deleted only at the instant when the simulation is paused; they are quickly nonlinearly driven back. In another simulation where the stable modes are untouched (retained), the spirals do not appear, however, and the time evolution immediately afterward remains almost identical to the initial stage [cf. Figs. 4(f) and Fig. 4(c)]. This exercise of instantaneous deletion of stable modes exposes the true nonlinear impact of the stable modes on magnetic field evolution, which remained hidden in Fig. 1(a).

This numerical experiment illustrates that the energy injected by the unstable modes into the large-scale fluctuations ($|k_{x}|<1$) via mean-fluctuation interactions, in the absence of stable modes, cascades in its entirety to smaller scales nonlinearly where it dissipates visco-resistively. This is a statement of the classical small-scale cascade by the advective nonlinearity, $(\mathbf{v}\cdot\bm{\nabla})\psi$, as commonly understood. When $\mathbf{v}$ is composed of comparable unstable and stable mode-amplitudes, though, the cascade is weakened by a factor of up to ten [Figs. 2, 3(c), and 4(b)]. Such a finding informs and can be used to improve reduced models of geo- and astrophysical instability-driven turbulence vandine2021 ; fuller2019 , where all the energy injected by the instability is assumed to pass through an inertial cascade, e.g., Refs. fuller2019 ; pessah2006 ; garaud2018 ; goodman1994 .

It may be noted that, although $2$D and $3$D turbulence have different conserved quantities and corresponding cascade processes, the stable modes of this study are the conjugate roots of the KH-instability, which exist in both $2$D and $3$D chandrashekhar1961 . Our preliminary studies of $3$D turbulence have indeed shown that virtually identical stable mode excitation and its impact in turbulence occur in $3$D as well. Details of such an investigation will be left for a future publication tripathi2022_testingQL .

This Letter demonstrates, for the first time, that stable modes via their counteracting straining motion of the flow, sequester magnetic energy at large scales, despite the conventional association of shear-flow with intensification and generation of small scales in MHD turbulence. Our analyses show that the stable modes greatly weaken the magnetic energy cascade to small scales and hence can thwart the magnetic-energy accumulation at sub-viscous scales in the high-magnetic-Prandtl-number regime, typically encountered in diffuse astrophysical plasmas; thus possibly allowing operation of a dynamo with magnetic fields concentrated at large scales, in accord with astrophysical observations schekochihin2002 ; tobias2013 . This work opens a new direction towards such studies.