Properties of Hall-MHD Turbulence at Sub-Ion Scales: Spectral Transfer Analysis

: We present results of a multiscale study of Hall-magnetohydrodynamic (MHD) turbulence, carried out on a dataset of compressible nonlinear 2D Hall-MHD numerical simulations of decaying Alfvénic turbulence. For the ﬁrst time, we identify two distinct regimes of fully developed turbulence. In the ﬁrst one, the power spectrum of the turbulent magnetic ﬂuctuations at sub-ion scales exhibits a power law with a slope of ∼− 2.9, typically observed both in solar wind and in magnetosheath turbulence. The second regime, instead, shows a slope of − 7/3, in agreement with classical theoretical models of Hall-MHD turbulence. A spectral-transfer analysis reveals that the latter regime occurs when the energy transfer rate at sub-ion scales is dominated by the Hall term, whereas in the former regime, the governing process is the dissipation (and the system exhibits large intermittency). Results of this work are relevant to the space plasma community, as they may potentially reconcile predictions from theoretical models with results from numerical simulations and spacecraft observations.


Introduction
Space and astrophysical plasmas are often found in a highly turbulent state. In particular, the solar wind is one of the few natural plasma environments (together with planetary magnetospheres) that is accessible to direct investigation, thanks to space missions and references therein, e.g., [1]. Among others, the turbulent nature of the solar wind is revealed by in situ measurements of the magnetic field. At large scales, where the solar wind behaves as a fluid and can be modeled within the framework of magnetohydrodynamics (MHD), the power spectrum of the magnetic field exhibits a power-law behavior with a spectral index close to −5/3 (the so-called inertial range), analogous to the Kolmogorov cascade in hydrodynamic turbulence [2]. When approaching the ion characteristic scales, the spectrum steepens, with a slope ranging from −2 to −4. The origin of such transition (possibly related to the onset of Hall currents and ion-kinetic dissipation, e.g., [3][4][5]) and the nature of sub-ion scales turbulence are still unknown.
Turbulence properties at sub-ion scales have been the main focus of several studies in the last decades, e.g., [6][7][8][9][10][11][12][13][14][15][16][17]. Early theoretical and numerical studies [7,8,[18][19][20] had predicted and/or measured a steepening of the power spectrum to a spectral index of −7/3, caused by either the onset of Hall physics or by the change in the wave-like nature of the nonlinear interactions (e.g., kinetic Alfvén waves) that mediate the turbulent energy cascade. Although such slope has been measured in early 3D Hall-MHD numerical simulations, e.g., [21], it was rarely observed by spacecraft [22][23][24]. Instead, a steeper slope with an average value around −2.8 or −3 is routinely reported, e.g., [25,26], as well as in later numerical studies [4,13,14,16,[27][28][29][30]. Refined models incorporating intermittency effects (e.g., assuming that energy dissipation mainly takes place in localized coherent structures such as current sheets [10,31,32] and vorticity structures that form in their vicinity [33,34]) are able to reproduce the ranges of observed slopes. However, there is no definite answer on the nature of the underlying physical mechanisms, nor on the conditions that determine the onset (if any) of different sub-ion scales turbulent regimes.
In this work, we observe the establishment of two different regimes of turbulence at sub-ion scales in 2D compressible viscous-resistive Hall-MHD numerical simulations: (i) the "classic" regime, where turbulence at sub-ion scales is highly intermittent and the magnetic field power spectrum shows a power law with a slope ∼−2.9; (ii) a pristine Hall-MHD regime, where the spectral index is −7/3. Furthermore, and for the first time, by carrying out a spectral energy transfer analysis [35], we link the appearance of either of the two regimes to specific terms in the equation for the evolution of the energy cascade, one related to the Hall currents, and the other related to dissipation.
The paper is organized as follows. In Section 2, we describe the numerical setup employed and the spectral transfer analysis method. In Section 3, we highlight the general evolution of the Hall-MHD simulations and describe their turbulent spectral properties. Finally, in Section 4, we summarize our results, also in the context of previous works.

Hall-MHD Simulations of Plasma Turbulence: Numerical Setup
The Hall-MHD model takes into account two-fluid effects generated by the ions decoupling from the electrons and from the magnetic field at scales smaller than the ion inertial length d i . Here, we use a simplified framework: the nonlinear viscous-resistive MHD equations are modified by only adding the Hall term in the induction equation. This is done by substituting the fluid velocity u with the electron velocity u e = u − J/en e , J being the current density, e the unsigned fundamental charge, and n e the electron numerical density. In their adimensionalized form, the Hall-MHD equations take the form where is the viscous stress tensor, Γ = 5/3 is the adiabatic index, and {ρ, u, B, T, P} denote the usual physical quantities. The pressure evolves according to the equation of state P = ρT. The above equations are normalized with respect to a characteristic length L = d i , a magnetic field amplitude B 0 , a plasma density ρ 0 , an Alfvén velocity c A = B 0 / 4πρ 0 = Ω i d i , a pressure P 0 = ρ 0 c 2 A , and a plasma temperature T 0 = k B c 2 A /m i . Ω i = eB 0 /m i c denotes the ion-cyclotron frequency and m i is the mass of the ions. With this normalization, the (adimensional) dynamic viscosity µ is in units of d i c A ρ 0 , and the magnetic resistivity η is in units of d i c A . Furthermore, we have a plasma β = 2P 0 /B 2 0 = 2Γ −1 (c s /c A ) 2 , where c s = ΓP 0 /ρ 0 is the adiabatic sound speed. Equations (1)-(4) are numerically solved by means of a pseudospectral 2.5D code [4,[36][37][38][39], that evolves the scalar fields {ρ, T} and the 3D vector fields {u, B} in a 2D (x, y) periodic domain, using Fourier decomposition to calculate spatial derivatives and a third-order Runge-Kutta scheme [40] for temporal integration. In Fourier space, a filter is applied according to the 2/3 Orszag rule [41], to avoid aliasing of the nonlinear quadratic terms. The aliasing of the nonlinear cubic terms is mitigated by the presence of a finite dissipation [42].
We carry out a series of three numerical simulations of decaying Alfvénic plasma turbulence. For all the simulations, we consider a 2D box of size L x × L y = 32 d i × 32 d i and a grid resolution of ∆x = ∆y = d i /32, corresponding to N x × N y = 1024 2 points. The system is initialized with a constant ambient magnetic field B 0 = B 0 e z out of the simulation plane, along the z direction (which we name as the parallel direction). The xyplane (i.e., the perpendicular plane) is filled with freely-decaying random (incompressible) Alfvénic-like sinusoidal fluctuations of same amplitude, characterized by a root-meansquare (rms) amplitude b rms = u rms , a zero mean cross helicity, and wavenumbers spanning from the smallest nonzero value contained in the box up to the injection wavenumber k ⊥,inj = 2π/ inj , where inj is the injection scale and k ⊥ = k 2 x + k 2 y (for more details on the initialization, see [14]). Thanks to such initialization, the energy reservoir located above the injection scale is enough to allow the simulations to reach and mantain a fully-developed turbulent state.
In the solar wind, the value of the ion plasma beta ranges from ∼0.01 to ∼100, with an average value of 1. If we assume that the electron beta is the same, this translates to a total plasma beta of β = 2. We set this value for the initial plasma beta in all the simulations and we vary the rms of the fluctuations as well as the injection scale and the value of viscosity µ and resistivity η (with µ = η). Table 1 reports the simulation dataset used in this work. Table 1. Physical and numerical parameters of the simulations used in this work (in code adimensional units). Box size L 0 , injection scale L inj , rms amplitude of the initial fluctuations (b rms = u rms ), nonlinear time at the injection scale τ nl (in units of Ω −1 i ), viscosity µ and resistivity η, (magnetic) Reynolds number (Rm)Re = L inj u rms ρ 0 /µ (since µ = η and ρ 0 = 1), and initial Mach number The nonlinear time is defined as τ nl = L inj /u rms . Q 0 is the total heating rate (in units of ρ 0 c 2 A Ω i ) measured in the time interval of fully developed turbulence considered (see Section 3.2).

Spectral Transfer Equations
The main property of fully developed turbulence is the onset of an energy cascade. The energy contained in the reservoir at the largest scales is transfered (i.e., "cascades") to smaller and smaller scales through nonlinear interactions of the fields between nearby spatial scales. In viscoresistive models, such cascade ends at the so-called Kolmogorov scale, where the plasma dynamics is dominated by the dissipation. The total energy transfer rate thoughout all scales results from different terms of different origin. Some terms may be, for example, dominant at large scales and subdominant at small scales, thus determining changes/transitions in the observed plasma regime. In this work, we use the spectral transfer approach [35,43,44] to measure and decompose the energy transfer rate into different contributions. Following [35], the temporal evolution ∂ t (ρu 2 + B 2 )/2 of kinetic and magnetic energy can be obtained by summing the two inner products of Equations (2) and (4) by u and B, respectively. Furthermore, by introducing w = √ ρu [45] and using Parseval's theorem to transform to Fourier space, we may write where the star denotes the complex conjugate, {...} is the real part, wide bars denote the discrete Fourier transform (e.g., w(k x , k y , t) = ∑ (x,y) w(x, y, t)e −i(k x x+k y y) /N x N y ), and each barred quantity is a function of the wave vector k ≡ k ⊥ = (k x , k y ) and time. The MHD and Hall transfer terms T MHD and T H are The magnetic+kinetic energy transfer rate from scales with wavenumber smaller or equal to k = |k| is given by In the case of fully developed turbulence (and in general, during the whole evolution [35]), as k increases, we expect to observe a negative energy transfer rate (∂ t E k < 0). Here, S MHD,k and S H,k denote the MHD and the Hall energy transfer term, while Ψ k and D k represent the pressure-dilatation and the dissipation term, respectively. The total heating rate Q 0 , which measures the dissipation of the energy density in the whole domain, correspond to the dissipation rate D k at the maximum wavenumber. In the following, we will focus on the Hall term and on the dissipation term.
We note that the energy transfer rates introduced in Equation (7) measure a global nonlinear transfer of energy from the sphere of wavectors with radius k. In the case of a turbulent cascade that is local in Fourier space, such energy transfer is the flux through the sphere's boundary and corresponds to the cascade rate. We also remark that when turbulence is well developed and a cascade is established, the spectral transfer equation is equivalent to the Von Karman-Howarth-Monin equation based on third-order structure functions [35].

General Evolution
The overall evolution in the three simulations listed in Table 1 is qualitatively the same. The initial large-scale Alfvénic-like fluctuations quickly evolve to form large vortices and other coherent structures, such as current sheets which then undergo magnetic reconnection. The system then evolves until turbulence has fully developed, corresponding to the maximum of the rms of the current density (for more details on the evolution of a similar setup, see [4]). The time of fully developed turbulence is different for each simulation (as shown in Figure 1). Run d16 takes roughly twice the time of Run d08a(b) to reach the maximum of J rms . This is because the eddy turnover time (also known as the nonlinear time) τ nl = L inj /u rms at the injection scale is roughly twice the nonlinear time of the other two runs. Indeed, once we normalize the time with respect to τ nl (see Figure 1b) and J rms with respect to its maximum, the three curves nicely overlap, with the maximum roughly located at the same time t = 1.25τ nl (marked with a vertical dotted line), followed by a plateau, indicating that a statistically quasi-stationary turbulent state (that lasts for ∼ 0.4τ nl ) has been reached and, finally, a slowly decaying phase. We note that the definition of the nonlinear time varies in the literature usually by a constant multiplicative factor of 2π or by using the size L of the full box instead of the injection scale L inj . However, such constant does not change qualitatively Figure 1b, as it only changes the normalization of the x-axis. Figure 2 reports, for all three runs, the amplitude of the magnetic field fluctuations (panels a, b, and c) right after the maximum of the turbulent activity has been reached, together with the amplitude of the current density (panels d to i). At large scales, vortices and filamentary structures in both Run d08a and d08b are roughly twice smaller than those in Run d16 (because of the different injection scale). At small scales, the morphological differences between the three runs at sub-ion scales are evident. The filamentary current structures (panels d, e, and g) are more space filling in Run d08a and b and with richer small-scales features (panels g, h, and i), while Run d16 shows more localized and thicker filaments. Finally, compared to Run d08a, the current density in Run d08b is more homogeneously distributed (panels e and f) and with a richer multiscale structure (see panels h and i). This picture was confirmed by measuring the intermittency in the three runs (not shown here). Indeed, high (small) values of intermittency were observed in Run d16 (d08a), while almost no intermittency was observed in Run d08b. The scenario just presented in Figure 2 is consistent with previous works which explored the properties of the dissipative structures that form in 2D, e.g., [17,[46][47][48], and 3D [49][50][51] numerical simulations of plasma turbulence.

Spectral Properties and Cross-Scales Energy Transfer
Once turbulence is fully developed, magnetic power spectra of all three runs show well-defined power laws, which last for the entire duration of the plateau and also in the following decaying phase (although the energy decreases). Figures 3 and 4 show the isotropized power spectra [13] (top panels) and the spectral energy transfer (bottom panels) averaged over a time interval of roughly the same length around the plateau, to improve the statistics and decrease the effects of the oscillations in the pressure-dilatation term Ψ k [35]. The spectral properties of Run d16 are qualitatively the same as those found in our previous works [4,39,48]. We observe the onset of an inertial range cascade of power-law −5/3 in the interval k ⊥ ∈ [k ⊥,inj , 2], followed by a transition to a steeper slope of −2.9 in the sub-ion range, shown in the top-left panel of Figure 3. In the inset in the same panel, we report the spectra multiplied by k 2.9 ⊥ , where we can appreciate the length of the interval of the −2.9 slope (6 k ⊥ d i 15). The −5/3 and the −2.9 power laws are in agreement with those found in our previous works, as well as in other studies [4,13,14,16,28,30]. However, compared to, e.g., Papini et al. [4], here, the transition from the inertial to the sub-ion range is smooth and more extended. Different to that observed in the MHD range in our previous works, the velocity power spectrum in Run d16 (as well as in the other two runs) almost overlaps with the magnetic spectrum, which hints that there is no residual energy present. for Run d16, d08a, and d08b, respectively. To better highlight the differences at small scales between the current structures in the three simulations, panels (g-i) show a zoom of the amplitude of the current density in a particular region.
The bottom-left panel of Figure 3 shows the spectral energy transfer of Run d16. The behavior is qualitatively similar to that observed by [3,35], but with some important differences. In the inertial range, the MHD term S MHD,k gives the main contribution where we observe the −5/3 power law, followed by a small interval in which the MHD, the Hall, and the dissipation term are comparable, with the Hall term S H,k being the largest. This interval (2 k ⊥ d i 6) is the same interval where we observe the smooth transition between the two power laws. In the sub-ion range, the dissipation term D k dominates, although S H,k still gives a non-negligible contribution. The resulting total energy transfer rate ∂ t E k is more or less constant at all wavenumbers above the injection wavenumber.
The top and bottom right panels of Figure 3 show the same plots, except for Run d08a, which differs from Run d16 only in the location of the injection scale, moved further down to 8d i (thus the Reynolds number of Run d08a is half that one of Run d16). In wavenumber, this translates to k ⊥,inj d i ∼ 0.8, i.e., at the edge of the ion characteristic scales. As a result, the size of the −5/3 MHD inertial range is further shortened. More importantly, the size of the interval where S H,k is the largest term increases, as well as the interval of the transition between the MHD and the sub-ion regimes. Finally, the sub-ion range now shows a more extended −2.9 power-law spectrum, almost up to the 2/3-cutoff wavenumber.  (6)) averaged in the interval from 70 to 90 Ω −1 i (∼0.31τ nl ) for Run d16 (left) and in the interval from 40 to 50 Ω −1 i (∼0.31τ nl ) for Run d08a (right). Each energy transfer rate term (reported in adimensional units [Q]) is normalized with respect to the total heating rate Q 0 (see Table 1). Vertical dash-dotted and dashed gray lines denote the location of the injection and of the 2/3-filter cutoff wavenumber, respectively. The dark gray area denotes the range in which the Hall term S H,k is bigger than the other spectral energy transfer terms.  Finally, we consider Run d08b, which has half the resistivity and viscosity of Run d08a and a higher rms of the initial pseudo-Alfvénic fluctuations (see Table 1). This is the most interesting run of this work. Figure 4 shows the power spectrum and the spectral energy transfer of Run d08b. The magnetic power spectrum only shows two power-law regimes. The usual −5/3 MHD range, which further extends up to k ⊥ d i ∼ 4.5, and a sub-ion range with a clear slope of −7/3 (see compensated spectra in the inset) that extends for almost a decade in the interval k ⊥ d i ∈ [4.5, 27.4] and nicely matches the interval where the Hall energy transfer term is dominant. In such range, the −7/3 slope is very stable in time, and lasts up to the end of the simulation, also in the decaying phase and for more than 4τ nl . Such slope is predicted by theoretical models of Hall-MHD turbulence, but it is rarely measured in spacecraft observations [22][23][24] and it has been only reported in early numerical studies of turbulence, e.g., [20,21].
The concurrent presence of a −7/3 slope and of a dominant Hall energy transfer term in Run d08b suggests to reexamine Run d16 and d08a. Indeed, in both runs and in the transition interval between the MHD and the sub-ion range, where S H,k is comparable to S MHD,k and D k , the power spectrum of magnetic fluctuations shows a hint of a −7/3 power law (see compensated spectra in the insets in Figure 3).

Discussion
The simulation dataset analyzed in this work clearly revealed a strict link between the spectral properties and the spectral energy transfer in Hall-MHD turbulence. Moreover, we have shown that the turbulent cascade at sub-ion scales can enter different regimes depending on which term gives the dominant contribution to the spectral energy transfer rate.
Our interpretation of the different scenarios we reproduced (see Figures 3 and 4) is the following. In the typical case, exemplified by Run d16, the energy injection scales are much larger than the ion characteristic scales; thus, a proper MHD inertial range cascade is established and vortices and coherent structures can form at scales larger than the ion characteristic ones. As a consequence, at sub-ion scales, the cross-scale energy transfer takes place in localized regions of enhanced dissipation. This is in agreement with our previous findings [48] and is coherent with the fact that in Run d16 (and also in other simulations previously analyzed [35]) the dissipation term D k gives the dominant contribution to the total energy transfer ∂ t E k . Indeed, in Run d16 we measured high levels of intermittency (not shown here) at sub-ion scales. This does not mean that the Hall term S H,k is unimportant, since the −2.9 slope would not form in MHD simulations. Indeed, such slope is present until D k 3S H,k in both Run d16 and d08a. Instead, in Run d08b the injection scale is so close to the ion characteristic scales that the formation of coherent structures is less efficient, the level of intermittency is very low, and the sub-ion scales turbulence is in a regime where the Hall term S H,k gives the main contribution to the total energy transfer. In this pristine Hall-MHD regime, we observe the −7/3 slope, as predicted and observed in previous Hall-MHD studies [18,21]. In particular, we believe that such slope was observed by [21] in their 3D Hall-MHD simulation because the box size was too small to allow for an extended MHD inertial range, which is necessary to develop the intermittency needed for a sub-ion scale regime with D k > S H,k . Other simulations (not shown in this work), performed by varying the parameters (rms of the fluctuations, resistivity, injection scale, box size), confirm this picture.
Our results are constrained by the 2D geometry, which prevents the onset of 3D intrinsic features of plasma turbulence (e.g., critically balanced turbulence [52]). Indeed, a complete description of plasma turbulence requires the use of 3D modeling, e.g., [53,54]. However, we expect that our results are, to some extent, relevant to a more realistic full 3D case. Indeed, reduced MHD models (that share similarities with the 2.5D MHD description, e.g., [55]) are able to reproduce the turbulent energy dissipation rate at fluid scales [56]. At sub-ion scales, 2.5D models qualitatively capture several properties of turbulent dissipation [46,57]. For instance, spectral turbulent properties in both 2.5D and 3D hybrid-kinetic simulations nicely match, and show remarkable agreement with solar wind observations [58,59]. Finally, our findings obtained with a 2.5D model provide the first explanation for the unfrequent appearence of the −7/3 power law in both observations and simulations. Future investigation will assess whether our findings hold in a 3D numerical setup.
In all the simulations presented in this work, we employed a value for the plasma β = 2. The effect of varying the plasma beta on the spectral index of the magnetic power spectra at sub-ion scales has been investigated by [28], who found that the spectral index increases with decreasing plasma beta (as well as the location of the spectral break), with a functional dependence on both the ion inertial length d i and the ion gyroradius ρ i . Unfortunately, the Hall-MHD model does not contain ρ i as a characteristic scale, therefore we are uncertain whether Hall-MHD is able to model the sub-ion scales turbulence at low betas. We will explore this aspect in future works.
We remark the limits of Hall-MHD in properly modeling collisionless plasmas. In particular, it has been shown that the dissipation term in the Von Kármán-Howarth equation (which is equivalent to the spectral transfer approach [35]) cannot properly model hybrid-kinetic simulations of plasma turbulence [3]. However, it is reasonable to assume that up to some extent (as in the case of the GEM challenge [60]), the dissipation term is able to mimick the nonthermal (kinetic) features or other dissipation mechanisms not included in the Hall-MHD description.
Results of this work are particularly relevant to the space plasma community and may reconcile predictions from theoretical models with both results from numerical simulations and spacecraft observations. In particular, the huge collection of solar wind data available, encompassing a much wider parameter range (in terms of, e.g., plasma beta, amplitude of the magnetic fluctuations, temperature anisotropy, cross-helicity, etc.) than the one that can be probed by simulations, allows us to test our findings and to make some predictions. Based on our interpretation, the −7/3 pristine Hall-MHD regime will be observed when the energy transfer at sub-ion scales is dominated by the Hall term, which implies that localized intermittent coherent structures of enhanced dissipation would be almost absent. This scenario could be easily tested in spacecraft observations. Indeed, the fact that the solar wind often shows high levels of intermittency is likely the reason why the −7/3 slope is rarely observed [24,61]. Finally, multispacecraft observations (that allow direct measurement of the Hall energy transfer rate, e.g., [62]) of the pristine Hall-MHD regime would be the ideal test to confirm our findings.