Evolution of Turbulence in the Kelvin-Helmholtz Instability in the Terrestrial Magnetopause

: The dynamics occurring at the terrestrial magnetopause are investigated by using Geotail and THEMIS spacecraft data of magnetopause crossings during ongoing Kelvin–Helmholtz instability. Properties of plasma turbulence and intermittency are presented, with the aim of understanding the evolution of the turbulence as a result of the development of Kelvin–Helmholtz instability. The data have been tested against standard diagnostics for intermittent turbulence, such as the autocorrelation function, the spectral analysis and the scale-dependent statistics of the magnetic ﬁeld increments. A quasi-periodic modulation of different scaling exponents may exist along the direction of propagation of the Kelvin–Helmholtz waves along the Geocentric Solar Magnetosphere coordinate system (GSM), and it is visible as a quasi-periodic modulation of the scaling exponents we have studied. The wave period associated with such oscillation was estimated to be approximately 6.4 Earth Radii ( R E ). Furthermore, the amplitude of such modulation seems to decrease as the measurements are taken further away from the Earth along the magnetopause, in particular after X ( GSM ) (cid:46) − 15 R E . The observed modulation seems to persist for most of the parameters considered in this analysis. This suggests that a kind of signature related to the development of the Kelvin–Helmholtz instabilities could be present in the statistical properties of the magnetic turbulence.


Introduction
Turbulence at the interface between the solar wind and the magnetosphere is an important subject in space physics. Properties of turbulence can help with understanding the transport of mass, momentum, and energy from the solar wind to the magnetosphere [1,2] and, more generally, the mechanisms of interaction between these two regions of the near-Earth space.
The Kelvin-Helmholtz Instability (KHI) can drive waves at the magnetopause. These waves can grow to form rolled-up vortices and facilitate transfer of plasma into the magnetosphere. This mechanism is considered one of the most important in the low-latitude boundary layer (LLBL) during

