We provide a tutorial on the paradigms and tools of magnetohydrodynamic (MHD) turbulence. The principal paradigm is that of a turbulent cascade from large scales to small, resulting in power law behavior for the frequency power spectrum for magnetic fluctuations $EB(f)$. We will describe five useful statistical tools for MHD turbulence in the time domain: the temporal autocorrelation function, the frequency power spectrum, the probability distribution function of temporal increments, the temporal structure function, and the permutation entropy. Each of these tools will be illustrated with an example taken from MHD fluctuations in the solar wind. A single dataset from the Wind satellite will be used to illustrate all five temporal statistical tools.

## I. INTRODUCTION

We present a tutorial of MHD turbulence. Our focus will be on magnetohydrodynamic (MHD) turbulence as measured in the solar wind, but the techniques described apply equally well to laboratory, and even simulated plasma turbulence. In Sec. I, we open with an introduction to turbulence, including the celebrated 1941 theory by Kolmogorov. The role of plasma turbulence in the solar wind is introduced. In Sec. II, we provide a pedagogical overview of the statistical analysis tools used in turbulence research. Analysis in the time domain is emphasized but analogs in the space domain are noted. Examples of five different statistical metrics calculated from a single solar wind dataset will be provided. Additional examples from other solar wind observations will also be discussed. Finally, we close in Sec. III with some conclusions.

^{1,2}More formally, a turbulent fluid has more energy in convective motion than is ultimately dissipated as heat. The ratio of those energies is the Reynolds number of the flow given by

*v*is a characteristic flow speed,

*L*is a characteristic length of the large scale flow, and

*ν*is the kinematic viscosity. A turbulent fluid has a large Reynolds number, which is to say convection dominates momentum diffusion, and there is a large separation of scales between the convective motion of the fluid and the scale at which the kinetic energy is dissipated. Turbulent flows are characterized by an energy cascade in which energy contained in the largest convective motion of the fluid is transferred via nonlinearities to ever smaller scales until it is dissipated. The range between the energy injection scale and the dissipation scale is known as the inertial range.

Fluctuation energy at different spatial scales is represented in Fourier space as a wavenumber power spectrum, *E*(*k*). The picture of the energy cascade begins with energy injected at the largest scales (smallest *k*) by stirring or interaction with boundaries. Initially, energy is added to the system at the largest scale *L*, or the smallest wavenumber $k=2\pi /L$. Energy accumulates at the largest scale (Figures 1(a) and 1(d)), but eventually nonlinearities transfer energy to smaller scales (larger spatial frequency *k*) without dissipation (Figures 1(b) and 1(e)), and a “cascade” develops in *k*-space (Figures 1(c) and 1(f)). It is only at the smallest scales (truly microscopic scales in the case of conventional fluid turbulence) where energy is ultimately dissipated. If the dissipation is small and the Reynolds number is large, then there is large separation of scales between injection and dissipation called the inertial range.

*k*to the next at a rate

*ϵ*. According to a hypothesis by Kolmogorov

^{3}(see below), the form of the wavenumber spectrum in the inertial range is

*η*is a resistivity (typically Spitzer

^{4}). Energy in turbulent magnetoplasmas can be dissipated by either viscosity or resistivity, though in the collisionless solar wind, the specific dissipation mechanism is complicated. Places in the fluid where there are sheared flows (i.e., vorticity given by $\u2207\xd7v$) give rise to viscous dissipation. Places in the fluid where there are sheared magnetic fields (i.e., currents given by $\u2207\xd7B$) give rise to resistive dissipation. The competition of these two effects is given by the magnetic Prandtl number

