A Maximum Entropy Principle for inferring the Distribution of 3D Plasmoids

TL;DR

This paper applies the maximum entropy principle to derive the distribution functions of 3D plasmoids in resistive MHD systems, revealing power-law and exponential behaviors for various physical quantities with minimal parameter dependence.

Contribution

It introduces a maximum entropy framework for modeling 3D plasmoid distributions, providing analytical power-law and exponential forms for key variables, and compares results with observational data.

Findings

Mass, flux, and helicity distributions follow power laws with exponential cutoffs.

Velocity distribution is flat at low speeds and follows a power law at high speeds.

Results are largely independent of free parameters.

Abstract

The Principle of Maximum Entropy, a powerful and general method for inferring the distribution function given a set of constraints, is applied to deduce the overall distribution of 3D plasmoids (flux ropes/tubes) for systems where resistive MHD is applicable and large numbers of plasmoids are produced. The analysis is undertaken for the 3D case, with mass, total flux and velocity serving as the variables of interest, on account of their physical and observational relevance. The distribution functions for the mass, width, total flux and helicity exhibit a power-law behavior with exponents of $-4/3$, $-2$, $-3$ and $-2$ respectively for small values, whilst all of them display an exponential falloff for large values. In contrast, the velocity distribution, as a function of $v = |{\bf v}|$, is shown to be flat for $v \rightarrow 0$, and becomes a power law with an exponent of $-7/3$ for $v \rightarrow \infty$. Most of these results are nearly independent of the free parameters involved in this specific problem. A preliminary comparison of our results with the observational evidence is presented, and some of the ensuing space and astrophysical implications are briefly discussed.

Figures

