Stochastic Flux-Freezing and Magnetic Dynamo

# Stochastic Flux-Freezing and Magnetic Dynamo

Gregory L. Eyink Department of Applied Mathematics & Statistics
and Department of Physics & Astronomy
The Johns Hopkins University, USA
July 8, 2019
###### Abstract

We argue that magnetic flux-conservation in turbulent plasmas at high magnetic Reynolds numbers neither holds in the conventional sense nor is entirely broken, but instead is valid in a novel statistical sense associated to the “spontaneous stochasticity” of Lagrangian particle trajectories. The latter phenomenon is due to the explosive separation of particles undergoing turbulent Richardson diffusion, which leads to a breakdown of Laplacian determinism for classical dynamics. We discuss empirical evidence for spontaneous stochasticity, including our own new numerical results. We then use a Lagrangian path-integral approach to establish stochastic flux-freezing for resistive hydromagnetic equations and to argue, based on the properties of Richardson diffusion, that flux-conservation must remain stochastic at infinite magnetic Reynolds number. As an important application of these results we consider the kinematic, fluctuation dynamo in non-helical, incompressible turbulence at unit magnetic Prandtl number. We present results on the Lagrangian dynamo mechanisms by a stochastic particle method which demonstrate a strong similarity between the and dynamos. Stochasticity of field-line motion is an essential ingredient of both. We finally consider briefly some consequences for nonlinear MHD turbulence, dynamo and reconnection.

###### pacs:
52.30.Cv, 52.35.Ra, 91.25.Cw, 52.35.Vd, 95.30.Qd

## I Introduction

Hannes Alfvén in a seminal paper in 1942 introduced the notion of flux-freezing in magnetohydrodynamic plasmas at infinite conductivity, noting that “every motion (perpendicular to the field) of the liquid in relation to the lines of force is forbidden because it would give infinite eddy currents” Alfvén (1942). In the years since, the property of flux-conservation has become a powerful tool in the analysis of many near-ideal plasma phenomena. For example, in his excellent monograph Kulsrud (2005), Kulsrud states that “The most important property of an ideal plasma is flux-freezing”, before proceeding to illustrate its many applications. Of course, physical plasmas in the the laboratory and in astrophysics are subject to various forms of non-ideality, including Spitzer resistivity, ambipolar diffusion, etc. The general assumption in the field of plasma physics, however, is that as long as such non-ideality is sufficiently “small”, then flux-freezing will hold in an approximate sense. This idea dominates the discussion of turbulent magnetic dynamo at high kinetic and magnetic Reynolds numbers. For example, it is commonplace to find statements in the literature such as the following: “The small-scale turbulent dynamo is caused by the random stretching of the (nearly) frozen-in field lines by the ambient random flow” Schekochihin et al. (2004), or “It is well established both analytically and numerically that a weak magnetic field can be amplified by the random motions of a highly conducting fluid [1-3]. This occurs because magnetic-field lines are generically stretched by the random motions of the fluid in which they are (almost) ‘frozen.’ ” Boldyrev et al. (2005). The physical idea which underlies these statements is that a tiny resistivity should diffuse field-lines only a short distance through the plasma. For example, the usual quantitative estimate is well expressed in this quote from Kulsrud’s monograph, Ch.13, on magnetic reconnection:

“Flux freezing is a very strong constraint on the behavior of magnetic fields in astrophysics. As we show in chapter 3, this implies that lines do not break and their topology is preserved. The condition for flux freezing can be formulated as follows: In a time , a line of force can slip through the plasma a distance

 ℓ=√ηct4π.(1)

If this distance is small compared to , the scale of interest, then flux freezing holds to a good degree of approximation.”

We shall argue that these commonplace ideas on flux-freezing are wrong. They contain an implicit assumption that the plasma fluid remains smooth and laminar for very small non-ideality. The quantitative estimate that field-lines slip through a resistive plasma only a diffusive distance in time is incorrect—by many orders of magnitude—in a turbulent plasma. Since laminar flow at very high kinetic and magnetic Reynolds is unstable to development of turbulence, flux-conservation in the conventional sense must be the exception rather than the rule in astrophysical plasmas. Indeed, we shall show that the standard views on flux-freezing must be incorrect, because the very notion of a Lagrangian fluid particle trajectory breaks down in turbulent flow with a spatially “rough” velocity field (i.e. with a power-law kinetic energy spectrum similar to that of Kolmogorov). Recent research has discovered a novel phenomenon of “spontaneous stochasticity”, according to which fluid particle trajectories are intrinsically random in high-Reynolds-number turbulence Bernard et al. (1998); Gawȩdzki and Vergassola (2000); Chaves et al. (2003); Kupiainen (2003); E and Vanden-Eijnden (2000, 2001); Jan and Raimond (2002, 2004). This surprising phenomenon is a long overlooked consequence of the fluid-dynamical effect of Richardson two-particle turbulent dispersion Richardson (1926). Because of spontaneous stochasticity, it makes no sense to assume that a field line follows “the” plasma fluid element, because there are infinitely many distinct fluid trajectories starting from the same point! But it is also not true that flux-freezing is completely broken. We shall argue below that magnetic flux-conservation remains valid in the ideal limit of high Reynolds numbers, but in a novel stochastic sense associated with the intrinsic stochasticity of the Lagrangian particle trajectories.

A correct formulation of flux-freezing is fundamental to understand a number of important astrophysical processes, such as turbulent dynamo and reconnection. In previous work Eyink and Neto (2010); Eyink (2010) we have shown how stochastic flux-freezing is involved in the small-scale “fluctuation dynamo” for a soluble model problem: magnetic fields advected by the Kazantsev-Kraichnan ensemble of velocity fields that are spatially rough and white-noise in time Kazantsev (1968); Kraichnan and Nagarajan (1967); Kraichnan (1968). It is worth noting, by the way, that “spontaneous stochasticity” is a rigorously established phenomenon for this model Bernard et al. (1998); Gawȩdzki and Vergassola (2000); Chaves et al. (2003); Kupiainen (2003); E and Vanden-Eijnden (2000, 2001); Jan and Raimond (2002, 2004). It was shown that presence of fluctuation dynamo effect at zero magnetic Prandtl number depends crucially on the degree of angular correlation between the infinite number of magnetic field vectors that are simultaneously advected by turbulence to the same spatial point Eyink and Neto (2010); Eyink (2010).