Three symmetries are important in descriptions of turbulence. Turbulence is called *stationary* if mean values, say, $\u3008b2\u3009$ are independent of the origin in time. Turbulence is called *homogeneous* if mean values are independent of position (e.g., if $\u3008b(r)2\u3009=\u3008b(r+\Delta r)2\u3009$, this is “weak homogeneity” strictly speaking). Finally, turbulence is called *isotropic* if mean values are independent of direction (e.g., if $\u3008bx2\u3009=\u3008by2\u3009$, and in fact all multipoint correlations). Turbulence never absolutely reflects these symmetries, for example, the flow direction away from the Sun is special in the solar wind. We find that for the SSX plasma wind tunnel,^{5–8} there are extended periods during which the turbulence is approximately stationary, homogeneous, and isotropic. The as yet unproven general ergodic theorem states that time averages are the same as ensemble averages, assuming the fluctuations are stationary (the ergodic theorem has been proven under certain conditions^{9}). We typically perform time averages over short epochs in the plasma wind tunnel,^{8} but also then perform averages over an ensemble of perhaps 80 realizations. For the solar wind data discussed here, we will need to invoke the ergodic hypothesis to compute averages, or else simply do calculations on the entire time series.

### A. Kolmogorov 1941 theory

^{2}Nonetheless, it is clear that turbulence is described by the Navier-Stokes equation (by numerical simulation, if nothing else) and the key criterion for turbulence, the Reynolds number, is derived from the ratio of the convective term (second on the left) to the dissipative term (second on the right).

*L*and the Kolmogorov scale is a function of the Reynolds number: $L/\u2113K=Re3/4$.

The essence of the Kolmogorov 1941 scaling argument for the omni-directional wavenumber spectrum for fully developed turbulence is that *E*(*k*) in the inertial range depends only on *k* (via a power-law) and the energy transfer rate *ϵ*. This is another of Kolmogorov's hypotheses. Kolmogorov thought about an energy transfer rate per unit mass: $\u03f5\u223cv2/\tau $. For MHD turbulence, it is total energy, $b2+v2$, that is transferred from one scale $\u2113\u221d1/k$ to the next in a characteristic time *τ*. So the energy transfer rate is $\u03f5\u223cb2/\tau $, where *b* is the fluctuating part of the magnetic field and *τ* is the time scale over which the energy is transferred.

*E*(

*k*) are such that

*τ*in the energy transfer rate depends on the physics of the transfer. For MHD, we consider an Alfvén crossing time at the scale $\u2113$

*v*is evaluated at the field at wavenumber

_{A}*k*. So now we do dimensional analysis:

*k*changes the scaling for

*E*(

*k*) at scales smaller than

*δ*

_{i}*C*, and where

_{H}*δ*appears for dimensional reasons)

_{i}### B. Solar wind

The solar wind is often referred to as the best studied turbulence laboratory (see Refs. 10–13 for a set of excellent reviews). Indeed, there are extended periods (hours to days) in which the solar wind is highly stationary. Aside from boundaries at planetary magnetospheres and the heliopause, there are no walls constraining the solar wind. Measurements in the solar wind require spacecraft and seldom are there more than a few spacecraft present to coordinate measurements (the Cluster group uses four satellites in a tetrahedral arrangement). It is known that the solar wind is anisotropic with different statistical character parallel and perpendicular to the local mean field. Magnetic field fluctuations in the solar wind tend to have minimum variance in the direction of the mean magnetic field. Much is known about the turbulence properties of the solar wind. Only a brief overview is presented here.^{10–13}

The solar wind is a high velocity ( $\u2265400\u2009km/s$), low density ( $10\u2009cm\u22123$) hydrogen plasma with embedded dynamical magnetic field (typically $B\u224510\u2009nT=100\u2009\mu G$ at 1 AU).^{10–13} The turbulent properties of the solar wind have been studied in great detail near Earth (1 AU) but some satellites (notably Voyagers 1 and 2) have made plasma measurements out to the heliopause, about 120 AU from the Sun. The flow is supersonic and super-Alfvénic (with M about 10), and there are periods of fast wind with velocities over $600\u2009km/s$. Temperatures in the solar wind plasma are typically about $10\u2009eV$ with $Ti\u2265Te$. Plasma beta ( $\beta =Wth/WB$), again at 1 AU, is approximately unity, indicating that neither magnetic pressure nor kinetic pressure dominates the dynamics. Interestingly, the solar wind is essentially collisionless; the mean free path for inter-particle collisions is approximately 1 AU. Nonetheless, the solar wind behaves in many ways like a collisional conventional fluid, with interactions mediated via processes other than collisions.

