1. Introduction
In the absence of any initial velocity field and without any type of forcing, an initially random magnetic field can only decay. This decay can be sped up by turbulent gas motions driven through the Lorentz force that is being exerted by the magnetic field itself. The decay of such a random field obeys power law behaviour where the magnetic energy density
$\mathcal{E}_{\mathrm{M}}$
decays with time
$t$
as
$\mathcal{E}_{\mathrm{M}}(t)\propto t^{-p}$
, and the magnetic correlation length
$\xi _{\mathrm{M}}$
increases as
$\xi _{\mathrm{M}}\propto t^q$
. For a helical magnetic field, we have
$p=q=2/3$
(Hatori Reference Hatori1984; Biskamp & Müller Reference Biskamp and Müller1999), while for a non-helical magnetic field, we have
$p=10/9$
and
$q=4/9$
(Hosking & Schekochihin Reference Hosking and Schekochihin2021; Zhou, Sharma & Brandenburg Reference Zhou, Sharma and Brandenburg2022). Such a decay has been seen in many hydromagnetic numerical simulations (Brandenburg, Kahniashvili & Tevzadze Reference Brandenburg, Kahniashvili and Tevzadze2015; Hosking & Schekochihin Reference Hosking and Schekochihin2021; Armua, Berera & Calderón-Figueroa Reference Armua, Berera and Calderón-Figueroa2023; Brandenburg et al. Reference Brandenburg, Sharma and Vachaspati2023), but not yet in plasma experiments. With the advance of high-powered lasers, it is already possible to amplify magnetic fields in the laboratory (Tzeferacos et al. Reference Tzeferacos2018) and similar advances may also allow us to achieve sufficient scale separation to perform meaningful inverse cascade experiments. However, such magnetic fields may be strongly anisotropic, so the question arises to what extent this affects the otherwise familiar decay dynamics.
Our goal here is to study the decay of an array of magnetic flux tubes with an electric current that is aligned with the magnetic field (Jiang, Pukhov & Zhou Reference Jiang, Pukhov and Zhou2021). Such a field is indeed highly anisotropic such that the correlation length in the direction along the tubes is much larger than that perpendicular to it. A simple numerical realisation of such a magnetic field is what is called the Roberts field I, which is more commonly also known as Roberts flow I. It is one of four flow fields studied by Roberts (Reference Roberts1972) in the context of dynamo theory. The field is fully helical, but with a slight modification, it can become a pointwise non-helical field, which is then called the Roberts field II. Both fields are here of interest. They are defined in § 2, along with a proper measure of anisotropy, the relevant evolution equations, and relevant input and output parameters. In § 3, we present numerical results for both flows using different magnetic diffusivities and scale separation ratios. Inverse cascading during the turbulent decay of helical and non-helical magnetic fields has applications to primordial magnetic fields in the radiation dominated era of the Universe, which are discussed in § 4. We conclude in § 5.
2. Our model
2.1. Roberts fields
To fix our geometry, we assume magnetic flux tubes to extend in the
$z$
-direction and being perpendicular to the
$xy$
-plane. Such a field can be realised by the so-called Roberts field I, i.e. the magnetic field
$\boldsymbol{B}$
is given by
\begin{equation} {\boldsymbol{B}}={\boldsymbol{B}}_{\mathrm{I}}\equiv {\boldsymbol{\nabla }}\times \phi \hat {{\boldsymbol{z}}} {}+\sqrt {2}k_0\phi \hat {{\boldsymbol{z}}} {}, \quad \textrm{where}\; \phi =k_0^{-1} B_0 \sin k_0 x\sin k_0 y \end{equation}
is an
$xy$
periodic field. Such a magnetic field has a component in the
$z$
-direction, but no variation along that direction, so it is highly anisotropic. This may change with time as the magnetic field undergoes a turbulent decay. The Roberts field I is maximally helical with
$\boldsymbol{A}\boldsymbol{\cdot} {\boldsymbol{B}}=\sqrt {2} k_0^{-1} B_0^2 (\sin ^2 k_0x+\sin ^2 k_0y)$
, so
$\langle \boldsymbol{A}\boldsymbol{\cdot} {\boldsymbol{B}}\rangle =\sqrt {2} k_0^{-1} B_0^2$
. Here,
$\boldsymbol{A}$
is the magnetic vector potential and
${\boldsymbol{B}}={\boldsymbol{\nabla }}\times \boldsymbol{A}$
. The Roberts field II, by contrast, is given by
\begin{equation} {\boldsymbol{B}}={\boldsymbol{B}}_{\mathrm{II}}\equiv {\boldsymbol{\nabla }}\times \phi \hat {{\boldsymbol{z}}} {}+k_{\textrm{f}}\tilde {\phi }\hat {{\boldsymbol{z}}} {}, \quad \textrm{where}\ \tilde {\phi }=k_0^{-1} B_0 \cos k_0 x\cos k_0 y, \end{equation}
where
$\tilde {\phi }$
is
$90^\circ$
phase shifted in the
$x$
- and
$y$
-directions relative to
$\phi (x,y)$
, and
$k_{\textrm{f}}=\sqrt {2}k_0$
is the eigenvalue of the curl operator for field I, i.e.
${\boldsymbol{\nabla }}\times {\boldsymbol{B}}_{\mathrm{I}}=k_{\textrm{f}}{\boldsymbol{B}}_{\mathrm{I}}$
, so
${\boldsymbol{B}}_{\mathrm{I}}\boldsymbol{\cdot} {\boldsymbol{\nabla }}\times {\boldsymbol{B}}_{\mathrm{I}}=k_{\textrm{f}}{\boldsymbol{B}}_{\mathrm{I}}^2$
, while
${\boldsymbol{B}}_{\mathrm{II}}\boldsymbol{\cdot} {\boldsymbol{\nabla }}\times {\boldsymbol{B}}_{\mathrm{II}}=0$
pointwise. In the Coulomb gauge, we have, for field II,
${\boldsymbol{B}}_{\mathrm{II}}={\boldsymbol{\nabla }}\times \boldsymbol{A}_{\mathrm{II}}$
, where
\begin{equation} \boldsymbol{A}_{\mathrm{II}}=k_{\textrm{f}}^{-1}\big({\boldsymbol{\nabla }}\times \tilde {\phi }\hat {{\boldsymbol{z}}} {}+k_{\textrm{f}}\phi \hat {{\boldsymbol{z}}} {}\big), \end{equation}
and therefore also
$\boldsymbol{A}_{\mathrm{II}}\boldsymbol{\cdot} {\boldsymbol{B}}_{\mathrm{II}}=0$
. Thus, for field II, not just the current helicity density vanishes pointwise, but in the Coulomb gauge, also the magnetic helicity density vanishes pointwise. Both for fields I and II, we have
$\langle {\boldsymbol{B}}^2\rangle =2B_0^2$
.
2.2. Quantifying the emerging anisotropy
To quantify the degree of anisotropy, we must separate the derivatives of the magnetic field along the
$z$
-direction (
${\boldsymbol{\nabla }}_\|$
) from those perpendicular to it (
${\boldsymbol{\nabla }}_\perp$
), so
${\boldsymbol{\nabla }}={\boldsymbol{\nabla }}_\|+{\boldsymbol{\nabla }}_\perp$
. We also decompose the magnetic field analogously, i.e.
${\boldsymbol{B}}={\boldsymbol{B}}_\|+{\boldsymbol{B}}_\perp$
. The mean current density can be decomposed similarly, i.e.
$\boldsymbol{J}=\boldsymbol{J}_\|+\boldsymbol{J}_\perp$
, but this decomposition mixes the underlying derivatives. We see this by computing
$\boldsymbol{J}\equiv {\boldsymbol{\nabla }}\times {\boldsymbol{B}}$
(where the permeability has been set to unity). Using this decomposition, we find
\begin{equation} \boldsymbol{J}={\boldsymbol{\nabla }}_\|\times {\boldsymbol{B}}_\perp +{\boldsymbol{\nabla }}_\perp \times {\boldsymbol{B}}_\|+{\boldsymbol{\nabla }}_\perp \times {\boldsymbol{B}}_\perp , \end{equation}
noting that
${\boldsymbol{\nabla }}_\|\times {\boldsymbol{B}}_\|=0$
. The term of interest for characterising the emergent isotropisation is the first one,
${\boldsymbol{\nabla }}_\|\times {\boldsymbol{B}}_\perp$
, because it involves only parallel derivatives (
$z$
-derivatives), which vanish initially. We monitor the ratio of its mean squared value to
$\langle \boldsymbol{J}^2\rangle$
.
The last term in (2.4) is just
$\boldsymbol{J}_\|={\boldsymbol{\nabla }}_\perp \times {\boldsymbol{B}}_\perp$
, but the first and second terms cannot simply be expressed in terms of
$\boldsymbol{J}_\perp$
, although
${\boldsymbol{\nabla }}_\|\times {\boldsymbol{B}}_\perp$
would be
$\boldsymbol{J}_\perp$
if the magnetic field only had a component in the plane and
${\boldsymbol{\nabla }}_\perp \times {\boldsymbol{B}}_\|$
would be
$\boldsymbol{J}_\perp$
if the magnetic field only had a component out of the plane. We therefore denote those two contributions in what follows by
$\boldsymbol{J}_{\perp \perp }$
and
$\boldsymbol{J}_{\perp \|}$
, respectively, so that
$\boldsymbol{J}_{\perp \perp }+\boldsymbol{J}_{\perp \|}=\boldsymbol{J}_\perp$
.
Thus, with the abovementioned motivation, to monitor the emergent isotropisation, we determine
$\langle \boldsymbol{J}_{\perp \perp }^2\rangle /\langle \boldsymbol{J}^2\rangle$
. For isotropic turbulence, we find that this ratio is approximately
$4/15\approx 0.27$
, and this is also true for
$\langle \boldsymbol{J}_{\perp \|}^2\rangle /\langle \boldsymbol{J}^2\rangle$
; see Appendix A for an empirical demonstration. In the expression for
$\langle \boldsymbol{J}^2\rangle$
, there is also a mixed term,
$\boldsymbol{J}_{\perp {m}}^2=-2\langle B_{x,z}B_{z,x}+B_{y,z}B_{z,y}\rangle$
, which turns out to be positive in practice. Here, commas denote partial differentiation. Thus, we have
\begin{equation} \langle \boldsymbol{J}^2\rangle = \big\langle \boldsymbol{J}_{\perp \perp }^2\big\rangle + \big\langle \boldsymbol{J}_{\perp \|}^2\big\rangle + \big\langle \boldsymbol{J}_{\perp \textrm{m}}^2\big\rangle + \big\langle \boldsymbol{J}_{\|}^2\big\rangle . \end{equation}
In the isotropic case, we find
$\langle \boldsymbol{J}_{\|}^2\rangle /\langle \boldsymbol{J}^2\rangle =1/3$
and for the mixed term, we then have
$\langle \boldsymbol{J}_{\perp {\mathrm{m}}}^2\rangle /\langle \boldsymbol{J}^2\rangle =2/15\approx 0.13$
.
2.3. Evolution equations
To study the decay of the magnetic field, we solve the evolution equations of magnetohydrodynamics (MHD) for an isotropic compressible gas with constant sound speed
$c_{\mathrm{s}}$
, so the gas density
$\rho$
is proportional to the pressure
$p=\rho c_s^2$
. In that case,
$\ln \rho$
and the velocity
$\boldsymbol{u}$
obey
\begin{equation} \frac {{\textrm{D} {}}\ln \rho }{{\textrm{D} {}} t}=-{\boldsymbol{\nabla }}\boldsymbol{\cdot} {\boldsymbol{u}}, \end{equation}
\begin{equation} \frac {{\textrm{D} {}}{\boldsymbol{u}}}{{\textrm{D} {}} t}=-c_{\mathrm{s}}^2{\boldsymbol{\nabla }}\ln \rho +\frac {1}{\rho }\left [\boldsymbol{J}\times {\boldsymbol{B}}+{\boldsymbol{\nabla }}\boldsymbol{\cdot} (2\rho \nu {\boldsymbol{\sf S}} {})\right ],\\[-12pt] \end{equation}
where
${\textrm{D} {}}/{\textrm{D} {}} t=\partial /\partial t+{\boldsymbol{u}}\boldsymbol{\cdot} {\boldsymbol{\nabla }}$
is the advective derivative,
$\nu$
is the kinematic viscosity and
${\boldsymbol{\sf S}} {}$
is the rate-of-strain tensor with components
$\mathsf{S}_{ij}=(u_{i,j}+u_{j,i})/2-\delta _{ij}{\boldsymbol{\nabla }}\boldsymbol{\cdot} {\boldsymbol{u}}/3$
. To ensure that the condition
${\boldsymbol{\nabla }}\boldsymbol{\cdot} {\boldsymbol{B}}=0$
is obeyed at all times, we also solve the uncurled induction equation for
$\boldsymbol{A}$
, i.e.
\begin{equation} \frac {\partial \boldsymbol{A}}{\partial t}={\boldsymbol{u}}\times {\boldsymbol{B}}-\eta \boldsymbol{J}. \end{equation}
As before, the permeability is set to unity, so
$\boldsymbol{J}={\boldsymbol{\nabla }}\times {\boldsymbol{B}}$
is the current density.
We use the Pencil Code (Pencil Code Collaboration et al. 2021), which is well suited for our MHD simulations. It uses sixth-order accurate spatial discretisations and a third-order time-stepping scheme. We adopt periodic boundary conditions in all three directions, so the mass is conserved and the mean density
$\langle \rho \rangle \equiv \rho _0$
is constant. The size of the domain is
$L_\perp \times L_\perp \times L_\|$
and the lowest wavenumber in the plane is
$k_1=2\pi /L_\perp$
. By default, we choose
$\rho _0=k_1={c_{\mathrm{s}}}=\mu _0=1$
, which fixes all dimensions in the code.
2.4. Input and output parameters
In the following, we study cases with different values of
$k_0$
. We specify the amplitude of the vector potential to be
$A_0=0.02$
for most of the runs with Roberts field I and
$A_0=0.05$
for Roberts field II. We use
$k_0=16$
, so
$B_0=k_0A_0=0.32$
for field I and
$0.8$
for field II. For other values of
$k_0$
, we adjust
$A_0$
such that
$B_0$
is unchanged in all cases. This implies
$\langle {\boldsymbol{B}}^2\rangle =2B_0^2=0.2$
and
$1.28$
, and therefore
$B_{\mathrm{rms}}=0.45$
and
$1.13$
, respectively. The initial values of the Alfvén speed,
$v_{\mathrm{A0}}=B_{\mathrm{rms}}/\sqrt {\mu _0\rho _0}$
, are therefore transonic. We often give the time in code units,
$({c_{\mathrm{s}}} k_1)^{-1}$
, but sometimes we also give it in units of
$(v_{\mathrm{A0}} k_0)^{-1}$
, which is physically more meaningful. However, we must remember that the actual magnetic field and therefore the actual Alfvén speed are of course decaying.
In addition to the Roberts field, we add to the initial condition Gaussian-distributed noise of a relative amplitude of
$10^{-6}$
. This allows us to study the stability of the field to small perturbations. To measure the growth rate, we compute the semilogarithmic derivative of
$\langle \boldsymbol{J}_{\perp \perp }^2\rangle /\langle \boldsymbol{J}^2\rangle$
for a suitable time interval.
The number of eddies in the plane is characterised by the ratio
$k_0/k_1$
. The aspect ratio of the domain is quantified by
$L_\|/L_\perp$
. The electric conductivity is quantified by the Lundquist number
${\textrm{Lu}}=v_{\mathrm{A0}}/\eta k_0$
, and the kinematic viscosity is related to
$\eta$
through the magnetic Prandtl number,
${\textrm{Pr}_{\mathrm{M}}}=\nu /\eta$
. In all our cases, we take
${\textrm{Pr}_{\mathrm{M}}}=5$
. This is an arbitrary choice, just like
${\textrm{Pr}_{\mathrm{M}}}=1$
would be arbitrary. The value of
$\textrm{Pr}_{\mathrm{M}}$
affects the ratio of kinetic to magnetic energy dissipation (Brandenburg Reference Brandenburg2014; Brandenburg & Rempel Reference Brandenburg and Rempel2019). While this topic is interesting and important, it is not the focus of our present study. Laboratory plasmas tend to have large values of
$\textrm{Pr}_{\mathrm{M}}$
, so the choice
${\textrm{Pr}_{\mathrm{M}}}=5$
instead of unity is at least qualitatively appropriate. Much larger values of
$\textrm{Pr}_{\mathrm{M}}$
would become computationally prohibitive. Furthermore, the choice
${\textrm{Pr}_{\mathrm{M}}}=1$
can lead to exceptional behaviour, particularly when the cross-helicity is finite; see figure 1 of Rädler & Brandenburg (Reference Rädler and Brandenburg2010).
Evolution of
$\langle \boldsymbol{J}_{\perp \perp }^2\rangle /\langle \boldsymbol{J}^2\rangle$
for (a) Roberts field I with
$k_0=4$
(blue),
$8$
(green),
$16$
(orange),
$32$
(red) and
$64$
(black dashed), and for (b) Roberts field II with
$k_0=2$
(black),
$4$
(blue),
$8$
(green),
$16$
(orange),
$32$
(red) and
$64$
(black dashed). The short thick line on the upper right indicates the value of 4/15, which is reached only at much later times outside this plot. The insets demonstrate that
$\langle \boldsymbol{J}_{\perp \perp }^2\rangle /\langle \boldsymbol{J}^2\rangle \to 4/15$
much later.
Important output parameters are the growth rate
$\lambda ={\textrm{d} {}}\ln (\langle \boldsymbol{J}_{\perp \perp }^2\rangle /\langle \boldsymbol{J}^2\rangle )/{\textrm{d} {}} t$
, evaluated in the regime where it is non-vanishing and approximately constant. It is made non-dimensional through the combination
$\lambda /v_{\mathrm{A0}} k_0$
. We also present magnetic energy and magnetic helicity variance spectra,
$\textrm{Sp}({\boldsymbol{B}})$
and
$\textrm{Sp}(h)$
, respectively. These spectra depend on
$k$
and
$t$
, so we denote the spectra sometimes also as
$\textrm{Sp}({\boldsymbol{B}};k,t)$
and
$\textrm{Sp}(h;k,t)$
, respectively.
Since
$\rho \approx \rho _0=1$
, the value of
$B_{\mathrm{rms}}$
is also equal to the instantaneous Alfvén speed,
$v_{\mathrm{A}}$
, and its square is the mean magnetic energy density,
$\mathcal{E}_{\mathrm{M}}=\langle {\boldsymbol{B}}^2\rangle /2$
. The latter can also be computed from the magnetic energy spectrum
${E_{\mathrm{M}}}(k,t)=\textrm{Sp}({\boldsymbol{B}})$
through
$\mathcal{E}_{\mathrm{M}}=\int\!\!{E_{\textrm{M}}}(k,t)\,{\textrm{d} {}} k$
. The integral scale of the magnetic field is given by
\begin{equation} \xi _{\mathrm{M}}(t)=\int k^{-1}{E_{\mathrm{M}}}(k,t)\,{\textrm{d} {}} k/\mathcal{E}_{\mathrm{M}}. \end{equation}
It is of interest to compare its evolution with the magnetic Taylor microscale,
$\xi _{\mathrm{T}}=B_{\mathrm{rms}}/J_{\mathrm{rms}}$
, where
$J_{\mathrm{rms}}$
is the root-mean-squared current density, i.e.
$({\boldsymbol{\nabla }}\times {\boldsymbol{B}})_{\mathrm{rms}}$
. (We recall that the permeability was set to unity; otherwise, there would have been an extra
$\mu _0$
factor in front of
$J_{\mathrm{rms}}$
.) Both in experiments and in simulations,
$\xi _{\mathrm{T}}$
may be more easily accessible than
$\xi _{\mathrm{M}}$
, so it is important to find out whether the two obey similar scaling relations.
During the decay,
$\mathcal{E}_{\mathrm{M}}=v_{\mathrm{A}}^2/2$
decreases and
$\xi _{\mathrm{M}}$
increases. The Alfvén time, i.e. the ratio
$\tau _{\mathrm{A}}\equiv \xi _{\mathrm{M}}/v_{\mathrm{A}}$
, therefore also increases; see Banerjee & Jedamzik (Reference Banerjee and Jedamzik2004) and Hosking & Schekochihin (Reference Hosking and Schekochihin2023a) for early considerations of this point. Both for standard (isotropic) helical decay with
$v_{\mathrm{A}}\propto t^{-1/3}$
and
$\xi _{\mathrm{M}}\propto t^{2/3}$
, as well as for non-helical decay with
$v_{\mathrm{A}}\propto t^{-5/9}$
and
$\xi _{\mathrm{M}}\propto t^{4/9}$
, the value of
$\tau _{\mathrm{A}}$
increases linearly with
$t$
, i.e.
\begin{equation} t\propto \tau _{\mathrm{A}}(t). \end{equation}
This is also consistent with the idea that the turbulent decay is self-similar (Brandenburg & Kahniashvili Reference Brandenburg and Kahniashvili2017). It was found that the ratio
$t/\tau _{\mathrm{A}}(t)$
approaches a constant that increases with the Lundquist number (Brandenburg et al. Reference Brandenburg, Neronov and Vazza2024). The difference between the quantity
$t/\tau _{\mathrm{A}}(t)$
and the factor
$C_{\mathrm{M}}$
defined by Brandenburg et al. Reference Brandenburg, Neronov and Vazza(2024) is the exponent
$p=10/9$
in the relation
$\mathcal{E}_{\mathrm{M}}\propto t^{-p}$
for non-helical and
$p=2/3$
for helical turbulence with
$t/\tau _{\mathrm{A}}=C_{\mathrm{M}}/p$
.
To compute the Hosking integral, we need the function
$\mathcal{I}_{\mathrm{H}}(R,t)$
, which is a weighted integral over
$\textrm{Sp}(h)$
, given by
\begin{align} \mathcal{I}_{\mathrm{H}}(R,t)&=\int _0^\infty w(k,R)\, \textrm{Sp}(h;k,t) \, \mathrm{d}k, \quad \textrm{where} \; w(k,R)=\frac {4\pi R^3}{3} \left [\frac {6j_1(kR)}{kR}\right ]^2,\end{align}
and
$j_1(x)=(\sin x-x\cos x)/x^2$
is the spherical Bessel function of order one. As shown by Zhou et al. (Reference Zhou, Sharma and Brandenburg2022), the function
$\mathcal{I}_{\mathrm{H}}(R,t)$
yields the Hosking integral in the limit of large radii
$R$
, although
$R$
must still be small compared with the size of the domain. They referred to this as the box-counting method for a spherical volume with radius
$R$
.
3. Results
3.1. Isotropisation
In figure 1, we show the evolution of
$\langle \boldsymbol{J}_{\perp \perp }^2\rangle /\langle \boldsymbol{J}^2\rangle$
for Roberts fields I and II. We see that, after a short decay phase, exponential growth commences followed by a saturation of this ratio. We expect the ratio
$\langle \boldsymbol{J}_{\perp \perp }^2\rangle /\langle \boldsymbol{J}^2\rangle$
to reach the value
$4/15$
at late times; see Appendix A. The insets of figure 1 show the degree to which this is achieved at late times. Especially in the helical case, when inverse cascading is strong, the peak of the spectrum has already reached the lowest wavenumber of the domain. This is probably the reason why the value of 4/15 has not been reached by the end of the simulation. However, also for the non-helical case, the system retains memory of the initial state for a very long time; see the insets of both panels.
The early growth of
$\langle \boldsymbol{J}_{\perp \perp }^2\rangle /\langle \boldsymbol{J}^2\rangle$
shows that both the Roberts fields I and II are unstable to perturbations and develop an approximately isotropic state. The normalised growth rates are given in table 1 along with the times
$t_{\mathrm{p}}$
of maximum growth. The normalised values are in the range 0.7–6, but mostly around unity for intermediate values with
$k_0=16$
. The normalised times,
$t_{\mathrm{p}}v_{\mathrm{A0}} k_0$
, tend to decrease with increasing values of
$k_0$
and are approximately 10–20 times larger for field I than for field II. This difference was also found in another set of simulations in which
$B_0$
was the same for fields I and II; see Appendix B.
Normalised growth rates
$\lambda$
and peak times
$t_{\textrm{p}}$
for different values of
$k_0/k_1$
. The hyphen indicates that no growth occurred.
Visualisations of
$B_z$
on the periphery of the computational domain are shown in figure 2 for Roberts fields I and II. The initially tube-like structures are seen to decay much faster for Roberts field II. At time
$t=100$
, the magnetic field has much larger structures for Roberts field I than at time
$t=1000$
for Roberts field II.
Visualisations of
$B_z$
on the periphery of the computational domain at times
$t=1$
, 10, 30 and 100 for Roberts field I (top) and at times
$t=1$
, 10, 100 and 1000 for Roberts field II (bottom).
3.2. Spectral evolution
In figure 3, we plot magnetic energy and magnetic helicity variance spectra for the Roberts field I. Note that the spectra are normalised by
$v_{\mathrm{A}}^2 k_0^{-1}$
and
$v_{\mathrm{A}}^4 k_0^{-3}$
, respectively. At early times, the spectra show spikes at
$k\approx k_{\mathrm{f}}$
and
$2k_0$
, respectively, along with higher harmonics. We also show the time evolution of the normalised values of these spectra at the lowest wavenumber
$k=k_1$
. For
$\textrm{Sp}(h)$
, we also scale by
$2\pi ^2/k^2$
, which then gives an approximation to the value of the Hosking integral (Hosking & Schekochihin Reference Hosking and Schekochihin2021). Again, we see a sharp rise in both time series when the fields becomes unstable.
Evolution of magnetic energy and magnetic helicity variance spectra,
$\textrm{Sp}({\boldsymbol{B}})$
and
$\textrm{Sp}(h)$
, respectively, for Roberts field I with
$k_0=16$
at different times
$t_i$
indicated by different colours and line types as seen in the time traces on the right. The open black symbols in panels (b) and (d) correspond to the dotted lines in panels (a) and (c).
We also see that at late times, a bump appears in the spectrum near the Nyquist wavenumber. This shows that the Lundquist number was somewhat too large for the resolution of
$1024^3$
. However, comparing with simulations at lower Lundquist numbers shows that the large-scale evolution has not been adversely affected by this.
In figure 4, we show the same spectra for the case of Roberts fields II. Again, we see spikes in the spectra at early times. Those of
$\textrm{Sp}({\boldsymbol{B}})$
are again at
$\sqrt {2}k_0$
, along with overtones, but those of
$\textrm{Sp}(h)$
are now at
$2\sqrt {2}k_0$
instead of
$2\,k_0$
, and there are no spikes of
$\textrm{Sp}(h)$
at
$t=0$
. This is a consequence of the fact that the field has zero initial helicity pointwise, and helicity is quickly being produced owing to the growth of the initial perturbations. The plot of
$\textrm{Sp}(h;k_1,t)$
shows nearly perfectly a constant level for
$tv_{\mathrm{A}} k_0=100$
. This indicates that the Hosking integral is well conserved by that time.
Same as figure 3, but for the Roberts field II at different times
$t_i$
as seen in the time traces on the right.
3.3. Spontaneous production of magnetic helicity variance
As we have seen from figure 4, the case of zero magnetic helicity variance is unstable and there is a rapid growth of
$\textrm{Sp}(h)$
also at small wavenumbers. This was already anticipated by Hosking & Schekochihin (Reference Hosking and Schekochihin2021) and the present experiments with the Roberts field II show this explicitly.
We now discuss the function
$\mathcal{I}_{\mathrm{H}}(R,t)$
; see Hosking & Schekochihin (Reference Hosking and Schekochihin2021) and Zhou et al. (Reference Zhou, Sharma and Brandenburg2022). The result is shown in figure 5. For small values of
$R$
,
$\mathcal{I}_{\mathrm{H}}(R)$
increases
$\propto R^3$
. This indicates that the mean squared magnetic helicity density is not randomly distributed on those scales. In the present case, the actual scaling is slightly shallower than
$R^3$
, which is probably due to the finite scale separation. For
$R\approx 1$
, corresponding to scales compatible to the size of the computational domain, we see that
$\mathcal{I}_{\mathrm{H}}(R)$
has a plateau. It is at those scales,
$R=R_\ast$
, that we must determine the Hosking integral
$I_{\mathrm{H}}(t)=\mathcal{I}_{\mathrm{H}}(t,R_\ast )$
. In figure 6, we show the time dependence of
$I_{\mathrm{H}}(t)$
for Roberts field II with
$k_0=16$
normalised both by
$v_{\mathrm{A0}}^4/k_0^5$
(which is constant) and by
$\mathcal{E}_{\mathrm{M}}^2\xi _{\mathrm{M}}^5$
(which is time-dependent). Note that the time axis is here also logarithmic. We see an early rapid growth of
$I_{\mathrm{H}}(t)$
proportional to
$t^8$
by over eight orders of magnitude. The detailed mechanism behind this initial generation of magnetic helicity variance is not clear. A comparison with a 20 times more resistive run shows the same initial growth
$\propto t^8$
. This suggests that it is not a resistive effect. We are therefore tempted to associate the magnetic helicity variance generation with the scrambling of the initially perfectly pointwise non-helical magnetic field. In figure 6, we have indicated this with a question mark to say that this is tentative.
$\mathcal{I}_{\mathrm{H}}(R)$
for Roberts field II with (a)
$k_0=4$
at
$t=1$
(black), 1.5 (blue), 2.2 (green), 3.2 (orange) and 4.6 (red). and (b)
$k_0=16$
at
$t=46$
(black), 147 (blue), 316 (green), 570 (orange) and 824 (red). The arrow indicates the sense of time.
Time dependence of (a)
$I_{\mathrm{H}}(t)$
(black solid line) along with
$\mathcal{E}_{\mathrm{M}}^2\xi _{\mathrm{M}}^5$
(red solid line) in units of
$v_{\mathrm{A}}^4k_0^{-5}$
as well as
$\mathcal{E}_{\mathrm{M}}^2/v_{\mathrm{A0}}^4$
(blue dashed line) and
$\xi _{\mathrm{M}}^5 k_0^5$
(orange dashed line) and (b) the ratio
$I_{\mathrm{H}}/\mathcal{E}_{\mathrm{M}}^2\xi _{\mathrm{M}}^5$
for Roberts field II with
$k_0=16$
. The plateaus at 0.03 and 3000 are marked by dotted lines. In panel (a), the dash-dotted straight lines indicate the slopes
$\propto t^8$
(black),
$\propto t^{3}$
(orange) and
$\propto t^{-3}$
(blue). The inset in panel (a) shows the growth of
$I_{\mathrm{H}}(t)$
in a semilogarithmic representation along with a line
$\propto e^{30 t}$
.
Previous work showed that the value of
$I_{\mathrm{H}}(t)$
can greatly exceed the dimensional estimate
$\mathcal{E}_{\mathrm{M}}^2\xi _{\mathrm{M}}^5$
(Zhou et al. Reference Zhou, Sharma and Brandenburg2022). Figure 6 shows that at late times,
$tv_{\mathrm{A0}} k_0\gt 100$
, this is also the case here. After the initial rapid growth phase, however, the normalised value of
$I_{\mathrm{H}}(t)$
is still well below unity (approximately 0.03). The growth of
$I_{\mathrm{H}}/\mathcal{E}_{\mathrm{M}}^2\xi _{\mathrm{M}}^5$
after
$tv_{\mathrm{A0}} k_0\gt 100$
is mostly due to the decay of
$\mathcal{E}_{\mathrm{M}}$
and it is later counteracted by a growth of
$\xi _{\mathrm{M}}$
. The dashed blue and orange lines in figure 6(a) show separately the evolutions for
$\mathcal{E}_{\mathrm{M}}^2/v_{\mathrm{A0}}^4$
and
$\xi _{\mathrm{M}}^5k_0^5$
, respectively.
If the Hosking scaling applies to the present case of initially anisotropic MHD turbulence, we expect
$\xi _{\mathrm{M}}\propto t^{4/9}$
and therefore
$\xi _{\mathrm{M}}^5\propto t^{20/9}$
. The actual slope seen in figure 6 is however approximately 3 at late times. For
$\mathcal{E}_{\mathrm{M}}$
, we expect a
$t^{-10/9}$
scaling and therefore
$\mathcal{E}_{\mathrm{M}}^2\propto t^{-20/9}$
, i.e. the reciprocal one of
$\xi _{\mathrm{M}}^5$
. Again, the numerical data suggest a larger value of approximately 3. In § 4.1, we analyse in more detail the anticipated scaling of
$\mathcal{E}_{\mathrm{M}}(t)\propto t^{-p}$
and
$\xi _{\mathrm{M}}\propto t^q$
. We find that the two instantaneous scaling exponents
$p$
and
$q$
are indeed larger than what is expected based on the Hosking phenomenology. However, the instantaneous scaling exponents also show a clear evolution towards the expected values.
It is interesting to observe that the evolution of
$I_{\mathrm{H}}$
proceeds in two distinct phases. In the first one, when
$tv_{\mathrm{A0}} k_0\lt 2$
,
$I_{\mathrm{H}}$
shows a rapid growth that is not exponential; see the inset of figure 6, where the growth of
$I_{\mathrm{H}}$
is shown on a semilogarithmic representation. The growth is closer to that of a power law, and the approximate exponent would be approximately eight, which is rather large. During this phase, the turbulent cascade has not yet developed, but a non-vanishing and nearly constant value of
$I_{\mathrm{H}}$
has been established. However, in units of
$\mathcal{E}_{\mathrm{M}}^2\xi _{\mathrm{M}}^5$
, its value is rather small (approximately 0.03).
In the second phase, when
$tv_{\mathrm{A0}} k_0\gt 100$
, turbulence has developed and a turbulent decay is established. It is during this time that the ratio
$I_{\mathrm{H}}(t)/\mathcal{E}_{\mathrm{M}}^2\xi _{\mathrm{M}}^5$
approaches larger values (approximately 3000) that were previously seen in isotropic decaying turbulence simulations (Zhou et al. Reference Zhou, Sharma and Brandenburg2022). The reason for this large value was argued to be due to the development of non-Gaussian statistics in the magnetic field. However, Brandenburg & Banerjee (Reference Brandenburg and Banerjee2025) presented an estimate in which the value of this ratio is equal to
$C_{\mathrm{M}}^2$
. With
$C_{\mathrm{M}}\approx 50$
, this would agree with the numerical findings.
4. Cosmological applications
4.1. Evolution in the diagnostic diagram
In view of the cosmological applications of decaying MHD turbulence, it is of interest to consider the evolution of the actual Alfvén speed
$v_{\mathrm{A}}(t)=\sqrt {2\mathcal{E}_{\mathrm{M}}/\rho }$
in an evolutionary diagram as a parametric representation versus
$\xi _{\mathrm{M}}(t)$
; see figure 7(a). With
$v_{\mathrm{A}}\propto t^{-p/2}$
and
$\xi _{\mathrm{M}}\propto t^{q}$
, we expect that
$v_{\mathrm{A}}\propto \xi _{\mathrm{M}}^{-\kappa }$
with
$\kappa =p/2q=1/2$
for the fully helical case of Roberts field I. This is in agreement with early work showing that
$v_{\mathrm{A}}\propto t^{1/3}$
and
$\xi _{\mathrm{M}}\propto t^{2/3}$
(Hatori Reference Hatori1984; Biskamp & Müller Reference Biskamp and Müller1999).
(a) Parametric representation of
$v_{\mathrm{A}}$
versus
$\xi _{\mathrm{M}}$
for Roberts fields I (red) and II (blue). The solid (dotted) curves are for
$\eta =2\times 10^{-7}$
(
$\eta =4\times 10^{-6}$
). Note that the red dotted line for
$\eta =4\times 10^{-6}$
starts at the same value
$v_{\mathrm{A}}=\sqrt {1.28}$
as the non-helical runs (blue lines). The similarity between the dotted and solid red lines shows that the initial amplitude does not matter much. The open (filled) symbols indicate the times
$t=10$
(
$t=100$
). The dash-dotted lines give the slopes
$\kappa =1/2$
and 5/4 for Roberts fields I (red) and II (blue), respectively. (b)
$pq$
diagram field fields I (red) and II (blue) with
$\eta =2\times 10^{-7}$
. Larger symbols indicate later times.
In figure 7(a), we have also marked the times
$t=10$
(open symbols) and
$t=100$
(filled symbols). The points of constant times depart significantly from the lines of constant Alfvén time,
$\tau _{\mathrm{A}}$
, for which
$v_{\mathrm{A}}=\xi _{\mathrm{M}}/\tau _{\mathrm{A}}$
grows linearly with
$\xi _{\mathrm{M}}$
. We expect the times to be larger by a factor
$C_{\mathrm{M}}$
than the corresponding values of
$\tau _{\mathrm{A}}(t)$
. This is indeed the case: the point
$t=100$
lies on the line
$\tau _{\mathrm{A}}=1$
, i.e.
$t/\tau _{\mathrm{A}}=100$
. This is twice as much as our nominal value of approximately 50.
There is an interesting difference between the cases of Roberts fields I and II in that for field II, there is an extended period during which
$\xi _{\mathrm{M}}$
shows a rapid decrease before the expected increase emerges. The fact that such an initial decrease of the characteristic length scale is not seen for Roberts field I is remarkable. The rapid development of smaller length scales is probably related to the breakup of the initially organised tube-like structures into smaller scales. In the helical case, however, the nonlinear interaction among helical modes can only result in the production of modes with smaller wavenumbers, i.e. larger length scales; see Frisch et al. Reference Frisch, Pouquet, Leorat and Mazure(1975) and Brandenburg & Subramanian (Reference Brandenburg and Subramanian2005) for a review. Such a constraint does not exist for the non-helical modes, where this can then reduce the effective starting values of
$\xi _{\mathrm{M}}$
and therefore also of the effective Alfvén time,
$\tau _{\mathrm{A}}=\xi _{\mathrm{M}}/v_{\mathrm{A}}$
, early in the evolution. In Appendix B, we present similar diagrams for different values of
$k_0$
, but with a drag term included that could be motivated by cosmological applications.
We inspect the time-dependences of
$t/\tau _{\mathrm{A}}=v_{\mathrm{A}} t/\xi _{\mathrm{M}}$
and
${\textrm{Lu}}=v_{\mathrm{A}}\xi _{\mathrm{M}}/\eta$
for Roberts fields I and II in figure 8. We see that
$t/\tau _{\mathrm{A}}(t)$
reaches values in excess of 100 for
$t=100$
in both cases. This is more than what has been seen before, but it also shows significant temporal variations.
(a)
$t/\tau _{\textrm{A}}$
and (b)
$\textrm{Lu}$
versus time for Roberts fields I (red) and II (blue).
4.2. Universality of prefactors in the decay laws?
The decay of a turbulent magnetic field is constrained by certain conservation laws: the conservation of mean magnetic helicity density
$I_{\mathrm{M}}=\langle h\rangle$
, where
$h=\boldsymbol{A}\boldsymbol{\cdot} {\boldsymbol{B}}$
is the local magnetic helicity density, and the Hosking integral,
$I_{\mathrm{H}}=\int h(\boldsymbol{x})h(\boldsymbol{x}+{\boldsymbol{r}})\,{\textrm{d} {}}^3r$
. When the magnetic field is fully helical, the decay is governed by the conservation of
$I_{\mathrm{M}}$
, and when it is non-helical, it is governed by the conservation of
$I_{\mathrm{H}}$
. The time of cross-over depends on the ratio
$t_\ast \equiv I_{\mathrm{H}}^{1/2}/I_{\mathrm{M}}^{3/2}$
(Brandenburg & Banerjee Reference Brandenburg and Banerjee2025). Specifically, the correlation length
$\xi _{\mathrm{M}}(t)$
, the mean magnetic energy density
$\mathcal{E}_{\mathrm{M}}(t)$
and the envelope of the peaks of the magnetic energy spectrum
${E_{\mathrm{M}}}(k,t)$
depend on the values of the conserved quantities with (Brandenburg & Larsson Reference Brandenburg and Larsson2023)
\begin{equation} \xi _{\mathrm{M}}(t)=C_i^{(\xi )} I_i^\sigma t^q,\quad \mathcal{E}_{\mathrm{M}}(t)=C_i^{(\mathcal{E})} I_i^{2\sigma } t^{-p},\quad {E_{\mathrm{M}}}(k)\leqslant C_i^{(E)} I_i^{(3+\beta )\,q} k^\beta , \end{equation}
where
$\sigma$
is the exponent with which the length enters in
$I_{i}$
:
$\sigma =1/3$
when the mean magnetic helicity density governs the decay (
$i=\textrm{M}$
) and
$\sigma =1/9$
for the Hosking integral (
$i=\textrm{H}$
). In figure 9, we show the appropriately compensated evolutions of
$\xi _{\mathrm{M}}$
and
$\mathcal{E}_{\mathrm{M}}$
such that we can read off the values of
$C_i^{(\xi )}$
and
$C_i^{(\mathcal{E})}$
for the helical and non-helical cases.
Compensated evolutions of
$\xi _{\mathrm{M}}$
and
$\mathcal{E}_{\mathrm{M}}$
allowing the non-dimensional prefactors in (4.1) to be estimated.
In table 2, we summarise the values for the six coefficients reported previously in the literature and compare with those determined here. The fact that the coefficients are now somewhat different under the present circumstances suggests that they might not be universal, although the anisotropy of the present set-up as well as the limited scale separation may have contributed to the new results. For the purpose of providing relevant information for future studies of anisotropic magnetic decay, we present in Appendix C the temporal evolution of the length scales and field strengths in the parallel and perpendicular directions.
Comparison of the dimensionless prefactors with values in earlier papers.
The question of universality is significant, however, because universality would mean that the decay laws of the form (e.g. Vachaspati Reference Vachaspati2021)
\begin{equation} \xi _{\mathrm{M}}(t)=\xi _{\mathrm{M}}(t_0)\,\left (t/t_0\right )^q,\quad \mathcal{E}_{\mathrm{M}}(t)=\mathcal{E}_{\mathrm{M}}(t_0)\,\left (t/t_0\right )^{-p} \end{equation}
could be misleading in that they suggest some freedom in the choice of the values of
$\xi _{\mathrm{M}}(t_0)$
and
$\mathcal{E}_{\mathrm{M}}(t_0)$
at the time
$t_0$
. Comparing with (4.1), we see that
\begin{equation} \xi _{\mathrm{M}}(t_0)/t_0^q=C_i^{(\xi )} I_i^\sigma \quad \textrm{and}\quad \mathcal{E}_{\mathrm{M}}(t_0)\,t_0^p=C_i^{(\mathcal{E})} I_i^{2\sigma }, \end{equation}
so they cannot be chosen arbitrarily, but they must obey a constraint that depends on the relevant conservation law.
5. Conclusions
We have seen that a tube-like arrangement of an initial magnetic field becomes unstable to small perturbations. The resulting magnetic field becomes turbulent and tends to isotropise over time. This means that tube-like initial conditions that could be expected in plasma experiments would allow us to study the turbulent MHD decay dynamics – even for moderate but finite scale separation of 4:1 or more. In other words, the number of tubes per side length should be at least four.
We have also seen that a pointwise non-helical magnetic field, as in the case of the Roberts field II, is unstable and develops magnetic helicity fluctuations. After approximately one Alfvén time, the Hosking integral reaches a finite value, but a fully turbulent decay commences only after approximately one hundred Alfvén times. From that time onwards, the value of the Hosking integral relative to that expected on dimensional grounds reaches a value of several thousand, a value that was also found earlier (Zhou et al. Reference Zhou, Sharma and Brandenburg2022).
Our present results have confirmed the existence of a resistively prolonged turbulent decay time whose value exceeds the Alfvén time by a factor
$C_{\mathrm{M}}\approx \tau /\tau _{\mathrm{A}}$
. As emphasised previously, the fact that this ratio depends on the microphysical magnetic diffusivity is in principle surprising, because one of the hallmarks of turbulence is that its macroscopic properties should not depend on the microphysics of the turbulence. It would mean that it is not possible to predict this behaviour of MHD turbulence by ignoring the microphysical magnetic diffusivity, as is usually done in so-called large eddy simulations.
The present results have shown that the decay time can exceed the Alfvén time by a factor of approximately 50–100, which is similar to what was found previously (Brandenburg et al. Reference Brandenburg, Neronov and Vazza2024). During intermediate times, however, the decay time can even be a hundred times longer than the Alfvén time. The dimensionless prefactors in the dimensionally motivated powerlaw expressions for length scale and mean magnetic energy density are also roughly similar to what was previously obtained from fully isotropic turbulence simulations.
Acknowledgements
We thank the two referees for detailed suggestions. In particular, we acknowledge Dr D. N. Hosking for suggestions regarding an anisotropic generalisation of the Hosking scaling.
Editor Alex Schekochihin thanks the referees for their advice in evaluating this article.
Funding
This work was supported in part by the Swedish Research Council (Vetenskapsrådet, 2019-04234), the National Science Foundation under grant no. NSF AST-2307698 and a NASA ATP Award 80NSSC22K0825. National Key R&D Program of China (No. 2021YFA1601700), and the National Natural Science Foundation of China (No. 12475246). We acknowledge the allocation of computing resources provided by the Swedish National Allocations Committee at the Center for Parallel Computers at the Royal Institute of Technology in Stockholm and Linköping.
Declaration of interests
The authors report no conflict of interest.
Data availability statement
The data that support the findings of this study are openly available on Zenodo at doi: https://doi.org/10.5281/zenodo.15739684 (v2025.06.25) or, for easier access, at http://norlx65.nordita.org/∼brandenb/projects/Roberts-Decay/. All calculations have been performed with the Pencil Code (Pencil Code Collaboration et al. 2021); DOI: https://doi.org/10.5281/zenodo.3961647.
Appendix A.
$\langle \boldsymbol{J}_{\perp \perp }^2\rangle /\langle \boldsymbol{J}^2\rangle$
for isotropic turbulence
We have examined the evolution of
$\langle \boldsymbol{J}_{\perp \perp }^2\rangle /\langle \boldsymbol{J}^2\rangle$
for isotropic turbulence using a set-up similar to that of Brandenburg et al. (Reference Brandenburg, Sharma and Vachaspati2023); see figure 10. The scale separation, i.e. the ratio of the peak wavenumber to the lowest wavenumber in the domain is 8 for this simulation and the Lundquist number, which is the r.m.s. Alfvén speed times the correlation length divided by the magnetic diffusivity, is approximately
$10^4$
. The other parameters are as in the earlier work of Brandenburg et al. (Reference Brandenburg, Sharma and Vachaspati2023); see the data availability statement of the present paper.
Evolution of
$\langle \boldsymbol{J}_{\perp \mathrm{m}}^2\rangle /\langle \boldsymbol{J}^2\rangle$
(green),
$\langle \boldsymbol{J}_{\perp \perp }^2\rangle /\langle \boldsymbol{J}^2\rangle$
(blue),
$\langle \boldsymbol{J}_{\perp\|}^2\rangle /\langle \boldsymbol{J}^2\rangle$
(orange),
$\langle \boldsymbol{J}_{\perp}^2\rangle /2\langle \boldsymbol{J}^2\rangle$
(red), and
$\langle \boldsymbol{J}_{\|}^2\rangle /\langle \boldsymbol{J}^2\rangle$
(black) for decaying isotropic turbulence with an initial peak wavenumber
$k_0/k_1=8$
using
$1024^3$
meshpoints (a) with helicity and (b) without helicity.
Appendix B. Diagnostic diagrams for different
$k_0$
In figure 7, we did already present a diagnostic diagrams of
$v_{\mathrm{A}}$
versus
$\xi _{\mathrm{M}}$
for
$k_{\mathrm{p}}=16$
. We also performed runs for different values of
$k_{\mathrm{p}}$
to compute the growth rates and the times
$t_{\mathrm{p}}$
of maximum growth in table 1, but not all the runs were long enough to compute similar tracks in the diagnostic diagram. In figure 11, we show such a diagram for a case in which a drag term of the form
$-\alpha {\boldsymbol{u}}$
is included on the right-hand side of (2.7). Here, we choose a drag coefficient that automatically changes in time so as to allow for a nearly self-similar decay. Using a multiple of
$1/t$
is an obvious possibility, but it would always be the same at all locations and for different types of flows. The local vorticity might be one possible option for a coefficient that varies in space and time, and has the right dimension. Another possibility, which is also the one chosen here, is to take
$\alpha$
to be a multiple of
$\sqrt {\mu _0/\rho _0}|\boldsymbol{J}|$
and write
$\alpha =c_\alpha \sqrt {\mu _0/\rho _0}|\boldsymbol{J}|$
, where
$c_\alpha$
is a dimensionless prefactor and
$\mu _0=\rho _0=1$
has been set. Again, as was already clear from figure 7, the tracks without helicity show a marked excursion to smaller values of
$\xi _{\mathrm{M}}$
before displaying a decay of the form
$v_{\mathrm{A}}\propto \xi _{\mathrm{M}}^{-\kappa }$
. The corresponding values of
$\lambda /v_{\mathrm{A0}} k_0$
and
$t_{\textrm{p}}v_{\mathrm{A0}} k_0$
are given in table 3.
Same as figure 7(a), but for
$c_\alpha =3$
, showing a parametric representation of
$B_{\mathrm{rms}}$
versus
$B_{\mathrm{rms}}/J_{\mathrm{rms}}$
and
$\xi _{\mathrm{M}}$
for Roberts field I (a) and II (b) with
$k_0=2$
(black),
$4$
(blue),
$8$
(green),
$16$
(orange),
$32$
(red),
$64$
(black) and 128 (blue). The open (filled) symbols in both plots indicate the times
$t=10$
(
$t=100$
).
Similar to table 1, showing normalised growth rates
$\lambda$
and peak times
$t_{\mathrm{p}}$
for different values of
$k_0$
, but with the photon drag term included. Here, unlike the case of table 1, the values of
$B_0$
are the same for Roberts fields I and II. The hyphen indicates that no growth occurred. Note that we used here what we called the rotated Roberts field.
Our definition of the Roberts fields follows the earlier work by Rheinhardt et al. (Reference Rheinhardt, Devlen, Rädler and Brandenburg2014). In the original paper by Roberts (Reference Roberts1972), however, the field was rotated by
$45\hbox{$^\circ $}$
. In that case,
$\phi =\cos k_0 x \mp \cos k_0 y$
, where the upper and lower signs refer to Roberts fields I and II. For this field, a lower eigenvalue of the curl operator, namely
$k_{\mathrm{f}}=k_0$
, can be accessed. In that case, we can accommodated one pair of flux tubes instead of four. This can be done both for fields I and II. They are given by
\begin{equation} {\boldsymbol{B}}_{\mathrm{I}}=\begin{pmatrix} \sin k_0 y\\ \sin k_0 x\\ \cos k_0 x -\cos k_0 y\end{pmatrix}, \quad {\boldsymbol{B}}_{\mathrm{II}}=\begin{pmatrix} \sin k_0 y\\ \sin k_0 x\\ \cos k_0 x +\cos k_0 y\end{pmatrix}, \end{equation}
which satisfies
${\boldsymbol{B}}_{\mathrm{I}}\boldsymbol{\cdot} {\boldsymbol{\nabla }}\times {\boldsymbol{B}}_{\mathrm{I}}=k_{\textrm{f}}{\boldsymbol{B}}_{\mathrm{I}}^2$
and
${\boldsymbol{B}}_{\mathrm{II}}\boldsymbol{\cdot} {\boldsymbol{\nabla }}\times {\boldsymbol{B}}_{\mathrm{II}}=0$
, just like the non-rotated field. However, here,
$k_{\mathrm{f}}=k_0$
is the eigenvalue of the curl operator.
Appendix C. Anisotropy
Given that the magnetic field remains anisotropic for a long time, it is useful to consider the possible effects of anisotropy. For this purpose, we define the length scales
\begin{equation} \left .\xi _\perp (t)=\int k_\perp ^{-1}{E_{\mathrm{M}}}(k_\perp ,t)\,{\textrm{d} {}} k_\perp \right / \int {E_{\mathrm{M}}}(k_\perp ,t)\,{\textrm{d} {}} k_\perp , \end{equation}
\begin{equation} \left .\xi _\|(t)=\int k_\|^{-1}{E_{\mathrm{M}}}(k_\|,t)\,{\textrm{d} {}} k_\| \right / \int {E_{\mathrm{M}}}(k_\|,t)\,{\textrm{d} {}} k_\|,\\[-10pt] \end{equation}
which represent the typical length scales in the directions perpendicular and parallel to the magnetic flux tubes, respectively. In figure 12, we plot the evolution of
$\xi _\perp (t)$
and
$\xi _\|(t)$
along with that of
$B_\perp (t)$
and
$B_\|(t)$
for the non-helical case of Roberts field II. We see that there are no clear power laws. During limited time intervals, however, the curves have the slopes
$\propto t^{4/9}$
and
$\propto t^{-5/9}$
for the length scales and field strengths, respectively, as expected from an isotropic evolution.
Scalings of (a)
$\xi _\perp (t)$
and
$B_\perp (t)$
, and (b)
$\xi _\|(t)$
and
$B_\|(t)$
in red and blue, respectively, both for the non-helical case. The expected slopes
$\propto t^{4/9}$
and
$\propto t^{-5/9}$
are indicated for reference.
We demonstrated already that the three-dimensional magnetic energy spectrum increases
$\propto k^4$
; see figure 4. This shows that there are no long-range correlations; see Hosking & Schekochihin (Reference Hosking and Schekochihin2023b) for a corresponding demonstration in the hydrodynamic case and Zhou et al. (Reference Zhou, Sharma and Brandenburg2022) for the application to magnetic fields. However, our two-dimensional spectra (see figure 13), and especially that of
${\boldsymbol{B}}_\|$
, as a function of
$k_\perp$
, increases
$\propto k_\perp ^3$
; see figure 13(b). This shows that there are no long-range correlations of the flux of
${\boldsymbol{B}}_\|$
over the
$xy$
-plane. Thus, even if the flux of
${\boldsymbol{B}}_\|$
over
$xy$
-planes might constitute an additional corresponding conserved quantity, it could not constrain the dynamics in the present case, because such a quantity vanishes in our case.
Spectra of (a)
${\boldsymbol{B}}_\perp$
and (b)
${\boldsymbol{B}}_\|$
as a function of
$k_\perp$
in both panels. The last time is shown as a thick line. The sense of time is also shown by the arrows in both panels.
Open access