Here we shall make a similar study for kinematic dynamo in non-helical, hydrodynamic turbulence, by a Lagrangian numerical method that employs data from a high Reynolds-number turbulent flow archived online Li et al. (2008); Perlman et al. (2007). We present results for unit magnetic Prandtl number which demonstrate—and quantify—the effect of Richardson diffusion and stochastic flux-freezing on the small-scale turbulent dynamo. We find, in fact, remarkable similarities between the Lagrangian mechanisms of small-scale dynamo in hydrodynamic turbulence at unit Prandtl number and in the Kazantsev model at zero Prandtl number. Previous numerical studies Haugen et al. (2004); Schekochihin et al. (2004) have instead found a close similarity of the unit Prandtl-number fluctuation dynamo with the solution of the Kazantsev model at infinite Prandtl-number, observing, in particular, the “Kazantsev spectrum” of magnetic energy at high wavenumbers. One study Schekochihin et al. (2004) went so far as to claim that “ It is clear that the kinematic dynamo in such [] runs is of the large- kind.” This interpretation is ruled out by our new results at much higher Reynolds numbers, which show that the inertial-range phenomenon of Richardson diffusion strongly affects the exponential growth rate of magnetic energy in the kinematic regime. However, stochastic flux-freezing is not a property of kinematic dynamo only but will hold also for fully nonlinear MHD turbulence and have important implications there for magnetic dynamo and reconnection.

The detailed contents of this paper are as follows: In the following section II we shall briefly review the phenomenon of “spontaneous stochasticity”, both its theoretical bases and its present confirmation from simulations and experiments. We aso present new numerical results of our own which support the essential predictions. In Section III we discuss stochastic flux-freezing. We begin with a demonstration of the stochastic flux-conservation properties of resistive MHD, which are novel results themselves. We employ Lagrangian path-integral methods that provide good physical insight. We then discuss the ideal case, via zero-resistivity and other limits. In section IV we employ the new results to discuss the turbulent kinematic dynamo. After reviewing the Lagrangian theory of dynamo, we present our numerical results and their comparison with analytical results for the Kazantsev model at zero Prandtl number and with previous numerical studies. In section V we briefly discuss some implications and open problems for nonlinear MHD turbulence and section VI contains our final discussion. An Appendix sketches the derivation of the path-integral formulas used in the main text.

## Ii Richardson Diffusion and Spontaneous Stochasticity

### ii.1 Richardson 2-Particle Dispersion

We briefly review Richardson’s theory Richardson (1926) of 2-particle or relative turbulent dispersion, emphasizing perspectives of recent research. See also more complete reviews Sawford (2001); Salazar and Collins (2009); Toschi and Bodenschatz (2009).

The object of Richardson’s study was the separation between a pair of passive Lagrangian tracer particles in a turbulent flow, such as ash particles in a volcanic plume. Richardson’s approach was semi-empirical. By estimating the “effective diffusivity” as a function of rms separation he inferred from data that there is a scale-dependent diffusivity coefficient

 K(ℓ)∼K0ℓ4/3. (1)

This is essentially a “running coupling constant” in the modern sense of renormalization group theory. Richardson proposed further that the probability density function of the separation vector would satisfy a diffusion equation

 ∂tP(BellSyst. Tech. J. ,t)=∂∂ℓi(K(ℓ)∂P∂ℓi(BellSyst. Tech. J. ,t)). (2)

Richardson observed that there is an exact similarity solution of his equation given by a stretched-exponential PDF, which he wrote explicitly in one space dimenson. Here we note its form

 P∗(BellSyst. Tech. J. ,t)=A(K0t)9/2exp(−9ℓ2/34K0t) (3)

in the more physically relevant case of three space dimensions. All solutions of (2) approach this self-similar form asymptotically at long times Eyink and Xin (2000). Averaging with respect to the self-similar density (3) yields

 ⟨ℓ2(t)⟩=γ0t3 (4)

with This is the famous Richardson -law.

Richardson’s work preceded the Kolmogorov 1941 (K41) theory of turbulence, but it was shown by Obukhov Obukhov (1941) to be fully consistent with that theory. This can be seen by a toy calculation in one space dimension. Assume that satisfies the initial-value problem

 ddtℓ(t)=δu(ℓ)=32(g0εℓ)1/3,ℓ(0)=ℓ0,

with velocity increment scaling as in K41 theory, where is the mean energy dissipation per unit mass. Separation of variables gives the exact solution

 ℓ(t)=[ℓ2/30+(g0ε)1/3t]3/2. (5)

If one defines a time which characterizes the initial separation then, for

 ℓ2(t)∼g0εt3. (6)

For sufficiently long times, the particles “forget” their initial separation and the Richardson law is obtained. The dimensionless parameter which appears in this form of the -law is usually called the Richardson-Obukhov constant. The physical mechanism of the explosive separation of particles, even faster than ballistic, is the relative advection of the pairs by larger, more energetic eddies as their separation distance increases.

This physics seems relatively simple and benign, but it has extraordinary consequences. As first pointed out in a seminal paper of Bernard, Gawȩdzki, and Kupiainen Bernard et al. (1998), Richardson’s theory implies a breakdown in the usual notion of Laplacian determinism for classical dynamics! This may already be seen in our toy calculation above. If we set then the solution (5) becomes the Richardson law for all positive times Thus, two particles started at the same point at time separate to a finite distance at any time The same oddity may be seen in Richardson’s similarity solution (3), which satisfies at initial time

 P∗(BellSyst. Tech. J. ,0)=δ3(BellSyst. Tech. J. ).