The solar wind is the only astrophysical collisionless plasma that we can study *in situ*. It is clear that the solar wind exhibits fully developed turbulence in the sense that an active cascade is present in all dynamical MHD quantities (**B**, **v**, n) as we will show in Sec. II (at least for **B**). The fluctuations tend to be Alfvénic insofar as **B** and **v** are either aligned or anti-aligned. Solar wind turbulence contains coherent structures that reveal themselves as temporal intermittencies in the flow. Since the solar wind plasma is effectively collisionless (the collisional mean free path is on the order of one AU), the wave processes that mediate the turbulent evolution create anisotropies at the smallest scales.

As we will discuss in Sec. II, turbulent dynamics can be studied in both the time domain and space domain. Analysis of solar wind plasma has focused on the time domain since typically only single spacecraft are available (with some exceptions) and the Taylor hypothesis is invoked ( $f=VSW/\lambda $, discussed below). The connection with cascade theory is through the wavenumber spectrum, which is approximated by Taylor hypothesis applied to frequency spectrum. Coherent structures and intermittency is manifest in non-Gaussian features in the probability distribution function of temporal increments (again, with the caveat of the Taylor hypothesis). Stochasticity of a waveform is revealed by a calculation of the permutation entropy.

## II. OVERVIEW OF STATISTICAL TOOLS FOR TURBULENCE

Here, we review the types of statistical and analytical tools used in the study of fully developed turbulence. These are described in much more detail in Batchelor^{1} and Frisch^{2} and have been largely developed in the study of turbulence in air in conventional wind tunnels and water in tidal channels.^{14} These tools can be applied to scalar fields such as density *n* and temperature *T*, as well as vector fields such as velocity $v$ and magnetic field $b$. We will focus on time variations of the magnetic field here.

In Figure 2, we show a sample magnetic field dataset from the Wind satellite. The Wind spacecraft provides high-cadence magnetic field observations of the solar wind using the MFI^{15} from the L1 Lagrangian point between the Earth and the Sun. Measurements are made 11 times per second using a flux gate magnetometer and then averaged to 3 s to remove the spacecraft spin signal from the data. We focus on one component of the magnetic field *B _{x}* over an eight-day period of so-called “fast wind,” with flow speed of about 600 km/s.

*B*is the sunward direction. We find for the tutorial purposes here, results are similar for all three components. The time cadence is 3 s, so there are 230 000 data points in the full record (Figure 2(a)). At 600 km/s, the sun-earth transit time is about 3 days, so the full dataset constitutes several AU's of turbulence. Here, we use multi-day long intervals of a fast wind stream (January $14\u221221$ 2008) with large scale magnetic fluctuations on the order of 10 nT.

_{x}In Figures 2(b)–2(d), we show shorter snippets from the full record (one day, 1 h, and 4 min, respectively). We will sometimes select an epoch of interest for study during a suitable period rather than analyzing the entire record. In addition, the epoch of interest should be during a period of otherwise stationary turbulence. By stationary we mean that average values are independent of the choice of time origin. Finally, the time series should persist for several dynamical times, in this case, several Alfvén times $tA=L/vA$, where the Alfvén speed is the characteristic velocity in a magnetized plasma: $vA=B/\mu 0Mn$, where *B* is the local magnetic field, *M* is the ion mass, and *n* is the number density. The Alfvén speed for this dataset is about 0.1 of the flow speed, i.e., $vA\u223c60$ km/s.

Long records should be available for the computation of higher order statistics such as structure functions and permutation entropy. The types of tools naturally divide into time domain and space domain, but we will focus on time domain analysis since the Wind spacecraft samples only at one spatial location. The extension of these tools to the space domain is straightforward, but in the case of an experimental measurement in the lab or in space (as opposed to a simulation), multiple single-time spatial diagnostics are required. In the solar wind, an armada of four synchronized spacecraft is the current state of the art. As such, the focus in the solar wind is on single spatial point time series. Five important temporal statistical tools of turbulence for time series are described below. We will illustrate each case with solar wind data. All analysis and results come from this single waveform from the Wind satellite (Figure 2).

### A. Autocorrelation function