b of Briggs *et al.* ([2011](#bib.bib143)). However, we wish to emphasize that the above statements must be viewed with due caution for the following reasons: (i) the number of plasmoids observed is typically low, (ii) the interactions between the plasmoids are likely to be relatively insignificant, and (iii) resistive MHD is not applicable to these systems. Hence, it is more instructive to consider the results from solar observations on account of the reasons that were identified in Secs. [I](#S1) and [II](#S2).

Equations

$$\frac{m}{\rho }=V=\pi R^{2}H=\lambda \pi R^3$$

$$M={\int }_{\Omega }\mathcal{F}mdmd\psi d\mathbf{v},$$

$$\Psi ={\int }_{\Omega }\mathcal{F}\psi dmd\psi d\mathbf{v}$$

$$E=\frac{1}{2}{\int }_{\Omega }\mathcal{F}\left[m{v}^2+\mathcal{C}{\psi }^2{\left(\frac{m}{\rho }\right)}^{-1/3}\right]dmd\psi d\mathbf{v},$$

$${\mathcal{E}}_B\sim {\ell }^{-1}\mathcal{H}\propto {\left(\frac{m}{\rho }\right)}^{-1/3}\mathcal{H},$$

$${\mathcal{E}}_B\ge \alpha \mathcal{H},$$

$$K=\frac{1}{2}{\int }_{\Omega }\mathcal{F}{\psi }^{2}dmd\psi d\mathbf{v},$$

$$S=-{\int }_{\Omega }\mathcal{F}\ln \mathcal{F}dmd\psi d\mathbf{v},$$

$$\delta A=0,A:=S-\beta E-\gamma M-\delta \Psi -\epsilon K.$$

$$\mathcal{F}=\exp \left(-\beta \mathcal{E}-\gamma m-\delta \psi -\frac{\epsilon {\psi }^2}{2}-1\right),$$

$$\mathcal{E}=\frac{m{v}^2}{2}+\frac{\mathcal{C}}{2}{\psi }^2{\left(\frac{m}{\rho }\right)}^{-1/3},$$

$${\mathcal{F}}_{m,\psi }=\int \mathcal{F}d\mathbf{v},$$

$${\mathcal{F}}_{m,\psi }$$

$${\mathcal{F}}_m={\int }_{\psi =0}^{\psi =\infty }{\mathcal{F}}_{m,\psi }d\psi ,$$

$${\mathcal{F}}_m$$

$$\sigma =\epsilon +\beta \mathcal{C}{\left(\frac{\rho }{m}\right)}^{1/3}.$$

$$\underset{m\to 0}{\lim }{\mathcal{F}}_m\propto m^{-4/3},$$

$$\underset{m\to \infty }{\lim }{\mathcal{F}}_m\propto m^{-3/2}{e}^{-\gamma m},$$

$${\mathcal{F}}_{m}dm={\mathcal{F}}_{R}dR,$$

$${\mathcal{F}}_R=3\pi \lambda \rho R^2{\mathcal{F}}_m\left(\pi \lambda \rho R^3\right),$$

$$\underset{R\to 0}{\lim }{\mathcal{F}}_R\propto R^{-2},$$

$${\mathcal{F}}_R^{(2D)}\propto R^{-2\Gamma +1},$$

$${\mathcal{F}}_R^{(3D)}\propto R^{-3\Gamma +2}.$$

$$\underset{R\to \infty }{\lim }{\mathcal{F}}_R\propto R^{-5/2}{e}^{-\pi \gamma \lambda \rho R^3}.$$

$${\mathcal{F}}_{\psi }={\int }_{m=0}^{m=\infty }{\mathcal{F}}_{m,\psi }dm.$$

$$\underset{\psi \to 0}{\lim }{\mathcal{F}}_{\psi }\propto {\psi }^{-3},$$

$${\mathcal{F}}_{\psi }$$

$$I\approx \exp \left(-\frac{4\gamma }{y_{\star }^2}\right){\int }_{y=0}^{y=\infty }\exp \left[-\frac{8\gamma }{3}{y}_{\star }^{-4}{\left(y-y_{\star }\right)}^2\right]dy,$$

$$y_{\star }={\left(\frac{\beta \mathcal{C}{\rho }^{1/3}{\psi }^2}{6\gamma }\right)}^{-3/8}.$$

$$I\approx \sqrt{\frac{3\pi }{8\gamma }}{y}_{\star }^2\exp \left(-\frac{4\gamma }{y_{\star }^2}\right),$$

Full text

I Introduction

It is now widely acknowledged that magnetic reconnection, which entails the conversion of magnetic energy into charged particle energy, drives many astrophysical phenomena in the Universe. Some well-known examples include magnetospheric substorms, gamma ray bursts, stellar and solar flares (Tajima and Shibata, 1997; Kulsrud, 2005; Birn and Priest, 2007; Zweibel and Yamada, 2009; Ji and Daughton, 2011; Kumar and Zhang, 2015). Under normal circumstances, this fast conversion would not be possible for highly elongated current sheets, as the conventional estimates for the reconnection rate are too low. Fortunately, it has been known since the 1980s that such current sheets are subject to a violent instability that was dubbed the plasmoid instability later (Bhattacharjee et al., 2009), as it results in the formation of plasmoids (Biskamp, 1986, 2000).

The plasmoid instability breaks up the elongated current sheets, and leads to very high reconnection rates that are nearly independent of the Lundquist number $S$ in the nonlinear regime (Bhattacharjee et al., 2009; Cassak et al., 2009; Daughton et al., 2009; Huang and Bhattacharjee, 2010; Loureiro et al., 2012; Huang and Bhattacharjee, 2013; Comisso et al., 2015). The emergence of this fast reconnection process, mediated by the plasmoid instability, has been shown to have several important consequences. Plasmoid reconnection can play an important role in regulating solar flares (Kliem et al., 2000; Shibata and Tanuma, 2001), coronal mass ejections (Lin et al., 2015), chromospheric jets (Ni et al., 2015), particle acceleration (Øieroset et al., 2002; Drake et al., 2006; Chen et al., 2008; Birn et al., 2012), pulsar winds (Sironi and Spitkovsky, 2014; Guo et al., 2014), blazar emissions and jets (Sironi et al., 2015; Petropoulou et al., 2016), as well as a wide range of laboratory and fusion experiments (Ebrahimi and Raman, 2015; Jara-Almonte et al., 2016; Ebrahimi, 2016).

In plasmoid-mediated reconnection, a great deal of analytic and numerical work has been undertaken for the linear regime. In Sweet-Parker current sheets, scaling relations for the growth rate and number of plasmoids produced were predicted by Tajima and Shibata (1997) and Loureiro et al. (2007) as a function of $S$, and were generalized to encompass arbitrary magnetic Prandtl numbers in Comisso and Grasso (2016). However, in dynamically evolving current sheets, the final Sweet-Parker aspect ratio cannot be achieved, as they disrupt prior to this stage (Biskamp, 1982; Comisso et al., 2016; Tenerani et al., 2016; Comisso et al., 2017; Huang et al., 2017). Recently, Comisso et al. (2016, 2017) utilized a principle of least time to formulate a general theory for the onset and development of the linear phase of the plasmoid instability. They demonstrated that the scalings were no longer true power laws, and depended on other parameters, namely the initial noise inherent in the system, the characteristic rate of current sheet evolution, and the nature of the thinning process.

From the observational standpoint, it is not very feasible to verify the linear stage of the plasmoid instability. As a result, a great deal of interest has been centred on the nonlinear phase of this instability. One of the most notable amongst them is the distribution of plasmoids as a function of their intrinsic properties, such as their flux and width. From the theoretical standpoint, there have been three major works in this subject, albeit only for plasmoids emergent from a 2D reconnection process. Uzdensky et al. (2010) and Huang and Bhattacharjee (2012) proposed phenomenological kinetic models for the flux distribution function, and arrived at scaling laws with exponents of $-2$ and $-1$ respectively. On the other hand, Fermo et al. (2010) prescribed rules for the generation, growth and merging of plasmoids, and arrived at an integro-differential equation for the distribution function that depended on the flux and the area.

There are, arguably, two important factors that need to be taken into account when paving the way for a more complete theory. Firstly, the distribution function is obviously statistical, as a result of which it would clearly be ideal to rely on a more statistical-mechanical approach based on first principles. Secondly, many of the theoretical and numerical studies in the area have relied on heuristic models with simplified geometries to arrive at their results. Hence, except in certain circumstances, the resultant parameters and the distribution function may not constitute true observables.

It is, therefore, the goal of this work to construct a 3D plasmoid distribution function in the highly nonlinear regime, characterized by a large number of plasmoids, that depends on three genuinely physical observables, namely the plasmoid mass, total flux and velocity. During the course of the paper, we shall use the nomenclature “plasmoids” to describe magnetic islands in three dimensions, i.e. flux ropes/tubes (Moldwin and Hughes, 1991). In 3D geometry, the magnetic field lines could become stochastic, and the flux surfaces may not be well-defined (Borgogno et al., 2005; Huang and Bhattacharjee, 2016). Yet, one can still identify coherent structures which can be treated on the same footing as flux ropes. Hence, in this paper, we shall operate with the implicit assumption that well-defined coherent structures (or flux ropes) are existent.

For the purpose of this work, we adopt the standard resistive magnetohydrodynamics (resistive MHD) as the underlying physical model. Our results are therefore relatively valid for astrophysical plasmas where resistive MHD constitutes a reasonable description of the plasma dynamics. Hence, for plasma environments such as planetary magnetospheres (Burch et al., 2016) where kinetic effects become important (Ji and Daughton, 2011), our results are not likely to be accurate. In contrast, resistive MHD constitutes a useful model for certain astrophysical systems like solar and stellar coronae (Priest, 2014) and the interstellar medium (Zweibel, 1999), and therefore we anticipate that our findings have a higher degree of applicability to these systems.

Our paper is centered around the Principle of Maximum Entropy (Jaynes, 2003), which has been successfully used in many fields of science. We shall derive the overall 3D distribution function, examine the various asymptotic limits, and eventually discuss its potential connections with astrophysics. The outline of the paper is as follows. A short introduction to the Principle of Maximum Entropy is offered in Sec. II. The physical parameters of interest, and the relevant constraints are introduced in Sec. III, and the overall plasmoid distribution function is derived. In Sec. IV, the single-variable distribution functions for the mass, width, total flux, helicity and velocity are obtained. Sec. V contains a discussion of the relevance of our work in space and astrophysical plasmas, and we finally summarize our results in Sec. VI.

II On the Principle of Maximum Entropy

In this Section, we shall briefly elucidate some of the salient features of the Principle of Maximum Entropy (MaxEnt), which has been used quite extensively in many areas of astronomy and astrophysics (Lynden-Bell, 1967; Narayan and Nityananda, 1986; Binney and Tremaine, 1987). We also discuss its relevance in the context of determining the statistical distribution of 3D plasmoids.

II.1 What is MaxEnt and where is it used?

In equilibrium statistical mechanics, it is known that the micro-, macro- and grand canonical distribution functions can be derived through a common principle. The central idea relies upon the maximization of the entropy $S=-\sum_{i}p_{i}\ln p_{i}$ with respect to the probabilities $p_{i}$, subject to holding the appropriate constraints on the system fixed (Tolman, 1938; Haar, 1955; Balescu, 1975; Callen, 1985). By drawing upon the classical ideas of Boltzmann, Gibbs and Shannon, E.T. Jaynes explored the connections between information theory and statistical mechanics and argued that statistical mechanics could be interpreted as a method for inferring the probability distributions based on the (limited) data available (Jaynes, 1957a, b, 2003).

It is instructive to consider Jaynes’ own words (Jaynes, 1957a) in this regard: “In the resulting “subjective statistical mechanics,” the usual rules are thus justified independently of any physical argument, and in particular independently of experimental verification; whether or not the results agree with experiment, they still represent the best estimates that could have been made on the basis of the information available.”

MaxEnt has enjoyed a considerable degree of success ranging from fields as disparate as ecology (Harte, 2011) and climate science (Paltridge, 1979; Ozawa et al., 2003) to gravitational dynamics (Lynden-Bell, 1967; Padmanabhan, 1990; Spergel and Hernquist, 1992; Nakamura, 2000; Herpich et al., 2017) and plasma physics (Hameiri and Bhattacharjee, 1987; Yoshida and Mahajan, 2008; Dewar et al., 2008). We shall not go into further details of MaxEnt as there exist several excellent reviews and summaries of the subject (Shore and Johnson, 1980; Grandy, 1980; Jaynes, 1982, 2003; Dewar, 2003, 2005; Martyushev and Seleznev, 2006; Pressé et al., 2013), and detailed analyses of the arguments for and against MaxEnt can also be found therein.

II.2 MaxEnt and the distribution of plasmoids

In dealing with the plasmoid distribution function, it is clear that the system cannot necessarily be regarded as being in equilibrium throughout, that the plasmoids themselves are complex entities (unlike, for e.g., ideal gas molecules), etc. Fortunately, as per the discussion in Sec. II.1, many of these points are rendered moot.

Firstly, we recall that MaxEnt must be understood as a fairly general statement regarding statistical inference 111This fact served as one of the primary motivations behind our use of the word “infer” in the title., and its generalization, Maximum Caliber (MaxCal), is valid when dealing with non-equilibrium systems (Jaynes, 1980; Pressé et al., 2013; Hazoglou et al., 2015). Since MaxEnt is applicable to complex systems close to equilibrium (Haken, 2006), it could potentially be used for deriving the plasmoid distribution function.

A comment regarding the choice of the entropy functional is in order. At times, there is a tendency to associate the Boltzmann-Gibbs-Shannon (BGS) entropy (Jaynes, 1965) only with the presence of an associated kinetic equation 222Although we shall use the nomenclature BGS entropy henceforth, it is important to recognize that there exist subtle differences between the Boltzmann, Gibbs and Shannon formulations (Jaynes, 1965).. In actuality, the BGS entropy has been “postulated” in a more general and axiomatic manner (Penrose, 1979; Callen, 1985). Moreover, the BGS entropy is endowed with the necessary mathematical and physical requirements, such as invariance, concavity, and additivity to name a few (Wehrl, 1978; Mackey, 1989).

We would like to mention some important caveats at this stage. Since we are dealing with a statistical treatment, there is an implicit assumption that there exist a sufficiently high number of plasmoids that interact amongst each other; as the theory is “thermodynamic” in a certain sense, a large number of interacting entities is ostensibly necessary. In the magnetospheres of planets (and satellites), plasmoids are produced in low numbers, as noted in Sec. V, implying that the MaxEnt formulation may not yield accurate results. On the other hand, since a large number of plasmoids are produced in the solar corona (Janvier, 2017), and given the applicability of resistive MHD in this regime, it could be more reasonable to employ MaxEnt in this context. The presence of a large number of plasmoids is also important from the standpoint of observations (or simulations), as otherwise proper statistics cannot be deduced.

Hence, bearing the above properties and caveats in mind, we are now in a position to adopt the standard approach of Jaynes (1957a) and construct a variational principle wherein the BGS entropy functional can be extremized whilst holding the necessary constraints fixed.

III The derivation of the Plasmoid Distribution Function

Here, we shall delineate some general properties of the 3D plasmoid distribution function $\mathcal{F}$, motivate the appropriate constraints, set up the variational principle, and thereby arrive at the plasmoid distribution function.

III.1 Properties of the plasmoid distribution function

Let us begin by considering the variables that the plasmoid distribution function $\mathcal{F}$ depends upon. We begin by noting that our “system” is now taken to comprise of a sufficiently large collection of plasmoids, and not the entire current sheet as a whole.

There are many properties that can characterize a particular plasmoid. These include:

• •

The mass $m$ enclosed within a given plasmoid.

• •

The density $\rho$ of a particular plasmoid. Notice that a knowledge of $m$ and $\rho$ suffices to determine the volume $V$ of the plasmoid.

• •

The total flux $\psi$ enclosed in the plasmoid. The use of the word “total” is deliberate, as we are dealing with the 3D setup. It also constitutes an observationally relevant physical variable, particularly in the solar context (Wang et al., 2005; Démoulin, 2008; Vieira and Solanki, 2010).

• •

The velocity ${\bf v}$ associated with the motion of the plasmoid.

• •

The position ${\bf x}$ describing the spatial location of the plasmoid.

In addition to the above quantities, the plasmoid distribution function can be expected to evolve over time as well. $\mathcal{F}$ could also depend on other properties, but we believe that we have listed all of the most salient candidates. Clearly, retaining all of these variables would make our procedure very onerous and rob it of physical transparency. Hence, we shall introduce the following simplifications:

• •

We are ultimately interested, by means of MaxEnt, to determine the final $\mathcal{F}$ that the system relaxes to. In this case, an explicit dependence on time is not expected to be present.

• •

We shall assume that the plasma inside the plasmoids is incompressible, i.e. we shall take $\rho=\mathrm{const}$, and hence it does not enter our model.

• •

Handling a spatial dependence is more difficult, and we shall consider the cases where $\mathcal{F}$ does not depend on ${\bf x}$, which may be valid when the system has a weak spatial dependence 333Alternatively, perhaps, one could interpret $\mathcal{F}$ as the coarse-grained distribution function where the ${\bf x}$ dependence has already been eliminated..

Thus, based on the above set of postulates, the distribution function has the form $\mathcal{F}\left({\bf v},\psi,m\right)$. As with the previous theories, the distribution function $\mathcal{F}$ and its constituent variables are continuous.

Although we know that $m\propto V$, we still do not know how the volume depends on the half-width $R$ and on the area. To resolve this matter, we must introduce an additional simplification. In 3D, it is natural to speak of flux ropes or flux tubes (Russell et al., 1990; Birn and Priest, 2007), and our choice of the word “plasmoids” should actually be taken to refer to such entities. These structures can be envisioned as cylinders with radius $R$ and height $H$, and we shall introduce the assumption $\lambda=H/R=\mathrm{const}$. Our ansatz is supported by recent numerical simulations undertaken by Huang and Bhattacharjee (2016) with a finite guide field, wherein cylindrical-shaped flux ropes were approximately characterized by a scale-independent anisotropy.

This does not significantly reduce the scope of our analysis, since one can easily make $\lambda$ a function of $R$. The current simplification amounts to the statement that all plasmoids have a self-similar structure, i.e. their volumes and areas may be different, but they have the same value of height to radius. Such ansatzen have already been used quite often in other areas of physics, where the basic “building blocks” are extended objects, such as liquid crystals and polymers (Reiss et al., 1959; Stephen and Straley, 1974; Vroege and Lekkerkerker, 1992). We finally end up with

[TABLE]

III.2 The choice of constraints

Before exploring the constraints that need to be included in the variational principle, a brief discussion of the dynamics is necessary (Bárta et al., 2008). There are three primary processes involved:

• •

Process I: New plasmoids are “born”, and they “grow” by accreting mass and flux, growing bigger in size.

• •

Process II: Two plasmoids can “collide” and thereby “merge”. Some of the merging rules have already been explored (in the 2D setting) by Fermo et al. (2010) and Huang and Bhattacharjee (2012).

• •

Process III: Plasmoids can be effectively said to have “died” when they have been ejected out of the current sheet.

Let us now take a look at some of the potential constraints that are feasible.

Total number: Although Processes I and III create and destroy plasmoids, there is no guarantee that they balance each other. In addition, Process II does not conserve the number of plasmoids, since the total reduces by unity. Therefore, taken collectively, there is no reason to believe that the total number $N=\int_{\Omega}\mathcal{F}\,dm\,d\psi\,d{\bf v}$ is conserved. Even if this constraint were incorporated, the overall effect would be qualitatively unimportant since the distribution function would just end up being multiplied by a constant factor.

Total mass: It is clear that the mass must be conserved during Process II. However, there are two other processes, one of which depletes the mass and the other that replenishes it. If we assume that, even though each of the processes are dynamical, the whole system can exist in a state of detailed balance, enabling the total mass to be conserved. Thus, the total mass is

[TABLE]

where $\Omega$ is the volume of the total phase space.

Total momentum: It is tempting to argue that a similar set of arguments should hold true for momentum conservation. In the 2D case, at least, one does not expect the total momentum (along each direction) to be conserved, owing to the fact that there exist clear inflow and outflow directions. Some processes can occur only along the former, and others along the latter alone. Thus, there is a fundamental asymmetry with respect to the axes for all the three processes, implying that momentum conservation along all the directions ought not be possible.

Total flux: When we consider flux, it is known that flux is not conserved during the merging of two plasmoids, viz. in Process II (Huang and Bhattacharjee, 2012). To compensate for this absence of conservation, Processes I and II would need to self-adjust in a very precise manner for flux conservation to hold true. If this condition were indeed valid, the total flux

[TABLE]

would be conserved. In (3), the distinction between the total flux $\Psi$ of the complete system and the total flux $\psi$ of a given plasmoid must be duly borne in mind; the latter has the standard definition, $\psi=\int{\bf B}\cdot d{\bf S}$.

Total energy: It is well known that the conservation of total energy is closely linked with time-translation symmetry (Landau and Lifshitz, 1960; Arnold, 1989). Taken individually, none of the processes can conserve the energy, but together, there is a possibility for doing so. If energy conservation is indeed valid for our system, then

[TABLE]

and the first term is clearly the kinetic energy, whilst the second term is the magnetic energy that is discussed further below. An important point worth noting here is that we have used the MHD energy, since we are interested in systems where collisionless effects are not important. For instance, if we consider a domain where extended MHD is applicable, a third term must be included in the above formula that will be proportional to the square of the electron skin depth in normalized units (Lüst, 1959; Keramidas Charidakos et al., 2014; Comisso and Asenjo, 2014; Lingam et al., 2015). Similarly, when Finite Larmor Radius (FLR) are important, they can be incorporated into the above formula by adding gyroaveraged contributions (Snyder and Hammett, 2001; Scott, 2010; Comisso et al., 2013). We have neglected thermal energy in our analysis, as this introduces a further level of complexity by means of an extra variable. Moreover, many plasmas are characterized by $\beta\ll 1$ (Goedbloed and Poedts, 2004; Kulsrud, 2005), which is consistent with the neglect of the thermal energy. Furthermore, even with the introduction of a finite thermal energy that is proportional to the temperature, we have verified that the overall scalings remain unaltered.

In the energy, $\mathcal{C}$ is a dimensionless constant that depends upon $\lambda$ and other numerical factors (such as $\pi$). We are not interested in the exact expression for $\mathcal{C}$ as we seek instead to evaluate the algebraic dependence of $\mathcal{F}$. We have also assumed that our system has a weak guide field. Our procedure can also be generalized to include an arbitrary guide field, and will result in an extra term that is linearly proportional to $m$ in (4). It can therefore be easily shown that the scaling relations for the mass (also the volume and radius) distribution function derived in Sec. IV.1 for the low-mass limit remain unchanged, as well as the overall behavior. Moreover, the flux and velocity distributions in both asymptotic regimes are also unaffected. Hence, the comparison of our theoretical results with the observational data undertaken in Sec. V is valid even when a strong guide field is present.

It is easy to verify that the second term represents the magnetic energy by carrying out a simple dimensional analysis. A second method for arriving at this expression is via the magnetic helicity. Suppose that the magnetic energy and helicity of a plasmoid are represented by $\mathcal{E}_{B}$ and $\mathcal{H}$ respectively. From dimensional considerations, it is once again evident that

[TABLE]

where the second equality follows by observing that $\ell$, the characteristic length scale, must correspond to either $R$ or $H$. Either of these two quantities can be determined in terms of $m$ via (1). Thus, if we show that $\mathcal{H}\propto\psi^{2}$, we will recover the magnetic energy in (4). This is explored below, when the total helicity is introduced as the next constraint.

The first equality in (5) can also be justified on more rigorous grounds than dimensional analysis. We have already stated that our plasmoids are now modelled as flux tubes. For some special configurations, wherein the sum of the twist and the writhe is zero, it is known that

[TABLE]

where $\alpha\propto V^{-1/3}$ (Ricca, 2008), and the $m$-dependence can be determined via (1). The above inequality, which also goes by the name of the ‘realizability condition’, was first proven by Arnold in the 1970s (Arnold, 1989; Arnold and Khesin, 1998).

The above approach relied upon topological fluid mechanics, and the reader may consult Lingam et al. (2016) for a recent overview of this subject for a diverse array of plasma fluid models. There is also a more direct means of arriving at $\mathcal{E}_{B}$ in (4) via this area of research. It was shown in Moffatt (1990) that the minimum value of the magnetic energy attained through a relaxation process is $\mu\,\psi^{2}V^{-1/3}$, where $\mu$ is an invariant that solely depends on the magnetic topology and can be computed exactly in certain scenarios (Chui and Moffatt, 1995). Clearly, upon using (1), this expression has the same form as the magnetic energy delineated in (4).

Total helicity: Most of the interest in the total magnetic helicity stems from its unique topological properties (Berger and Field, 1984; Berger, 1999). There is considerable experimental and theoretical evidence that supports the robustness of helicity conservation even in the presence of magnetic reconnection (Berger, 1984; Pfister and Gekelman, 1991; Kleckner and Irvine, 2013; Scheeler et al., 2014). Hence, it is quite plausible that the total magnetic helicity of the system is also conserved. Then, we end up with

[TABLE]

and the factor of $\psi^{2}$ could be justified, as earlier, through dimensional analysis. It can also be obtained on mathematical grounds - it was shown by Moffatt and Ricca (1992) that the helicity, modulo topological invariants (such as the Gauss and Calugareanu-White linking numbers), was quadratic in the flux; we also refer the reader to Moffatt (1990); Chui and Moffatt (1995); Berger (1999); Ricca (2008). The relation $\mathcal{H}\propto\psi^{2}$ is also backed by a fairly high degree of empirical evidence from solar observational data (Démoulin, 2007, 2008; Balogh et al., 2014). In reality, (7) would also involve an extra numerical factor, which can be easily eliminated by simply rescaling the definition of $K$.

III.3 The formulation of the MaxEnt variational principle

The BGS entropy functional for our model corresponds to

[TABLE]

with $\mathcal{F}$ denoting the plasmoid distribution function. We are now in a position to write down our variational principle as all the relevant pieces have been assembled. Let us commence our analysis with an important observation. Of all the constraints discussed in Sec. III.2, only the energy exhibits a complex dependence on all three variables, namely $m$, $\psi$ and ${\bf v}$. In contrast, all of the other functionals are dependent on just one variable.

To summarize, the mass (2), flux (3), energy (4) and helicity (7) have been chosen as our constraints. Thus, our variational principle is given by

[TABLE]

We have chosen the notation $A$ and $\beta$ for the target functional and the Lagrange multiplier preceding the energy. Our choices are deliberate since these quantities can be rightfully seen as the generalizations of the Helmholtz free energy and the inverse temperature (from thermodynamics) respectively (Tolman, 1938; Callen, 1985).

Upon variation, we arrive at

[TABLE]

and if we take $\{\gamma,\delta,\varepsilon\}=0$, the similarity with the canonical ensemble distribution function (Balescu, 1975) is not very surprising, given that only the energy constraint has been included. Here, it must be noted that

[TABLE]

is the total plasmoid energy in the MHD regime.

IV Analysis of the plasmoid distribution function

We are now in a position to start analyzing the properties of the plasmoid distribution function (10). We shall present the derivation of the distribution functions for the mass, width, total flux, helicity and (3D) velocity in this Section.

IV.1 The mass and width distributions

We commence our derivation by introducing

[TABLE]

which represents the joint distribution function of the plasmoids as a function of $\psi$ and $m$. An inspection of (10) reveals that it is dependent only on the modulus of ${\bf v}$. Hence, we are free to use the transformation $d{\bf v}=v^{2}\,\sin\theta\,dv\,d\theta\,d\phi$. In Sec. III.2, we stated that we shall consider a system where the guide field is absent, implying that one cannot define motion parallel and perpendicular to the guide field. Hence, it seems reasonable to suppose that the conventional ranges for $\theta$ and $\phi$, namely $\left[0,\pi\right]$ and $\left[0,2\pi\right]$ respectively, are valid. This leads us to $d{\bf v}=4\pi v^{2}\,dv$, and substituting it into (12) results in

[TABLE]

where the limits of integration are from $v=0$ to $v=\infty$. Next, one can analytically calculate the mass distribution as follows

[TABLE]

and the integration limits are taken to be positive as per the standard convention. The final answer is given by

[TABLE]

where we have introduced the auxiliary function

[TABLE]

We are now in a position to extract the two limits of the mass distribution. We find that

[TABLE]

but it is necessary to recognize that the immediate higher-order contribution is $m^{-7/6}$, which is smaller than the leading-order contribution (17) only by a factor of $m^{1/6}$ which exhibits a weak dependence on $m$. Hence, the $-4/3$ exponent will not prevail for long, as it will be subsumed by the higher-order contributions at a fairly early stage. In the other limit, we have

[TABLE]

and the exponential cutoff is clearly dominant in this region, subsuming the power-law behavior. We have plotted (15) in Fig. 1 - the power-law limit for small $m$, namely (17), has also been included. It is clear that, for larger values of $m$, the exponential falloff becomes more important. The overall structure of Fig. 1 is quite similar to the 2D numerical results of Lynch et al. (2016) (see their Fig. 10), although our (3D) power-law slope is slightly steeper.

Under a change of variables, the probability contained in a differential area is known to be constant (Kolmogorov, 1950). Let us denote the width distribution of the plasmoids by $\mathcal{F}_{R}$. Using this principle, we find

[TABLE]

and we can now make use of (1). After simplification, we arrive at

[TABLE]

where $\mathcal{F}_{m}\left(\pi\lambda\rho R^{3}\right)$ implies that the function $\mathcal{F}_{m}$ is evaluated at $m=\pi\lambda\rho R^{3}$. For the sake of completeness, we shall present the two limits for the width distribution below. We begin with

[TABLE]

and it is important to reiterate that our analysis is for the 3D case. A similar result was also proposed in the 2D case by Uzdensky et al. (2010) via their phenomenological arguments and numerical results. We also refer the reader to the 2D simulations by Shen et al. (2013) and Sironi et al. (2016) where the same scaling law was obtained. The slopes (for the width distribution) can be evaluated in the 3D case when compared to the 2D case given the same mass distribution. Of course, one expects the 2D and 3D cases to have different mass distributions as well, but it is still instructive to consider the above argument. Let us assume that $\mathcal{F}_{m}\propto m^{-\Gamma}$ in the limit $m\rightarrow 0$, where $\alpha$ is positive. Then, in the 2D case, we find that

[TABLE]

since $m\propto R^{2}$ in the 2D case. For the 3D case, we use $m\propto R^{3}$ that leads us to

[TABLE]

From (22) and (23), we see that the 3D width distribution is steeper than its 2D counterpart only when $\Gamma>1$ is valid. It is clear from (17) that this condition is indeed met since $\Gamma=4/3$, leading to a steeper 2D distribution. A consequence of our analysis is that the width distribution of $-2$ would not be possible in 2D, if the mass distribution were still the same as (17). The other limit for the width distribution is found from (18), and is of the form

[TABLE]

We have plotted (20) in Fig. 2, and the power-law limit for small $R$, namely (21), is also clearly illustrated. In a manner similar to the mass distribution, the exponential cutoff dominates when larger values of $R$ are considered.

We conclude this section by observing that a pure power law distribution for (15) follows if only the energy constraint is retained, i.e. if we let $\{\gamma,\delta,\varepsilon\}=0$. For this simplified case, we find that the $-4/3$ power law holds true for all values of $m$. Hence, the result coincides with the more general distribution function in the limit $m\rightarrow 0$ as seen from (17).

IV.2 The flux and helicity distributions

In order to obtain the total flux distribution function, we must use (13), thereby obtaining

[TABLE]

The integral can be performed analytically, but the result is tremendously complex since (13) has a fairly complicated dependence on $m$. Hence, we shall not present it here, but we have duly carried out the series expansion for $\psi\rightarrow 0$. The leading order contribution is

[TABLE]

indicating that the slope of the distribution function is steeper than the 2D analyses conducted by Uzdensky et al. (2010) and Huang and Bhattacharjee (2012), who arrived at slopes of $-2$ and $-1$ respectively. But, we wish to point out two notable differences: (i) our work assumes a 3D geometry, and (ii) the flux $\psi$ that we have used is the total flux and is therefore not the same variable as the one that they had used. The fact that the slope for the 3D case is steeper than its 2D counterpart is quite consistent with the general observations for the width distribution that were presented in Sec. IV.1.

There is also another means of arriving at (26), which we describe below. If we introduce $y=m^{-1/2}$, we have

[TABLE]

If we let $\psi\rightarrow 0$ in the above expression, by dropping the first term, it is evident that the integral does not converge. The transformation $\tilde{y}=\psi\,y^{1/3}$ yields the factor of $\psi^{-3}$ that can be pulled outside the integral, and the latter converges to a finite value upon setting $\psi=0$. Consequently, we have illustrated how the limit (26) can be obtained.

In Sec. IV.1, we proved that the inclusion of the energy constraint alone or the general case led to the same result for $\mathcal{F}_{m}$ in the limit $m\rightarrow 0$. Hence, it is instructive to study the two cases for $\psi\rightarrow 0$ as well. The general case leads to (26), whilst the case with only the energy constraint is easily evaluated by setting $\{\gamma,\delta,\varepsilon\}=0$. The integral can be evaluated analytically by using the fact that $\int_{r=0}^{r=\infty}\exp\left(-Ar^{2/3}\right)\,dr\propto A^{-3/2}$. We end up with (26), as before, thereby confirming that the two cases yield identical results when $\psi\rightarrow 0$.

Let us now proceed to evaluate the other limit. We shall formally consider the case where the coefficients preceding the two terms inside the integral (IV.2) are very “large”. In general, note that this condition will be fully satisfied when $\psi\rightarrow\infty$. With this assumption, the first term in the expression is unity for small $y$ and rapidly falls off at large $y$, whilst the opposite behavior is seen for the second term.

This enables us to construct a heuristic approximation along the lines described in de Bruijn (1958); Bleistein and Handelsman (1975); Hinch (1991) for Laplace’s method 444Alternatively, we could introduce the change-of-variables $\bar{y}=\psi^{3/4}\,{y}$, thereby enabling us to directly apply the conventional form of Laplace’s method. The asymptotic expansion of the integral for $\psi\rightarrow\infty$ follows as a natural consequence.. The trick is to Taylor expand the function (inside the exponential) around the maximum, and carry out the integration. Denoting the integral by $I$, we find that it simplifies to

[TABLE]

where we have introduced the variable

[TABLE]

(28) can now be further simplified to yield

[TABLE]

where we have used the fact that $\psi\rightarrow\infty$ to eliminate an error function $\left(\mathrm{Erf}\right)$ that approaches unity. Hence, we can conclude that

[TABLE]

where we have chosen to deliberately retain only the highest power of $\psi$ in the exponential factor arising from (IV.2) and (30). In general, there are also terms involving $\psi$ and $\psi^{3/2}$ inside the exponential function.

We have plotted the numerical solution of (IV.2) in Fig. 3. Akin to the mass (or width) distribution, the existence of an exponential falloff at large values of $\psi$ is quite apparent. The veracity of the power law analytical limit (26) for small $\psi$ has also been verified through comparison with the numerical solution of (IV.2).

In Sec. III.2, it was pointed out several times that the total flux and the helicity are proportional to one another. The constant of proportionality is taken in some instances to be $0.5$ (Balogh et al., 2014), and this was also adopted in (7). We use the property of probability conservation in a differential area to write

[TABLE]

where $\mathcal{H}=\psi^{2}/2$ and $\mathcal{F}_{\mathcal{H}}$ denotes the helicity distribution of plasmoids. Thus, we find that

[TABLE]

and $\mathcal{F}_{\psi}$ is evaluated at $\psi=\sqrt{2\mathcal{H}}$. The two limiting cases of $\mathcal{F}_{\mathcal{H}}$ are therefore found to be

[TABLE]

[TABLE]

and we have retained only the highest power in the exponential factor in the above expression.

We have plotted the numerical solution of the helicity distribution function in Fig. 4. This has been obtained by using the numerical solution for $\mathcal{F}_{\psi}$ in conjunction with (33). The limiting behavior for small $\mathcal{H}$, which is seen to be a power law from the figure, also matches the analytical limit obtained in (34). For higher values of $\mathcal{H}$, the exponential cutoff is increasingly dominant.

IV.3 The velocity distribution

To compute the velocity distribution, we first introduce

[TABLE]

which leads to the analytical expression

[TABLE]

where $\sigma$ was defined in (16). The velocity distribution follows from

[TABLE]

The integral cannot be performed analytically, and even the asymptotic expansions are somewhat intricate. Hence, we shall resort to an alternate strategy to deduce the two limits. From (16), we see that $0<\sigma^{-1}<1/\varepsilon$, and the two limits are attained at $m=0$ and $m=\infty$ respectively. This implies that

[TABLE]

and, similarly, we also find

[TABLE]

The above relations indicate that the two functions are bounded by finite values from above and below and their product is also expected to have finite lower and upper bounds. Hence, we do not expect the functional dependence of (38) on ${\bf v}$ to be significantly altered by these functions, especially for certain choices of $\delta$ and $\varepsilon$. In contrast, $\sigma^{-1}$ and $\exp\left(-\zeta m\right)$ do not have finite lower bounds since they vanish for $m\rightarrow 0$ and $m\rightarrow\infty$ respectively, where

[TABLE]

Hence, we shall center our attention on the integral

[TABLE]

which will help us elucidate the general behavior of $\mathcal{F}_{\bf v}$. Our approach is motivated partly by classical asymptotic analysis (de Bruijn, 1958; Bleistein and Handelsman, 1975; Bender and Orszag, 1978), but it has a greater range of validity since there are no small or large parameters involved. We have verified that the full numerical solution of (38) is well represented by the integral (42), which has the analytical solution (43) discussed below. A comparison of these two plots is provided in Fig. 5, and the accuracy of (42) stands out. Observe that the two curves are not exactly the same, since we had neglected some constant numerical factors such as $e^{-1}$ whilst defining (42).

After a careful evaluation of the integral, the result is analytical, expressible in terms of the Meijer G-function (Bateman et al., 1953; Abramowitz and Stegun, 1965), given by

[TABLE]

If we take the limit ${\bf v}\rightarrow 0$, it follows that $\zeta\rightarrow\gamma$. This leads us to

[TABLE]

In the other case, we consider ${\bf v}\rightarrow\infty$, which is equivalent to letting $\zeta\propto v^{2}$ to leading order, evident from the definition (41). The asymptotic expansion of (43) for large values of $\zeta$ enables us to conclude that $\mathcal{I}\propto\zeta^{-7/6}$. Rewriting this expression in terms of $v$, we end up with

[TABLE]

Hence, we have successfully arrived at the two asymptotic limits for $\mathcal{F}_{\bf v}$, and they are both power laws. Such behavior is rather unusual since the other distribution functions exhibited an exponential falloff for larger values of their respective variables.

Suppose that we consider the case where all the constraints apart from the energy are dropped, which amounts to setting $\{\gamma,\delta,\varepsilon\}=0$ in (37). The integral (38) can then be performed analytically via the use of the Gamma function, and a power law behavior in $v^{2}$ is observed. The exponent turns out to be exactly the same as the one obtained in (45). In this regard as well, the velocity distribution function is rather unusual. For the mass and flux distribution functions, retaining only the energy gave the same power-law behavior for small values. In contrast, we have seen that the same power-law behavior for the velocity distribution is obtained for large values of ${\bf v}$.

Finally, a clarification regarding $\mathcal{F}_{\bf v}$ is necessary. Although we have determined, and plotted, it as a function of $v$, the distribution function $\mathcal{F}_{\bf v}$ is really a function of ${\bf v}$. This stems from the fact that it actually depends on $v^{2}\equiv{\bf v}\cdot{\bf v}$. In other words, the distribution function is isotropic in ${\bf v}$, but it encompasses the full 3D velocity dependence. If our system was endowed with momentum or angular momentum dependence, the resultant distribution function would have been anisotropic - a similar situation also exists in gravitational dynamics (Binney and Tremaine, 1987).

In Fig. 6, the numerical solution of the 3D velocity distribution function (38) is provided. Note that we have plotted $\mathcal{F}_{\bf v}$ as a function of $v$ (instead of $v^{2}$) since this variable is more physically transparent for visualization purposes. An inspection of the figure reveals that the analytical asymptotes (44) and (45), valid for small and large values of ${\bf v}$ respectively, capture the correct behavior of $\mathcal{F}_{\bf v}$ in these regimes. Unlike all of the previous distribution functions, there is no exponential cutoff since we observe a power-law behavior in both limits.

One can also compute the distribution function for the modulus of the velocity alone, viz. effectively a 1D distribution, denoted by $\mathcal{F}_{v}$ hereafter. Using the probability conservation law described in Sec. IV.1, we have

[TABLE]

and we can now use $d{\bf v}=4\pi\,v^{2}\,dv$ as per the discussion in Sec. IV.1. Note that $\mathcal{F}_{\bf v}$ must be expressed in terms of $v$, which has already been done since $\mathcal{F}_{\bf v}$ depended on $v^{2}$. Hence, we end up with

[TABLE]

and the LHS and RHS are solely functions of $v$. We are now in a position to derive the two limiting cases of $F_{v}$, which are presented below.

[TABLE]

[TABLE]

and this implies that the number of particles increases with $v$, attains a maximum, and then decreases for higher values of $v$.

V The connections with space and astrophysical plasmas

Since the 1970s and 1980s, detailed observations of flux ropes/tubes (that we have referred to as “plasmoids”) in the Earth’s magnetotail have been available (Russell and Elphic, 1978; Paschmann et al., 1982; Saunders et al., 1984). It was also during this period that a great deal of progress, both from a theoretical and observational standpoint, was achieved concerning the formation and properties of plasmoids (Hones, 1977, 1979). By the 1980s, it was well understood that plasmoids were inherently three-dimensional in nature, corresponding to flux rope/tube structures (Hones et al., 1982; Hughes and Sibeck, 1987). In contrast, many theoretical and numerical studies carry out 2D analyses.

There have been many statistical analyses of 3D plasmoids (flux ropes) over the past few decades in the terrestrial, solar and planetary contexts (Richardson et al., 1987; Moldwin and Hughes, 1992; Ieda et al., 1998; Hu and Sonnerup, 2001; Vignes et al., 2004; Jackman et al., 2011; Vogt et al., 2014). From these studies, a few broad inferences can be drawn, and we shall begin by considering the data from the magnetotail.

• •

It has been shown that the width distribution is subject to an exponential falloff for large values of $R$ (Fig. 7 of Fermo et al. (2011)), which is loosely consistent with (24).

• •

Earlier papers by Richardson et al. (1987) (see Fig. 16) and Moldwin and Hughes (1992) (see Fig. 13) also appear to show a somewhat similar trend.

• •

In Fig. 9 of Moldwin and Hughes (1992), the velocity distribution declines for larger $v$, thereby proving to be qualitatively similar to the trend predicted in Sec. IV.3.

Furthermore, in the Martian magnetosphere, the width distribution may also display an exponential falloff for large values, as predicted by (24) - see Fig. 4b of Briggs et al. (2011). However, we wish to emphasize that the above statements must be viewed with due caution for the following reasons: (i) the number of plasmoids observed is typically low, (ii) the interactions between the plasmoids are likely to be relatively insignificant, and (iii) resistive MHD is not applicable to these systems. Hence, it is more instructive to consider the results from solar observations on account of the reasons that were identified in Secs. I and II.

There exist some commonalities between the flux rope structures observed in the magnetotail and in coronal mass ejections (CMEs), as noted in Lin et al. (2008) and Linton and Moldwin (2009). In the latter context, several statistical analyses of solar flux rope properties (such as the width) have been conducted in recent times. In Janvier et al. (2014), small flux ropes (Moldwin et al., 2000) were shown to obey a power-law width distribution with an exponent of $-2.4$ (other cases such as $-1.8$ and $-2.1$ were also considered), whilst magnetic clouds (larger flux ropes) were characterized by a Gaussian distribution.

These results are mostly consistent with our theoretical predictions (21) and (24), especially the former where the power-law exponent of $-2$ was derived. We also refer the reader to Feng et al. (2007) and Feng et al. (2008) for the size distributions of solar flux ropes, which were also analyzed by Janvier et al. (2014) in their work. A summary of the current status of these statistical studies can be found in Section 4.3 (see also Fig. 17) of Janvier (2017). The combination of a power law at small scales, and the exponential cutoff at larger scales was also obtained in Guo et al. (2013) for post-CME sheets through the use of topological methods. Clearly, this result is qualitatively consistent with our methodology which predicts the absence of a universal power law across all scales.

Although Guo et al. (2013) suggest that the power-law exponent for small values may be approximately $-1$, this prediction is partially motivated by the relation $\psi\propto R$ that was proposed in Uzdensky et al. (2010). Moreover, our result does not appear to directly contradict their findings since the limit (21) is formally valid only in the limit $R\rightarrow 0$. The analytical prediction in this regime cannot be easily probed by either simulations or data analysis (since there are an insufficient number of plasmoids).

Hitherto, we have focused on comparing the general trends (and specific results) predicted by our model with some of the existent observations. We close this Section by pointing out a couple of avenues where our scalings may prove to be useful in understanding astrophysical systems. The first is the capacity of plasmoid motions and plasmoid-mediated reconnection for driving particle acceleration, which has been studied extensively (e.g Øieroset et al., 2002; Drake et al., 2006; Chen et al., 2008; Oka et al., 2010; Birn et al., 2012; Liu et al., 2013; Nishizuka and Shibata, 2013; Nishizuka et al., 2015; Nalewajko et al., 2015; Innes et al., 2015). The velocity distribution of the plasmoids computed herein may therefore enable a better understanding of how particle acceleration can occur via interactions with plasmoids.

Lastly, we observe that plasmoid reconnection has been studied extensively in the context of solar phenomena such as flares and CMEs (Kliem et al., 2000; Shibata and Tanuma, 2001; Bárta et al., 2011; Shibata et al., 2013; Lin et al., 2015; Shibata and Takasao, 2016; Dudík et al., 2016). Most of these processes are likely to be operational on other stars as well (Maehara et al., 2012; Shibayama et al., 2013; Lingam and Loeb, 2017a). A better understanding of the statistical properties of the ejected plasmoids could contribute to our understanding of the distribution of superflares and CMEs on M-dwarfs (Candelaresi et al., 2014; Shibata and Takasao, 2016). In turn, we expect that gaining such knowledge would enable us to make some progress towards understanding the prebiotic chemistry (Segura et al., 2010; Airapetian et al., 2016; Lingam and Loeb, 2017b) and atmosphere escape rates of exoplanets, especially around M-dwarfs (Shields et al., 2016; Dong et al., 2017; Lingam and Loeb, 2017c).

VI Conclusion

In statistical physics, the Principle of Maximum Entropy has proven to be very successful in deducing the probability distribution function given a set of invariants. The advantage of this approach has been argued to lie in its simplicity, physical transparency and versatility; the latter arises from the fact that this method can be interpreted as a form of statistical inference that can be applied to near-equilibrium systems (Jaynes, 2003).

After adopting this approach, we arrived at the final distribution function for 3D plasmoids, which depended on the mass, velocity and total flux, each of which constitutes physically meaningful variables. Our treatment is primarily applicable to systems with large numbers of interacting plasmoids, and where resistive MHD constitutes a reasonable description of the underlying dynamics. To the best of our knowledge, a 3D analysis along these lines is novel, and the distributions for these variables have never been derived before. Our final expression for the three-variable distribution function, viz. (10), depended on four free parameters, as well as another one that quantified the geometric shape of the plasmoids. At first glimpse, one may therefore argue that there is a high degree of arbitrariness, owing to the large number of unspecified parameters.

However, we showed that a certain degree of generality (within the resistive MHD model) still emerged when the mass, width, total flux, and helicity distribution functions were determined. They all possessed a power-law behavior for small values, and were subject to an exponential falloff for large values. Moreover, the power-law exponents at small values were very robust, and were nearly independent of the choices of free parameters. The power-law exponents were $-4/3$, $-2$, $-3$ and $-2$ for the mass, width, flux, and helicity distributions respectively. Even for larger values, the exponential functions were shown to have some commonality, regardless of the actual choices of the parameters.

The emergence of this behavior was even more striking when the distribution function for the (3D) velocity was considered. We showed that the small and large limits of the velocity distribution both led to power laws, with exponents of [math] and $-7/3$, when computed as a function of the variable $v=|{\bf v}|$. The two exponents remained mostly unaffected by the choices of the free parameters prevalent in our overall distribution function, which lends credence to the hypothesis that these asymptotes may be quite robust insofar our model is concerned.

Finally, we have undertaken a preliminary comparison of our model with the observational data, and concluded that the former yields results that are mostly consistent with the latter. It must, however, be recalled that this study was somewhat rudimentary, and in-depth observational (and numerical) data analyses (Schanche et al., 2016; Ni et al., 2017) are undoubtedly necessary for obtaining a clear picture of the strengths and weaknesses of the theoretical model. Along the way, we also delved briefly into the attendant astrophysical implications for phenomena ranging from particle acceleration to solar flares and CMEs.

As the formalism proposed in this paper might be applicable to other scenarios, several avenues open up for subsequent research. The same analysis can be redone by introducing a new set of constraints for the 2D system that are clearly distinct from the 3D case (Kraichnan and Montgomery, 1980; Morrison and Hazeltine, 1984), and the ensuing scaling laws are also expected to be different compared to their 3D versions.

Another possibility is to obtain the distribution of relativistic plasmoids, which could then be compared against recent simulations (Sironi et al., 2016; Petropoulou et al., 2017). The inclusion of turbulent kinetic and magnetic energies might also pave the way for analyzing the statistical properties of turbulent plasmoid-mediated reconnection (Huang and Bhattacharjee, 2016); see also Lazarian et al. (2015). Based on the simplicity and elegance of the MaxEnt approach, we suggest that further theoretical analyses along these lines may be warranted.

Acknowledgements.

The authors are very grateful to Amitava Bhattacharjee and Yi-Min Huang for the enlightening discussions and the perceptive remarks. The beneficial conversations with Eero Hirvijoki, Russell M. Kulsrud, Takuya Shibayama and Roscoe B. White are also acknowledged. The authors were supported by the Department of Energy Contract No. DE-AC02-09CH11466 and the National Science Foundation Grant Nos. AGS-1338944 and AGS-1552142 during the course of this work.