All particles start with separation However, is a smooth density for so that with probability one at later times. Richardson’s theory thus implies that two particles advected by the fluid velocity which start at the same initial point

 ddtx(t)=u(x(t),t),x(0)=x0

can follow different trajectories. This seems to violate the theorem on uniqueness of solutions of initial-value problems for ODE’s. However, such theorems assume that the advecting velocity is Hölder-Lipschitz continuous in the space variable

 |u(x1,t)−u(x2,t)|≤C|x1−x2|h (7)

with exponent A turbulent velocity field in a Kolmogorov inertial range has instead Hölder exponent and the uniqueness theorem need not apply. Our toy calculation earlier is just the standard textbook example for failure of uniqueness (see HartmanHartman (2002), p.2) . In that example, for any non-negative “waiting time”

 ℓ(t)=(g0ε)1/2(t−τ)3/2+

is a solution of the initial-value problem with . (Here for and otherwise.)

The above considerations may seem fairly technical and mathematical. It has been shown, however, that this breakdown in uniqueness of trajectories can appear in various physical limits for turbulent advection. Even more remarkably, the solutions of the deterministic classical dynamics become intrinsically stochastic! See the important series of papers Bernard et al. (1998); Gawȩdzki and Vergassola (2000); Chaves et al. (2003); E and Vanden-Eijnden (2000, 2001). Advanced probablistic techniques have obtained the most refined results Jan and Raimond (2002, 2004). A very clear and concise review of the subject Kupiainen (2003) is available, which the reader may consult for further details. Below we briefly describe the most basic theory and results.

### ii.2 High-Reynolds-Number Limit and Spontaneous Stochasticity

The easiest way to understand the phenomenon is via the problem of stochastic particle advection,

 ddt˜x(t)=uν(˜x(t),t)+√2κ˜\boldmathη(t),x(t0)=x0 (8)

with advecting velocity perturbed by a Gaussian white-noise multiplied by and with velocity assumed spatially smooth at subviscous length-scales , for a finite viscosity The transition probability for a single particle in a fixed (non-random) velocity realization can be written using a “sum-over-histories” approach as a path-integral Drummond (1982); Shraiman and Siggia (1994); Bernard et al. (1998):

 Gν,κu(xf,tf|x0,t0)=∫x(t0)=x0Dxδ3(xf−x(tf)) (9) ×exp(−14κ∫tt0dτ|˙x(τ)−uν(x(τ),τ)|2). (10)

Since this formula plays an important role in our analysis, we provide a self-contained derivation in the Appendix. A physical motivation to study such random advection is the problem of the evolution of a passive scalar, such as a temperature field or dye concentration. These fields solve the scalar advection-diffusion equation

 ∂tθ+(uν\boldmath⋅\boldmath∇)θ=κ△θ, (11)

with the molecular diffusivity. The exact solution of (11) is given by the Feynman-Kac formula Shraiman and Siggia (1994); Bernard et al. (1998); Drummond (1982); Sawford (2001):

 θ(x,t)=∫d3x0θ(x0,t0)Gν,κu(x0,t0|x,t) (12) =∫a(t)=xDaθ(a(t0),t0) (13) ×exp(−14κ∫tt0dτ|˙a(τ)−uν(a(τ),τ)|2) (14)

for This corresponds to solving backward in time the stochastic equation

 ddτ˜a(τ)=uν(˜a(τ),τ)+√2κ˜\boldmathη(τ)

from to with the condition The present value of the scalar field is thus the average, along stochastic Lagrangian paths, of its earlier values.

It naively appears by an application of the Laplace asymptotic method to (10) that the transition probability collapses to a delta-function

 Gν,κu(xf,tf|x0,t0)→δ3(xf−x(tf)) (15)

as with the solution of the ODE for initial condition Only for such time-histories is the action vanishing in the exponent of the path-integral. However, it was shown Bernard et al. (1998); Gawȩdzki and Vergassola (2000); Chaves et al. (2003); Kupiainen (2003); E and Vanden-Eijnden (2000, 2001); Jan and Raimond (2002, 2004) that (15) may not hold if simultaneously (or the Reynolds-number ) and if in that limit the velocity field for a rough (non-smooth, singular) . In that case, as

 Gν,κu(xf,tf|x0,t0)→Gu(xf,tf|x0,t0) (16)

for a nontrivial probability density The Lagrangian trajectories can remain random as ! This phenomenon has been called “spontaneous stochasticity” Chaves et al. (2003), because of the analogy with spontaneous symmetry-breaking in condensed matter physics and quantum field-theory, where, for example, a ferromagnet may retain a non-vanishing magnetization even in the limit of vanishing external magnetic field. It is important to appreciate, however, that “spontaneous stochasticity” is a very different type of randomness than is usual in turbulence theory, associated to a random ensemble of velocity fields. Instead, the randomness in (16) is for a fixed (non-random) ensemble member . The limiting distribution consists of time-histories that are all solutions of the same deterministic initial-value problem

 ˙x=u(x,t),x(t0)=x0. (17)

As is clear from (10), the limiting probability measure is in a certain sense the uniform or equal-weight distribution over all such solutions. We note in passing that Kneser’s theorem in the mathematical theory of ODE’s implies that, whenever there is more than one such solution, then there is in fact a continuous infinity of solutions (see Hartman Hartman (2002), section II.4).

The above results have been rigorously established Bernard et al. (1998); Gawȩdzki and Vergassola (2000); Chaves et al. (2003); Kupiainen (2003); E and Vanden-Eijnden (2000, 2001); Jan and Raimond (2002, 2004) for some model problems of turbulent advection, primarily the Kraichnan model of advection by a Gaussian random velocity field with zero mean and covariance

 ⟨uνi(x,t)uνj(x′,t′)⟩=[D0δij−Sνij(x−x′)]δ(t−t′).