*R*is a tensor if we consider correlations of different components of $b$, but we will focus on single components (the diagonal elements of the tensor), and just

_{ij}*R*here.

_{xx}If we had many solar wind realizations in an ensemble, we would compute the correlation coefficient for a particular *τ* by averaging over a time interval during a stationary phase of the turbulence, then averaging this result over several realizations of an ensemble. If the turbulence is truly stationary, then the function $R(\tau )$ should be independent of the choice of the origin of t. This is a good functional check of stationarity. For a single waveform such as in Figure 2, we are forced to invoke the ergodic theorem and employ time averages in place of true ensemble averages. Here, we compute the correlation coefficients for the entire eight-day waveform of Figure 2. Correlation coefficients are computed in this way for a range of *τ*'s in order to construct $R(\tau )$. $R(\tau )$ is an even function, i.e., $R(+\tau )=R(\u2212\tau )$. We define the de-correlation time *τ _{C}* as the time at which $R(\tau )$ drops by some factor: 1/2 or $e\u22121$. A more general definition involving the function (normalized by

*R*(0)) is $\tau C=\u222bR(\tau )d\tau $. If we perform time averages, we demand that averages are taken over many de-correlation times. The stationary phase of the turbulence should persist for many de-correlation times.

In Figure 3, we show an example of a temporal autocorrelation function from the Wind solar wind dataset for *b _{x}*. Note that as we see from visual inspection of Figure 2, the autocorrelation time is about 1/2 h, and fluctuations rapidly de-correlate for times larger than one hour. We calculate a similar range of autocorrelation times: $\tau 1/2=0.32\u2009h,\tau e=0.54\u2009h,\tau C=0.59\u2009h$.

The temporal autocorrelation function in the solar wind has been measured several times before. The notion of stationarity in the solar wind (i.e., that average properties of $B(t)$ do not depend on the origin of time) has also been tested. In a classic set of papers, Matthaeus and Goldstein^{16,17} analyzed magnetometer data from Voyager, ISEE 3, and IMP satellites and found correlation times as high as 50 000 s, but can be an order of magnitude smaller depending on solar wind speed and other parameters. Our sample of solar wind magnetic fluctuations has a *τ* over an order of magnitude smaller than they report. A high degree of variability of estimates of correlation time is expected, particularly if solar rotation effects are included in a long sample.

Matthaeus and Goldstein also found that the solar wind magnetic field is statistically time stationary, at least in the “weak” sense. Weak stationarity suggests that the simple two-time $R(\tau )$ defined above (*N* = 2) should be independent of the choice of the origin of t, while strict stationarity requires that all higher order correlations ( $N\u22652$) are independent of time origin. The first proper two-point single time measurements of the spatial correlation function in the solar wind plasma were performed by Matthaeus *et al*.^{18} The spatial correlation coefficients are computed exactly analogously to the temporal correlation coefficients discussed above. They used simultaneous magnetic field data from several spacecraft, including the four Cluster spacecraft flying in tetrahedral formation. Simultaneous measurements were performed with separations as small as 150 km (using pairs of Cluster satellites) to as large as $350\u2009RE$ ( $2.2\xd7106\u2009km$). Matthaeus *et al.*^{18} find a spatial correlation length of $193RE$ or $1.2\xd7106\u2009km$. This is comparable to our measured autocorrelation time $\tau C=0.59\u2009h$ multiplied by the wind speed of 600 km/s.

### B. Frequency power spectrum $EB(f)$

*b*(

*t*) can be obtained with a Fourier transform or wavelet transform. Typically, we deal with the purely real power spectrum $EB(f)$ or $EB(\omega )$

*b*(

*t*) anywhere in the plasma. In a turbulent flow, the frequency power spectrum is most useful if spatial structures are frozen into the flow. This is the Taylor hypothesis,

^{19}meaning that if a structure of size

*δ*is convected by a spacecraft at velocity

*V*, then a frequency of order $f=V/\delta $ is registered in the power spectrum. In this way, information on spatial fluctuations is encoded in the time series (i.e., time derivatives can be converted to spatial derivatives). The hypothesis pertains as long as the magnetic field of the structure changes slowly during the time the structure is advected across the spacecraft. Another way to state it is that the fluctuation velocity

