Multi-Beam Energy Moments of Compound Measured Ion Velocity Distributions
Martin V. Goldman, David. L. Newman, Jonathan P. Eastwood, Giovanni Lapenta
1 Multi-Beam Energy Moments of Complex Measured Ion Velocity Distributions
M. V. Goldman, D. L. Newman, J. P. Eastwood and G. Lapenta Physics Department, University of Colorado at Boulder. MG Orcid 0000-0001-8039-9143; DN Orcid 0000-0003-0810-1204 The Blackett Laboratory, Imperial College London, London SW7 2AZ, UK. Orcid 0000-0003-4733-8319 Department of Mathematics, KU Leuven, University of Leuven, Celestijnenlaan 200B, 2001, Leuven, Belgium. Orcid 0002-3123-4024
Martin V. Goldman ([email protected])
Key Points:
Multi-beam ion distributions, f i ( v ), measured by spaceraft and in PIC simulations are analyzed in terms of appropriate velocity moments. Three methods are developed for approximating f i ( v ) as a sum of beams: a "visual" method, a "k-means" method and a least-squares-fit method. Multi-beam energy moments of f i ( v ) found by taking standard moments of each beam and summing over beams are free of pseudothermal contributions. 2 Abstract
Complex ion distributions, f i ( v ), have been measured with high-time resolution by NASA's Magnetospheric Multi-Scale Mission (MMS) and have been found in reconnection simulations. A complex distribution, f i ( v ), consisting, for example, of essentially disjoint pieces will be called a multi-beam distribution and modeled as a sum of "beams," f i ( v ) = f ( v ) + ... + f N ( v ). Velocity moments of f i ( v ) are taken beam by beam and summed. Such multi-beam moments of f i ( v ) have advantages over the customary standard velocity moments of f i ( v ), for which there is only one mean flow velocity . For example, the standard thermal energy moment of a pair of equal and opposite cold particle beams is non-zero even though each beam has zero thermal energy. We therefore call this thermal energy pseudothermal. By contrast, a multi-beam moment of two or more beams has no pseudothermal energy. We develop three different ways of decomposing into a sum and finding multi-beam moments for both a multi-beam f i ( v ) measured by MMS in the dayside magnetosphere during reconnection and a multi-beam f i ( v ) found in a PIC simulation of magnetotail reconnection. The three methods are: • A visual method in which the velocity centroid of each beam is estimated and its density determined self-consistently • A k-means method in which particles in a particle-representation of f i ( v ) are sorted into a minimum energy configuration of N (= k) clusters • A nonlinear least squares method based on a fit to a sum of N kappa functions. Multi-beam energy moments are calculated and compared with standard moments for the thermal energy density, pressure tensor, thermal energy flux (heat plus enthalpy fluxes), bulk kinetic energy density, RAM pressure and bulk kinetic energy flux. Applying this new formalism to real data demonstrates in detail how multi-beam techniques may provide significant insight into the properties of observed space plasmas. 3
Particle velocity distributions, f( v ), in space plasmas have been measured with increasing resolution over the past two decades by means of electrostatic analyzer instruments (Fazakerley, et al., 1998). The launch of the Fast Plasma Instrument (FPI) (Pollock et al., 2016) on board the Magnetosphere Multi-Scale (MMS) Mission (Burch, et al, 2016) now enables the full 3D electron and ion distribution to be measured over unprecedented time-intervals of 30 ms and 150 ms cadence respectively (or 7.5 ms and 37.5 ms in certain operating modes (Rager, et al, 2018). A common and useful method of analyzing particle velocity distributions has been to take standard velocity moments of f( v ). Moments enable the study of f( v ) at a given spatial point and time (or over a space-time average) in terms of fluid variables, such as flow velocity, energy density, temperature, pressure, etc. The meaning of fluid moments is well understood and not controversial for contiguous and effectively single-peaked velocity distributions, f( v ). However, recently measured electron and ion velocity distributions (Burch, et al, 2016, Eastwood, et al, 2014) can be effectively disjoint and are not easily interpreted in terms of standard velocity moments. We shall refer to such f( v ) as " multi-beam" velocity distributions. An alternative formal way of taking multi-beam moments of multi-beam velocity distributions, of form f( v ) = f ( v ) + ... + f N ( v ), has been developed recently (Goldman, et al 2020). A multi-beam moment of f( v ) results from taking a standard moment of each of the "beams," f ( v ), ..., f N ( v ), and then summing over the beams. Multi-beam moments are often more easily interpreted than standard moments; they don’t contain pseudothermal energy, such as that present in the standard thermal moment of an f( v ) consisting of two equal and opposite cold beams. The pseudothermal part of any standard thermal moment is defined as the difference between the standard thermal moment and its multi-beam counterpart. In the present paper we examine how this formalism can be applied to real data, using a multi-beam ion velocity distribution, f i ( v ), measured during dayside magnetopause reconnection (Burch, et al, 2016) and an f i ( v ) found in a Particle-In-Cell (PIC) tail reconnection simulation (Eastwood, et al, 2014; Goldman et al, 2016). For a given f( v ), we calculate both the standard and multi-beam moments, thereby identifying and eliminating the pseudothermal parts of standard thermal moments. Although each of the three different methods developed in this paper gives (slightly different) multibeam moments, we find they provide a reliable means to calculate multi-beam moments of complex distribution functions. The key step in taking multi-beam moments is to find a way to approximate a measured ion multi-beam f i ( v ) as a sum of N effective "beam" contributions, f i ( v ) = f ( v ) + ... + f N ( v ). Three different methods are developed for doing this and applied both to the multibeam f i ( v ) measured during magnetospheric reconnection and to the multibeam f i ( v ) from a PIC simulation of mag-netotail reconnection. The three methods are • A visual method in which the velocity centroid of each beam is estimated • A k-means method in which particles in a particle-representation of f i ( v ) are sorted into a minimum energy configuration of N (= k) clusters • A nonlinear least squares method based on a sum of N kappa functions . Review of standard and multi-beam moments of particle velocity distributions, f( v ) Kinetic theory definitions of standard velocity moments are reviewed and distinguished from their multi-beam counterparts. Moments treated are: the standard bulk kinetic energy density, the thermal energy density, the bulk kinetic energy density flux, the thermal energy density flux (= enthalpy flux + heat flux) and the standard pressure tensor and pressure ellipsoid.
Pseudo-thermal moments are defined as differences between standard and multi-beam thermal moments.
They can be expressed in terms of beam densities and velocities. Pseudothermal parts of the standard thermal energy moments are shown to be equivalent to bulk kinetic energy deficiencies of the same moments. Visual method of finding multi-beam moments of a multi-beam f( v ) Visually estimating the centroid velocities of four assumed beams (with some constraints) enables one to find the beam densities and multibeam moments of a velocity distribution f( v ) = f ( v ) + f ( v ) + f ( v ) + f ( v ) whose standard moments are known. Fewer than four beams may be used, with additional constraints on the centroid velocities. k-means and least squares methods of approximating f( v ) and finding multi-beam moments In addition to the visual method , there are two methods for decomposing a given velocity-space distribution, f( v ), into a sum of N beams without any constraint on N. The first method, typically referred to as “ k -means,” divides the given f( v ) into k disjoint “clusters” in such a way that the thermal energy summed over the k clusters is a minimum. (The thermal energy is defined relative to the centroid of each cluster) The second method uses a nonlinear “least squares” algorithm to minimize the cumulative square of the difference between a chosen distribution, f( v ), and the corresponding value specified by a parameterized model comprising a superposition of populations. A description of each method is provided together with an explanation of how the measured velocity distribution, f( v ), is preprocessed. Multibeam moments of ion distribution from a PIC simulation of magnetotail reconnection
The three methods of finding multibeam moments are applied to an ion velocity distribution, f( v ), taken from a PIC simulation of magnetotail reconnection. Three-dimensional isosurface contour plots of f( v ) are generated and centroid velocities of the N beams displayed. In most cases N is assumed to be 4, and in one case, 2. For the k-means and least-squares methods, the N beam distributions are shown explicitly. The multibeam and standard energy densities are compared, as are the energy density flux vectors. Ion heat fluxes are found to be negligible. Multibeam moments of ion distribution, f( v ) measured by MMS during dayside reconnection Multibeam moments are found by each method for an ion velocity distribution, f( v ), measured by MMS during magnetopause reconnection. Three-dimensional isosurface contour plots are 5 used to illustrate f( v ) and the beams. Multibeam and standard energy densities are compared, as are the energy density flux vectors. Ion heat fluxes are found to be non-negligible. Conclusions
A summary is presented of the three methods for approximating an f( v ) measured in a PIC simulation of magnetotail reconnection and an f( v ) measured by MMS in the dayside magnetosphere by a sum of N beams. The quantitative findings regarding velocity moments are presented concisely. The practical significance of multibeam moment analysis for magnetic reconnection events is explained. 6 Review of standard and multi-beam moments of particle velocity distributions, f( v ) In Goldman, et al (2020) it is shown how to define and take multi-beam moments of a multi-beam f(v) from the densities and velocities of N beams in the approximation f(v) = f (v) + ... + f N (v). 2.1 REVIEW OF STANDARD MOMENTS IN TERMS OF KINETIC THEORY First, we review definitions and notation used here for standard moments of f( v ): Table 1
Standard velocity moments of f(v)
Name of standard moment Symbol Definition (dimensional units)
Number density n ò d v f( v ) Number density flux n u ò d v f( v ) v , u = flow (centroid) veloc. Bulk kinetic energy density U bulk nmu /2 Thermal energy density Frame invariant U thermal (m/2) ò d v f( v )( d v) , d v = v - u Undecomposed energy density U U bulk + U thermal = (m/2) ò d v f( v )v Bulk kinetic energy flux vector Q bulk n u mu /2 Heat flux vector Frame invariant Q heatflux (m/2) ò d v f( v )( d v)· ( d v) , d v = v - u Enthalpy flux vector Q enthalpy (m/2) ò d v f( v )[ u ( d v) +2 u· ( d v) ( d v) ] Thermal flux vector Q therm Q heatflux + Q enthalpy Undecomposed energy flux
Q Q bulk + Q thermal = (m/2) ò d v f( v ) v v Pressure tensor
Frame invariant P m ò d v f( v ( d v) ( d v) RAM pressure tensor P RAM nm uu Stress tensor (undecomposed moment) T P RAM + P = m ò d v f( v ) vv Thermal energy density in terms of pressure U therm Tr P /2 Enthalpy flux in terms of pressure Q enthalpy u Tr P /2 + u · P It is evident from Table (1) that each standard undecomposed moment is the sum of the corresponding (coherent) bulk and (incoherent) thermal moments. 7 2.2 MULTI-BEAM MOMENTS There are at least two ways to take multi-beam moments of a multi-beam f( v ). One, already mentioned in the Introduction , is to take a standard moment of each of the terms (beams) in f( v ) = f ( v ) + ... + f N ( v ) and then sum. The other is to use only the densities, n , ... n N and the centroid velocities of the N beams u , ... u N as described in Goldman, et al (2020). From that paper we reproduce, in Table (2), expressions for the multi-beam thermal moments in terms of the “pseudothermal” parts of the thermal moments, D U therm , D Q t herm , and D P using standard moments and only the N partial beam densities, h j = n j /n and centroid velocities, u j . Note the dimensionless units introduced in the middle column in Table 2. In Goldman, et al (2020) it is proven that an undecomposed multi-beam moment is always equal to the corresponding undecomposed standard moment. Thus, if a standard thermal moment is reduced by a pseudothermal part, D , when becoming multi-beam moment then the corresponding multi-beam bulk moment is increased from the standard bulk moment by the same quantity, D . Table 2: Formulas for determining multi-beam moments from standard thermal moments, beam densities and beam centroid velocities (dimensionless unit scaling factor in 3rd column ) f( v ) moment Pseudothermal part of standard moment Units Multibeam moment Energy Density
Pressure tensor
Energy Flux vector Δ U ≡ u n η j u j − u j = N ∑ ⎡⎣⎢⎢ ⎤⎦⎥⎥ mnu n U thermMB = U therm − Δ U Δ P ≡ u n η j u j u j − uu j = N ∑ ⎡⎣⎢⎢ ⎤⎦⎥⎥ mnu n P MB = P − Δ P Tr P = U therm , Tr P MB = U thermMB Δ Q = u n η j u j u j − u u j = N ∑ ⎡⎣⎢ ⎤⎦⎥ mnu n Q thermMB = Q therm − Δ Q Q therm ≡ Q enthalpy + Q htflux , Q enthalpy ≡ ˆ u Tr P + u · P "pseudothermal" parts of standard moments are equivalent to differences between multi-beam bulk moments and the corresponding standard bulk moments. For example, D U = [U
MBbulk - U bulk ] is also equal to [U therm - U MBtherm ] because of the equivalence of multi-beam and standard undecomposed moments, U MB = U. Hence, we may express, (1a) U MBtherm = U - U
MBbullk
One can think of the pseudothermal part of a standard thermal moment as a deficiency in the bulk kinetic part of the corresponding standard bulk moment. The pseudothermal part and the bulk kinetic energy deficiency are equal to each other whether they are scalars (as for energy density moments), tensors (as for pressure moments, or vectors (as for energy density flux moments). In the construction of multi-beam moments using any of the methods described in this paper, one must first decide on the number of beams , N, in the representation f( v ) = f ( v ) + ... + f N ( v ). The visual method requires N ≤ 4, whereas there is no restriction on N for the k-means and least squares methods. A recent study based on a sum of Maxwellian f j ( v ) suggests that the most statistically likely values of N are low single digits (Dupuis, et al, 2020). For the standard thermal energy density, a choice of larger N tends to give a smaller U MBtherm and a larger pseudothermal energy, D U. This is consistent with the fact that, in the limit of large N, we recover the kinetic result that there are only discrete particles: If N is equal to the number of particles there can be no temperature or thermal energy density. These are strictly fluid concepts. In this paper we will mainly use N = 4 (i.e., four beams) for all methods, to provide a clear illustration of their implementation and to enable comparisons. We will also include one N = 2 two-beam example (for the visual method). The two-beam case gives the largest value of U MBtherm and the smallest pseudothermal energy (or, equivalently, the smallest deficit of bulk kinetic energy). The pseudothermal parts of the standard pressure moment is a tensor and the pseudothermal part of the standard thermal energy flux is a vector, so one must consider direction as well as magnitude in the pseudothermal parts. In other words, Q MBbulk can change in magnitude and direction as beams are added to N. Once again, (1b). In the examples treated in this paper | D Q therm | << | Q therm |. Q thermMB = Q − Q bulkMB Visual method of finding multi-beam moments of a multi-beam f( v ) In this Section it is explained how to use Table (2) to find all multi-beam moments by the visual method . In Sections (5) and (6) the visual method and other methods are used to find the multi-beam moments of multi-beam ion velocity distributions found in a simulation and measured in the magnetosphere 3.1 STRATEGIES IN THE
VISUAL METHOD
The visual method begins with a visual representation of the measured f( v ) at a given place and time or over a narrow space-time average near one place and time. This may consist of plots of the three mutually orthogonal reduced velocity distributions of f( v ) or consist of a 3-D contour map of f( v ). By inspection, one then estimates the velocity centroids, u j of each of the putative N "beams," f ( v ), ..., f N ( v ), as well as the number of beams, N. As will be shown below, if N = 4, the densities of each of the beams, n , n , ..... n N are fully determined by Eqs. 2, below. The density and centroid velocity of beam j are explicitly defined via standard zero and first-order velocity moments of f j ( v ), n j = ò d v f j ( v ) and n j u j = ò d v f j ( v ) v . However, the full form of the distributions f j (v) is neither determined nor needed. The velocity moment integrals are never actually performed on the beams since, after estimating the locations of the beam centroid velocities u j , the values of the beam densities, n j follow from Eqs. 2. Secs. (3.2) and (3.3) show how to find the densities, n j for N = 4 and for N ≠ 4. 3.2 HOW TO DERIVE BEAM DENSITIES FROM FOUR BEAM VELOCITIES We first demonstrate that, for N = 4, the beam densities , n j (or fractional densities , h j =n j /n) are completely determined in terms of the four beam velocities found, for example, by visual inspection of the measured f( v ). This is because conservation of particles and particle flux constitute exactly four scalar eqns. Other values of N may be employed but are less useful in association with this method, as we show later. The zero-order-moment of f( v ) = f ( v ) + f ( v ) + ... + f N ( v ), yields conservation of particles (when integrated over volume), Here, n j is the density of beam j, n is the cumulative density of the ensemble of beams and h j = n j /n is the fractional density of each beam. The first-order moment of f( v ) yields conservation of particle flux. When divided by n, this becomes (2b) , where u j is the centroid velocity of the j th beam and u is the (weighted) mean velocity of the ensemble of beams (i.e. the standard flow velocity of f( v )). Each beam velocity, u j , in the weighted mean, u , is multiplied by the fractional beam density, h j . (2 a ) d v f v ( ) / n ∫ = = η + η + ... + η N , η j ≡ n j / n d vv f v ( ) / n ∫ = u = η u + η u + ... + η N u N
10 It is evident from Eqs. (2) that the fractional densities, h j and beam velocities, u j are not independent. Consider fitting a measured multi-beam f( v ) with exactly four beams (N = 4) with beam velocities, a , b , c , and d , where, the components of each beam are given by In practical applications, the implicit {v x , v y , v z } coordinate system may be simulation coordinates, or, for satellite data, GSE, or field-aligned, {v || , v ⊥ , v ⊥ }. Equations (2) may be expressed in terms of 4-dimensional matrix eqns. A matrix, M , acts on the four-vector, h = { h , h , h , h } to yield a four-vector composed of the three components of u and a fourth component equal to 1. Assuming M has a non-zero determinant, the inverse gives the four-vector, h : (3) Assuming the centroids locations allow the determinant of M to be non-zero, the four fractional beam densities are then found in terms of the inverse matrix, M -1 by. Eqn, (3) Thus, the four fractional beam densities, { h , h , h , h } can be evaluated by choosing ion beam velocities, a , b , c , d and a mean ion flow velocity u. The values of the four fractional beam densities, { h , h , h , h } will therefore depend on the values of the four fitted beam velocities and the given value of u . In applying this method, we must be very careful to choose beam velocities which guarantee that all four h j > 0. This will be guaranteed if u lies within the velocity-tetrahedron whose four vertices are at the four points a , b , c , d . The beam velocities chosen must also be chosen so as to guarantee that the eigenvalues of the multi-beam pressure tensor are all positive, as will be discussed later. a ≡ u x , u y , u z { } , b ≡ u x , u y , u z { } , c ≡ u x , u y , u z { } , d ≡ u x , u y , u z { } η η η η ⎡⎣⎢⎢⎢⎢⎢ ⎤⎦⎥⎥⎥⎥⎥ = M − · u x u y u z ⎡⎣⎢⎢⎢⎢⎢ ⎤⎦⎥⎥⎥⎥⎥ , M · η η η η ⎡⎣⎢⎢⎢⎢⎢ ⎤⎦⎥⎥⎥⎥⎥ = u x u y u z ⎡⎣⎢⎢⎢⎢⎢ ⎤⎦⎥⎥⎥⎥⎥ , where, M ≡ a x b x c x d x a y b y c y d y a z b z c z d z ⎡⎣⎢⎢⎢⎢⎢ ⎤⎦⎥⎥⎥⎥⎥
11 3.3 FINDING DENSITIES WITH FEWER/MORE THAN FOUR BEAMS. Suppose we try to fit the multi-beam f( v ) with two beams instead of four. The two beams must satisfy the vector Eqs. (2), whose three components give three expressions for h : (4a) In Fig. (1), a ¢ , b ¢ , and u ¢ are the projections of a , b , and u onto the x-y plane, whose origin has been taken at b ¢ = 0. Since b x = b y = 0, the x and y components of Eqn (4a) yield, (4b) Let j u be the angle u ¢ makes to the x-axis and let j a be the angle a ¢ makes to the x-axis. Then Eqn. (4b) becomes tan f u = tan f a , or f u = f a , which requires the {x-y}-projected flow velocity, u ¢ , to lie on the line segment joining a ¢ , and b ¢ . (In order for h to be positive, u must lie between a and b) . These constraints can make it harder to find even an approximately self-consistent estimate of the centroid locations, a and b , of the two beams. Three beams also present constraints in fitting to f( v ), which are absent for four beams. Three beams define a triangle and hence a plane in 3D velocity space. The constraint is that the mean flow velocity must lie in this plane. For five or more beams, the fractional density of the fifth and higher order beams is not determined by particle and particle flux conservation. For the five-beam case, if the fractional density, h , and velocity, u of the fifth beam are determined by other means, then the same matrix, M , used above in the four-beam case, may be employed, together with the source term now containing h and u to find the other four fractional densities: (5) η a + −η ( ) b = u ⇒ η = u x − b x a x − b x = u y − b y a y − b y = u z − b z a z − b z > u y u x = a y a x M × η η η η ⎡⎣⎢⎢⎢⎢⎢ ⎤⎦⎥⎥⎥⎥⎥ = u x − η u x u y − η u y u z − η u z − η ⎡⎣⎢⎢⎢⎢⎢⎢ ⎤⎦⎥⎥⎥⎥⎥⎥ , η η η η ⎡⎣⎢⎢⎢⎢⎢ ⎤⎦⎥⎥⎥⎥⎥ = M − × u x − η u x u y − η u y u z − η u z − η ⎡⎣⎢⎢⎢⎢⎢⎢ ⎤⎦⎥⎥⎥⎥⎥⎥ Figure 1:
Visual method constraints for two-beam fit
4. K-means and least square methods of approximating f(v) and finding multibeam moments
In addition to the visual method , we have found two complementary methods for decomposing a given velocity-space distribution, f( v ), into a specified number of populations (e.g., beams). • The first method, typically referred to as “ k -means ,” divides the given f( v ) into N = k disjoint “clusters” in such a way that the thermal energy summed over the N clusters is a minimum. The k-means method is a widely used technique, well-known in image analysis. • The second method uses a nonlinear “ least squares ” algorithm to minimize the cumulative square of the difference between f( v ) (or a function of f( v )) and the corresponding value specified by a parameterized model comprising a superposition of populations. In each case, we employ algorithms from the MATLAB R analysis package; with comparable algorithms implemented in other software (e.g., Mathematica R ). The methodology is described more fully in the subsections below. The two "given" velocity distributions, f( v ), to be used later in this paper will be • from a PIC simulation (Eastwood, et al, 2014) (Sect. (5), (7)) • from the MMS mission (Sect. (6), (7)) One important way in which the k -means and least-squares algorithms differ is in how the data is represented. While the least-squares algorithm operates directly on f( v ) in 3D velocity space, the k -means algorithm requires a collection of velocity-space coordinates of discrete equally weighted particles distributed according to f( v ). Our discussion of the two decomposition models will assume the satellite and simulation data are in the necessary form, with a description of how the data is brought into the required form deferred until the end of this section. 4.1 K-MEANS DECOMPOSITION The starting point for the k-means decomposition is a set of N p particles of uniform weight, each with velocity v i , for i = 1 , ..., N p , distributed according to the velocity distribution function f ( v ) as described below. The object is to divide the N p particles into N = k clusters . In the continuum picture this is equivalent to f( v ) = f ( v ) + ... + f k ( v ), where none of the f j ( v )'s are overlapping. To be consistent with the visual method of beam determination discussed above, we will apply the method with N = 4, although it works for N equal to any integer. The k-means algorithm finds the 'center' velocity, u j , of cluster j, with j = 1 , .. , k N , labelling each of the N clusters. It also yields the number, N j , of particles associated with each cluster. The set of N j particles in cluster j is associated with one (and only one) center velocity, u j . The N p particles are divided into k disjoint sets (clusters) of N j particles each, with, In addition, the k-means algorithm determines the k center velocities, u j , so as to satisfy the following conditions: (6 a ) N jj = N ∑ = N p
13 Eqn. (6b) states that each cluster center, u j , is the centroid of particle velocities in that cluster. The condition in Eqn. (6c) is equivalent to the requirement that the total thermal energy of particles relative to their respective cluster centers be a local minimum. The left side of Eqn. (6c) is essentially the multi-beam thermal energy, based on the standard U therm in Table (3). From the derived expression, it is clear that the k-means procedure will simultaneously maximize the pseudothermal energy, D U, as defined in Table 2. The k-means algorithm is nonlinear and performs an iterative search for the k velocities u j , yielding a local minimum of the total multibeam thermal energy while meeting the requirement that the cluster centers also be the cluster centroids. There is therefore no guarantee that the solution is a global minimum. By repeating the search with different initial guesses, one can ascertain whether the local minimum is likely close to a global one. For the examples presented in this study, this multibeam thermal energy appear to be close to the global minimum with high likelihood. The fractional density in each cluster is simply the fraction N j /N of all particles in cluster j . Furthermore, because each center is the centroid of the respective cluster it immediately follows that the standard bulk flow velocity, u , satisfies particle density flux conservation, Because the k -means algorithm identifies the cluster center for each of the original N particles the full pressure tensor P j and heat-flux, Q j-heatflux, can be computed for each cluster, j. Therefore, terms such as Q j − heatflux can be determined for each cluster as well as for the sum over all clusters Q MBheatflux . This is an advantage over the visual method which does not yield the pressure tensors and heat fluxes of individual beams. Another advantage to this method is that it should be easy to set up automatically. In the Conclusion (Section 8) we compare the three methods of calculating multibeam moments 4.2 NONLINEAR LEAST-SQUARES DECOMPOSITION The nonlinear least-squares (NLSQ) decomposition attempts to find an optimal representation of given velocity distribution function f( v ) as a superposition of simple analytic distributions functions, g j ( v ). Here, we will use “kappa” distributions for our base set, g j ( v ), j = 1, ... , N. (6 b ) u j = N j v jii = N j ∑ , and (6c) v ji − u j ( ) = local minimum, i = N j ∑ j = k ∑ where, v ji = i th partice in j th cluster U thermMB = m d v f j v ( ) v − u j ( ) ∫ j = k ∑ U thermMB = U therm − Δ U (6 d ) u = η j u jj = k ∑
14 The density, bulk velocity, temperature tensor, and index κ of each g j ( v ) are taken as parameters to be determined through the fitting algorithm. The “data” to be fit is, in general, a scalar function of the 3D f( v ) on a cubic Cartesian grid of dimension N vx × N vy × N vz along the three orthogonal basis directions. The MATLAB algorithm we use assumes a 1D array of data values. This entails simply concatenating the 3D data into a one-dimensional data vector. Since the computation time and memory demand increase with the number of data points (as well as the number of model parameters), we subsample our distribution function before performing the fit. A grid of N g = 93 ≈ 5 . data points is found to afford a reasonable trade off between velocity-space resolution and computational demand. The model distribution we will use to approximate the satellite or simulation data in this study will be a superposition of triaxial kappa distribution, for which the tri-Maxwellian distribution is a limiting form as . The essential feature of such distributions (in their CM frame) is that each g j ( v ) is constant on nested ellipsoidal surfaces of uniform shape and orientation. Specifically, each g j ( v ) depends on velocity only through the dimensionless scalar function, h( v ), defined as follows: , where ( v , v , v ) are the components of velocity, v , in the coordinate system for which the pressure tensor is diagonal (it can always be diagonalized, due to the symmetry property P ij ≡ P ji ). The mass of the species being considered is denoted by M. The scalars, T i , for i = 1 , P /n (i.e., the diagonal elements of P /n in that coordinate system). Note that v T and v have the components of the velocity in the original coordinate system and are expressed as row and column vectors, respectively. Each kappa distribution is a function of velocity, v , characterized by 11 parameters: the density, n (scalar); the local mean velocity vector, u , (three components); the symmetric temperature tensor, P , (which has six independent components); and the scalar index κ . The explicit form of each parameterized kappa distribution, f kappa takes the following form: , Note f kappa has the usual physical units, n/v . The fitting distribution g( v ) = consists of a superposition of N independent kappa functions, each labeled by subscript j. Since g j ( v ) = , the fitting distribution, g( v ) has 11 N free parameters. The MATLAB NLSQ algorithm allows an upper bound, lower bound, or both to be set on any parameter. We constrain n to be positive for each component and restrict κ to the range 2 ≤ κ ≤ κ max , but do not impose any constraints on the other parameters. κ → ∞ (7) h v ( ) ≡ nM v T P − v ( ) = M v T + v T + v T ⎛⎝⎜ ⎞⎠⎟ (8) f kappa v ; n , u , P , κ ( ) ≡ n Γ κ + ( ) / Γ κ − ( ) π κ − ( ) T T T [ ] / M + h v − u ( ) κ − ⎛⎝⎜ ⎞⎠⎟ − κ + ( ) g j v ( ) j = N ∑ f kappa v ; n j , u j , P , κ j ( )
15 The upper bound ( κ max = 50 for the PIC simulation and κ max = 10 for the FPI data) differs relatively little from a Maxwellian (especially for κ max = 50) except well out on the power-law ( g ∼ v −2 κ max ) tail. The lower bound κ = 2 is larger than the minimum value κ = 3 / P ) of the distribution breaks down. Never-theless, with g ( v ) ∼ v −4 for larger v higher moments of the distribution such as the heat flux, may diverge. This, however, does not affect the validity of fitting with a small index such as κ = 2, because the fit is restricted to the domain of the finite velocity-space data grid. Outside of this domain, the tails of any physical distribution will likely decrease more rapidly with v . Note, however that the heat flux can remain finite even with a v -4 tail on the velocity distribution, provided the tails become sufficiently symmetric at high velocity (see Section (4.2) on tri-Maxwellians) In its simplest form, the NLSQ algorithm varies the model parameters to minimize the cumulative square error where the integer subscript n labels grid points ranging from 1 to N g . Note that any scalar function of the distribution can be used instead. This is particularly useful for modeling the observed FPI distributions where the cumulative square error due to differences in the tail of the distribution are insignificant and therefore cannot effectively influence the value of κ in the model. However, by first taking the logarithm of the distribution [both the observed f ( v ) and the model g ( v )], the controlling feature of the distribution shifts from the local maxima to the power-law tails. A hybrid approach in which both f ( v ) and log [ f ( v )] are used in the NLSQ fitting process is described in connection with the modeling of the actual FPI data. As described in the preceding sections, the k -means algorithm requires the particle distribution be represented by an ensemble of equally weighted discrete particles while the NLSQ algorithm requires the distribution function be specified on a 3D Cartesian grid. While the output of a PIC simulation might be directly usable as input to the k -means algorithm, the iPic3D code used in the reconnection simulation described here uses particles of unequal weight and therefore the simulation output requires additional processing. Since we need the velocity distribution on a Cartesian velocity-space grid for the NLSQ method, we first generate such a distribution from either the PIC simulation output or from the FPI distribution skymap and then generate an ensemble of equally weighted particles from the resulting distribution grid. A brief description of the methodology follows: 4.3.1 FROM IPIC3D PARTICLE DATA TO A CARTESIAN f( v ) The ensemble of weighted particles culled from the iPic3D reconnection simulation comprise the totality of particles located in a square 0 . d i on a side centered at the location where the two trajectory bundles in Fig. 3e and 3f of Goldman, et al, (2020) cross ( x = 120 d i and y = 14 d i ). There are ∼ macro-particles in the sample region. The velocity distribution function is created by adding the weight of each particle to the cumulative value for the cubic cell containing the particle’s velocity in a pre-defined Cartesian velocity-space grid. Since only ∼ . particles have a significant weight, the resulting f ( v n ) − g ( v n ) n = N g ∑ ,
16 distribution has relatively large statistical fluctuations. We apply five iterations of nearest-neighbor smoothing to reduce the fluctuation level as a visualization aid. The reduced distributions in the figures are based on the smoothed distributions. 4.3.2 FROM FPI SKYMAPS TO A CARTESIAN f( v ) The FPI skymaps take the form of a distribution function over a 16×32×32 grid in θ , ϕ , E space where E is the logarithmically spaced energy in eV. While the value of the distribution function can be interpolated directly onto a Cartesian velocity-space grid, we employ an intermediate step whereby the logarithm of the distribution on a shell of constant E is fit to a superposition of spherical harmonics Y lm ( θ,ϕ ) up to a predetermined order l max . In the present study we use l max = 13 for which there are ( l max + 1) = 196 parameters (i.e., Y `m amplitude coefficients) to fit the 32 × 16 = 512 angular data point (subject to a floor value in those angular bins where there are zero counts). The reconstruction of the distribution on each energy shell from the corresponding Y lm ( θ, ϕ ) is then used for the interpolation onto the predetermined Cartesian velocity-space grid. 4.3.3 FROM CARTESIAN f( v ) TO AN ENSEMBLE OF EQUALLY WEIGHTED PARTICLES The procedure for generating the desired ensemble of N p equally weighted particles distributed according to the 3D velocity distribution function on a Cartesian grid consists of the following steps: 1. Map the 3D Cartesian velocity-space grid of dimension ( N vx , N vy , N vz ) into a 1D vector of length N g = N vx × N vy × N vz . The gridded velocity distribution f ( v x ,v y ,v z ) then maps onto the data vector f j , j =1..., N g . 2. Construct the cumulative distribution vector, with F tot = F Ng being the sum over all grid points of f . 3. Generate a set of N p random numbers uniformly distributed between 0 and F tot . These values will be used to identify the velocities of the N p particles in the target ensemble. 4. Interpolate each random number onto the vector F j defined in step 2, rounding up to the next higher integer j . 5. Invert the mapping defined in step 1 from j back to the gridded 3D velocity ( v x , v y , v z ). 6. Add a random velocity perturbation to each particle so that its final velocity is randomly distributed with the corresponding cubic grid cell. Note that we are using v x , v y , and v z as generic velocity components. For the 2D PIC magnetotail reconnection simulation studied in Sec. (6), the initial state is a Harris current sheet of finite thickness in y, with current in the (out-of-plane) negative-z-direction, producing reversed magnetic field lines in the ±x directions. For this PIC simulation, the velocity components, { v x , v y , v z ,} are in these simulation coordinates. (This differs from the GSE convention, in which x is the same but y and -z are exchanged) For the case of the dayside MMS FPI-measured velocity distribution, v x , v y , and v z will be interpreted as magnetic field-aligned velocities { v || , v ⊥ , v ⊥ }. It is essential that the cubic cells in the Cartesian velocity-space grid all have the same volume for either case to work without modification. F j = f ii = j ∑ ,
5. Multibeam moments of f( v ) from PIC simulation of magnetotail reconnection A visual representation of the ion velocity distribution whose moments are to be taken is the starting point for the visual method. Although not intrinsic to the methods, a visual representation is also helpful in visualizing the results of the k-means and the least squares methods of finding multibeam moments. In both this Section, which addresses velocity distributions from PIC simulation, and in the next (dealing with MMS-FPI ion velocity measurements during dayside reconnection) two visual representations are employed • three orthogonal reduced
2D ion velocity distributions. •
3D fly-around surface-contour plots of the ion velocity distributions The underlying ion velocity distribution is of necessity somewhat different for each of the three methods of finding multibeam moments. In the present section the visual method determines the effective beam velocities by inspection of the ion velocity distributions obtained by processing and smoothing the PIC results (see Section (4.3.1)) to obtain the simulation "box" ion distribution, f box ( v ). The k-means method uses f box ( v ) to create the distribution of ion particles necessary to employ the k-means method (Sec. (4)). Finally, the least squares method is based on the least-squares velocity distribution, f fit ( v ) fitted to f box ( v ). For the PIC case addressed in this section the differences in the multibeam moments based on different underlying velocity distributions are essentially negligible. The multibeam moments found in this section by each of the three methods assumes that f box or f fit can be interpreted as a sum of beams, f( v ) = f ( v ) + ... + f N ( v ), with N = 4 or 2. The pseudothermal parts of the standard moments of f( v ) are found by comparing the standard moments with the multibeam moments in each case. The multibeam moments U MBtherm and Q MB therm will differ from the standard moments U therm and Q therm by amounts equal to the pseudothermal part of each standard moment (Table 2). The standard moments will differ for each of the three methods for finding multibeam moments since, in general, f( v ) ≠ f box ( v ) ≠ f fit ( v ). The pseudothermal parts of the standard moments will also differ for each of the three methods for finding multibeam moments. For the FPI case to be addressed in Section (6), the visual method is based on the FPI standard moments file, with beams guessed-at visually by inspecting the velocity distribution, f box ( v ), processed from the FPI skymap. The Cartesian velocity coordinate f box ( v ) is based on the FPI-measured spherical velocity-coordinate skymap, as described in Section (4). 5.1 REDUCED ION VELOCITY DISTRIBUTION, f( v ), FROM PIC
SIMULATION The reduced ion velocity distributions of f box ( v ) for the magnetotail reconnection PIC simulation are shown in Fig. (5), in terms of simulation coordinates, {v x , v y , v z }, and reduced velocity distibutions, F a , F b , F c . The 2D spatial simulation coordinates are {x, y, z}, where z is out-of-plane. The x-velocity, v x , is the component of v colinear with the initial in-plane magntic field; the y-velocity, v y , is in the vertical (other in-plane) direction, and v z is out of plane. Field aligned coordinates are not necessary here since the ions are not strongly magnetized in this region.
18 The reduced velocity distibutions, F a , F b , F c are defined in terms of {v x , v y , v z } as, (9a) F a (v x , v y ) º (9b) F b (v x , v z ) º (9c) F c (v y , v z ) º The three reduced distributions are depicted in three orthogonal planes in Fig. (2), below. The reduced distributions may not be best-suited for visualizing the dominant beam-like portions of f( v ). This is because there may be multiple beams along one direction in a reduced velocity distribution compressed by the velocity integration along that direction. A more comprehensive "fly-around" 3D representation of f( v ) will be used shortly, resulting in different beam centroid velocities and multi-beam moments found from the visual method. Two principal beams are especially apparent in both F a (v x , v y ) and F b (v x , v z ) in Fig. 2. The third reduced distribution, F c (v y , v z ), peaks in a crescent shape, which is harder to interpret (hence, the advantage of a "fly-around" 3D representation of f( v )). Around the pair of peaks is a high-velocity cloud or "halo." Clouds are present and much more extended to high velocities in the MMS f( v ) to be considered in Section (6). This has consequences for the heat flux moment which will be discussed later. dv z f ( v x , v y , v z ) −∞∞ ∫ dv y f ( v x , v y , v z ) −∞∞ ∫ dv x f ( v x , v y , v z ) −∞∞ ∫ Figure 2:
Contour plots of logarithm of reduced velocity distributions of the PIC ion f( v ): ( a) F a (v x , v y ), (b) F b (v x , v z ), and (c) F c (v y , v z ) defined in Eqns. (9). Units of velocities and of velocitiy distributions are arbitrary. Crossing points of horizontal and vertical white lines is origin, {0,0}. Components of the mean flow, u , are indicated by the + sign in each of the three planes. STANDARD MOMENTS
OF f( v ) FROM PIC SIMULATION In standard energy-transport analyses one generally uses the standard decomposed moments of f( v ). The visual method determines the multibeam moments by subtracting the pseudothermal moments from the corresponding standard moments , using the equations in Table (2). The standard moments of f box ( v ) which can be used to find multibeam moments for the PIC simulation are given below in Table (3). Table 3:
Standard moments of f box ( v ) from PIC magnetotail reconnection simulation (Fig. (1)) . The standard mean flow velocity, u is in arbitrary units. The standard pressure tensor, P , standard energy density, U, and standard energy flux vectors, Q bulk , Q enthalpy and Q are in the dimensionless units defined in Table (2), with u n = u . Hence, the standard energy flux vectors are all in units of Q bulk and | Q bulk | = 1. The standard moments u , P, U and Q are all calculated from the same velocity distribution, f box ( v ). Table (4) gives the standard mean flow velocity, u , in arbitrary units as well as the components of the pressure tensor and energy flux vectors (bulk, enthalpy and heat fluxes) in appropriate dimensionless units (Table (2), with u n = u ). Note that the heat flux is negligible in magnitude compared to the enthalpy flux in this PIC tail simulation. This is associated with the relatively narrow velocity "haloes" in f( v ), as seen in Fig. (1). In a later section, MMS-measured velocity distributions in the magnetopause with more extended haloes will be seen to be associated with large heat fluxes. This association arises because heat fluxes are third order veocity moments which are strongly influenced by higher velocities (i.e., on the tail of f( v )). The standard pressure tensor, P , is symmetric, so it possesses only six independent matrix elements instead of nine. The trace of the pressure tensor, Tr P , is about 2 in these units, and is identically equivalent to the dimensionless standard thermal energy density , U therm . (Table (1)). Hence there is no need to calculate U therm separately. Moreover, since the bulk energy density is 1 in these units, the total standard ion energy density moment is U = 1 + U therm = 3.06 Standard moments of f box ( v ) from TAIL PIC SIMULATION x y z Norm, Trace or Eigenvalues Mean flow, u
10 226 -123 | u | =257.5 Presure, P x,i P = U therm = 2.06 Eigenvalues( P ) = {1.09, 0.50, 0.45} Pressure, P y,i -0.208 0.705 -0.220 Pressure, P z,i Q bulk Q bulk | =1.00 Q enthalpy -0.40 3.23 -1.96 | Q enthalpy | = 3.80 Q heatflux -0.01 0.15 0.12 | Q heatflux | =0.19 Q therm = Q enthalpy + Q heatflux -0.41 3.38 -1.85 | Q therm | = 3.87 Q = Q bulk + Q therm -0.37 4.26 -2.33 |Q| = P (both pseudo and true thermal parts), specify the thermal asymmetry of f( v ). The three eigenvalues of P give three generally different pressures, P , P and P associated with each of the three principal axes defined by the directions of the three orthogonal eigenvectors. The nature of the pressure anisotropy is illustrated in Fig. (3). The red, green and blue vectors point in the directions, , of the three pressure eigenvectors, but have magnitudes equal to the square roots of the corresponding eigenvalues, √P , √P , and √P . These define the semi-major and semi-minor axes of a pressure ellipsoid in velocity space. With velocity v in units of the magnitude, u, of the flow velocity, u , The orientation of the ellipsoid in simulation coordinates is oblique to the x-axis. From Table (3), √P = 1.04, √P = 0.71 and √P = 0.67 are the magnitudes of the R, G and B vectors. The magnitudes of the G and B velocity vectors are almost equal, while the magnitude of the R vector is larger, so the pressure ellipsoid in this case is essentially a prolate spheroid . Note there is a related temperature ellipsoid which can be expressed in terms of the three orthogonal temperatures, T = P /n, T = P /n, and T = P /n. This is simply a rescaled version of the pressure ellipsoid. 5.3 PIC f( v ) MOMENTS FROM 3-D VISUAL METHOD AND K-MEANS METHOD We next construct a of the PIC ion velocity distribution, f( v ), whose reduced distributions are displayed in Fig. (2). This and other 3-D velocity contour maps are used to visualize the results of the different methods for approximating f( v ) as a sum of ion beams. The centroid velocities required by the visual method are easier to estimate from the 3-D visualization than from the reduced distributions and the centroid velocities for the other two methods (k-means and least squares) can be displayed for comparison Multi-beam moments are found by each method and organized in Tables for comparison each other and with the standard moments given in Table (3). Multibeam and standard pressure ellipsoids associated with each of the three methods are displayed in the Figures to give a clearer picture of the nature and origin of the pseudothermal nature of the standard pressure tensor for each f( v ). In Section (5.4) we consider multi-beam moments found by the least squares method. In Sect. (5.5) we compare multi-beam energy and energy flux vector moments . The green and yellow surfaces in Fig. (4a) comprise a 3-D view of the ion distribution, f( v ), displayed in Fig. (2) in terms of reduced distributions. In both plots, f( v ) is the distribution ˆ x ', ˆ y ', ˆ z ' { } v x '2 P + v y '2 P + v z '2 P = Figure 3:
Orthogonal eigenvectors of the standard pressure tensor for the tail PIC simulation are in the directions of the red, blue and green vectors. The colored vectors each have a magnitude equal to the square root of the corresponding eigenvalue (Table (3)). Velocities are dimensionless, scaled by the flow speed, u. The R, G, and B vectors define semi-major and semi-minor axes of this pressure ellipsoid. ���� (cid:1) ��
21 in the PIC simulation box. The multi-connectedness of f( v ) is much easier to see in Fig. (4) than in Fig. (1). The bold magenta plus-mark indicates the location of the standard moment flow velo-city, u . The other four plus-marks show the (guessed-at) centroid velocities of the beams in a four-beam approximation to f( v ). Table (4) gives the magnitudes of the five marked velocities. In Fig (4c) the dark green, blue and red clusters at 25% f max are visible, but the black cluster is blocked. In the almost opposite view in (d), the black cluster is visible but the green one is hidden. The light-colored green, blue, red and grey cluster surface contours are at 10% f max . Note the similarity of the cluster shapes to beam-like features in f(v) visible in (a). The standard pressure pink ellipsoid is again shown and seen to contain most of f(v). The multi- Figure 4: 3D visualizations for PIC tail simulation: (a)
Surface contour map of PIC velocity distribution, f( v ). Green surfaces are at 25% of global maximum of f( v ) and yellow surfaces are at 10%. Magenta star marks the standard flow velocity u . Colored plus-marks indicate the four visually-selected beam velocities in Table (4), below. (b ) Pressure anisotropy ellipsoids based on eigenvalues and eigenvectors of pressure tensors for f( v ). Standard pressure ellipsoid is pink; multibeam pressure ellipsoid is black. Standard and individual beam velocities are again shown. False pressure is contained in space between ellipsoids. ( c ), ( d ) Two surface-contour views (almost 180° apart) showing partition of f( v ) into four clusters obtained by k-means method with four beams (colored black, red, green and blue). All properties of the four beams are completely determined by k-means .
22 beam pressure ellipsoid determined by the k-means method is black and oriented almost per-pendicular to the multi-beam pressure ellipsoid determined by the visual method seen in (b)
In Table (5) the four-beam velocity moments found by the visual method are displayed in the dimensionless units of Table (2). The four-beam moments should be compared with the standard moments shown in Table (3). The four-beam thermal energy density in Table (5) is approximately 1, compared to a a normal thermal energy density of approximately 2 in Table (3). Thus, about half of the standard thermal energy density is pseudothermal (and shows up as an increase in the four-beam bulk kinetic energy density, which is not shown). In like fashion, the magnitude of the normal thermal energy flux vector is reduced by about a factor of two, with the difference showing up as an increased magnitude of the four-beam bulk kinetic energy density flux vector. The eigenvalues and eigenvectors of the four-beam pressure tensor enable one to construct a four-beam pressure ellipsoid and compare with the standard pressure ellipsoid in Fig. (3). The two ellipsoids are compared in Fig. (4b). The four-beam ellipsoid (in black) is smaller than the normal ellipsoid (in pink) due to pseudothermal contributions to the standard pressure tensor.
Table 5 : Multi-beam pressure tensor, eigenvalues and trace and energy density flux vectors from visual method applied to 3D f box (v) from PIC simulation. See Fig. (4) above
Figs. (4c) and (4d) are visualizations relevant to the k-means method, which here divides f( v ) into four non-overlapping "clusters." The clusters are colored green, blue, red and black in the PIC simulation 3D f box (v) fitted visually with four beams
Beam 1 Beam 2
Beam 3 Beam 4 x-comp of centroid velocity
140 -240 y-comp of centroid velocity
90 350
125 250 z-comp of centroid velocity -35 -290
300 120
Inferred
Density
Fraction, h j MBx,j P MB = U MBtherm = 0.960 P
MBy,j P MB = {0.502, 0.284, 0.174} P MB z,j -0.042 0.026 0.202 Q MBbulk -0.606 2.420 -1.593 |Q
MBbulk | = 2.960 Q MBtherm = Q MBenth + Q MBhtflux
MBtherm | = 1.993
Q = Q
MBbulk + Q MBtherm --0.373 4.260 -2.323 | Q | = 4.866 Table 4:
Four beam visually-found centroid velocities and etas from 3D f box ( v ) for PIC simulation. See Fig. (3) above.
23 front and back views of Figs. (4c) and (4d). Also shown in these two figures are the standard pressure ellipsoid (in pink), and the smaller four-beam pressure ellipsoid (in black). The four plus marks show the locations of the centroids of the four clusters, which are completely determined by the k-means method rather than the visual method. Shown below in Table (6) are the four-beam velocity moments of f( v ), as found by the k-means method. Comparing with the standard moments in Table (3), once again the four-beam thermal energy density is smaller than the standard thermal energy density. Here they are in the ratio 2.06 to 0.803 so this method yields even more pseudothermal energy than the visual method. Table 6:
Multibeam moments found by k-means partitioning of f box (v) from PIC simulation: four-beam pressure tensor moment, P MB , and four-beam energy fluxes Q MB bulk , Q MB enthalpy , Q MB heatflux , Q MB thermal , and Q MB . PIC k-means 4-beam mo-ments of tail simultion f(v) x comp. y comp. z comp. Related scalar/vector P MBx,i Tr P MB = U MBtherm = MB y,i -0.020 Eigenvalues of P MB = {0.318, 0.269, 0.216} P MB z,i Q MBbulk -0.375 3.079 -1.725 | Q MBbulk | = 3.55 Q MBtherm = Q MBenth + Q MBhtflux -0.004 1.253 -0.630 | Q MBtherm | = 1.40 Q MBhtflux -0.002 -0.013 -0.006 | Q MBhtflux | = 0.014 Q MB = Q MBbulk + Q MBtherm -0.380 4.332 -2.359 | Q MB | = 4.947 24 5.4 PIC f( v ) MOMENTS FROM LEAST SQUARES METHOD Finally, we address the least-squares method of approximating the PIC ion distribution f( v ) as a sum of four kappa-function ion beams, in accordance with Section 4.2, above. Figure (5a) shows the results of the four kappa-function beam least-squares fit to the PIC f( v ) in Fig. (4a). The standard flow velocity, u , is marked by the magenta star and the other colored plus-signs mark the velocity centroids of the four beams whose free parameters have been determined from the least-squares fit. Once again, the pink ellipsoid is based on the standard pressure tensor, and the black ellipsoid is based on the four-beam pressure tensor, with pseudo-thermal pressure in between the two ellipsoids. Figure (5b) exhibits in a different color each of the four kappa-function beams which are summed in Fig. (5a). The velocity centroid of each beam is indicated by a plus-sign of the same color as the beam. Figure
5: 3D visualization of least squares fit to f( v ) from PIC tail simulation: (a) Surface contour map of least squares four-beam fit, f ls ( v ) to PIC velocity distribution, f( v ) shown in image (a) of previous Figure (4). Green surfaces are at 25% of global maximum of f ls ( v ) and yellow surfaces are at 10%. Magenta star marks standard flow velocity u . Colored plus-marks indicate least squares beam velocity centers. Also shown in (a) are pressure anisotropy ellipsoids based on eigenvalues and eigenvectors of pressure tensors, P , for f ls ( v ). Standard pressure ellipsoid is pink; multibeam pressure ellipsoid is black. Pseudothermal pressure is contained in the space between ellipsoids. (b) the four different kappa-functions (each in a different color) from the least squares fit shown in (a). All free parameters have been determined by the least-squares method.
Table (7) gives the velocity moments of the fitted distribution in Fig. (4a). Note the four-beam energy density is 1.16, which is slightly more than half of the standard energy density of the PIC f( v ) given in Table (3). Also note that the multi-beam heat flux is zero, while the normal heat flux is non-vanishing. This is because the kappa-function-beams used in the fit have no heat flux, due to symmetry. This is a minor disadvantage of the assumed fitting functions; it does not signify that the normal heat flux is entirely pseudothermal. Table 7:
Four-beam moments from least squares fit, f fit (v) , to f(v) from PIC simulation: pressure tensor moment, P MB , multibeam bulk K.E. flux vector, Q MBbulk and thermal energy flux vector, Q MBtherm . Tr P MB = Thermal part of U
MBthermal , in dimensionless units.
Another comparison worth making is to contrast the standard moments of the PIC f( v ) with the standard moments of the fitted distribution in Fig. (5a). This is a measure of the suitability of the fit for finding moments of f( v ). The standard moments of the fitted distribution are given in Table (8), below. Consulting Table (3), we see that the standard moments of f fit ( v ) are acceptably close to the standard moments of f( v ). Note that the standard heat flux of f fit ( v ) is not zero. Table 8:
Standard moments of least squares fitted distribution , f fit ( v ) to PIC f box ( v ) PIC least squares 4-beam moments of distribution , f fit ( v ), fitted to PIC f box ( v ) x comp. y comp. z comp. Related scalar/vector P MBx,i P MB = U MBtherm = P MB y,i Eigenvalues( P MB ) = {0.509, 0.395, 0.255} P MBz,i -0.066 Q MBbulk -0.290 2.107 -1.290 | Q MBbulk | = 2.487 Q MBhtflux
0 0 0 Q MBhtflux = 0 Q MBtherm = Q MBenth + Q MBhtflux Q MBtherm | = 2.273 Q MB = Q MBbulk + Q MBtherm -0.245 4.192 -2.193 | Q MB | = 4.737 PIC standard moments of least sqs. fitted distribution , f fit ( v ), for x y z Norm, Trace or Eigenvalues Presure P x,i P = U therm = 2.186 Eigenvalues( P ) = {1.085, 0.615, 0.486} Pressure P y,i -0.162 0.691 -0.219 Pressure P z,i Q bulk Q bulk | = 0.875 Q heatflux -0.090 0.210 0.153 | Q heatflux | = 0.275 Q therm = Q enthalpy + Q heatflux -0.2952 3.415 --1.793 | Q therm | = 3.868 Q = Q bulk + Q therm -0.245 4.191 -2.191 |Q| = v ) FOUND BY DIFFERENT METHODS A visualization of comparative moments is desirable in order to better understand the similarities and differences among the numerical values of moments of the PIC f( v ) in Tables 3-6 and in the Supplement. The pressure ellipsoids in Figs. 3,4 and 5 already constitute a visual comparison of the pressure tensors. The multi-beam pressure ellipsoids found by the three different methods are generally smaller than the standard pressure ellipsoid of Fig. (3) and embedded within it. In this Section we display the scalar values of the bulk and thermal energy densities in terms of a bar graph (Fig. (6)) and the bulk and thermal energy density fluxes in terms of vector diagrams (Fig. (7)). The bottom bars in purple at the bottom of the bar graph show the standard energy densities of the PIC f( v ) given in Table 3. The right purple bar is U thermal and the left purple bar is U bulk . The sum of two is the undecomposed energy density U, which is also equal to the sum of the multi-beam bulk and thermal energy densities shown in the five bars above the purple standard energy densities. The light orange bars are from the visual four-beam moments in Table (5) and the pink bars are from a different visual four-beam. The yellow bars are from a two-beam visual method fit. The grey and light blue bars at the top are from the k-means method (Table 6) and the least squares method (Table 7). In all cases the pseudothermal energy density is equal to the difference between the standard thermal energy density bar and the multibeam energy density bar. The pseudothermal energy density for all cases is roughly 50% except for the two-beam case, where it is closer to 40%. In all cases the lost pseudothermal ���� (cid:1) ��
012 0 1 2
4B Least sqrs4B k - means2B Visual4B Visual - U MBthermal
Standard U thermal U M B bu l k U bulk Figure 6: Comparison of PIC ion energy density moments:
Purple:
Standard bulk and thermal energy density moments from PIC velocity distribution, f( v ); Five colors above purple: multibeam bulk (left) and thermal (right) energy density moments found by visual, k-means and least squares methods. Among these three are different examples of the visual method: the four-beam and two-beam cases based on 2D planar reduced distribution visualization in the Appendix and a four-beam case based on the 3D visualization in Fig. (4a). All (dimensionless) energies are in units of standard bulk energy from PIC data.
27 energy shows up as a gain in bulk energy density. The sum of the bulk and thermal energy densities (standard or multibeam) is 3.1 in all cases. Next, we summarize our results for the standard and multi-beam energy flux vectors associated with the PIC f( v ). Fig, (7) shows the thermal energy flux vector found by using the different methods. Recall that the thermal energy flux vector is the sum of the enthalpy flux vector and the heat flux vector, which is much smaller in magnitude than the enthalpy flux for this PIC f( v ). The black vector is the standard thermal moment, Q therm, from the PIC simulation, f( v ). It has a greater magnitude than the other vectors, whch are all multi-beam thermal flux moment vectors. All of the multibeam vectors, Q MBtherm are bunched together fairly close in angle to the standard Q therm . However, their magnitudes cover a larger range of values. The 2-beam visual method Q is the largest in magnitude of the multi-beam Q MBtherm and the 4-beam k-means method Q is the smallest in magnitude. This is ordering is consistent with the ordering of the thermal energy densities in Fig. (6), in which the PIC thermal energy density, U , found visually with two beams is the largest of the U MBtherm and the U found from k-means is among the smallest. It is tempting to conclude from the ordering of the magnitudes, of the multi-beam thermal energy fluxes, | Q MBtherm |, that | Q therm | - | Q MBtherm |, or | (Q therm - Q MBtherm |, is a reasonable measure of the magnitude of the pseudothermal energy flux. The large magnitude of the 2-beam thermal flux, | Q | leads to the smallest magnitude of the pseudothermal energy flux contained in | Q therm |. Figure 7
Comparison of standard and multi-beam thermal energy flux vectors for PIC f( v ). Standard Q therm is black. Three multi-beam thermal energy flux vectors, Q MBtherm , found by the visual method are displayed in warm colors, orange, magneta and red, representing, respectivey, N = 4 beams and N = 2 beams, from reduced distributions (see Supplement), and N = 4 beams from PIC 3D distribution in Fig. (4a). The vector moment, Q , found from k-means is colored cyan and the Q found from least squares is colored blue.
28 The sum of thermal and bulk energy density fluxes of the PIC f( v ) are shown in the vector diagrams of Fig. (8). The thick arrows show the sum of standard energy flux vectors: the bulk kinetic flux (thick blue) plus the thermal energy density flux (thick red). The heat flux is negligibe compared to enthalpy flux in all of these thermal energy flux vectors. The thin arrows are all four-beam energy fluxes. The thin red, green and tan arrows are the thermal energy flux vectors found, respectively, by the four-beam visual, k-means and least-squares methods. In comparison to the standard thermal energy flux vector they are all smaller in magnitude and make a small angle to it. Each is vectorially added to a different thin blue vector (multi-beam bulk K.E. flux) to equal an undecomposed energy density flux vector (not shown explicitly). These undecomposed energy density flux vectors are almost all the same length, although for the least squares method there is a noticeable difference (gap). Their magnitudes may be found in the Tables. Figure 8: Multibeam total energy flux moment vectors (thin arrows) compared to standard total energy flux moment vectors taken directly from tail PIC simulation (thick arrows) . The flux vectors have all been made dimensionless by dividing by the magnitude of the PIC bulk energy flux vector, mnu /2. The thick blue vector is the standard PIC bulk energy flux , of magnitude one. The standard PIC enthalpy flux is the thick red vector, and the much smaller PIC heat flux is not visible. The thin blue arrows are the multibeam bulk energy flux vectors in the three multibeam approximations. Each thin blue vector is joined to a different-colored thin vector, which is the multibeam thermal energy flux moment (mostly enthalpy flux ), found by a different method: visual method indicated by thin red vector; k-means method indicated by thin green vector and least-squares method indicated by thin tan vector. ���� (cid:1) ��
6. Moments of ion distribution, f( v ) measured by MMS-4 during dayside reconnection In this Section, four-beam moments are calculated for the velocity distribution measured by the MMS Fast Plasma Instrument in the dayside magnetosphere during the MMS-4 reconnection event of October 2015, at 13:07:06 (Burch, et al, 2015). This distribution and the corresponding moments differ in significant ways from the PIC simulation case treated above in Section (5) (see Conclusions). The multi-beam moments calculated by the visual method is based on inspection of the"box"-processed ion velocity distribution, f box ( v ), derived from the FPI skymap (Sec. 4.3.2). However, the multi-beam moments will be compared both to the standard moments of f box ( v ) and to the FPI standard moments file The multi-beam moments calculated by the k-means method is based on the ensemble of equally weighted particles derived from f box ( v ) (Sect. 4.3.3) and compared to the standard moments of f box ( v ) The multi-beam moments calculated by the least-squares method is based on the least-squares fitted ion velocity distribution, f fit ( v ) and compared to the standard moments of f box ( v ). The ion velocity distributions appropriate for each method are summarized in Table (9) Table 9 : Ion velocity distribution used by each method to find multi-beam moments and standard moments for comparison. Dayside reconnection MMS4 Oct. 2015, 13:07:06 (Burch, et al, 2015). Once again
3D fly-around surface-contour plots are employed to visualize both f( v ) and the effective ion " beams" corresponding to each of the three different methods for finding multi-beam moments. Reduced
2D velocity distributions are used below (in Sect. 6.1) to visualize the MMS f box ( v ) and to compare it to the PIC f box ( v ). They are also used in the Appendices to visualize the locations of beam centroid velocities in another way. Field-aligned velocity coordinates are used in both visualizations, with v || parallel to B , v ^ is in the direction of z x B , where z is the GSE coordinate perpendicular to the plane of the ecliptic and v ^ is in the direction of B x ( z x B) . Method for finding multi-beam moments Distribution used to find multi-beam moments Compare multi-beam moments with
Visual method f box ( v ) • FPI standard moments file • standard moments of f box ( v ) k-means method Particle distribution derived from f box ( v ) • standard moments of f box ( v ) Least-squares method
Fitted ion distribution, f fit ( v ) • standard moments of f fit ( v ) • standard moments of f box ( v ) 30 6.1 REDUCED
ION VELOCITY DISTRIBUTION, f box ( v ), FROM MMS-4 The reduced ion velocity distributions of f box ( v ) derived from the FPI Skymap on MMS-4 during reconnection at the magnetopause are shown in Fig. (9), below, in terms of field-aligned coordinates, {v || , v ^ v ^ }, and reduced velocity distributions, F a , F b , F c , defined as (10a) F a (v || , v ^ ) º (10b) F b (v || , v ^ ) º (10c) F c (v ^ , v ^ ) º Fig. (8), displays a logarithmic contour plot of the reduced distributions of the MMS ion f box ( v ) in field -aligned velocity coordinates. The standard moment flow velocity is indicated by the white plus signs. The two dark areas smaller than 200 km/s indicate that the MMS-measured ion distribution, f( v ), is multi-beamed. A high-velocity ion cloud is seen to extend out to 1000 km/s, although the magnitude of the MMS ion f( v ) at such high velocities is down by at least five orders of magnitude. The high-velocity cloud is an important feature not present in the PIC ion f( v ) (Fig. 2). By contrast, the PIC f( v ) falls to zero at around 500 km/s. The presence of a high-velocity ion cloud in the MMS f( v ) requires that special care be taken in applying the methods for finding moments (especially for the least-squares method). The high-velocity cloud also leads to a relatively greater ion heat flux for the MMS f( v ) than for the PIC f( v ). 6.2 STANDARD MOMENTS OF MMS ION VELOCITY DISTRIBUTION Table (9) indicates the three different kinds of standard moments required to compare with the multi-beam moments found by the three different methods for the dayside reconnection event measured onboard MMS-4, Oct. 2015, 13:07:06 (Burch, et al, 2015). These standard moments are given below in Tables 10-12 d v ⊥ f v || , v ⊥ , v ⊥ ( ) −∞∞ ∫ d v ⊥ f v || , v ⊥ , v ⊥ ( ) −∞∞ ∫ Figure 9: Log of MMS reduced ion velocity distributions defined in Eqns. 10 . Velocities are in units of km/s.
Note the high velocity cloud surrounding the beams in the MMS f( v ). It is not present in the PIC f( v ) reduced ion velocity distributions shown in Fig. (2) d v || f v || , v ⊥ , v ⊥ ( ) −∞∞ ∫ Table 10:
Standard moments from FPI moment file during magnetopause reconnection Mean flow velocity, u, in km/s. Dimensionless thermal energy, U thermal in units of U bulk =Mnu /2, pressure in units of mnu and fluxes in units of |Q bulk | = Mnu /2. Vector and matrix components are in field-aligned coordi-nate system in which parallel direction is along B and Perp1, Perp2 are in directions z x B and B x ( z x B ) All three sets of standard moments are in dimensionless units based on the standard U bulk and |Q bulk | from the FPI moments file in Table 10. The thermal energy densities, U thermal from Tables 10, 11 and 12, respectively are 5.0, 5.2 and 4.9. The magnitudes of the undecomposed energy density flux vectors, | Q |, are, respectively, 11.4, 12.8 and 13.24. Thus, there is a 13% difference between the values the sum, | Q |, from the FPI moments file and the value of, | Q | based on f box . Fig. (10) shows the three standard moments, Q bulk , Q enthalpy and Q heatflux which add up to | Q | in the FPI standards moments file and add up to a different | Q | for f box The three different standard moments each have slightly different flow velocities and may be regarded as moments in different frames. This affects U bulk and Q bulk the most (especially for the least squares case) but does not affect the frame invariant moments, U therm , U, Q heatflux and Q , which differ slightly from each other anyway in each of the three cases from FPI standard moments file ⊥1 ⊥2 || Norm, Trace or Eigenvalues Mean flow, u -72.6 126.4 73.4 u = | u | = 163 km/s Pressure P ⊥1,i P ) = {1.914, 1.663, 1.446} Pressure P ⊥2,i P ||,i -0.0973 0.0500 1.533 Tr P = U thermal = Q bulk -0.445 0.774 0.450 | Q bulk | =1.00 Q enthalpy -3.60 6.53 3.80 | Q enthalpy | =8.37 Q heatflux -1.06 2.22 -0.805 | Q heatflux | =2.59 Q thermal = Q enthalpy + Q heatflux Purple 8.75 3.00 | Q thermal | =10.36 Q = Q bulk + Q thermal -5.11 9.52 3.45 | Q | = 11.36 n i Mu = − J / m n i Mu = mW / m , n = 8.15 cm − Figure 10:
Comparison of standard Q bulk , Q enthalpy , Q heatflux and their sums, Q from FPI moments file (thick, Red, Green, Blue) and from f box moments (thin, Magenta, Cyan, Yellow) Table 11:
Standard moments of f box ( v ) (processed FPI) . Mean flow velocity, u , in km/ss. PRESSURE IN UNITS OF mnun AND FLUXES IN UNITS OF | Q BULK | = mnun /2, FROM FPI MOMENTS FILE (TABLE 10). Vector and matrix components in field-aligned coordinate system in which the parallel direction is along B , Perp1 is in direction z x B , and Perp2 is in direction B x (z x B)
Table 12 FPI standard moments of least square fit velocity distribution, f fit (v)
Mean flow velocity, u, in km/ss. Pressure in units of mnu n2 and fluxes in units of |Q bulk | = mnu n3 /2, from FPI moments file (Table 10). Vector and matrix components in field-aligned coordinate system in which the parallel direction is along B, Perp1 is in direction z x B , and
Perp2 is in direction
B x (z x B)
Standard moments of f box ( v )- processed FPI skymap ⊥1 ⊥2 || Norm, Trace or Eigenvalues Mean flow, u -83.40 130.8 81.6 u = | u | = 175 km/s Pressure P ⊥1,i P ) = {1.991, 1.795, 1.438} Pressure P ⊥2,i P ||,i -0.096 0.112 1.512 Tr P = 5.223 = U thermal Q bulk -0.592 0.928 0.579 | Q bulk | =1.24 Q enthalpy -4.42 7.22 4.41 | Q enthalpy | = 9.55 Q heatflux -1.73 2.28 -1.79 | Q heatflux | = 3.38 Q thermal = Q enthalpy + Q heatflux -6.15 9.50 2.62 | Q thermal | = 11.62 Q = Q bulk + Q thermal -6.74 10.43 3.20 | Q | = 12.82 FPI Standard moments of f fit ( v ) = least squares distribution fit to f box ( v ) ⊥1 ⊥2 || Norm, Trace or Eigenvalues Mean flow, u -85 144 104 u = | u | = 197 km/s Pressure P ⊥1,i P ) = {1.889, 1.723, 1.322} Pressure P ⊥2,i P ||,i -0.145 0.123 1.455 Tr P = U thermal = 4.933 Q bulk -0.761 1.289 0.933 | Q bulk | =1.764 Q enthalpy -4.265 7.530 5.372 | Q enthalpy | = 10.185 Q heatflux -1.235 2.031 -2.018 | Q heatflux | = 3.118 Q thermal = Q enthalpy + Q heatflux -5.500 9.560 3.355 | Q thermal | = 11.528 Q = Q bulk + Q thermal -6.261 10.850 4.287 | Q | = 13.24 33 MMS f( v ) VISUAL METHOD AND k-MEANS MULTI-BEAM MOMENTS As in Fig. 3 for the PIC f( v ), Fig. 11 shows a 3-D velocity contour map of the MMS ion velocity distribution , f( v ), whose reduced distributions are displayed in Fig. (9). Fig 11a shows a surface contour map of the MMS f box ( v ) with the standard flow velocity and four visually-chosen beam centroids marked. In Fig. 11b the standard and visual multi-beam pressure ellipsoids are displayed. They are closer together in this MMS case than in the PIC analysis (Fig. 3b). Figs. 11c and 11d give front and back views of the clusters (and their centroid velocities) found by the k-means method as well as pressure ellipsoids. Multi-beam moments found by each method are organized below in Tables (13-14) for comparison with each other and with the standard moments given in Tables (10-12). In Section Figure 11: MMS f(v), 3D visualizations. ( a): Surface contour map of MMS f box ( v ). Green (inner) surfaces are at 10% of global maximum of f( v ) and yellow (outer) surfaces are at 1%. Magenta star marks standard flow velocity u . Colored plus-marks indicate four visually-selected beam velocities. (b ): Pressure anisotropy: standard pressure ellipsoid is pink; multibeam pressure ellipsoid is black. ( c ) and ( d ): Two surface-contour views (almost 180° apart) for partition of f( v ) into four clusters obtained by k-means method. Inner contours are at 10% f max . Note that there is no 10% surface for two of the k-means clusters. Outer contours are at 1%. Pressure ellipsoids are also shown.
34 (6.4) we consider moments found by the least-squares method; in Sec. (6.5) we compare multi-beam and standard energy and energy flux moments found by the different methods. More detailed analysis based on reduced distributions can be found in the Supplement. The visual method is based on the FPI standard moments file. From Table (13) it is seen to yield a 4-beam thermal energy density, U = 4.42 , which is less than the FPI moments file stan-dard thermal energy density, U therm = 5.02. The difference of D U = 0.60, means that about 12% of the thermal energy in the FPI moments file is pseudo thermal energy.
This is a small correction to U therm compared to the PIC pseudothermal energies found in Section 5. If the 4-beam energy density U is compared to the moments of f box instead of to the moments in the FPI moments file, the pseudothermal energy is even greater (see Table (9)) Note that the total 4-beam energy flux found by the visual method is |Q MB | = 11.35, which compares very well with the the total standard energy density flux from the FPI moment file (Table 10), | Q | = 11.34. However, the magnitude | Q | of the undecomposed energy density flux vector is 12.82, indicating again that care must be taken in comparing with the moments of f box . Table 13. Four-beam moments of MMS f(v) from visual method (based on reduced distribution visualization in Appendix also displayed in Fig. 11a). Pressure in units of mnu n2 and fluxes in units of |Q bulk | = mnu n3 /2, from FPI moments file (Table 10). Table (14) shows multi-beam moments found by the k-means method with k =4. Two of the four clusters and their centroids are visible in Figs. (11c, d). The pseudothermal parts of the thermal moments are generally larger for this method than for the visual method. This can be seen by comparing the four-beam energy density of 3 from Table (14) with the four-beam energy density of 4.2 from Table (13) and by comparing the pressure eigenvalues in Table (14) with those in Table (13). ⊥1 comp. ⊥2 comp. || comp. Related scalars P MB ⊥1,i -0.03 Eigenvalues of P MB = {1.65, 1.44, 1.33} P MB ⊥2,i -0.05 1.56 0.135 P MB ||,i -0.03 0.14 1.41 Tr P MB = U MBtherm = 4.42 Q MBbulk -0.67 1.30 0.72 | Q MBbulk |= 1.63 Q MBtherm = Q MBenth + Q MBhtflux -4.44 8.22 2.74 | Q MBtherm | =9.74 Q MB = Q MBbulk + Q MBtherm -5.11 9.52 3.46 |Q MB | = 11.35 35 Table 14:
Multi-beam moments of MMS f box ( v ) from k-means method with four beams (k =4). Pressure in units of mnu n2 . Fluxes in units of |Q bulk | = mnu n3 /2, from FPI moments file (Table 10). MMS f( v ) MULTI-BEAM MOMENTS FOUND BY MODIFIED LEAST SQUARES METHOD Finally, we describe how multi-beam moments of the MMS f( v ) are found using a modified least squares method. Here we introduce for the first time a novel way of dealing with the high energy cloud in the MMS f( v) velocity distribution shown in Fig. 8. Modified least-squares fit to an MMS f( v ) with a high-energy cloud. • First the cloud, f cloud , is modeled from the MMS velocity distribution, f( v ) by fitting the logarithm of the distribution at high velocities to the logarithm of the sum of three kappa distributions (beam 1, beam 2 and beam 3). The sum of these three distributions is identified as f cloud ( v ). • Next, the cloud distribution is subtracted from f(v) and the difference [f(v)-f cloud ( v )] is fit separately with a sum of three different kappa functions (beam 4, beam 5 and beam 6), which are identified as core beams to differentiate them from the cloud . • Rather than adding together the standard moments of beams 1, 2, 3, 4, 5 and 6 to get (six-beam) multi-beam moments, a new strategy is employed to improve the cloud modeling: • One standard moment of the ensemble of the three cloud beams , 1, 2, and 3 is taken and added to the three standard moments of each of the core beams
4, 5, and 6 to produce (four-beam) multi-beam moments of f( v ). • This has the advantage that the cloud contribution to the multi-beam moments of f( v ) has a non-vanishing standard heat-flux moment, whereas there is no standard heat-flux for any of the single kappa-function beams, due to their symmetry. The beam centroid velocities and partial densities from the fits are given explicitly in Table (15), and in the number-density flux-vector diagram in Fig. (12). Note that n c u c = n u + n u + n u k-means four-beam moments of MMS f(v) ⊥1 comp. ⊥2 comp. || comp. Related scalar/vector P MB ⊥1,i Eigenvalues of P MB = {1.154, 1.017, 0.848} Tr P MB = U MBtherm = 3.02 P MB ⊥2,i MB ||,i -0.084 0.051 0.919 Q MBbulk -3.466 5.581 2.135 | Q MBbulk| = Q MBtherm = Q MBenth + Q MBhtflux -3.378 4.900 0.998 | Q MBtherm | = Q MBhtflux Q MBhtflux | = Q MB = Q MBbulk + Q MBtherm -6.844 10.478 3.133 |Q MB | = 12.902 36 is the standard number density flux of the set of three beams, 1, 2, and 3 representing the cloud. Also, n u is the standard number density flux of the six beams representing the modified least squares fit to the MMS f( v ). Table 15:
Centroid velocities and fractional densities of the beams in a (modified) least-squares fit to the MMS f(v). Beams 1, 2, and 3 belong to the cloud, the cloud flow is from the standard flux of the three cloud beams and beams 4, 5 and 6 are core beams.
Fig. (13) exhibits the shapes of beams 4, 5, 6 and of the cloud (sum of beams 1, 2, 3) in this modified least-squares fit and their velocity centroids, u , u , u and u c , together with the multibeam and standard pressure ellipsoids. The multi-beam (four-beam) moments of the modified least squares fit, f fit ( v) to the MMS f box ( v ) are given in Table. (16). The multi-beam energy density is 4.5, which is close to the value 4.4 from the visual method (Table 13) and above the value 3 from the k-means method (Table 14). The pressure tensor eigenvalues are closer to those obtained by the visual method than to those obtained by k-means. The magnitude of the heat-flux is 1.9, which is higher than the k-means heat flux, 1.3, but lower than the standard heat flux, 3.1, from the Standard Moment Table (12) of the fit distribution , f fit , or the standard heat flux moment magnitude of f box (Table 11), with value 3.4. Hence, there is approximately 40% pseudothermal heat flux magnitude by either of these measures! Although the standard heat flux from the MMS FPI moment file has a magnitude of only 2.6, comparing it with the magnitude of the least squares multi-beam heat flux of 1.9 there is a 27% pseudothermal heat flux magnitude. Note all energy moments of MMS ion velocity distribution are normalized in the same way (see Fig. captions). As it should be, the magnitude of the undecomposed energy flux Q MB is equal to Q = 13.24, from the Standard Moment Table (12). Least-squares fit to MMS f( v ) Mean flow, u Beam1 u Beam2 u Beam u Cloud flow , u c Beam u Beam u Beam u ⊥1 -85 -107 -141 -133 -114 52 39 -136 ⊥2
222 193 142 215 -18 || -37 98 157 90 92 NORM
261 245 218 236 165 h j Figure 12:
Density flux vectors for the six beams in Table 15 and standard moment flux vectors of cloud (n c u c ) and overall (n u ). n c u c n u Table 16:
Multi-beam moments of modified least-squares fit , f fit ( v ) to MMS f( v ). Pressure in units of mnu n2 and fluxes in units of |Q bulk | = mnu n3 /2, from FPI moments file (Table 10). Four-beam moments of least squares fit, f fit , to MMS f box (v) ⊥1 comp. ⊥2 comp. || comp. Related scalar/vector P MB ⊥1,i Eigenvalues of P MB = {1.663, 1.588, 1.220} P MB ⊥2,i MB ||,i -0.184 0.117 1.439 Tr P MB = U MBtherm = Q MBbulk -0.979 1.937 1.225 | Q MBbulk | = MBenth -5.092 8.576 4.915 | Q MBenth | = 11.119 Q MBhtflux -0.188 0.331 -1.851 | Q MBhtflux | =1.890 Q MBtherm = Q
MBenth + Q
MBhtflux -5.280 8.910 3.060 | Q MBtherm | = MB = Q MBbulk + Q
MBtherm -6.261 10.850 4.287 | Q MB | = 13.24 Figure 13: 3D visualizations of modified least-squares fit to MMS f(v) in different colors. ( a): Surface contour map of fitted distribution. Green surfaces (core beams) are at 10% of global maximum of f( v ) and yellow surface (high energy cloud) at 1%. Magenta star marks standard flow velocity u . Three colored +marks indicate three core beam velocities and black +sign indicates cloud standard flow velocity, u c . Standard pressure ellipsoid is pink; multibeam ellipsoid is black. ( b ) Each of four "beams" shown in different color, with 1% outer surface in lighter color for three core beams.
38 6.5 COMPARISON OF MULTI-BEAM MOMENTS OF MMS f(v) Here we summarize the results and significance of the three methods for finding multi-beam moments of the MMS ion velocity distribution, f( v ). There are three different standard moments: those from the MMS FPI Moments File (Table 10), those based on f box (Table 11) and those based on f fit (Table 12). The FPI Moments File is based on the FPI Skymap, which includes higher energy ions than in f box or f fit . (The higher energy ions must be arranged in velocity-space to make the magnitude of the standard heat flux from the FPI Moments File (Table (10)) lower than the magnitude of the standard heat flux for f box ( v ) (Table 11).) To carry out the visual method for determining multibeam moments values of U thermal , P and Q therm the FPI Moments File (Table 10) is used as for the standard moments in Table (2). This guarantees that the undecomposed multi-beam moments found by the visual method are identical to the corresponding undecomposed standard FPI Moment File moments. This is an advantage of the visual method, but it also has disadvantages which will be discussed later. In Fig. (14) scalar multi-beam energy densities, U MBtherm and U
MBbulk found by the three methods are compared with each other and with the standard energy densities from the FPI Moments File, from f box and from f fit (the least squares fit to f box ). Thermal energy densities are on the right, and bulk thermal energy densities are on the left. For each pair, the undecomposed
U is the sum of U
MBtherm and U
MBbulk (i.e., the sum of the lengths of the right and left bars of a given color). Comparing these energy densities with those in Fig. (6) for the PIC f( v ), we see that there is far more false thermal energy in the former. This is because there are fewer high energy ions (i.e., no halo) in the PIC f( v ) and more lower energy beam-like structures. Subtracting the visual U from the standard U thermal from the FPI Moments File, using Tables (10) and (13) gives a pseudothermal energy density, D U thermal = 5 - 4.4 = 0.6. This is 12% of U thermal . The pseudothermal energy lost from U thermal adds to U bulk . In that sense, D U thermal can be thought of as a (false) deficit in the kinetic energy in U bulk . Since U bulk is 1 in our units, U MBbulk becomes 1.6, a significant 60% correction to U bulk . The value of U found by the least squares method is similar to the value found by the 4-beam visual method. However, the value found by the k-means method is quite different, possibly due to the lack of veloocity-space overlap in the four clusters (beams).
Figure 14:
Standard and multibeam thermal and bulk energy densities. All (dimensionless) energies are in units of standard bulk energy from the FPI moments file. U therm FPI moment file U therm f box standard U therm l-sq. standard U visual 2beam U visual 4beam U least square U k-m U bulk FPI U bulk f box U bulk momU U U U
39 In in order to understand the dependence of multi-beam moments on the number, N, of beams chosen, we have included in Figure (14) the results of the visual method after removing two of the four beams beams marked in Fig. (11a) with + signs. We keep only the two beams with centroid velocities in the green peaks in the MMS f( v ). The visual method now yields values for U and U , shown by the brown bars in Fig. (14). U is now 4.6 instead of 4.4 and U is now 1.45 instead of 1.6. Hence, there is slightly less pseudothermal energy density and slightly more false bulk K.E. in the undecomposed energy density, U Next, we compare standard and multi-beam energy density flux moments of the MMS velocity distribution, f( v ). These moments are vectors rather than scalars. Emphasis will be on thermal energy density flux, Q MBtherm , and multi-beam bulk kinetic energy density fluxes, Q MBbulk . The significance of Q MBthermal is explained in Goldman, et al (2020). The energy density flux moments, Q MBtherm , are higher order in velocity than the energy density moments, U
MBtherm , so the effect of the high energy cloud should be more pronounced. In Fig. (15) we plot vector diagrams showing how the undecomposed multi-beam energy flux moments, Q MB found from the MMS f( v ) by each of the three methods might be decom-posed. The decomnposition is into a sum of a bulk (kinetic) energy flux, Q MBbulk and a thermal energy flux, Q MBtherm , so that Q MB = Q MBbulk + Q MBtherm . The least-squares and k-means Q MBtherm vectors are not shown explicitly in Fig. (15). However, Q MB enthalpy and Q MBheatflux , which add together to form Q MBtherm , are each shown explicitly. 40 The blue vectors are all bulk kinetic energy fluxes. The undecomposed Q MB found by the visual method is the sum of the green Q MBthermal vector and the blue bulk K.E. vector to which it is attached. This Q MB is identical to the undecomposed Q from the FPI Standard Moments File found by adding together the FPI standard moments, Q bulk , Q enthalpy , and Q heatflux (the sum of the smallest blue vector, the red vector and the purple vector). The least squares and k-means methods do not have undecomposed Q MB equal to the undecomposed Q from the FPI Moments File. Instead, each is equal to the undecomposed Q based on the standard moments of f fit ( v ) and f box ( v ), respectively. These Q vectors are not shown in Fig. (15); if shown, they would connect to the end of the least squares heat flux vector (red) or the end of the k-means heat flux vectors (magnenta), Figure 15:
Comparison of energy flux vectors to standard energy flux vectors from FPI moments file . Units of energy fluxes are | Q bulk | from standard FPI moments file. Blue vectors are bulk
K.E. fluxes from FPI Standard Moments File; smallest is Q bulk of magnitude 1, which adds to red Q enthalpy and purple Q heatflux to yield Q (not shown), all from standard FPI moments file. Next longest blue vector is multibeam Q from visual method, which adds to Q to yield Q . Remaining two blue Q bulk vectors are joined to Q and Q , found by different methods: the k-means fluxes in yellow and magenta and the least-squares in brown and pink. These add up to different Q -vectors (not shown) v || v ^ v ^ k - m ea n s Q B e n t h a l py l ea s t s qu a r e s Q B e n t h a l py Q h e a t f l u x Q e n t h a l py v i s u a l Q B t h e r m a l
41 The multi-beam heat flux moment, Q MBheatflux , cannot be found by the visual method. The Q MBheatflux vectors found by least squares (red) and by k-means (magenta) point in very different directions from each other and from the direction of the standard heat flux vector (purple) from the FPI Moment File and other standard moments (Table 11 and 12). The magnitude of the measured standard ion heat flux, | Q heatflux | in Fig. (15) is much less than the magnitude of the measure ion enthalpy flux, | Q enthalpy |, a result also found in other magnetospheric measurements (Eastwood, et al, 2013). For this and other reasons we will not dwell on Q MBheatflux , but will, instead, compare and interpret thermal flux vectors without decomposing them into enthalpy flux and heat flux vectors. Since thermal flux vectors are not frame-invariant, the amount of pseudothermal energy flux will be a frame-dependent vector. (In the frame with zero flow velocity, Q thermal consists entirely of heat flux and there is no bulk kinetic energy flux.) In Fig. (16), the three standard Q thermal vectors and the multibeam Q MB thermal vectors are displayed together. They are confined to a relatively narrow cone of angles, unlike the large angular spread of the heat flux vectors in Fig. (15). In addition to the four-beam thermal flux vectors, constructed by the visual, k-means and least squares methods, Fig. (16) shows another thermal flux vector, the green Q vector constructed by the visual method with two beams rather than four, as discussed in connection with the scalar energy densities shown in Fig. (14). It is close in magnitude and direction to the Q vector from the visual method. Once again, the multi-beam moment found by k-means (magenta) is an outlier —in angle as well as in mag-nitude. However, the Q found by the least-squares method (brown vector) is not very different from Q or Q found by the visual method. Among the standard Q therm moments, the moment from the FPI Moment File has the smallest magnitude. The Q thermal of the box distribution f box ( v ) and the Q thermal of the least-squares fitted distribution, f fit ( v ), both have larger magnitudes. Thermal flux vector moments, MMS f(v)
Figure 16:
Three standard Q therm fluxes (black, tan and yellow) from FPI moments file, f box ( v ), and f fit ( v ) are shown together with a four-beam and a two-beam Q MBtherm, found by the visual method (red and green) and two Q found by k-means and least squares (l-sq). methods (magenta and brown). SeeTables 10, 11, 12, 13, 16 and 14.
42 In Fig. (17) and Table (17) we see the pseudothermal flux vectors determined by sub-tracting the multi-beam thermal flux from standard thermal flux. The change in the corresponding bulk kinetic energy flux vectors is in Table (17).
Table 17. Magnitudes and angles for Fig. (17)
In Fig. (17) four pseudothermal fluxes are shown: The red, green and purple, Q th - Q , Q th - Q and Q t h -Ls Q vectors show the pseudothermal flux vector parts of the standard moment Q th from the MMS FPI Moments File, based respectively on a visual four-beam, a visual two-beam and a four-beam least squares construction. The blue vector, fit Q th - Ls Q , shows the pseudo thermal flux vector part of the standard moment Fit Q th of the least squares fit distribution, f fit ( v ) based on the least-squares four-beam fit thermal flux, lsq Q . Table (17) shows the magnitudes and orientations of these pseudothermal flux vectors. The magnitudes of the pseudothermal vectors are all around 0.6 or 0.7, more than an order of magnitude smaller than magnitudes of the standard thermal moments , | Q th | ~ 10 to 12 (from Tables 10, 11, 13). Thus, when an asymmetric high energy cloud is present there is not much pseudothermal energy flux; the cloud contains true thermal energy flux. The magnitude of the pseudothermal vector part is a much larger fraction of the magnitude of the standard bulk KE fluxes, as is evident from Table 17, where the fraction is seen to be around 60%. Several other conclusions can be drawn from Table 17: • The (green) pseudothermal part of Q therm based on the two-beam visual method makes an approximate 45° angle to the (red) pseudothermal part of Q therm based on the four-beam visual method. • The (purple) pseudothermal part of Q therm based on the least squares method makes a126° angle to the (red) pseudothermal part of the thermal moment, Q therm based on the four-beam visual method, indicating that this is not a good measure of the pseudothermal part. • The (purple) pseudothermal part of fit Q therm (the standard thermal flux moment of f fit ( v )) based on the least squares method makes a tiny 4° angle to the (red) pseudothermal part of Q therm based on the visual method and is around the same magnitude, meaning that this measure of the pseudothermal part is in good agreement with the visual method for finding the pseudothermal part of Q therm Pseudothermal flux Mag-nitude % stndrd therm flux Angle with Q % stndrd bulk KE flux Q th - Q Q th - Q Q th - Ls Q Q th -Ls Q Figure 17.
Pseudothermal flux parts of Q th , from FPI Moments File and of standard Fit Q th of least squares fit f fit ( v ) to f box ( v )
7. Overall Summary
After reviewing concepts and methods developed in this paper and in the earlier paper, Goldman, et al, (2020), we summarize key results. Multi-beam moments of complex velocity distributions, f( v ), taken both from simulations and from MMS reconnection measurements during reconnection at the magnetopause are described and interpreted. We then discuss implications for emerging new measurements of multi-peaked and other complex velocity distributions. 7.1 METHODS 1. Standard energy moments (Table 1) of complex ion velocity distributions, f( v ), can contain pseudothermal parts, which may be misinterpreted as true temperature, pressure, etc. Pseudothermal parts are equivalent to deficits in the standard moments of both bulk kinetic energy and kinetic energy flux. Complex velocity distributions include f( v ) with pronounced multiple peaks ("beams") and f( v ) with extended two and three dimensional-ly shaped peaks ("ridges"). An alternate way of taking (multi-beam) energy moments of f( v ), with no pseudothermal parts, is achieved by first approximating f( v ) as a sum of N "beams," f( v ) = f ( v ) + ... + f N ( v ). (Goldman, et al, 2020). Reasonable values of N are low single digits (Dupuis, et al, 2020) 2. Three methods developed for approximating a complex f( v ) as a sum of N beams (Goldman, et al, 2020) are applied to an f( v ) measured in a PIC-simulation and to an MMS-measured f( v ). They are o A visual method, in which one estimates the N beam centroid (flow) velocities o A k-means method in which f( v ) is divided into N non-overlapping clusters o A least-squares method in which each beam, f j ( v ), is expressed parameters found from a least-squares fit, f fit ( v ) = [f ( v ) + ... + f N ( v )] to f( v ) 3. In the visual method, multi-beam thermal moments of f( v ) are constructed as follows: o With N chosen to be four or less, one visually estimates the N beam velocity centroids, u j , by referring to a 2D or 3D graphic display of f( v ) (i.e., reduced f( v ) planes, 3D fly-around, etc.) o The N beam-densities, n j , are found from the N velocities, u j , the standard flow velocity, u , and the standard density, n (Section 3). o The densities and velocities are then used (via Table (2)), together with known standard moments, U therm and Q therm to find multi-beam thermal moments, U MBtherm and Q MBtherm . o The multi-beam bulk energy moments U
MBbulk and Q MBbulk are then found by subtracting U
MBtherm and Q MBtherm from the (known) standard undecomposed, 44 moments, U and Q . The multi-beam pressure tensor, P MB , is found in an analogous manner. o Advantages of the visual method : • Intuitive • Fast, once the centroids are specified (which could be time-consuming). • Computer usually not needed to find multi-beam moments of f( v ). • Method directly gives pseudothermal parts of a given set of standard thermal moments. • Standard flow velocity, u , and density, n, of f( v ) are pre-specified, so electron or ion beam fluxes are guaranteed to add up to n u. o Disadvantages of the visual method : • Must estimate beam centroid-velocities, u j from 2D or 3D images of f( v ) while avoiding choices leading to negative beam densities or negative multi-beam pressure eigenvalues. • Number, N, of beams must be be four or less. • Only beam densities and centroid velocities are determined explicitly. Beam shapes are never fully determined. • Multi-beam enthalpy flux , Q MBenthalpy , and heat flux , Q MBheatflux cannot be found separately; in the visual method, only their sum is found: the thermal energy flx, Q MBthermal = Q MBenthalpy + Q MBheatflux . 4.
In the k-means and least- squares methods multi-beam moments of f( v ) are determined by taking a standard moment of each of the N beams in f ( v ) + ... + f N ( v ) and adding. (Goldman, et al, 2020). N can be any integer. Beam shapes are fully determined. o Advantages of k-means method • Method is systematic once N (= k) is chosen, and may be readily automated • Method yields a lowest energy partition of f( v ) o Disadvantages of k-means method • The beams it yields are non-overlapping, non-interpenetrating clusters arising from the partitioning of f( v ) • The standard moment is that of the f box ( v ) used to generate particles but does not depend on the choice of N or on the resulting cluster properties • Method may not yield reasonable-looking beam centroid velocities, especially if an inappropriate value of N is chosen. (Note, this may be remedied with a special treatment of clouds for k-means now being developed) o Advantages of least-squares method • Method is systematic once N is chosen, and method may be automated • Beams can overlap and interpenetrate unlike in the k-means method • Method can be modified to treat high energy velocity space clouds in f( v ) 45 o Disadvantages of least squares method • Method requires choosing parametrized beam functions of velocity which may not be sufficiently general. • Standard moments are an output, not an input. They are moments of a fitted distribution, f fit ( v ), and may differ from the standard moments of the measured f( v ), especially when high energy clouds are present in f( v ) . • High energy clouds in f( v ) may require special treatment when specifying the beam decomposition. • In such cases the heat flux moments are of f fit ( v ), not of f( v ). 5. Method-specific modified versions of MMS f( v) (Table 9) o The FPI standard moments file is used in the visual method . o A smoothed Cartesian f box ( v ) generated by interpolation from the MMS Sky Map is used to find beam velocities in the visual method . o f box ( v ) is also used to generate the particle distribution in the k-means method. o In the least squares method , f fit ( v ) is the least squares fit to f box ( v ). 6. Broad high-energy clouds in velocity distributions such as the MMS f( v ). o The high-energy cloud in the MMS-measured f( v ) can lead to non-negligible differences among the above modified velocity distributions. o The cloud in f fit ( v ) may be better-represented by more than one beam in the least squares method of finding multibeam moments. (Table 15). This is especially true when the beams in the fit are taken to be symmetric functions of velocity and therefore lack heat flux moments. o High-energy clouds in f( v ) may contribute significantly to the heat flux moment. o There are significant differences in the relative magnitudes of the energy density and energy flux density moments associated with the MMS f(v) in comparison with those associated with the PIC f( v ) (Table 18 and next section). 46 Table 18:
Comparison of multibeam energy and energy flux moments of PIC f(v) and of MMS f( v ) Issues PIC simulation of tail reconnection MMS FPI measurement during dayside reconnection
Is there a significant high energy cloud ? No. Yes How large are the ion pseudothermal energies , D U therm as a fraction of U therm? visual method 1.3/2.1 = visual method .58/ 5 = k-means 1.3/2.1= k-means 2.2/5.2 = least squares 1.0/2.2 = least squares 0.45 / = 9.3% How large a fraction of the thermal flux is the heat flux, | Q heatflx | /| Q thermal |? PIC f box ( v ) standard moment 0.19/3.9 = FPI moment file 2.6/8.7 =
Magnitudes | Q |, (determined by different velocity distributions) (Table 9) visual method 4.86 visual method 11.34 k-means 4.95 k-means 12.82 least squares 4.73 least squares 13.24 47 7.2 QUANTITIVE FINDINGS REGARDING MOMENTS Standard energy moments of complex ion distributions often do not distinguish properly how much of an energy density moment is thermal (incoherent) and how much is bulk-kinetic (coherent). A classic example is the complex distribution, f( v ), consisting of a system of two equal and opposite cold beams of electrons or of ions. Standard energy moments of f( v ) indicate that the system has no kinetic energy density and all thermal energy density, whereas the opposite would appear to be more physical. Of course, this difficulty is easily resolved by taking the thermal and kinetic energy density moments of each beam and adding to obtain the “ multibeam ” energy density moments of this two-beam system. The multibeam thermal energy density moment, (U therm ) MB does not have any of the pseudothermal energy density present in the standard thermal energy density moment, U therm . In this paper we have taken generalized multibeam moments (Goldman, et al, 2020) of a complex ion velocity distribution, f( v ), in a PIC simulation of magnetotail reconnection and of another f( v ) taken from an MMS-satellite-measured f( v ) during reconnection in the dayside magnetosphere. The multibeam thermal energy density moments, (U therm ) MB and the multibeam thermal energy density flux moments, ( Q therm ) MB = ( Q enthalpy ) MB + ( Q heatflux ) MB, will lack the pseudothermal parts of the standard moments, , U therm and Q therm but will retain any truly incoherent thermal parts contained in the standard moments. The pseudothermal energy density is defined as, D U therm = U therm - (U therm ) MB and the pseudothermal energy density flux vector as D Q therm = Q therm - ( Q therm ) MB . The pseudothermal part of the standard thermal energy density of f( v ) which has been removed in the multibeam thermal energy density shows up as multibeam bulk-kinetic energy density, so that the sum of thermal and bulk-kinetic energy density moments is the same whether the moments are standard or multibeam. One implication of our finding is that the traditional standard way of taking moments of a complex ion distribution, f( v ), at a given location and time, can exaggerate the standard thermal energy density, U therm , of f( v ), and hence exaggerate the ion temperature , T i = U therm /n at that location and time. In addition, the multibeam bulk kinetic-energy density, (U bulk)MB will be increased from U bulk by the same amount, D U therm , that (U therm ) MB is decreased from U therm . The same holds for the pseudothermal energy density flux vector, D Q thermal . To obtain multibeam moments, D Q thermal is subtracted from Q thermal and added to Q bulk . Multibeam energy moments therefore help sort out how much of a fluid element’s energy at a given place and time is random (characterized by a true temperature) and how much of it is directed kinetic energy. This is significant because a set of multibeam energy moments taken at nearby places and times could help distinguish particle heating from particle acceleration , A key issue is: by how much do multibeam thermal moments differ quantitatively from standard thermal moments, both for the PIC f( v ) analyzed in Sec. 5 and for the MMS f( v ) analyzed in Sec. 6? In other words, how large are pseudothermal effects in the two cases and can they be com-pared? In Fig. 6 the standard and multibeam thermal energy densities for the PIC f( v ) found by the various methods were compared (see units in Table 2). In most cases four beams were assumed in taking multibeam moments. In all cases the multibeam thermal energy densities are roughy 48 half the standard thermal energy densities (also see Table 18). By contrast, the standard and multibeam thermal energy densities for the MMS f( v ) (Fig. 14), are fairly close to each other. From Table 18 we can get some idea of differences between the energy moments of the PIC simulation f( v ), which does not have a significant high-energy ion “ cloud ” and the MMS-measured f( v ) which does have a significant high-energy ion cloud. For the PIC f( v ), the ion pseudothermal energy density, D U therm , (based on four ion beams) was found to be around 50% - 60% of the standard energy density, U therm , whereas for the MMS velocity distribution the ion pseudothermal energy density, D U therm , it was around 10% of the standard energy density, U therm (assuming the outlying k-means result is less reliable). A probable explanation is that the cloud possesses a large amount of true thermal energy density that dominates both the standard and multibeam thermal energy density moments. In this scenario the (smaller) pseudothermal energy density would arise mainly from the lower-energy, colder beams. Figures 7 and 8 show the standard and multibeam thermal and bulk-kinetic energy density flux vectors found by the three different methods, as applied to the PIC f( v ). The directions of the vectors are not significantly different from each other but the magnitude of the PIC standard bulk kinetic energy density flux vector, Q bulk is significantly smaller than the magnitudes of the corresponding multibeam vectors, ( Q bulk ) MB , found by the various methods. In addition, the magnitude of the standard thermal energy density flux vector, Q therm (mainly enthalpy flux for the PIC case), is much larger than the corresponding multibeam vectors, ( Q therm ) MB . Figures 15,16,17 show the standard and multibeam thermal and bulk-kinetic energy density flux vectors found by the three different methods as applied to the MMS f( v ). The differences here are related to the high-energy ion velocity cloud in the MMS f( v ), which was absent from the PIC f( v ). The main difference is that the MMS standard heat flux, Q heatflux is no longer negligible in comparison to the standard enthalpy flux, Q enthalpy in the sum, Q therm = Q enthalpy + Q heatflux . This is evident in Figure 15. It is quantified in Table 18, which shows, for the PIC case, that, | Q heatflux | is 5% of | Q thermal |, so that | Q thermal | is dominated by | Q enthalpy |. However, for the MMS case, | Q heatflux | is significant: It is equal to 30% of the magnitude of | Q thermal |. However, the heat flux can only be enhanced by a cloud if the cloud is highly asymmetric or if its velocity centroid is sufficiently skewed from the standard flow velocity. In Figure 17 and Table 17 MMS pseudothermal flux vectors, D Q therm = Q therm – ( Q therm ) MB , found by the visual method (with 4-beams and with 2-beams) and by the least-squares method with 4-beams are compared. There are once again complications having to do with representing the high-energy cloud as a beam and with differences between the modified distributions. However, the magnitudes of all of the pseudothermal vectors, | D Q therm | turn out to be comparable (~0.6) and ~5% of the magnitude of the standard thermal flux, | Q therm |, (Table 10). However, all of the | D Q therm | are around 60% of the magnitude of the standard bulk kinetic energy density flux. As seen in Fig. 17, the D Q therm vectors found by the different methods are significantly skewed in direction. Finally, we turn to the pressure tensor moments, P , found by the three different methods for the PIC f( v ) and for the MMS f( v ). For the PIC f( v ) Table 3 gives the six independent components of the standard pressure tensor and its eigenvalues, {1, .5, .5}, which are represented by a 49 pressure ellipsoid, in Fig. 3. Tables 5, 6 and 7 give the multibeam pressures and their eigenvalues for the three methods. The visual method, for example, yields multibeam pressure eigenvalues, {.5, .3, .2}. The multibeam pressure eigenvalues are roughly half the standard pressure eigenvalues, as was the case for the energy densities which are, in fact, equal to the sum of the three eigenvalues, Tr( P ). In view of this equality, the comparative PIC energy densities in the bar graph of Fig. 6 give a pretty good measure of the smaller multibeam pressures in comparison to the standard pressure. For the MMS f( v ), Table 10 gives the six independent components of the standard pressure tensor from the FPI moments file and its eigenvalues, {1.91, 1.66, 1.45}. Tables 13, 14 and 15 show the multibeam pressure tensors and their eigenvalues for the three methods. The visual method yields multibeam pressure eigenvalues, {1.65, 1.44, 1.33}. Thus, for the MMS f(v), the multibeam pressure eigenvalues are roughly 15% smaller than the standard pressure eigenvalues. The standard and multibeam pressure ellipsoids are shown in velocity-space for the PIC f( v ) and MMS f( v ) cases in Figs. 4, 5, 11 and 13. For the PIC f( v ) the multibeam ellipsoids differ in both size and orientation from the standard ellipsoid, whereas for the MMS f( v ) the multibeam ellipsoids are fairly tightly embedded in the standard ellipsoid with the same orientation. 7.3 PRACTICAL SIGNIFICANCE OF MULTIBEAM MOMENT ANALYSIS How might multibeam moments of complex velocity distributions impact studies of magnetic reconnection in the magnetosphere and elsewhere? First, we point out that the potential impact is not limited to complex ion velocity distributions but extends to electron distributions as well. For example, during reconnection, electron crescents have been measured together with electron core distributions which may consist of both cold and hot components. We have presented multibeam electron moment studies of complex electron f( v ) containing crescents at Workshops and Meetings (e.g., AGU, 2019). Finding reliable moments of complex electron and ion distributions is germane to the identification of particle diffusion regions during reconnection as well as to particle heating and acceleration. Off-diagonal pressure tensor elements, for example may break the frozen-in condition of electrons or ions, since the pressure force enters into the particle momentum equation together with electric and magnetic forces (Aunai, et al 2011; Dai, et al, 2015). Our results show that pseudothermal effects are less pronounced for the MMS ion f( v ) than for the PIC f( v ) due to a high energy ion velocity-space cloud present in the former but missing from the latter. For example, the MMS f( v ) has a higher magnitude ion heat flux, | Q heatflux | relative to that of the ion enthalpy flux, | Q enthalpy |, when compared to the size of | Q heatflux | in the PIC f( v ). This might appear to suggest that the PIC simulation should have included higher-energy ion velocities in order to properly evaluate higher order moments such as the ion heat flux. A possible alternate explanation is the following: The PIC simulation is based on magnetotail reconnection measurements while the MMS f( v ) is from a dayside magnetopause reconnection event. Measurements by Eastwood, et al (2013) of energy flux moments during antiparallel 50 symmetric reconnection in the magnetotail suggest that the energy flux is dominated by ion enthalpy, with contributions from the electron enthalpy and heat flux and ion kinetic energy flux. The electron kinetic energy and ion heat flux appear to be negligible. Even if they are present during magnetotail reconnection high-energy ion clouds may be less relevant to heat flux moments if the clouds are symmetric and centered at the standard flow velocity. The PIC simulation studied here may not be a bad approximation. At any rate, the magnetotail should be an interesting future arena for carrying out multibeam analyses of complex ion distributions. The different methods (visual, k-means and least-squares-fit) used in this paper to identify the ion beams whose sum is an approximation to a measured complex ion f( v) have enabled multibeam moments of a measured f( v ) to be taken practically and perhaps eventually automated. By comparing the multibeam thermal energy density moments, (U thermal ) MB to the standard moment, U thermal , over a neighborhood of nearby space-time locations, one may ultimately be able to distinguish between whether a species has undergone heating or has undergone particle acceleration . Hopefully, the multibeam methods presented here will also be instructive for future applications to complex particle velocity distributions measured in the magnetosphere and to emerging measurements of complex distributions in the solar wind (e.g. from the Solar Orbiter and the Parker Solar Probe missions). Acknowledgments
We wish to thank Dr. Steven Schwartz for helpful conversations leading to this paper. JPE is supported by UKRI (STFC) grant ST/ S000364/1. M. Goldman, D. Newman and G. Lapenta are supported by NASA Grants NNX08AO84G, 80NSSC17K0607 and . We thank NASA for the computational resources used to generate the graphics. The MMS data is accessible through the public link provided by the MMS science working group teams: http://lasp.colorado.edu/mms/sdc/public/. We are grateful for the dedicated efforts of the entire MMS mission team, including development, science operations, and the Science Data Center at the Univ. of Colorado.
References
Aunai, N., A. Retin`o, G. Belmont, R. Smets, B. Lavraud, and A. Vaivads (2011), The proton pressure tensor as a new proxy of the proton decoupling region in collisionless magnetic reconnection,
Annales Geophysicae , 29, 1571–1579, doi: 10.5194/angeo-29-1571-2011 Burch, J.L., et al, (2016) Electron-scale measurements of magnetic reconnection in space,
Science , DOI: 10.1126/science.aaf 2939. Dai, L., C. Wang, V. Angelopoulos, and K. H. Glassmeier (2015), In situ evidence of breaking the ion frozen-in condition via the non-gyrotropic pressure effect in magnetic reconnection,
Annales Geophysicae , 33, 1147–1153, 2015, doi:10.5194/angeo-33-1147-2015 Dupuis, R., Goldman, M. V., Newman, D. L., Amaya, J., Lapenta, G. (2020) Characterizing magnetic reconnection regions using Gaussian mixture models on particle velocity distributions,
Astrophysical Journal,
Physical Review Letters, , https://doi.org/10.1103/PhysRevLett.110.225001 Eastwood, J.P., M. V. Goldman H. Hietala D. L. Newman R. Mistry G. Lapenta, (2014), Ion reflection and acceleration near magnetotail dipolarization fronts associated with magnetic reconnection, Journal of Geophysical Research , https://doi.org/10.1002/2014JA020516. Fazakerley, A. N., S. J. Schwartz, G. Paschmann (1998). Measurement of plasma velocity distributions, in analysis methods for multi-spacecraft data (eds. G. Paschmann and P. W. Daly),
International Space Science Institute/ESA , ISSI SR-001, p 91-124. Goldman, M.V., D.L. Newman J.P. Eastwood and G. Lapenta, (2020), Multi-beam Energy Moments of Multi-beam Particle Velocity Distributions,
Journal of Geophysical Research, https://doi.org/10.1029/2020JA028340 Goldman, M.V., D.L. Newman and G. Lapenta, (2016), What can we learn about magnetotail reconnection from 2d PIC Harris-sheet simulations?
Space Science Review,
Vol. 199 , Space Science Reviews , Vol. 199 (1-4), pp.331-406. Rager, A. C., Dorelli, J. C., Gershman, D. J., Uritsky, V., Avanov, L. A., Torbert, R. B., … Saito, Y. (2018). Electron crescent distributions as a manifestation of diamagnetic drift in an electron‐scale current sheet: Magnetospheric Multiscale observations using new 7.5 ms Fast Plasma Investigation moments.