The Data: Geotail and THEMIS Magnetopause Crossing
The description of the properties of turbulence is customarily performed through the high-order statistical analysis of the field increments. This requires that each interval should have a statistically significant number of measurements. As a rule of thumb, the correct estimate of the qth order moment requires 10 q measurements in order to ensure convergence. Based on this, in order to enable the computation of the fourth-order moments, we have selected a set of intervals each including a number of data points N 10 4 . Furthermore, we used intervals that are not affected by an excessive presence of data gaps, so that missing points are less than 10%. According to these requirements, we have selected 17 Geotail events from the list of 19 reported by Hasegawa et al. [15]. Only two THEMIS events from the set of 14 reported by Lin et al. [30] have been chosen as the only ones with an interesting position for our analysis, i.e., behind the dawn-dusk terminator on the left-side where we have fewer Geotail events. The observed condition of all events are listed in Table 1 and in Table A1 in Appendix A. Table 1 shows the following information (from Hasegawa et al. [15]): date; time interval; GSM position measured in units of the Earth radius R E ; IMF condition; ion mixing status; fluctuation period, related to the rolled-up vortices; magnetosheath mean bulk velocity v ms ; wavelength λ, obtained multiplying the magnetosheath flow speed by the fluctuation period, and dividing by 1 R E . The ion mixing status is defined "mixed" if a significant amount of cool magnetosheath-like ions was present on the magnetospheric side of the magnetopause, where the density, evaluated using the full phase space distribution, is n > 1/cm 3 ; or "weakly-mixed", when magnetosheath-like ions were found on magnetospheric side, with density lower than n < 1/cm 3 . The fluctuation period corresponds to the perturbations in the flow that are interpreted by Hasegawa et al. [15] as being due to vortical motions of plasma (e.g., Fujimoto et al. [14]), whereas those in the magnetic field are due to deformation of the field lines when those near the magnetopause are brought into rolled-up vortices [23,31]. The event recorded on 25 March 2002 (event F) is split into two different sets F a and F b because it presents a large gap in its central part. Table 1. Event List of rolled-up vortices detected by Geotail over nine years from 1995 to 2003, adapted from Hasegawa et al. [15]. The last two events were detected by THEMIS probe C, adapted by Lin et al. [30]. The locations of all the selected events are presented in Figure 1. As can be seen, most of the events are behind the dawn-dusk terminator, and, more precisely: nine on the dawn-side (left-side in Figure 1) and nine on the dusk-side (right-side in Figure 1). Events E and E2 present the same X and Y (GSM) coordinates, so they are over-plotted on top of each other. From now on, we refer to the GSM coordinates with X, Y and Z in capital letters and for the magnetic field component in the same direction with x, y, z in lower case letters. Before proceeding to the analysis of the turbulence properties of the data, it is necessary to test the validity of the Taylor hypothesis, which allows the switch between time and space measurements [32]. In particular, the time series of a field can be assumed to be an instantaneous spatial scan of the field if the typical velocities associated with the dynamics are slower than the probe speed inside the medium. To be more precise, for a turbulent distribution of modes in wavevector space, the plasma-frame frequency term ω and spatial advection term k · v ms , (where v ms is the mean speed of each sample and the subscript "ms" indicates "magnetosheath" region, where our dataset is located), both contribute to the spacecraft-frame frequency [33], according to the relation [32]: (1) Figure 1. Locations of rolled-up Kelvin-Helmholtz events identified from Geotail (light-blue/blue dot) and THEMIS (red diamond).
As stated above, assuming ω ∼ kv A , and for a super-Alfvénic flow, v ms v A , then Equation (1) would imply ω sc k · v ms , thereby relating the spacecraft-frame frequency directly to the wavenumber of spatial fluctuations. This would correspond to the Taylor hypothesis [32,34].
In order to test the validity of the Taylor hypothesis in the samples studied here, values of the magnetosheath mean bulk speed and of the Alfvén speed were estimated for each interval under study. Unfortunately, their values are of the same order for almost all of the samples, so that the Taylor hypothesis does not hold, as shown in Table 2. However, it is possible to invoke a phenomenological approximation as an alternative to the Taylor hypothesis, previously presented by Stawarz et al. [35], and utilized in other magnetospheric studies of turbulence [36,37]. These authors assume that fluctuations are mainly Alfvénic (at least in the large-scale domain), so that their frequency in the plasma frame can be estimated as ω ∼ k · v A = kv A cos θ, where θ is the angle between k and B 0 . Indicating by φ the angle between k and v ms , the advection term is k · v ms = kv ms cos φ. Since in our samples v ms ∼ v A , the advection term dominates over the frequency term (thus allowing space-time conversion in the time series) only when cos θ cos φ, it is well known that, in magnetohydrodynamic (MHD) turbulence, the energy cascade tends to develop in the directions perpendicular to the mean magnetic field B 0 , so that perpendicular wave-vectors k ⊥ dominate over parallel wave-vectors k || . Thus, we can expect that θ is close to π/2. Examining our data, we verified that the mean magnetic field B 0 is mainly oriented in the z-direction (in GSM coordinates), therefore k is essentially in the xy plane. Instead, the plasma bulk velocity v ms is in the x-direction, i.e., mainly perpendicular to B 0 . This is confirmed by calculating the angle α between B 0 and v ms as which is consistently close to 90 • for all samples ( Table 2). This implies that the condition cos θ cos φ is satisfied in our database, so that the argument given by Stawarz et al. [35] holds, and the time series can be interpreted in terms of spatial measurements. For the present analysis, we used the magnetic field time series sampled at dt = 0.0625 s obtained by the ASCII listings of Geotail MGF (Magnetic Field Measurement) high resolution magnetic field data (1/16 s sampling), where we got the data of magnetic field. Instead, the ASCII listings of Geotail LEP ion moment data (12 s sampling) was adopted to get the data of density and plasma velocity. The official website http://themis.ssl.berkeley.edu provides THEMIS data. We have worked with magnetic field data by an FGM (flux gate magnetometer) instrument at high resolution (1/128 s sampling), density and plasma velocity data taken by MOM (on-board moments) instrument at low resolution (3 s sampling). The frame of reference is the Geocentric Solar Magnetospheric (GSM) that has its X-axis towards the Sun and its Z-axis is the projection of the Earth's magnetic dipole axis (positive North) on to the plane perpendicular to the X-axis. The direction of the geomagnetic field near the nose of the magnetosphere is well ordered by this system. Thus, it is considered the best system to study the effects of interplanetary magnetic field components (e.g., B z ) on magnetospheric and ionospheric phenomena (see Hapgood [38], Russell [39]).