*v*in the moving plasma frame is small, $v/V\u226a1$. This is a good assumption for high flow speeds and small structures. It is an especially important assumption in the turbulent analysis of the solar wind ( $v/V\u223c0.1$ for the dataset studied here).

In Figure 4, we show the frequency power spectrum for our Wind dataset. The full eight-day record was used. Note that between $10\u22123$ and $10\u22121$ Hz, the frequency power spectrum exhibits a $\u22125/3$ power-law behavior. The standard tool for computing the frequency power spectrum is the fast Fourier transform or FFT. The FFT affords an improvement in computational speed of the discrete Fourier transform from $ON2$ to $ONlnN$ (hence “fast”). Typically, however, the FFT is taken over a long time duration, often the entire record. It is often useful to analyze the spectrum as a function of time. This can be accomplished with a windowed FFT.

Some groups have adopted the more flexible wavelet transformation.^{20} The idea of the wavelet transform is to isolate shorter portions of the waveform for analysis while providing some weight to the entire time series. The time localization and weighting are performed by selection of the “mother wavelet.” There are three typical mother wavelets (derivative of Gaussian or DoG, Paul, and Morlet), and we use the fourth order Morlet wavelet in Figure 4.^{8}

Another excellent example of a frequency power spectrum from the solar wind was discussed by Sahraoui *et al*.^{21} In this measurement, very high frequency ( $100\u2009Hz$) solar wind magnetic and electric data were analyzed from the Cluster spacecraft. Data from a three-hour epoch were studied. During this time, the solar wind speed was $640\u2009km/s$ so at $100\u2009Hz$, structures as small as $6.4\u2009km$ could be detected. The plasma density was $n\u223c3\u2009cm\u22123$ and the mean magnetic field was $B\u223c6\u2009nT$. At that density, $\delta i=c/\omega pi\u223c130\u2009km$. The plasma temperatures were $Tp\u223c50\u2009eV$ and $Te\u223c12\u2009eV$, so the proton gyro radius was $\rho i\u223c120\u2009km$ and the local proton gyro frequency is $fcp\u223c0.1\u2009Hz$.

Sahraoui *et al.*^{21} calculated the frequency power spectrum from low frequency data (below about $1\u2009Hz$) from the Flux Gate Magnetometer (FGM) and higher frequency data are from the Cluster STAFF Search-Coil (SC). The data were resolved into fluctuations parallel and perpendicular to the mean magnetic field. They found that there is more energy in the perpendicular fluctuations so the turbulence is anisotropic, but the slopes are similar indicating that during this three-hour epoch, the anisotropy seems to be independent of scale. The low frequency part of the spectrum is consistent with the Kolmogorov prediction of $k\u22125/3$, again assuming the Taylor hypothesis so $\omega =kV$. The key result of their paper was that there are two break points in the spectrum, and the breakpoints are associated with structures advected across the satellite at solar wind speed and not with characteristic frequencies in the plasma frame. In other words, because the wind velocity is so high ( $V\u226bvtp$), the frequency $f\rho p=V/\rho p$ is much higher than the proton gyro frequency $fcp\u223cvtp/\rho p$ and much more consistent with the measured break point in the spectrum.^{22}

### C. Temporal increment

^{23–26}If the fluctuations in a stationary time series are truly random, then after some delay

*τ*(beginning at any time

*t*in the time series), we expect as many upward changes in the signal as downward changes, and we expect large increments to be rare. This can be quantified by constructing a record of increments for some time lag

*τ*

^{25,26}that PDFs of increments are much broader than expected from a Gaussian or normal distribution. These “fat tails” can be quantified using the statistical metric flatness (or kurtosis) for each time lag

*τ*: $F(\tau )=\u3008\Delta b4\u3009/\u3008\Delta b2\u30092$. A Gaussian distribution of a scalar variable has

*F*= 3 but turbulent PDFs of increments can have a flatness an order of magnitude higher.