The velocity fields are temporal white-noise and spatially rough for with a Hölder exponent but smooth for We discuss the model only for incompressible flow, in which case

 Sνij(r)={D1[(1+h)δij−h^ri^rj]r2h,ℓν≪r≪Lu,D1ℓ2h−2ν[2δij−^ri^rj]r2,r≪ℓν. (18)

The velocity realizations are divergence-free and Hölder continuous with exponent for The kinetic energy spectra are power-laws for the “inertial range” with and thus The key feature which makes this model analytically tractable is the Markovian property in time, which leads to the exact validity of Richardson’s 2-particle diffusion equation in the form

 ∂tP(r,t)=∂∂ri(Sνij(r)∂P∂rj(r,t))+2κ△rP(r,t), (19)

with Note that Richardson’s original equation is obtained for rather than the Kolmogorov value a peculiarity of the white-noise approximation. The main features of the solutions of this equation can be inferred from a study of the dispersion. Multiplying (19) by and integrating over leads to

 ddt⟨rk(t)rℓ(t)⟩=2⟨Sνkℓ(r(t))⟩+4κδkℓ (20)

An analysis of (19) and (20) leads to the following key results for the ensemble of particles which are all started at the same point () .

For the dispersion is that of two independent Brownian motions in three space dimensions

 ⟨r2(t)⟩∼12κt (21)

for short times with but of Richardson-type

 ⟨r2(t)⟩∼gh(D1t)1/(1−h) (22)

for longer times with a constant independent of For the behavior is a bit more complex. The result (21) still holds at very short times but at large Prandtl numbers there is an intermediate range of exponential growth

 ⟨r2(t)⟩∝κtνe2λνt (23)

for Here is the leading Lyapunov exponent of the smooth advecting velocity field. The dispersion at longer times again follows the Richardson law (22). In either case, the Richardson law is valid once and thus holds for arbitrarily small times in the limit as

A non-vanishing dispersion implies that the Lagrangian trajectories must stay random in the limit. Although the diffusion equation (19) has been averaged over velocity realizations it must be the case that for and for a set of with nonzero probability, or otherwise the average over would also be a delta-function! The physical mechanism of spontaneous stochasticity is clearly the “forgetting” of the length-scales by Richardson diffusion for sufficiently long times, with those times also vanishing in the limit In the case of incompressible flow, the limiting distribution is completely independent of how the limit is taken. We note in passing that this is not true in general for compressible flows and that the possibility for Lagrangian particles to stick as well as to stochastically split allows there to be different limits, depending upon the Prandtl number in the limit as Gawȩdzki and Vergassola (2000); E and Vanden-Eijnden (2001); Jan and Raimond (2004). However, for incompressible flow the limiting distribution is very universal and robust.

The same limit is obtained for incompressible flow even with or if the randomness is introduced through the initial conditions rather than stochastic noise. Consider the solution of the initial-value problem for the Kraichnan ensemble of velocities

 ddt˜x(t)=uν(˜x(t),t),˜x(t0)=x0+ϵ˜% \boldmathρ

where is a zero-mean, unit-variance random vector with probability density One may interpret as the size of error in measuring the initial position of the particle. This corresponds to solving the Richardson diffusion equation (19) with and initial condition so that as However, if the limits are taken first and subsequently, then the solution of the Richardson equation does not degenerate to a delta-function for This may be seen by solving for the dispersion from eq. (20) with and When

 ⟨r2(t)⟩∼2ϵ2e2λνt (24)

for times but for longer times follows the Richardson law (22). When instead, then the short-time behavior is diffusive

 ⟨r2(t)⟩∼2ϵ2+(const.)D1ϵ2ht, (25)

for times Dispersion for such a “random cloud” of initial positions was first considered by Batchelor Batchelor (1950, 1952), who obtained instead ballistic growth for hydrodynamic turbulence. The diffusive result above is an artefact of the white-in-time velocity. At times the Richardson law (22) holds. As in the previous cases, the Richardson law holds for any time if with vanishing slower than

As this last discussion should make clear, the phenomenon of “spontaneous stochasticity” is not especially connected with the random perturbation of the motion equations in (8). Instead it is the advection of the particle by a rough velocity field and the “forgetting” of the initial separations which makes the Lagrangian particle motions intrinsically stochastic. Spontaneous stochasticity should not be confused with “chaos”, as that term is used in dynamical systems theory Chaves et al. (2003); Kupiainen (2003). For chaotic dynamical systems with a smooth velocity field, one sees only exponential growth of deviations as in eq.(24). Because this result is proportional to the initial separation is never “forgotten” and, for all times, as For chaotic dynamics, any imprecision in the initial data is exponentially magnified, leading to loss of predictability at long enough times. Spontaneous stochasticity corresponds instead to . The solution is unpredictable for all future times, even with infinitely precise knowledge of the initial conditions!

We remark finally that spontaneous stochasticity has very important implications for turbulent advection of a passive scalar Bernard et al. (1998). Note that

 ddt∫d3xθ2(x,t)=−2κ∫d3x|\boldmath∇θ(x,t)|2,

so that, naively, the integral is conserved for However, in the infinite Reynolds-number, fixed Prandtl-number limit ( with constant)

 θ(x,t)=∫d3x0θ(x0,t0)Gu(x0,t0|x,t),

using the results (14) and (16). The scalar at the present time is a nontrivial average along stochastic Lagrangian trajectories of its values at an earlier time, with molecular mixing replaced by turbulent mixing. In particular, the scalar intensity—for every velocity realization, without ensemble averaging— is still dissipated

 ∫d3xθ2(x,t)<∫d3xθ2(x,t0),t>t0

even as ! This is the scalar analogue of the dissipative anomaly of Onsager for fluid turbulence Onsager (1949); Eyink and Sreenivasan (2006); Eyink (2008). The Lagrangian mechanism is spontaneous stochasticity.

### ii.3 Experimental and Numerical Results

Our quantitative discussion above has been based upon the original Richardson theory and, in particular, his diffusion equation (19). This equation is exact for the Kraichnan white-in-time velocity ensemble but is only an approximation for hydrodynamic turbulence, where several of its quantitative predictions are known to be incorrect. We have already mentioned the diffusive short-time growth in dispersion for a particle cloud, eq.(25), which is ballistic for fluid turbulence. (Note that the ballistic regime is correctly predicted by Richardson’s diffusion equation for a suitably time-dependent eddy diffusion tensor Batchelor (1952); Kraichnan (1966).) The diffusion equation (19) holds in the Kraichnan model also for backward-in-time dispersion. However, actual dispersion rates are different forward and backward in time because of the negative skewness of turbulent velocity increments Kraichnan (1966); Sawford et al. (2005). There is presently no quantitative theory of turbulent dispersion which successfully accounts for all aspects of the phenomenon. It is necessary to stress that the prediction of “spontaneous stochasticity” has more general grounds in the mathematical theory of ODE’s and is not dependent upon the diffusion approximation (19). Nevertheless, in the absence of any fully successful, quantitative theory, it is important to develop understanding from numerical simulations and laboratory experiments. We shall here briefly review the empirical studies of turbulent dispersion and the status of Richardson’s theory. In particular, we shall present some new numerical results of our own on stochastic particle advection according to eq.(8) for a turbulent velocity field.

We confine our discussion to just some of the latest studies by experiments Ott and Mann (2000); Bourgoin et al. (2006); Ouellette et al. (2006); Xu et al. (2008) and simulations Ishihara and Kaneda (2002); Biferale et al. (2005); Sawford et al. (2008) at the highest Reynolds numbers. For more complete surveys of the literature, see those papers and also recent review articles Salazar and Collins (2009); Toschi and Bodenschatz (2009). Athough the -law (4) and the stretched-exponential PDF (3) are probably the most famous predictions of Richardson’s theory, even more important for our discussion is the “forgetting” of initial separations. If is the initial particle separation distance and is energy dissipation per mass, then for times much greater than both and should become independent of As we have seen, this is the crucial physical mechanism underlying spontaneous stochasticity. In general, it has proved rather difficult to observe in a completely consistent and convincing way all of these predictions of Richardson’s theory.

Experiments of Ott and Mann at maximum Taylor-scale Reynolds number observed both a law and the Richardson PDF, but varied only by a factor of around the value (for the Kolmogorov dissipation length). Thus, they provide no information on collapse independent of A series of experiments by Bodenschatz and collaborators Bourgoin et al. (2006); Ouellette et al. (2006); Xu et al. (2008) at substantially higher Reynolds numbers up to fail to see a law and instead produce results consistent with Batchelor’s ballistic range. However, their smallest achievable value of the initial separation was only about so that it is arguable that they need longer times and still higher Reynolds numbers. At their smallest value of they did observe Richardson’s PDF (3) and for they see results roughly consistent with Richardson’s predictions for the quantity and also some tendency to collapse independent of (see Xu et al. Xu et al. (2008), Fig. 1).

Numerical simulations of Ishihara and Kaneda Ishihara and Kaneda (2002) at claimed to observe an inertial-range law for in a range between but no tendency whatsoever for collapse at long times independent of However, Biferale et al. Biferale et al. (2005) in a simulation at nearly identical and with in a range of see exactly the opposite: only a slight indication of a -law for their smallest separation but a strong tendency to collapse at long times. They also observe Richardson’s PDF (3) for More recently, Sawford et al. Sawford et al. (2008) have performed simulations with maximum They find the best evidence yet of Richardson’s predictions for the dispersion (see their Fig.4), with a reasonable range for Other values of in the range of do not give a convincing -law but do verify a tendency for collapse at long times.

In view of the incomplete verification of Richardson’s predictions, we decided to undertake our own numerical investigation. Unlike previous works, however, which have nearly all studied deterministic fluid particles with variable initial separations we have instead studied the problem of stochastic particle advection according to the eq.(8). The velocity field is obtained from a pseudospectral numerical simulation of forced, statistically stationary turbulence at The flow data is available online at http://turbulence.pha.jhu.edu and fully documented there and in papers Li et al. (2008); Perlman et al. (2007). The entire flow history for about one large-eddy turnover time is archived at a time resolution suitable for particle-tracking experiments, with spatial and temporal interpolation implemented within the database. This is very convenient for our purposes, since it permits us to study particle dispersion backward in time as well as forward. As we have discussed above, it is backward dispersion that is most relevant for turbulent mixing Sawford (2001).

We have studied two values of the Prandtl number and . We solved (8) using the simplest Euler-Maruyama scheme and also, for convergence analysis, an explicit, 1.5th-order strong scheme (Kloeden and Platen Kloeden and Platen (1992), Section 11.2). We took time discretization which guaranteed that particles moved a fraction of under both turbulent advection and Brownian diffusion at each time-step. The velocity field between stored data-points was interpolated by 6th-order Lagrange polynomials in space and piecewise-cubic Hermite polynomials in time. The results were verified to be converged in both by comparison with the higher-order method and with the Euler scheme at a halved step size. For each initial location we evolved independent particle realizations, giving 523,776 particle pairs, over the whole time range of the database. For forward tracking we started particles at and for backward tracking at (the final time in the database). We then averaged all results over 256 initial locations obtained by choosing 4 independent, uniformly distributed points from each of subcubes of the whole flow domain.

Our results for particle dispersion are given in Fig.1. We present there a log-log plot of (normalized by versus time (normalized by ). For backward dispersion we take to facilitate comparison with the forward-in-time results. In the viscous units that we employ, the early-time diffusive separation (21) becomes This regime is clearly seen for both backward and forward cases and for both and Furthermore, for the unit Prandtl number cases, we see a convincing transition to a -law for This occurs slightly earlier for backward dispersion than for forward. Also, we find that the asymptotic Richardson-Obukhov constant is greater for backward dispersion than for forward, in agreement with earlier results Sawford et al. (2005). An average of the local constants in the -scaling range gives values of for backward dispersion and for forward dispersion. The latter agrees perfectly with a recent theoretical prediction Franzese and Cassiani (2007) and both are generally consistent with previous values Sawford et al. (2008). Pure cube scaling laws with these coefficients are plotted in Fig.1 for comparison with the numerically obtained mean dispersions. The agreement is obviously quite good at long times, especially for the backward case. Even more importantly, the dispersions show a very clear trend to approach the same cubic laws at sufficiently large times or, in viscous units, Our Fig.1 thus provides strong evidence of the “forgetting” of the molecular diffusion time-scale by turbulent Richardson diffusion at long times, which is the essential ingredient of spontaneous stochasticity.

To be completely conclusive, we would need to see collapse of the dispersion curves for different in a range where both show -scaling. We see no clear Richardson for the cases in the time ranges plotted. We cannot continue the time-integration further for two reasons. First, the velocity field from the turbulence archive contains no data for longer times. Second, the rms dispersion distance at the final time has reached a value , just slightly smaller than the velocity integral length-scale for the flow. (This is expected, since is about one large-eddy turnover time.) To integrate further to see a conclusive collapse, we would need a numerical simulation at higher Reynolds numbers and integrated to longer times. However, our confirmation of Richardson diffusion at using stochastic Lagrangian trajectories is comparable to, or even better than, the results of Sawford et al. Sawford et al. (2008) at using deterministic Lagrangian trajectories. (See Fig.4 in that paper). There are two plausible arguments why this should be so. In the first place, using stochastic trajectories, all the particles start at the same point. At the time when (for ), the particles are not randomly placed in the flow, but have already been experiencing relative advection by different-sized eddies at the onset of the inertial range. Thus, they begin to experience Richardson diffusion at that time. However, using the usual technique of seeding the flow with particles at initial separations one would still need to wait some additional time for the initial configuration to be “forgotten”. A second reason is that backward dispersion is faster than forward dispersion, so that the range of scaling occurs even earlier in that case. The technique of stochastic Lagrangian trajectories appears to be promising in the numerical study of Richardson diffusion.

In order to make a completely convincing case that we are observing Richardson diffusion, we have also numerically calculated the PDF of the particle-separations. Our results for are presented in Fig.2, with the normalization As has been previously observed Ott and Mann (2000); Biferale et al. (2005); Salazar and Collins (2009), Richardson’s analytical formula for the long-time PDF of separation distances implies that all the PDF’s at different times will collapse when scaled with In fact, equation (3) is equivalent to

 L3(t)P(r,t)=exp⎡⎢⎣−α(rL(t))2/3+β⎤⎥⎦

with numerical values

 α=(1287/8)1/3≐5.4387

and

 β=ln(335(143)3/2√2π)≐4.7617.

Thus, Richardson’s theory makes a parameter-free prediction that a log-linear plot of versus should give a straight line with slope and -intercept In Fig.2, therefore, we have plotted our PDF’s in this way, at three times all lying in the range of scaling. We have also plotted the straight line predicted by Richardson’s theory. We see that the PDF’s scaled in this way collapse very nicely. Furthermore, except for some deviation at small in the backward dispersion case, they very closely agree with the predictions of Richardson’s theory.

## Iii Stochastic Flux-Freezing

The standard views on flux-freezing in high-conductivity plasmas are inconsistent with the phenomenon of spontaneous stochasticity. It is nearly ubiquitously argued that flux-freezing should hold better as magnetic diffusivity However, high magnetic Reynolds numbers are usually associated also with high kinetic Reynolds numbers. If kinematic viscosity simultaneously with the resistivity and the plasma becomes turbulent, then Lagrangian trajectories will no longer be unique. Which fluid trajectory shall a magnetic field line follow if there are infinitely many such trajectories? This is the paradox of flux-freezing.

As we shall argue below, a form of flux-freezing does survive at small resistivities and viscosities, but in a novel stochastic sense. Before we make this argument, however, we shall first discuss the related subject of flux-freezing properties of resistive hydromagnetics.

### iii.1 Resistive Hydromagnetics

In this subsection we shall discuss magnetic fields that satisfy the resistive induction equation

 ∂tB=\boldmath∇\boldmath×(u\boldmath×B−λ\boldmath∇% \boldmath×B), (26)

with the magnetic diffusivity. It is important to stress that our analysis applies here applies to a very general velocity field It may be incompressible or compressible. It may be externally prescribed or it may satisfy a dynamical equation that contains itself. For example, may be the plasma velocity that obeys the standard magnetohydrodynamic momentum equation or it may be taken to be the electron fluid velocity in Hall magnetohydrodynamics 111We take this opportunity to correct a typo in our previous paper Eyink (2009) where it was incorrectly written that with the mass density and the electron mass. Instead the latter should have been the ion mass.. Our only assumption in this section shall be that is a smooth vector field.

A priori there is no obvious way how to describe magnetic-line motion for non-ideal plasmas. One approach that has been widely employed in discussions of magnetic reconnection Boozer (1990) (or Kuslrud Kulsrud (2005), Section 3.4) is to introduce a “slip velocity” . In that case, one may attempt to introduce an “effective velocity” of the field lines. Unfortunately, this approach is not generally successful because if and only if (or for Ohmic non-ideality). As has been emphasized Priest et al. (2003); Hornig and Priest (2003), no effective-velocity approach is satisfactory for discussions of three-dimensional magnetic reconnection. In fact, those authors show that, even if the non-ideality is spatially localized, there generally exists no smooth velocity field whatsoever such that for a non-ideal plasma.

For magnetic fields that obey (26), however, there is a natural and consistent way to describe line-motion as a process of stochastic advection. Such approaches have already been employed for some time in discussion of kinematic magnetic dynamos, at least for incompressible velocity fields Molchanov et al. (1984); Drummond and Horgan (1986). Recently, we gave a rigorous proof of stochastic flux-conservation properties for nonlinear hydromagnetic models using mathematical methods of stochastic analysis Eyink (2009). We shall present here a more physical demonstration of these results using path-integral methods which, also, extends their validity to compressible fluid models.

To begin, we note that the induction equation (26) may be rewritten as

 ∂tB+(u\boldmath⋅\boldmath∇)B=(B\boldmath⋅\boldmath∇)u−B(\boldmath∇\boldmath⋅u)+λ△B. (27)

In this form it is the same as the scalar advection equation (11), except for the additional two terms on the righthand side. The path-integral formula (14) for the scalar solution may thus be easily adapted to this situation. The solution of (27) with initial condition is given by the “sum-over-histories” formula

 B(x,t) = ∫a(t)=xDaB0(a(t0))\boldmath⋅\boldmathJ(a,t) (29) ×exp(−14λ∫tt0dτ|˙a(τ)−uν(a(τ),τ)|2)

where where is a matrix and is interpreted as a 3-dimensional row vector. satisfies the following ODE along the trajectory :

 ddτ\boldmathJ(a,τ) = \boldmathJ(a,τ)\boldmath∇xu(a(τ),τ) (31) −\boldmathJ(a,τ)(\boldmath∇x\boldmath⋅u)(a(τ),τ),

with initial condition It is easy to check by taking the time-derivative of (29) and using (31) to show that the induction eq.(27) is satisfied. Just as for the scalar problem, the condition on the path-integral trajectories implies that they correspond to solutions of the stochastic equation

 ddτ˜a(τ)=u(˜a,τ)+√2λ˜\boldmathη(τ),˜a(t)=x (32)

integrated backward in time from to

However, the stochastic equation (32) may also be integrated forward in time from to . In that case, the same ensemble of trajectories may be obtained by considering only those particles with initial locations carefully selected to arrive at at time , for a given realization of the white-noise With a slight change of notation, we may characterize this ensemble of time-histories as those which solve

 {ddτ˜x(a,τ)=u(˜x(a,τ),τ)+√2λ˜% \boldmathη(τ),τ>t0˜x(a,t0)=a (33)

such that the inverse map to specifies the starting point by Notice that (33) is a stochastic generalization of the usual equation for a Lagrangian flow map of a particle with initial “label” and that is the “back-to-labels” map. It is easy to show furthermore by applying to (33) that

 ˜\boldmathJ(a,t)≡1det(\boldmath∇a˜x(a,t))% \boldmath∇a˜x(a,t). (34)

solves equation (31) with initial condition It is therefore possible to re-express the path-integral formula (29) as

 B(x,t)=¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯1det(\boldmath∇a˜x(a,t))B0(a)% \boldmath⋅\boldmath∇a˜x(a,t)|˜a(x,t). (35)

The overbar represents the average over realizations of the random white noise process in (33).

We shall call the above result the stochastic Lundquist formula, since it is the stochastic generalization of the standard Lundquist formula Lundquist (1951) (or Kulsrud Kulsrud (2005), section 4.8). It may be cast into a more familiar form by noting that the determinant that appears there can be interpreted as the ratio of initial and final mass densities 222This interpretation requires some caution. The quantity satisfies the stochastic (Stratonovich) equation with initial condition It thus represents the mass density in the random flow with noise realization However, the noise-average is not the physical mass density! In fact, by converting the previous stochastic equation to Ito form and taking the average, it follows that which is not the correct continuity equation.:

 det(\boldmath∇a˜x(a,t))=ρ0(a)˜ρ(˜x(a,t),t)

It follows that the vector field is stochastically “frozen-in” and advected along stochastic Lagrangian trajectories, where is defined to be the quantity under the overbar in (35). Notice, therefore, that the average in (35) is not over the frozen-in field but rather over the magnetic field itself. This is necessary in order to reproduce the Laplacian term in (27), which has the form and not

The use of the stochastic Lunquist formula is illustrated in Fig. 3. The aim is to calculate the magnetic field at spacetime point The first step is to generate an ensemble of stochastic Lagrangian trajectories solving (32) backward in time from at time to random locations at the initial time The path-integral formula (29) sums over all such random time-histories. We show in Fig. 3 (top) three stochastic trajectories generated numerically from the turbulence database together with the starting magnetic field vectors indicated by arrows at the starting locations The next step is to transport each of the field vectors in the usual “frozen-in” fashion along the stochastic Lagrangian trajectories to the final spacetime point The result is an ensemble of field vectors at that point, stretched and rotated by the flow. These are illustrated in Fig. 3 (middle) by the collection of three arrows at obtained by transporting the three initial vectors. In the usual deterministic Lundquist formula there would be just one trajectory and one vector at the final point, which would give the desired magnetic field. Now however as the final step one must average over the ensemble of random vectors in order to obtain the resultant magnetic field This is illustrated by the black arrow in Fig. 3 (bottom). In contrast to the previous transport step which preserved line-topology (in each individual realization), the final averaging step resistively “glues” the transported lines together and changes the magnetic field-line topology.

There is an elegant reformulation of the stochastic Lundquist formula which must be mentioned here, both because of its conceptual simplicity and also because of its potentially greater generality (see next subsection). Consider any smooth, oriented surface at final time Then the formula (35) may integrated in over the surface with respect to the vector area element and the ensemble-average and surface-integration interchanged on the righthand side. Because the expression under the overbar is the one that appears in the usual Lundqust formula, the standard multi-variable calculus manipulations convert this into a surface integral over the surface randomly advected backward in time to the initial time As before is the “back-to-label map” for the stochastic forward flow. The result is the following stochastic Alfvén theorem

 ∫SB(x,t)\boldmath⋅dA(x)=¯∫˜a(S,t)B0(a)\boldmath⋅dA(a),t>t0. (36)

This result generalizes a previous theorem Eyink (2009) to compressible plasmas. Eq.(36) expresses the conservation of magnetic flux on average, as illustrated in Fig. 4. An initial loop boundary of the surface is shown there in black. This is stochastically advected backward in time to give an infinite ensemble of loops at the initial time These are represented by the three colored loops. The ensemble-average of the magnetic flux through the collection of loops at the initial time is equal to the magnetic flux through the loop at the final time

The stochastic Alfvén theorem is an example of what is called a “martingale property” in probability theory. The magnetic flux through each advected loop at the earlier time is unequal to the magnetic flux through at time Nevertheless, the mean flux remains the same. Note that this result implies an irreversibility or an “arrow of time” since it only holds for backward stochastic advection of loops. Backward-in-time is the causal direction, since the magnetic flux at the present must be obtained as an average of past values and not of future values. If we assumed a “forward martingale” property then we would obtain instead the magnetic induction equation (27) with a negative resistivity term Note, in fact, that the stochastic Alfvén theorem (backward in time) is mathematically equivalent to the usual resistive induction equation (26) or (27) Eyink (2009).

### iii.2 High-Reynolds-Number Limit

We now consider the limit of large kinematic and magnetic Reynolds numbers. For simplicity we shall assume that remains fixed as

Consider the Feynman-Kac formula (29). By a naive application of the Laplace method, one would assume that the path-integral collapses to a single deterministic trajectory as with rms fluctuations of order for small but nonzero This is precisely the heuristic estimate of line-slippage made by Kulsrud Kulsrud (2005) which was quoted in the Introduction. This estimate is rigorously correct if the velocity and magnetic fields are assumed to remain smooth in the limit Thus, the heuristic estimate is correct if the plasma flow remains laminar, but this will be the exception rather than the rule at high Reynolds numbers. In a turbulent flow, the behavior will be quite different. As we can see from our Fig. 1 for incompressible hydrodynamic turbulence, the heuristic estimate is only valid for very short times smaller than the resistive time At longer times, the rms slip distance of the field lines instead follows the Richardson law independent of and The quantitative behavior will be different in plasmas with strong magnetic fields, due to the effects of the Lorentz force, as discussed more in section V. However, the qualitative behavior must be the same whenever the advecting velocity is turbulent and spatially rough.

The Feynman-Kac formula (29) is not very well-suited to analyzing the limit of high Reynolds numbers, however, because the velocity-gradients that appear in the definition of the matrix diverge in that limit. Likewise, the gradients of the Lagrangian flow map that appear in the definition (34) of are expected to diverge. The integrated form of flux-conservation, the stochastic Alfvén theorem (36), is more likely to remain meaningful in the limit of infinite Reynolds number. The backward-advected loops at finite values of are expected to approach well-defined curves as which, however, are not rectifiable but fractal Mandelbrot (1976); Sreenivasan and Meneveau (1986); Fung and Vassilicos (1991); Villermaux and Gagne (1994); Kivotides (2003); Nicolleau and Elmaihy (2004). To make mathematical sense of magnetic flux through such fractal loops, we may introduce the vector potential and rewrite the flux through surface as a line-integral around its perimeter We may then further transform by change of variables to a line-integral around the original loop as:

 ∮˜a(C,t)A0(a)\boldmath⋅da=∮CA0(˜a(x,t))\boldmath⋅d˜a(x,t).

The integral on the right may be interpreted as a generalized Stieltjes integral, which is well-defined as long as the map is suitably Hölder continuous Young (1936, 1938).

It may seem from our arguments to this point that the validity at very high Reynolds numbers of the stochastic flux-freezing result (36) is dependent upon the particular stochastic representation of resistive effects employed in eqs.(33) and (35). However, the same “martingale property” can be obtained in the limit by a different argument that employs only the standard Lagrangian flow Eyink (2007). Define the deterministic flow, as usual, by

 {ddτxν(a,τ)=uν(xν(a,τ),τ),τ>t0xν(a,t0)=a (37)

where the superscript is a reminder that the dynamical equation for the advecting velocity field contains a certain viscosity Correspondingly one defines the inverse map Stochasticity can be introduced by assuming small random perturbations of the loop, taking where is a random loop from a well-behaved ensemble 333Note that the addition of loops is pointwise. That is, if two loops are given parametrically by periodic functions then the sum-loop is parameterized by . Thus, can be regarded as the spatial-resolution in determining the precise form of the loop We may then argue that, at least for incompressible flow, the ensemble of loops obtained from stochastic advection in the limit coincides with the ensemble of loops obtained from deterministic advection taking the limits first and then As discussed in the previous section, this is rigorously known to be true for point-particles advected by velocities selected from the Kraichnan white-in-time ensemble E and Vanden-Eijnden (2000); Jan and Raimond (2004). The physical mechanism is just turbulent Richardson diffusion. We therefore conjecture that the same result holds for loops. If this is so, then the double limit

 limϵ→0limν,λ→0∮aν(C+ϵ˜C,t)A0(a)\boldmath⋅da

gives precisely the same ensemble of fluxes with the same distribution as in the previous approach. In that case, (36) must again hold in the limit, or, more precisely,

 ∮CA(x,t)\boldmath⋅dx=limϵ→0limν,λ→0¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯¯∮aν(C+ϵ˜C,t)A0(a)\boldmath⋅da, (38)

where and the overbar now indicates average over the ensemble of loop perturbations This is a nontrivial result because for an individual loop

 ddt∮xν(C,t)Aλ(x,t)% \boldmath⋅dx=−λ∮