Analysis
In order to characterize the properties of the fluctuations, any event was analyzed using the standard diagnostics for intermittent turbulence. For each event, we have obtained: the autocorrelation function, which gives useful information about the correlation scale of the field; the associated energy power spectrum, whose power-law scaling exponent has to be compared with Kolmogorov-like spectrum observed at MHD scales, while a steeper power law is suggested below proton scales; the Probability Distribution Functions (PDFs) of the scale-dependent increments, whose deviation from Gaussian will qualitatively illustrate the presence of intermittency, and finally the kurtosis with its scaling exponent. Figure 2 shows examples of the autocorrelation function versus time scale, for two different samples of the whole collection, one chosen among the mixed status (E2) and the other one among the weakly-mixed status (O). The behaviour of the autocorrelation function is typical of turbulent fields, i.e., with roughly parabolic shape near the origin, indicating the field smoothness in the "dissipative" range, followed by a slower decay to zero, indicative of the inertial range of turbulence [40]. For the example of Figure 2, the values of the time-scale at which the autocorrelation function approaches zero, are the following. Event E2:

Magnetic Field Spectral Properties
The values obtained for the three magnetic field components for all the intervals under study vary between 13 s and 58 s, in agreement with the typical values in this region [35,41].  Figure 3 where we show two characteristic frequencies: the frequency related to the correlation time f corr = 1/τ corr and the frequency f d i , related to the ion inertial length of a thermal ion d i = c/ω pi , i.e., the ratio between light speed and ion plasma frequency. The latter has been estimated assuming the validity of the argument given by Stawarz et al. [35].
Specifically, using the Taylor approximation of relation ω sc k · v ms and replacing the angular frequency ω with the spatial frequency f → ω = 2π f , we obtain 2π f kv ms , and being the wave vector k = 1/ the inverse of a typical length (in this work, we use . At large scales, the correlation frequency very well represents the large-scale boundary of the spectral inertial range. Similarly, the inertial range clearly breaks around the frequency associated with the ion inertial scale f d i , where kinetic plasma effects start being non-negligible, and in agreement with the usual observations of solar-wind and magnetosheath turbulence [42].
In the MHD range of scales, i.e., above the ion inertial length, the spectrum is well represented by a power law. We report even for the spectra, the measurements of the two selected events, E2 chosen among the mixed status and O for the weakly-mixed status. For the example, in the left panel of  The observed spectra seem to indicate that the typical behaviour of space plasma is retrieved in these samples, and that the turbulence might be developing as a consequence of, or superimposed to, the Kelvin-Helmholtz instability. Figure 4 shows how the spectral exponent α kol and the exponent α ion are distributed around the typical expected values. The histograms, collecting the exponents for all magnetic components and for all events, show clearly that the inertial range exponent distribution is sharply peaked around the expected Kolmogorov value α kol −5/3. On the contrary, and in agreement with solar wind observations [43], the small-scale exponents are more broadly distributed around their mean α ion 2.44.

Intermittency
The statistical properties of turbulent fields cannot be fully described by spectra, and the intermittency in particular requires a more accurate statistical treatment of the fluctuations. For this reason, we have estimated the PDF for all components of the magnetic field increments at different scales, for all samples here considered. For each scale, the magnetic field increments were standardized by normalizing to their standard deviation to allow the comparison between the different scales. PDFs were thus computed as histograms normalized to the total number of field increments and to the equally spaced bin size, so that the integral of the whole PDF is one. Figure 5 shows two examples of the increment PDFs at five different scales, for the same two intervals and components already shown in previous figures. The black dashed line represents a reference Gaussian distribution. The probability distribution functions are characterized by the increasing deviation from Gaussian towards smaller scales [44][45][46], typical of intermittent turbulence. At small scales, the distributions show heavy tails indicating the presence of particularly intense magnetic field fluctuations, usually related to the presence of intermittent structures. Similar results were found for all intervals and for all magnetic field components. A way to quantify the deviation from a Gaussian distribution is the normalized fourth-order moment of the increments, the kurtosis K. In all the intervals under study, the kurtosis of all magnetic field components has the typical Gaussian value K = 3 at large scales, roughly down to the correlation scale, and then increases toward small scales as a power law K(l) ∼ l κ . For all cases, the power-law fitting range was approximately located between the correlation time (at large time scales) and the time associated with the ion-inertial scale, via the Taylor hypothesis (at small time scale). This range is mostly consistent with the observed spectral inertial range, although in some occasions a small shift towards small scales is present. The scaling exponent κ gives a quantitative estimate of the intermittency, i.e., of the anomalous scaling of the magnetic fluctuations [47].
In Navier-Stokes turbulence, it is often observed that κ 0.1 [48]. In Figure 6, the scaling exponent of kurtosis is larger, consistent with a more efficient intermittency, for both datasets E2 and O.
Note that at small scales a saturation and decrease of the kurtosis is observed in about half of the cases, as shown in Figure 6 for the E2 event on the left side. For the other cases, the kurtosis keeps its increasing trend or shows a saturation ( Figure 6 event O on the right side). The first behaviour, i.e., a saturation and decrease at small scales, is typical of the solar-wind magnetic fluctuations, while it was not usually observed in the magnetosheath (e.g., [35,49]). Although we can exclude instrumental noise effects (the observed scales are sufficiently far from the noise level, as also seen from the spectra), the reason for this behaviour is not fully understood, and it is outside of the scope of the present paper.

Evolution of the Turbulence along the Flank Magnetopause
Following the analysis of the properties of turbulence in each individual event, we now focus on the possibility to identify signatures of its evolution along the flank magnetopause by comparing the different samples under study. The position along the magnetopause likely corresponds to different stages of the KH evolution, whose possible influence on the turbulence characteristics will be discussed. Turbulence in the dayside magnetopause region is usually strong, and can have several sources besides the KHI [50,51]. These include for example magnetic reconnection poleward of the cusp [16,17,[52][53][54][55][56], or kinetic Alfvén waves mode converted in the non-uniform magnetopause region from compressional fluctuations present in the dayside magnetosheath [36,37]. The turbulence observed in our sample could thus be either generated by the identified KHI, or pre-existent, due to other drivers, and modulated by the KHI.
We are interested in understanding if the statistical properties of the field fluctuations, which may be used as an indicator of the state of turbulence, evolve following the development of the KHI along the magnetopause. The majority of the intervals selected for this work were identified as collected in the mixing region of the KHI. This makes our database roughly homogeneous in terms of the macroscopic KHI geometry and of the spacecraft position within the instability structures in the shear direction (approximately corresponding to X in our reference frame). The upstream wind conditions were relatively similar for most of the intervals, as specified in Table A1 in Appendix A. Furthermore, all samples present comparable roll-up period of 2-4 minutes. Using the tailward flow speed listed in Table 1 to transform the period in wavelength, the resulting average KH wavelength is ∼ 6.4 R E , consistent with that reported in the literature [6,12,13,15,23]. In this approximation, the different intervals are thus samples of turbulence at different distances from the instability region, along KHI structures of comparable scale, allowing to study the modulation of turbulence caused by the instability.
Therefore, we will show the variation of the turbulence and intermittency parameters estimated in previous sections as a function of the spacecraft position in the magnetosheath. In Figure 7,  (in the sub-ion range, i.e., for scales smaller than the typical ion scales) (bottom left), and the kurtosis scaling exponent κ are plotted as a function of the -X coordinate. The figure suggests a possible fluctuating behaviour of the parameters as the spacecraft position spans the X-coordinate, visible as a quasi-periodic modulation of the spectral exponents. Although such behaviour does not fit a periodic function, the approximate wavelength associated with such oscillation can be roughly estimated as The amplitude of such modulation seems to decrease as the measurements are taken further away from the dusk-dawn line, and for X −15 R E the parameters seem to become more stabilized. Note that the fluctuations of the power-law index are typically larger than the estimated fitting error on each value (as clearly visible from the figures).  [43,[57][58][59]. A modulation is suggested, during the departure along -X coordinate, consistently because the error that affect measures are significantly smaller than the α-index values. In the right-bottom panel, the fitted power-law index of the kurtosis is plotted as a function of -X coordinate. The reference value is κ = 0.101, which is typical value observed in Navier-Stokes turbulence [48].
This result suggests that, for the observed KHI events, the vortex roll-up periodicity provides a modulation of the associated turbulence at all times, possibly due to the quasi-periodic alternation of dominating plasma (i.e., on both sides of the KHI interface) along the flanks. Eventually, turbulence converges towards a fully developed state, i.e., with a stabilized spectral index, after about ∼ 15 R E , where the two regions are finally mixed. Note that this type of evolution was also recently observed in numerical simulations of turbulence in the KHI [60]. This interpretation is corroborated by the fact that the observed modulation is persistent for most of the parameters obtained in this analysis. Indeed, the scaling exponent of kurtosis κ (z) , also shown in the bottom panel of Figure 7, also presents the same modulation, with similar periodicity. For a more comprehensive overview of our observation, the fitted power-law indexes α kol , α ion and the kurtosis scaling index κ are plotted together as a function of the -X, Y, Z coordinate (Figure 8), in order to better visualize the overall behavior which characterizes the KH vortices.  α, κ α, κ Figure 8. All fitted power-law index, i.e., α kolm at MHD scales (blue symbols), α ion at ion scales (red symbols) and κ scaling exponent of the kurtosis (green symbols), as a function of -X, Y, Z coordinate. Different shades of color and symbol shape refer to the different component (see legend). The mean value of each sample is reported as grey or black square. The overall fluctuating behaviour is seen for all three indexes.
The spectral index α kol at MHD scales is indicated by different shades of blue symbols for his three field components; similarly, α ion is described by different red symbols for the three component, and the index κ of the kurtosis is indicated by green symbols. We have also plotted the mean of three component for each sample as grey squares. Moreover, reference lines at 5/3, −2.44 and 0.101 are plotted for each exponent, respectively. The fluctuating behaviour is seen for all three indexes. This confirms the possible signature of transition to turbulence in the region downstream of the KH instability, with spectral and intermittency properties evolving while the KHI vortices roll-up. This interpretation requires the assumption that all KHI events studied here might be originated at the same distance from the dusk-dawn line, e.g., at magnetosphere nose, and that the development of the KHI under similar conditions is characterized by a similar evolution of turbulence.
It is noteworthy that the emergence of a clear anti-correlation between the κ and the two spectral indexes, although it is more evident between κ and α kol . In order to look for correlations between the variations in the turbulence parameters, we directly compare the fluctuations around their means of the inertial range spectral exponent and of the kurtosis scaling exponent, as shown in Figure 9 as a function of -X for the z-component of the field. There is a hint of some correlation in the trend of the two indexes.
To quantitatively measure the amount of such correlations between all parameters, we computed the (linear) Pearson and (rank) Spearman correlation coefficients (respectively, C P and C S ) for all pairs of turbulent indexes. As also visualized in the scatter plots in Figure 10

Conclusions
We have studied and characterized the properties of plasma turbulence and intermittency, along the tail-flank magnetopause and its boundary layer, when Kelvin-Helmholtz instability was reported. We have surveyed the Geotail and THEMIS data, recognized as rolled-up vortices by Fairfield et al. [12], Hasegawa et al. [15], Fairfield et al. [21], Fujimoto et al. [28], Stenuit et al. [29], Lin et al. [30], taken during satellite magnetopause crossings. Firstly, we have applied time-series analysis techniques to the collection of 20 samples, in order to obtain the autocorrelation function, the power spectrum, the probability distribution functions of the field increments and their kurtosis. The behaviour of the autocorrelation functions is standard, with values of the correlation scales τ corr that vary between 13 and 58 s, in agreement with typical values observed in this region. In the MHD range of scales, the spectrum is well represented by a power law with exponent ∼ −1.69, not far from the Kolmogorov value −5/3. Below the typical proton scales, the spectrum is instead compatible with a steeper power law with exponent which we find in the range between −1.89 and −2.76, with a mean value α ion = −2.44. The inertial range clearly breaks around the frequency associated with the ion inertial scale f d i , where kinetic plasma effects start being non-negligible, and in agreement with the usual observation of solar-wind and magnetosheath turbulence [42]. Probability distribution functions are characterized by high tails and the deviation from Gaussian increases towards smaller scales [44][45][46]. The fat tails are due to particularly intense magnetic field fluctuations, usually related to the presence of structures. Finally, we have analysed the behaviour of the kurtosis. The range of scales where we showed the presence of a power-law is generally consistent with the spectral inertial range, and the scaling exponent κ gives a quantitative estimate of the intermittency [47]. In light of the results obtained, we have investigated the behaviour of several parameters as a function of the progressive departure along the Geocentric Solar Magnetosphere coordinates, which roughly represent the direction in which we expect the KHI vortices to evolve towards fully developed turbulence. It appears that a fluctuating behaviour of the parameters exist, visible as a decreasing, quasi-periodic modulation with an associated periodicity, estimated to correspond to approximately 6.4 R E . Such observed wavelength is consistent with the estimated vortices roll-up wavelength reported in the literature for these events [6,12,13,15,23]. The observed modulation seems robust, as it exists for most of the parameters considered in this analysis, which also present moderate correlations among each other. If the turbulence is pre-existent, it is possible that the KHI modulates its properties along the magnetosheath, as we observed.
On the other hand, if we assume that the KHI has been initiated near the magnetospheric nose and develops along the flanks, then the different intervals we study may be sampling the plasma at different stages of evolution of the KH-generated turbulence, after the instability has injected energy in a cascading process as large-scale structures.
Of course, a possible role of the KH boundaries (i.e., of the actual mixing conditions) cannot be excluded, although in that case more random fluctuations of the parameters would be expected. Instead, the presence of a modulation suggests that the KHI events are initiated at the nose, and their structure is relatively similar for all the events studied here. The initial injection scale is persistent as the scale of the modulation of the transition to turbulence. The evolving nature of turbulence is represented by the initially broad, then decreasing fluctuations of the spectral and intermittency indexes. Such regime could be associated with the transient competition between linear dynamics and the emergence of secondary instabilities, which later evolves to fully developed turbulence. In our sample, the observed parameters roughly converge to the typical values of fully developed turbulence at a distance ∼ −15 R E from the dusk-dawn line, corresponding to a factor of 2-3 in terms of the KHI vortex roll-up. This observation might be an indication of the typical distance for the observation of fully developed turbulence initiated by a KH instability. If this behaviour is general, it could be of relevance for KHI near interplanetary shocks, in the solar corona, and in the magnetotail, providing useful information for currently operating and future space missions, such as MMS, the Parker Solar Probe and Solar Orbiter. Table A1. Event List of rolled-up vortices detected by Geotail over nine years from 1995 to 2003, adapted from Hasegawa et al. [15]. The last two events were detected by THEMIS probe C, adapted by Lin et al. [30]. We collected the daily conditions of solar wind: the magnetic field magnitude |B sw |, the mean values of three components of the field B sw x , B sw y , B sw z , the average values of bulk speed |V sw |, the mean values of three components of the velocity V sw x , V sw y , V sw z , the values of the proton density n sw and proton temperature T sw .

Event
Interval