An example of this technique, comparing a solar wind measurement and simulation was done by Greco *et al*.^{26} A 27-day time series of magnetic data were studied from the Magnetic Field Experiment on the ACE spacecraft. The data were subdivided into 12-h subintervals, and increments $\Delta b$ were computed for $\tau =4\u2009min$, normalized to the standard deviation for each 12-h subinterval. PDFs of increments were constructed and compared to unit variance Gaussians.

The PDFs of normalized increments of one component of magnetic field from ACE data were plotted with PDFs from both 2D and 3D MHD simulation data. The key result was that the 2D simulation more closely matches the ACE solar wind data; both have non-Gaussian “fat tails” suggestive of a preponderance of small scale coherent structures. A close analysis of the 2D simulation shows what kind of structure contributes to the non-Gaussian tails and the culprit is a sea of small-scale current-sheet-like structures that form the sharp boundaries between magnetic flux tubes.

In Figure 5, we show the increment analysis for our Wind dataset. We also see that for shorter lag times (from 3 s up to 10 min), there are non-Gaussian “fat tails” in the PDF, suggesting a preponderance of large excursions or even discontinuities in the dataset. These discontinuities are associated with current sheets in the MHD turbulence and are most pronounced at the shortest lag times (see, for example, Figure 2(c)). Note that 3–300 s at a 600 km/s solar wind speed correspond to spatial scales of 1800 to 180 000 km, or 10's to 1000's of proton gyro-orbits. This represents the range of spatial scales of the structures.

Further studies have associated ion heating with the “spontaneous cellularization” of solar wind turbulence.^{27–31} The idea is that as MHD turbulence evolves, flux tubes, discontinuities, and thin current sheets form as part of the temporal evolution and relaxation processes. The discontinuities and current sheets, revealed by the PDF of increment technique described here, can become sites of local plasma heating if magnetic reconnection ensues.

### D. Temporal structure function

*τ*. These are all raised to the power

*p*(where

*p*is not necessarily an integer), and we compute the average. Again, for a single waveform, we are forced to invoke the ergodic theorem and employ time averages in place of true ensemble averages. Here, we pick a value of

*τ*and construct the table for the entire eight-day record. The process is repeated for a series of

*τ*'s. Alternatively, the PDF of increments can be constructed for a range of time lags

*τ*and the

*p*th moment can be taken. Note that structure functions have already been discussed above in computing the flatness. Flatness can be viewed as fourth order structure function suitably normalized.

Hydrodynamic turbulence theory predicts that if the turbulence is self-similar and fully developed, then higher order structure functions should scale linearly with the order of the structure function: $SBp(\Delta s)\u223c\Delta s\zeta $, where *s* is typically a spatial displacement but connected to a time series by the Taylor hypothesis. The Kolmogorov 1941 (K41) prediction for fluid turbulence is $\zeta =p/3$,^{2,3} while the Iroshnikov-Kraichnan (IK) prediction for MHD is $\zeta =p/4$.^{32,33} The extent to which there is intermittency and coherent structures in the flow is manifest in departures from a linear relationship of the scaling exponents. Dissipation is likely to occur in these localized coherent structures whether they are viscous vortex filaments or resistive current sheets. Indeed, the dissipation need not be collisional but in the case of magnetic dissipation, almost certainly involves collisionless dissipation mechanisms at electron scales.

In Figure 6, we construct the temporal structure functions from the solar wind time series of Figure 2. We plot the structure function vs. time lag *τ* for orders $p=1\u22126$. Note that the slope of the structure function increases with increasing order p. In Figure 7, we plot the slope of the structure function as determined in Figure 6 as a function of order, but for orders up to *p* = 10. In addition, we can compute the slope of the structure function for fractional orders p, so that we can display $\zeta (p)$ as a continuous function. Note that at low order, the slope of $\zeta (p)$ tracks the K41 theory well ( $\zeta =p/3$) but rapidly departs from that model for $p\u22652$.

Temporal structure functions have previously been studied in the solar wind (see the review by Marsch and Tu^{34} and references therein, as well as^{35}). Assuming again the Taylor hypothesis that the rate of evolution in the plasma frame is slow compared to the rate at which structures are advected by spacecraft, the prediction from Kolmogorov turbulence theory is $Sp(\tau )\u221d\tau p/3$.^{2} In other words, a log-log plot of $SBp(\tau )$ for some power *p* should have a linear region corresponding to the inertial range. Famously, the second order structure function grows like $\tau 2/3$ in numerous wind tunnel experiments.^{2} Growth of the second order structure function like $\tau 2/3$ is closely connected with the frequency power spectrum (another second order statistical quantity) dropping like $f\u22125/3$ in Kolmogorov turbulence. Data from long time series of the Voyager and Helios spacecraft show structure functions with slopes much flatter than the Kolmogorov prediction. In the Marsch and Tu review, data from structure functions up to $20th$ order are shown. It is only above sixth order that departures from the Kolmogorov prediction are observed in that dataset.

In a recent series of observations using the Cluster spacecraft and the FGM and STAFF-SC instruments discussed above, Kiyani *et al.*^{35} measured high order structure functions in a stationary interval of fast solar wind. They show in a log-log plot of $SBp(\tau )$ vs *τ* an increase in the slope as the order increases from *p* = 1 to 5. An increase is observed for both the inertial range ( $\tau >10\u2009s$) as well as the dissipation range ( $\tau <1\u2009s$). Just as in the Marsch and Tu review and our results here, the Kiyani *et al.* results show structure functions with slopes flatter than the Kolmogorov prediction at the higher orders due to intermittency.

### E. Permutation entropy

Finally, we consider permutation entropy and complexity of a turbulent waveform.^{36,37} The idea is to study the ordinal pattern of a sequence of values in a time series. If we consider *N* = 5 sequential points in a waveform, we ask in what order do they appear? One possibility is that they appear in ascending order $1\u22122\u22123\u22124\u22125$. Another is that the largest value of the five appears first, followed by the next highest value, then lowest, third, second. This ordinal pattern would be represented $5\u22124\u22121\u22123\u22122$. Some examples are shown in Figure 8. There are $N!=120$ such permutations if *N* = 5 (called the embedding dimension). We are interested in the frequency each ordinal pattern appears in a long time series.

*P*consisting of all 120 frequencies of occurrence

*P*of a given length 5 ordinal pattern in all 5-value segments of the time series. Following Bandt and Pompe,

_{i}^{36}we define the Shannon permutation entropy of the time series as

*et al.*

^{37}added another metric to the study of time series related to the embedded structure in the waveform and based on the notion of the so-called “disequilibrium.” They define a Jensen-Shannon complexity

*P*is the uniform probability distribution function that admits the maximal entropy discussed above. $QJ[P]$ is a normalized quantity, and it quantifies how different P is from the uniform distribution

_{e}*P*. We see immediately that $CJS[P]$ is zero when $H[P]=0$ (e.g., the linear ramp), and $CJS[P]$ is also zero when

_{e}*P*=

*P*(the maximal entropy case).

_{e}If we plot $CJS[P]$ versus $H[P]$ (in the so-called CH-plane), we see that for entropies between these extremes, there is a range of possible complexities $CJS[P]$, reflecting the fact that there are often differing degrees of structure which can exist in systems which appear equally random. The complexity can be interpreted as a measure of the “non-triviality” of the distribution of a systems ordinal patterns, reflecting correlational structures neglected in the calculation of the entropy. For example, while deterministic chaos is highly unpredictable, reflected by moderately high entropies, there are intricate structures embedded in chaotic dynamics, reflected by near-maximal complexities on the CH-plane. It can be shown that if *C _{JS}* is plotted as function of

*H*, the position of any time series on the CH plane is constrained to fall within a crescent-shaped area, outlined in black in Figure 9.

We can compute $CJS[P]$ and $H[P]$ for any time series and plot it on a CH-plane. What do we see? For the two extreme cases above, we find first that the linear ramp has zero permutation entropy and zero complexity (lower-left corner), while the uniform probability distribution *P _{e}* has normalized entropy of unity and zero complexity (lower-right corner). We can generate waveforms from deterministic chaotic systems such as the logistic map, the skew tent map, Henon map, and the Lorenz map (all defined in Rosso

*et al.*

^{37}), and we find that for these systems the normalized entropy is about 0.5 but the complexity is maximal (also about 0.5). Numerically generated fractional Brownian motion waveforms (fBm)

^{37}have close to maximal entropy (and therefore C near zero).

Finally, we need not select our *N* = 5 points sequentially in a time series, but rather, we can opt to select every other point, or every 10th. The idea here is that the physics of interest may well be slower than digitization rate of our instrument. We call the number of points skipped the embedding delay, and it allows us to study permutation entropy at longer time scales, and therefore larger spatial scales. An interesting future study would be to establish a connection between coherent structures manifest in the “fat tails” in the PDF of increments discussed above, and non-trivial ordinal patterns in the time series manifest by finite complexity *C _{JS}*.

For the solar wind data shown in Figure 2, we find that the entropy is nearly maximal (and C near zero) so that solar wind fluctuations are highly stochastic, even more stochastic than fractional Brownian noise.^{38} In Figure 9(a), we plot the entire CH plane, and show the CH locations for the solar wind data of Figure 2 with a range of embedding delay times from 3 s to 3 min. We also show the CH locations of chaotic skew tent, Henon, and logistic maps, as well as fractional Brownian motion. Figure 9(b) shows a zoom in of the same data from the high entropy corner of the CH plane. We see that as we increase the embedding delay, the solar wind dataset has even higher permutation entropy (*H* = 0.97), and less complexity (*C* = 0.05). It is remarkable that solar wind fluctuations have higher entropy and are therefore more stochastic than numerically generated fractional Brownian noise.

## III. CONCLUSIONS

We have presented a tutorial on the paradigms and tools of magnetohydrodynamic (MHD) turbulence. The principal paradigm is that of a turbulent cascade from large scales to small, resulting in power law behavior for the frequency power spectrum for magnetic fluctuations. A single dataset from the Wind satellite was used to illustrate five temporal statistical tools. The five statistical tools for MHD turbulence in the time domain include: the temporal autocorrelation function, the frequency power spectrum, the PDF of temporal increments, the temporal structure function, and the permutation entropy.

We have several findings that corroborate prior measurements in the solar wind. These include an autocorrelation time of about 1/2 h, and a power-law scaling of the frequency power spectrum with a spectral index of 5/3. In addition, we corroborate that for our Wind satellite data set, the PDF of temporal increments has fat tails at short lag times, indicating a population of discontinuities in the magnetic field time series. We note that for structure functions at low order, the slope of $\zeta (p)$ tracks the Kolmogorov theory well ( $\zeta =p/3$) but rapidly departs from that model for $p\u22652$. A new finding is that the solar wind time series has particularly large normalized permutation entropy, suggesting that the turbulent solar wind is nearly maximally stochastic. We suggest that applying different tools to the same dataset could illuminate new physics. For example, coherent structures manifest in the “fat tails” in the PDF of increments could be revealed by non-trivial ordinal patterns in the time series manifest by modest complexity *C _{JS}*. While entropy is a measure of missing information, complexity is a measure of correlational structure. So finite complexity

*C*in the solar wind could indicate a preponderance of correlational structure.

_{JS}## ACKNOWLEDGMENTS

This work was supported by grants from the Department of Energy (OFES), and the National Science Foundation (Physics Frontier Center for Magnetic Self Organization, CMSO). The authors gratefully acknowledge Robert Wicks (NASA) for providing the data, the assistance of technicians S. Palmer and P. Jacobs, and experimentalist students Adrian Wan'15, Peter Weck'15, Emily Hudson'17, Jeffrey Owusu-Boateng'16, Alexandra Werth'14, Darren Weinhold'12, Xingyu Zhang'12, Ken Flanagan'12, Bevan Gerber-Siff'10, Anna Phillips'10, Vernon Chaplin'07, and postdocs Tim Gray'01, Chris Cothran who all contributed to turbulence studies at SSX. The authors would also like to acknowledge the expert suggestions of our referee, who helped clarify several aspects of this tutorial.

## REFERENCES

*Theory of Homogeneous Turbulence*

*Turbulence: The Legacy of A. N. Kolmogorov*

*Random Functions and Turbulence*