Structure functions and intermittency in ionospheric plasma turbulence

Low frequency electrostatic turbulence in the ionospheric E-region is studied by means of numerical and experimental methods. We use the structure functions of the electrostatic potential as a diagnostics of the fluctuations. We demonstrate the inherently intermittent nature of the low level turbulence in the collisional ionospheric plasma by using results for the space-time varying electrostatic potential from two dimensional numerical simulations. An instrumented rocket can not directly detect the one-point potential variation, and most measurements rely on records of potential differences between two probes. With reference to the space observations we demonstrate that the results obtained by potential difference measurements can differ significantly from the one-point results. It was found, in particular, that the intermittency signatures become much weaker, when the proper rocket-probe configuration is implemented. We analyze also signals from an actual ionospheric rocket experiment, and find a reasonably good agreement with the appropriate simulation results, demonstrating again that rocket data, obtained as those analyzed here, are unlikely to give an adequate representation of intermittent features of the low frequency ionospheric plasma turbulence for the given conditions.


Introduction
The study of the structure functions associated with the fluctuating velocity is an important tool to characterize turbulence of neutral incompressible flows.It is well known Correspondence to: H. L. Pécseli (hans.pecseli@fys.uio.no)(Chandrasekhar, 1957) that the second order structure function, as a function of spatial separations, can be obtained by simple dimensional arguments, apart from a numerical constant.For the longitudinal second order velocity structure function in the universal Kolmogorov-Oubokhov range of homogeneous isotropic turbulence we thus find 2 (r)≡ (u (0)−u (r)) 2  =C 2 (r ) 2/3 , (1) in terms of the energy dissipation per unit mass and a universal Kolmogorov constant C 2 which is experimentally found to be in the range 2.1-2.5.In Eq. ( 1), the notation indicates the velocity component parallel to the separation vector r.The result in Eq. (1) has found extremely solid experimental support (Hinze, 1975).One could attempt to model higher order structure functions by similar arguments, finding trivially that n ≡ |u (0)−u (r)| n =C n (r ) n/3 .Experiments demonstrate, however, that for n>3, this analytical result no longer agrees with observations, the deviations becoming more and more pronounced with increasing n.The explanation is found in the intermittent nature of turbulence, implying that energy is dissipated in concentrated "spots" or localized regions of space (Hinze, 1975;Anselmet et al., 1984).A more specific definition is given by Rollefson (1978), stating that "a variable with zero mean will be called intermittent if it has a probability distribution such that extremely small and extremely large excursions are much more likely than in a normally distributed variable".
The universal scaling law given by Eq. ( 1) is reflected also in the turbulent power spectrum of the velocity fluctuations, as expressed in the Kolmogorov-Oubokhov spectrum, which is given as 2/3 k −5/3 apart from a universal Published by Copernicus Publications on behalf of the European Geosciences Union and the American Geophysical Union.
L. Dyrud et al.: Structure functions and intermittency in plasma turbulence constant (see also Appendix A).Since power spectra are easily obtained by spectrum analyzers, many studies prefer to use this representation for studying turbulence in fluids (Hinze, 1975) as well as plasmas (Chen, 1965;Pécseli et al., 1983;Krane et al., 2000).
The first observations and discussion of intermittency effects seemingly originate from studies of fluid turbulence.The basic ideas will apply also for plasma turbulence and many studies have been carried out, numerically as well as experimentally.Magneto-hydrodynamic (MHD) turbulence in the solar wind has been reported by Tu and Marsch (1995) and by Bruno and Carbone (2005).MHD turbulence is in a sense more complicated than its counterpart in incompressible flows since in plasmas generally two vector quantities are involved, the magnetic field and the plasma flow velocity.A plasma can however also support a simpler form of wave phenomena: electrostatic waves, which can be adequately described by the space-time variation of a scalar quantity, the electrostatic potential.Such waves are often spontaneously excited in nature by plasma instabilities and have been frequently observed also in the Earth's ionosphere.Intermittency effects have been studied in the ionospheric plasma by, for instance, Tam et al. (2005), where their work refers to ∼700 km altitudes.Other relevant studies of space plasma turbulence can be found in the work by Chang and Wu (2008).In fusion plasma studies it has been found that intermittency effects are often related to anomalous turbulent transport (Boedo et al., 2003;Xu et al., 2005), an observation also supported by earlier laboratory studies (Huld et al., 1991).Intermittency effects have been recognized in several different laboratory plasma devices also by e.g.Fredriksen et al. (2003a,b) and Kervalishvili et al. (2008).The analysis is not necessarily based on structure functions as discussed in the present work.Conditional sampling methods have been used, for example.
In the present paper we analyze turbulent fluctuations in magnetized partially ionized plasmas in the ionospheric Eregion, where collisions between charged particles and neutrals dominate the effects of ion-electron collisions.The fluctuations are electrostatic and we study the turbulent electrostatic potential φ(r, t) associated with the low frequency ionospheric plasma turbulence.We analyze the space and time evolutions of the structure functions in the form n (r, t) ≡ (φ(0, 0)−φ(r, t)) n ≡ n φ(r, t) for cases to be discussed in the following (Rose et al., 1992;Krane et al., 2000;Dyrud et al., 2006), assuming locally homogeneous and time stationary conditions.

Gaussian random processes
The second order structure function is directly related to the correlation functions of the signal, which for Gaussian random processes with zero mean contain all available information (Bendat, 1958).The joint probability function for two scalar variables with zero mean, say φ 1 and φ 2 , is in this case given by where σ 2 1, 2 ≡ φ 2 1, 2 and ρ(1, 2)≡ φ 1 φ 2 /(σ 1 σ 2 ) is the correlation function for the two times t 1 and t 2 or, if spatial variations are considered, the two positions r 1 and r 2 .For stationary and homogeneous conditions σ 1 =σ 2 ≡σ .The normalized structure function 2(1−ρ(1, 2)) depends in this case only on the separation of the two sampling positions (or times), and not on their absolute values.
Introducing the difference and sum variables, ≡φ 1 −φ 2 and ≡φ 1 +φ 2 , we can readily rewrite Eq. ( 3).After integration with respect to , we obtain the probability density for .For Gaussian processes we find , which is independent of ρ(1, 2) for all n.It is here perfectly feasible to let n be a continuous variable.For Gaussian random processes, the ratio | | n /( 2 ) n/2 is thus scale invariant, being independent of the separation of the two sampling positions, here labeled 1 and 2, or corresponding sampling time separations.
The power-law exponent for | | n is consequently directly proportional to n for this case.This property can be used as a characterization of Gaussian random processes.Upon division by n/2, we can introduce the compensated exponent α which is a constant for Gaussian random processes.In the appendix we discuss some relations between power-law spectra and structure functions.

Ionospheric turbulence
Low frequency electrostatic fluctuations are frequently observed in the lower parts of the Earth's ionospheric E-region, in the equatorial as well as in the polar ionospheres.
Several candidates for instabilities giving rise to these waves have been proposed (Rogister and D'Angelo, 1970).For the present analysis, we focus on the the Farley-Buneman instability that arises in a plasma with a large ion-neutral collision frequency, ν i > ci while at the same time ν e ce , when a dc-electric field is imposed in the direction perpendicular to the Earth's magnetic field (Farley, 1963;Buneman, 1963).The instability can have importance also in other environments, meteor tails, for instance (Dyrud et al., 2002).
We present here a simplified version of the linear dispersion relation as obtained by a fluid plasma model.The real and imaginary parts of the frequency are denoted ω r and ω i , respectively.We have (Fejer and Kelley, 1980) the approximate expressions with where ϕ= ν e ν i ce ci and L n denotes the scale length of a possible large scale plasma density gradient in the direction perpendicular to B, while V d is the difference between the electron and ion drift velocities, and θ is the angle between V d and k.Since quasi-neutrality is assumed, the results only apply for wavelengths much longer than the Debye length, λ D .The term −2β R n 0 accounts for the damping effect of recombinations, with β R being the recombination coefficient and n 0 the local plasma density (Fejer et al., 1984).Equations ( 5) and (6) are valid in the limit of very small growth rates, 0<ω i ω r , and almost B-perpendicular wave propagations, k k ⊥ .We note that a gradient in plasma density contributes to an instability at any drift velocity (second term in the parenthesis of Eq. 6).We will argue that for the relevant plasma conditions analyzed in the following, we can ignore large scale plasma density gradients perpendicular to B. The relative drift velocity V d between electrons and ions has to exceed the ion sound speed C s in order to give unstable waves, otherwise it has a damping effect.In this simple model, the first waves to become unstable are those where k⊥B.Since ω ce ν e and ci ≤ν i for the relevant ionospheric conditions, waves with large k give large ϕ and therefore small ω r , and will consequently remain linearly stable for realistic values of V d .
The enhanced non-thermal fluctuations were first discovered by radar scattering off the ionosphere, and later investigated by in-situ measurements by instrumented rockets.In a sense, the rocket and the radar represent complementary types of diagnostics: the radar selects a constant wavelength determined by the wavenumber matching condition, while the rocket data are evidently dominated by the largest amplitude signal, irrespective of its characteristic wavelength.
While radar scattering can be an important diagnostic in some respects, it can evidently not provide detailed information on the space-time evolution of the instability and its saturated turbulent state.More information can be gained by an instrumented rocket traversing the unstable region, but even here only a time-varying signal will be available, reflecting the properties of the fluctuations along the rocket trajectory.(Rose et al., 1992).
In addition, we have a practical problem with rocket data: since no absolute potential reference ("ground"-potential) is available, potential variations have to be detected by taking the potential difference between two probes.In principle, potential variations could be measured with respect to the rocket body, but experience has shown that this gives rise to very "noisy" signals, presumably because the probes are in general outside the Mach cone, and the rocket body inside.The standard configuration (Bahnsen et al., 1978;Rose et al., 1992), as addressed also in this work, consists of booms carrying the probes, as illustrated in Fig. 1, where the potential differences can be obtained between probes on the same boom or alternatively between probes on different booms.
It is not evident that the available probe difference signals are sufficient for recovering features of intermittency effects that may be present in the plasma turbulence.This question is addressed in the present study.We compare here data from numerical simulations with those obtained from the insitu rocket observations.In the numerical simulations, all information is available, in principle, and we can here unambiguously identify intermittency effects as represented by the structure functions.The input data for the simulations are chosen to be representative of the extreme values of the ionospheric parameters at the time of launch.The simulation data are analyzed by two methods: a simple analysis of one-point structure functions followed by an analysis where we "mimic" the signal as it would be obtained by a potential difference measurement carried out by an instrumented rocket.For the ionospheric rocket only the latter option is available.Our analysis includes structure functions up to 8th order, being aware that the accuracy of the estimate decreases for increasing order.Our analysis refers, as stated, to one particular plasma instability.The Farley-Buneman instability is driven by a current (i.e. the E 0 ×B-electron flow through unmagnetized or weakly magnetized ions), and is thus likely to have properties in common with other current driven instabilities.We therefore anticipate that our results are qualitatively relevant for other plasma instabilities.

Numerical simulations
The numerical simulations were conducted in two spatial dimensions in the plane perpendicular to the imposed magnetic field, using a Particle-in-Cell (PIC) code (Birdsall and Langdon, 1991) for the ion component and a fluid model for the electrons (Oppenheim et al., 1995;Oppenheim and Otani, 1996;Dyrud et al., 2006).In the present analysis, the electron inertia is ignored.Details of the simulation code are presented by Oppenheim et al. (2003).We have performed also smaller box simulations with the same parameters used in the simulations shown here, but with finite electron inertia, and found no substantial difference in the resulting evolution or spectral characteristics of the present parameters.The plasma parameters used in the present study are summarized in Table 1.For E 0 we use the largest value that is relevant for the rocket experiment discussed in Sect. 4. For our conditions, the recombination coefficient is β R ≈3×10 −7 cm −3 s −1 .Recombination effects are not included in the simulations, since they give only small corrections for the present strongly driven case (Fejer et al., 1984).We deal with low-β plasmas, and the magnetic field is assumed constant.We use a value for the electron temperature which is consistent with the present observations.Evidence can be found for anomalous electron temperature enhancements for increasing dc electric field in the ionospheric Eregion (St.-Maurice et al., 1999;Noël et al., 2005), where the effective electric field needs to be considered in case we have neutral winds.Unfortunately, we have no means for obtaining information concerning neutral winds for the ROSEexperiment.For E 0 ≈40 mV/m, i.e. for the up-leg conditions of the rocket experiment discussed in Sect.4, the increase in T e is expected to be minute, but for the somewhat larger down-leg fields, E 0 ≈60-70 mV/m, nontrivial enhancements of T e are anticipated, but not observed for the present conditions.The wave propagation velocities, for instance, as found by Iranpour et al. (1997), Krane et al. (2000) and Dyrud et al. (2006) are best explained by an electron temperature of approximately 400 K. Also other reports (Pfaff et al., 1992) noted the lack of electron temperature enhancements for conditions similar to ours.A previous study (Dyrud et al., 2006) attempted to explain the low electron temperatures by thermal conduction to the colder regions below the enhanced wave activity, but used too low numerical values for the electron energy loss per collision, by taking this energy loss to be at most an order of magnitude larger than for inelastic collisions.By far, the dominant cooling rate is due to inelastic collisions in the E-region.It may be that the analysis of Dyrud et al. (2006) applies for conditions in laboratory experiments, where the dominant collisional electron energy loss will usually be for elastic collisions with neutral inert gases.The collisional energy losses for elastic collisions are generally much smaller than for the inelastic collisions.
In the related rocket experiments over Greenland (Bahnsen et al., 1978;Pécseli et al., 1989) the electron temperature was determined.The propagation speed for the fluctuations that was found there agreed well with the sound speed obtained by an average ion mass and the experimentally obtained electron temperature (Pécseli et al., 1989).
An illustrative result from the simulations is shown in Fig. 2 for three times in physical units, t=5, 10, and 22 ms.The axes are in physical units as well.We note the evolution of small scale structures in the linear initial phase of the instability.Eventually, in the nonlinear phase, larger scale structures develop and a saturated turbulent stage of the instability is reached.Typically, the saturated potential fluctuations have a characteristic wavelength of ∼2 m, and a peak value of ∼0.3 V.A typical root-mean-square value of the potential fluctuations is ∼0.08 V, corresponding to acelectric fields ∼3×10 −2 E 0 , for the given conditions.The fluctuations in density are relatively modest, typically below 20%, even though we can observe larger spikes (Dyrud et al., 2006).
A sample of the time-series is shown in Fig. 3.We select 25 such series, taken at separations corresponding to 3 m in the ionosphere, as the basis for the structure function analysis.These separations are sufficiently large to let us assume The initial red part contains the non-stationary initial growth phase and is omitted from the analysis.We verified that our results are robust with respect to small variations in the lengths of the omitted time-sequences.
that the time-evolutions of the small scale structures are statistically independent.Each of these samples contains approximately 500 time steps.We omit the initial ∼200 time samples when analyzing the data, since they contain an initial non-stationary exponential growth phase.The amplitude probability density for the signals used in the analysis are given in Fig. 4. For this figure we used all data from the available points and all the time-series, except the initial omitted part.Figure 4

Data analysis
We study the structure functions associated with temporal and spatial variations of the signal.To study the time variations, we consider a set of 25 time-series for the fluctuating potential φ obtained in a 5×5 grid with 3 m separation.This grid should not be confused with the simulation grid, which is much finer.To study the fixed-time spatial variations, we consider 25 samples with full spatial resolution, taken at different times in the saturated stage.

One-point statistics for temporal variations
We obtain first the temporal structure functions |φ(t 1 )−φ(t 2 )| n for n=1,. .., 8, with t 1 and t 2 being two times in the same record.The averaging is performed over the individual time samples and then over the 25 sets of the data.Results are shown in Fig. 5 in a double logarithmic presentation.The structure functions are normalized to the first time sample.Note that for n≥7 we have a non-trivial uncertainty in the estimate of the corresponding structure function.We perform a power-law fit t α n to these structure functions in the interval 3-10, showing in Fig. 6 the exponent α n for different values of n.By varying the length of the time interval used for obtaining the structure functions we find that the values of α n up to n=6 are robust, while they become increasingly uncertain for larger n.For n>8, we do not consider the estimates for α n to be reliable.The power-law index α n of the structure functions shown in Fig. 6  expected for Gaussian random signals.This deviation is conspicuous already for n=4.In Appendix B we give a more detailed discussion of the uncertainty of the estimators of the structure functions due to finite record lengths.

Two-point potential difference statistics, temporal variations
We note that the structure functions obtained by the foregoing analysis can not be directly compared to rocket observations as obtained by many instrumented rockets, for reasons outlined in the introduction.In order to make the analysis more directly relevant for comparisons with rocket data, we consider the potential difference between two positions separated by 3 m, which is representative for many rocket geometries, in particular also for those to be discussed later in this paper (Rose et al., 1992).This difference can be taken in basically two directions in the available two dimensional geometry.The corresponding values of the exponent are denoted by α ⊥ and α , respectively, where the subscripts ⊥ and refer to the E 0 ×B-direction.Figure 7 shows the variation of α ⊥ and α with n.We find these results to be significantly different from those summarized in Fig. 6.Thus, with relevant separations, the two-point difference signal has statistical properties significantly different from those found when analyzing the one point signal.
The results summarized in Fig. 7 corresponds to a rocket at rest in the ionosphere.To make the analysis more realistic we should in principle analyze the simulated signal corresponding to a spinning and coning rocket moving along a prescribed trajectory.A complete analysis taking into account Nonlin.Processes Geophys., 15, 847-862, 2008 possible values of all the parameters entering the problem will make the analysis extremely lengthy.We argue in the following that the modifications are unlikely to be significant for the comparison with the available rocket data, to be discussed later.

Structure functions for spatial separations at fixed times
The analysis of Sect.3.1 and 3.2 refers to temporal separations when calculating the structure functions.Similar results can be obtained by a spatial sampling of the potential at a fixed time, and then varying the separation.
We now obtain the structure functions of the potential difference between two spatial positions separated along x and along z, respectively.Based on a set of 5 samples at different times distributed over the available time-interval for the fluctuating potential φ in the entire available plane, we obtain the structure functions |φ(x 1 , z)−φ(x 2 , z)| n and |φ(x, z 1 )−φ(x, z 2 )| n for n=1, . .., 8.The averaging is performed over the spatial samples and then over the 5 sets of data, taken at different times.Results are shown in Fig. 8 in a double logarithmic presentation.The exponents corresponding to the spatial structure functions are shown in Fig. 9 for different n.
The relation between the spatio-temporal variation of the ionospheric signal and the time varying signal obtained from the rocket have already been discussed by Pécseli et al. (1989).Basically, we find a Doppler shift due to the rocket motion and a frequency and amplitude modulation due to the rocket spin.If the rocket spin frequency is small as compared to relevant wave-frequencies (which  8. Fixed-time structure functions as a function of spatial separations, in units of numerical spatial samples.On top we have spatial separations in the z-direction (perpendicular to E 0 ×B) and in the bottom for the x-direction (parallel to E 0 ×B).The order parameter n=1, . .., 8 increases from bottom to top.Note that for n≥7 we have also here a non-trivial uncertainty in the estimate of the corresponding structure function.
is often the case), we may ignore the latter effects.Concerning a rapidly moving rocket, with velocity U , we may argue that the Taylor hypothesis (or the "frozen turbulence approximation") can be applied (Shkarofsky, 1969).Physically, the Taylor hypothesis assumes the transit time of the turbulent eddies to be much less than the characteristic evolution time, implying that the observed frequencies can be approximated by the Doppler shifts.Under relatively mild assumptions (Shkarofsky, 1969) 9. Variation with n of the two exponents α x and α z , for the potential structure functions for spatial separations taken at a fixed time, see also Fig. 8.The power-law variation is fitted in the interval 3-15 spatial lags, except for n=8 where it is 3-12.The circles refer to the difference in the z-direction, asterixes to differences in the x-direction.conditions, we find the following relation with φ(t)≡φ(t 1 )−φ(t 2 ) and φ(r)≡φ(r 1 )−φ(r 2 ), with t≡t 1 −t 2 and r≡r 1 −r 2 .Alternatively, we can assume the observer to be fixed and the wave field to be propagating with a large velocity U , and the Taylor hypothesis can again be applied.The plots in Fig. 2 indicate that all frequency components propagate uni-directionally, at least to a good approximation.The propagation velocity is approximately the sound speed.In either case, the spatial and temporal correlation functions, and consequently also the corresponding structure functions, will be related as 2 φ(r=0, t) ≈ 2 φ(r=U t, t=0) .If we assume that both 2 φ(t) and 2 φ(r) have a power-law variation in a significant interval, we can argue by Eq. ( 7) that the characteristic exponents α t and α z for the two cases are directly related, α t ≈α z , taking the z-coordinate to be along U .By comparing Figs. 6 and 9 we find a reasonable agreement up to n=3 for the results argued on basis of the Taylor hypothesis.However, the latter is not applicable on the structure functions for spatial separations perpendicular to U .

Discussion of the analysis of the simulation data
We find that the first few values of α n are close to be directly proportional to n, while for n>3 we find a pronounced deviation from a linear variation.It thus seems that we have found a clear indication for intermittency effects in the turbulence.For n=1, 2 the values of α-exponents in Fig. 6 and 9 are close, indicating that the frozen turbulence approximation has a certain region of applicability, but also note that the differences increase rapidly with n as soon as n≥3, indicating that the higher order structure functions are very sensitive to deviations from the assumption of frozen turbulence.It is also interesting that for n=1 and n=2 we find no difference between the α-values in Fig. 9, which could be taken as a sign of spatial isotropy, i.e. lack of distinction between the x and z-directions for the smallest scales.This fact is consistent with an intermittency model based on secondary instabilities associated with gradients of larger scale structures, which are in turn again influenced by even larger scale structures.When we approach the smallest scales being resolved, we have an approximate local isotropy of the spatial potential variations.Similar observations have been made for electrostatic drift wave turbulence (Okabayashi and Arunasalam, 1977;Pécseli, 1982), where similar arguments apply (Hallatschek and Diamond, 2003).

Analysis of rocket data
We analyze data from the ROSE4 rocket (Rose et al., 1992).The ionospheric conditions and details of the instrumentation relevant for the present dataset were discussed in a special issue of Journal of Atmospheric and Terrestrial Physics (54,1992).Here we present only a short summary: dc-electric field values of typically 40 and 70 mV/m were measured by the rocket instruments on up-leg and down-leg passages of the E-region, respectively.The corresponding E 0 ×B/B 2 velocities are approximately 800 and 1400 m/s.These values are of a sufficient magnitude to excite the Farley-Buneman instability.The threshold value for the electric field is approximately 20 mV m −1 .
The extremely low frequency (ELF) signals analyzed here were obtained by means of gold-plated spherical probes of 5 cm diameter, mounted on two pairs of booms, one near the top of the payload (labeled 1 and 2) and the other 185 cm lower (labeled 3 and 4), oriented at an angle of 90 • with respect to the first pair, as illustrated in Fig. 1 (Rinnert, 1992).The length of each boom was 180 cm.We analyzed the following fluctuating signals where φ j (t) for j =1, 2, 3, 4 is the potential on the j -th probe with respect to a suitably defined common ground.There is an evident redundancy in the available signals, which can be used to check the performance of individual probes.For wavelengths much larger than the probe separations, it is evident that the potential difference signals can be used to estimate the fluctuating electric fields, while the interpretation becomes more complicated when the spectra contain wavelengths comparable to or smaller than the probe separation.In general, the difference signal can be interpreted in terms of a filtering operation of the spatial potential variation.
The space-time varying electric field fluctuations of the electrojet were originally sampled with a 4 kHz sampling frequency.By averaging sampling points two-by-two, we here increased the sampling interval to 0.5 ms, giving a Nyquist frequency of 1000 Hz.The electric circuits give an effective frequency limitation being noticeable for frequencies larger than 600 Hz.The signals were digitized with 12 bit resolution.Amplitudes of the potential differences were typically in the range 15-30 mV.The amplitudes of the relative density fluctuations, n/n 0 , were in the range 1-3%.The amplitude probability density of the detected potential difference fluctuations is non-Gaussian (Larsen et al., 2002), but it is important to emphasize that this conclusion refers to the filtered signal.The probability amplitude of the non-filtered signal is, on the other hand, significantly affected by the rocket spin.We have analyzed all probe-combinations, but show here only results for U 6 , the differences between e.g.U 6 and U 2 being small, these two difference signals referring to probe-sets perpendicular to the rocket axis.For the small time separations relevant for the present analysis we do not find significant differences between signals such as U 6 and U 5 either, but note that the correlation times for these signals are somewhat different (Iranpour et al., 1997;Krane et al., 2000).
In Fig. 10 we show the structure functions n φ(t) for n=1,. .., 8, as a function of temporal separations.These data are obtained from the up-leg part of the flight, where the fluctuation amplitude level is somewhat smaller than for the down-leg part.We note an overall similarity with the results obtained from numerical simulations (see Fig. 5).Also in the present case we can fit a power-law t α n and show in Fig. 11 the variation of α n with n.The points lie on a curve which is close to a straight line, so we also show the compensated variation α n /n.In this representation, the deviations from the Gaussian results become more noticeable.
In order to test the significance of the result, we also carry out a similar analysis for randomized (or "surrogate") data (Schreiber and Schmitz, 2000;Wernik, 1996;Pécseli and Trulsen, 1993), with results shown in Fig. 12, including also here the variation of the compensated exponents.randomizing the phase information by a standard random number generator.The power-spectrum of the resulting dataset is then the same as the one obtained for the original data, but the relative phases of the various Fourier components are unrelated to the previous ones.Coherent structures that may be present in a signal will be characterized by distinct phase relations of their Fourier components.The randomization of the phases will consequently destroy such possible coherent structures.These surrogate data have approximately Gaussian properties and the compensated exponents α n /n should be approximately constant.The properties of the surrogate data depend also on the seed of the random number generators.In order to demonstrate the effects of this dependence, we generated many sets of synthetic data, and show the range of variability of the resulting values of α n by a green-shaded area in Fig. 12.The thin red line gives the average curve.The results in Fig. 11 fall outside the shaded region for large n, but only marginally so.
We extended the analysis to include also a sequence from the down-leg part of the flight in a time-interval 250.0-254.1 s, where the overall fluctuation level is somewhat increased as compared to the up-leg part (Krane et al., 2000;Dyrud et al., 2006).In Fig. 13 we show the structure functions for varying order as functions of temporal separations for this down-leg time interval.The analysis of the exponents have been carried out also for these data, as shown in Fig. 14.In particular, Fig. 15 displays the analysis of surrogate data, as in Fig. 12.Also here, we show the statistical scatter by a green shading.The values of α n found in the original dataset fall somewhat outside the shaded region, in particular for large n.The slight narrowing of the green areas in Figs. 15 and 12 are an artifact due to the finite number of realizations of the random number generator seeds.For the surrogate data, the statistical spread in the estimate on α n /n increases nonlinearly with n.
The statistical significance of the results in Fig. 13 is better for the down-leg part of the flight as compared to the results shown in Fig. 10.This observation is consistent with the increase in the fluctuation level, which is expected to enhance the phase couplings in the turbulent spectrum.

Discussion of the analysis of the rocket data
Comparing the results from analyzing the data from the numerical simulations with the corresponding analysis from the rocket data we find somewhat similar results.The part of the analysis of the simulation data where a comparison is appropriate (i.e.what concerns the potential differences) the structure functions have a subrange characterized by a clear power law variation.The exponent varies systematically with the order n of the structure function, but the compensated exponent is not a constant.The electrojet turbulence is intermittent in the sense that it has measurable differences from the results expected for Gaussian random processes.This conclusion is supported by the analysis of the surrogate data, emphasizing the statistical Fig. 15.Results for the synthetic data, showing the variation of the exponent with the order of the structure function together with the corresponding variation of the compensated exponent for varying order of the structure function.The results correspond to Fig. 14.The green shaded areas represent also here the statistical scatter obtained by varying the seed of the random number generators for the surrogate data.
significance of the results.It is however clear that the intermittency effects found in the numerical simulations are more evident than in the ionospheric data.
We emphasize one of the basic differences between the data from the numerical simulations and those originating from the ROSE4-rocket: Apart from an initial growth phase, the simulation data represent approximately time-stationary and spatially homogeneous (but anisotropic) plasma turbulence.The rocket data, on the other hand, represent long time-records obtained by instruments traversing changing plasma conditions, with nontrivial differences in the driving dc-electric field in the upleg and downleg conditions.

L. Dyrud et al.: Structure functions and intermittency in plasma turbulence
The approximate local isotropy of the smallest scales discussed in Sect.3.4 makes the rocket spin immaterial for the present analysis.

Conclusions
We have demonstrated the presence of a power-law subrange for the structure functions associated with the electrostatic potential in turbulent plasma fluctuations for conditions appropriate for the ionospheric E-region.We found clear indications for intermittent fluctuations in the sense that the power-law index for the structure functions of order n deviates from a simple n-proportionality.These features are summarized in Figs. 6, 7 and 9.In order to explain the physical reason for the intermittency, we give particular attention to secondary instabilities developing on the gradients of larger scale structures (Sudan, 1983).We note that such secondary instabilities can be found also for other instabilities (Hallatschek and Diamond, 2003).These secondary instabilities are, by nature, of a "bursty" appearance, requiring the presence of local large scale gradients, associated with a long wavelength component.The presence of secondary instabilities could be anticipated already by inspections of Eq. ( 6), where the local gradients can be considered as being associated with large scale waves.As evident from the analysis, our diagnostic is based on structure functions of the electrostatic potential.Other related works (Tam et al., 2005) are based on wavelet transforms.They studied the degree of intermittency on different scales and found electric field fluctuations to be more intermittent on smaller scales.
Similarly, we can argue 2 φ(t) ≈ (∂φ(t)/∂t| t=0 ) 2 t 2 , see also the discussion in Appendix B. The origin of time (as well as of position) variables is arbitrary because of the stationarity and homogeneity of the turbulence.We observe neither any t 2 nor an r 2 dependence of the structure function, implying that the range of validity of the previous approximation is very limited, and most likely constrained by collisions.This collisional timescale is not resolved by the simulation, nor by the sampling period of the rocket instruments.The relevant smallest length scales are not resolved by the simulations.
We found several interesting features of the ionospheric plasma turbulence.First of all, intermittency, as evidenced by a lack of proportionality between the exponents α n and the order n of the structure function, is much more evident for ionospheric turbulence as compared to turbulence in neutral incompressible flows (Anselmet et al., 1984).For n≥4 there is a significant difference between the fixed-position temporal intermittency (see Fig. 6) and the one associated with the fixed time spatial-difference variable (see Fig. 9).On the other hand, we note some overall similarity between the variation with the order parameter n of the structure functions obtained for the time varying potential difference between two fixed positions (see Fig. 7), and the structure function taken at a fixed time with varying spatial separation (see Fig. 9), both cases referring to numerical simulations.
We find that the low frequency electrostatic turbulence in the ionospheric E-region is likely to be strongly intermittent for dc-electric field values that are common (i.e. in excess of 50 mV/m), but on the other hand we also find that standard rocket probe set-ups, as illustrated in Fig. 1, are not well suited for recovering such features.Evidence for intermittency can be found, but only by detailed investigations of the data, where the use of surrogate data can be an important tool for assessing the statistical significance of the results (Wernik, 1996).
The power-law exponents α n with n=1 found in the simulations are somewhat smaller than those for the rocket data, even when we consider the down-leg part which is most unstable, and the difference becomes more conspicuous for n≥2.The simulations show slightly stronger intermittency effects than the rocket data even when potential difference signals are considered.Part of the explanation deals with the sampling rate of the rocket data, which is too small to resolve the smallest time-scales.Similarly, the grid-resolution and the finite time-step in the numerical simulations prohibit the finest details of the space-time variations of the physical instability to be resolved completely.Considering these shortcomings, we might argue that the magnitudes of the exponents α n for the simulations and the rocket observations for n=1 and 2 agree quite well.
The parameters chosen for the simulations are representative for the most unstable conditions on the down-leg part of the rocket flight.If we average over the entire up-leg and down-leg parts, we find average electric fields smaller than the 70 mV/m used here.In spite of the strong fluctuation levels, we find that the detection method based on potential differences between two probes with a large separation gives signals that are close to exhibiting characteristics of Gaussian signals.We recover the strongly intermittent features only by making a one point analysis of the data.Rockets equipped as for the Rose campaign (Rose et al., 1992) are very useful for detecting the bulk features of plasma conditions and fluctuations, but inadequate as soon as finer details, such as intermittency effects, of the E-region turbulence are studied.Some rockets, the TOPAZ II and TOPAZ III rockets for instance, has a somewhat different set-up (Vago et al., 1992) and it may be worthwhile to investigate the signals from these probes for studying intermittency effects.
Indirectly, our results thus emphasize the importance of detailed numerical simulations and laboratory experiments for the understanding of these instabilities.We note for instance that for waves propagating exactly perpendicular to the magnetic field (at zero aspect angle) we have one threshold E 0 /B-velocity and at larger aspect angles we have a slower one (Pécseli et al., 1989).The transition is due to the change in electron dynamics, which is adiabatic for small aspect angles and isothermal for larger angles (Pécseli et al., 1989).The zero aspect angle is likely to be a part of the earlier evolution of the amplitude.The full evolution of the structures may be three dimensional and such that once the amplitude has increased enough for the growth rate to slow down through the nonlinear effects (St.-Maurice and Hamza, 2001), then shears and rotations can introduce a fast evolving aspect angle that destroys the structures while heating the electrons (J.-P.St.-Maurice, private communications, 2008).The largest amplitudes may be met when the phase speed has slowed down to be the threshold speed, i.e. isothermal ion-acoustic speeds at large aspect angles.Also the altitude dependence of the collision frequency can introduce an important aspect angle effect on the properties of the non-linear wave structures as they approach saturation.We find it unlikely that these details can be recognized by an instrumentation as the one shown in Fig. 1, and foresee that numerical simulations can have an important role in this discussion.The bulk of the rocket observations outlined here (although not in all detail, as discussed by Dyrud et al., 2006) can be accounted for by a two-dimensional numerical simulation as the one discussed in the present work.
We emphasize that the structure functions as obtained in the present study refer to relatively small spatial and short temporal scales.We might add a large amplitude, low frequency, long wavelength component which will make any signal significantly non-Gaussian, but such a wave will have negligible consequences on the present structure functions, by adding a slowly varying bias to our data.

Appendix A
The correlation function ρ is related to the power spectrum S of the fluctuations by the Wiener-Khinchin theorem (Bendat, 1958).Considering the case with temporal variables we have ρ=ρ(τ ) where τ =t 1 −t 2 .The frequency power spectrum is then obtained by the cosine transform of ρ(τ ).Assuming that we have a range {τ a :τ b } of τ -values where ρ≈1−Aτ α we have For ω-intervals where the three integrals in Eq. (A1) are slowly varying with ω, we have a power spectrum S(ω)∼1/ω α+1 in that interval, relating the exponent in the power-spectrum to the exponent in the structure function.
For spatial separations, we have similar expressions in terms of wavenumbers (Hinze, 1975).If we, as an illustration, consider again the universal range of the second order structure function in fully developed incompressible turbulence, we have an ∼( r) 2/3 variation in terms of the separation r and the specific energy dissipation rate , while the wave-number power spectrum varies as ∼ 2/3 k −5/3 , consistent also with the foregoing estimates.
A spectral representation can be convenient from an experimental point of view and several studies of plasma turbulence analyzed turbulent spectra.Results from laboratory experiments that were particularly relevant for the E-region fluctuations (Mikkelsen and Pécseli, 1980;Pécseli et al., 1983) have been compared to spectra obtained from rocket experiments (Krane et al., 2000).

Appendix B
The present appendix deals with the consequences of finite time sequences for the estimates of the structure functions.The analysis is limited by considering only a time-varying signal.We consider this as an insignificant restriction.
Taking one record, we can obtain an estimate for the n-th order structure function where φ(t a , t b )≡|φ(t a )−φ(t b )|, noting the implied simplifying assumption that the integration interval T is independent of τ .The estimate E n =E n (T , τ ) is statistically varying over the ensemble of realizations (Bendat, 1958;Pécseli, 2000).It has an average value We now assume that the process is time stationary.We then have n φ(t 1 , t 1 +τ ) = n φ(0, τ ) ≡ n (τ ) independent of t 1 , giving E n = n (τ ), since the integral in Eq. (B1) becomes trivial.This result for E n will be used in the following.
First we assume that τ is small, so we can make the approximation φ(t a , t b )≡|φ(t a )−φ(t b )|≈|φ a | |t a −t b |, where φ a ≡dφ/dt| t=t a .In the limit of small τ we have n (τ ) ≈ |φ | n τ n , where we here can omit the subscript on φ because of the assumed time-stationarity of the process.Similarly, we have ℵ n (τ, ν)≈ |φ 1 | n |φ 2 | n τ 2n for small τ .Consequently, we have for this limiting case where |φ 1 φ 2 | n implicitly depends on ν≡t 2 −t 1 .For large ν where φ 1 and φ 2 can be assumed to be statistically independent, we have , so the integrand vanishes in this limit.At ν≈0 we have ≥ |φ | n 2 by the Schwartz-Cauchy inequality.
The result Eq. (B7) demonstrates that the root-meansquare error σ n increases with time separation as τ n for constant T .This observation can be used for fixed τ and varying n, or vice versa.
In order to illustrate the variation with T , we postulate a simplified model of the correlation function for |φ | n in the form ρ n (ν)≈ |φ | 2n exp(−|ν|/τ c )+ |φ | n 2 , with τ c being a correlation time, and use this model in Eq. (B7).The proposed correlation function will be accurate for the case where the potential derivative |φ | is a random Gaussian Markov process with non-zero mean (Bendat, 1958;Pécseli, 2000).In Fig. B1 we show the variation of the normalized variance σ n τ n |φ | 2n 1/2 with varying normalized record length, T /τ c .It is interesting that the result shown in Fig. B1 is independent of n for this model.We find that the signal to noise ratio is significantly improved when T is increased from zero, but an increase from 10T /τ c to 20T /τ c gives a comparatively much smaller improvement.As far as the rocket data are concerned, we have approximately 10 s of data for upleg and a similar time interval for the downleg conditions, to be compared with typically 30-50 ms correlation times (Krane et al., 2000;Dyrud et al., 2006), which is much shorter than the available record length.For the numerical simulations, we have shorter time durations in comparison, but have here 25 records available for averaging.
The present analysis assumes small τ .For arbitrary τ we have to model the entire variation of ℵ n (τ, ν).This will not be discussed here.A discussion of the uncertainty of the estimate of structure functions for the spatial variations of the potential at a given fixed time can be carried out as shown before and need not be discussed here.
The error estimates discussed in this appendix assumes that we have one time record available.In the numerical simulations discussed in Sect. 2 we had 25 such records of equal lengths.The estimate E n can be generalized as 1 T T 0 n φ m (t 1 , t 1 +τ )dt 1 , with m being the label of the record.In our case we have N=25.If we can assume the N records to be statistically independent (in practice by assuming that they are obtained at positions separated by more than a correlation distance), it is relatively straightforward to generalize the foregoing analysis.We find that the error σ n scales as ∼1/ √ N. The present appendix assumes continuous functions, where we in our foregoing analysis had sampled space-time varying functions.It is however evident even from the present appendix that the number of samples alone can not determine the accuracy of an estimate.If we have a very dense sampling within a time sequence shorter than the correlation time, our estimate will be inaccurate under all circumstances.It is important to distinguish the number of samples of φ(t a , t b )≡|φ(t a )−φ(t b )| obtained from different and statistically independent realizations, and the number of time samples in one record.

Fig. 1 .
Fig.1.Schematic diagram for the positioning of the probes on some instrumented rockets, with length scales relevant for the present analysis of data from the ROSE experiment(Rose et al., 1992).

Fig. 2 .Fig. 3 .
Fig. 2.Summary plots illustrating the electrostatic potential for three times as obtained from the numerical simulations.The magnetic field is perpendicular to the plane of the paper, and the E 0 ×B-drift is in the vertical direction, with E 0 =70 mV/m in the positive x-direction.

Fig. 5 .
Fig.5.Structure functions as functions of time separations ("time lags"), in units of numerical time samples.The electrostatic potential is sampled in one fixed spatial position.The order parameter n=1. ..8 increases from bottom to top.

Fig. 6 .
Fig.6.Exponents in the fit t α n of the fixed-point temporal structure functions for different values of n, see also Fig.5.The power-law variation is fitted in the interval 3-10 time lags.

Fig. 7 .
Fig. 7. Variation with n of the two exponents α ⊥ and α , for the temporally varying potential difference structure functions, taken at two spatial positions separated by 3 m, in the x (asterixes) and z (circles) directions, respectively.
Fig.8.Fixed-time structure functions as a function of spatial separations, in units of numerical spatial samples.On top we have spatial separations in the z-direction (perpendicular to E 0 ×B) and in the bottom for the x-direction (parallel to E 0 ×B).The order parameter n=1, . .., 8 increases from bottom to top.Note that for n≥7 we have also here a non-trivial uncertainty in the estimate of the corresponding structure function.
Fig.9.Variation with n of the two exponents α x and α z , for the potential structure functions for spatial separations taken at a fixed time, see also Fig.8.The power-law variation is fitted in the interval 3-15 spatial lags, except for n=8 where it is 3-12.The circles refer to the difference in the z-direction, asterixes to differences in the x-direction.

Fig. 11 .
Fig. 11.Variation of the exponent with the order of the structure function, shown together with the corresponding variation of the compensated exponent for varying orders of the structure function, corresponding to Fig. 10.The time interval is 112.0-116.1 s, during the up-leg part of the flight.The dashed line on the figure for α n for varying n is determined by the spectral index of the frequency power spectrum as discussed in the appendix.

Fig. 12 .
Fig. 12.Results for the synthetic data, corresponding to Fig.11, showing the variation of the exponent with the order of the structure function together with the corresponding variation of the compensated exponent for varying order of the structure function, corresponding to Fig.11.The green shaded areas represent a statistical scatter obtained by varying the seed of the random number generators for the surrogate data.

Fig. 13 .
Fig. 13.Structure functions for the time interval 250.0-254.1 s on the down-leg part of the flight.

Fig. 14 .
Fig. 14.Variation of the exponent with the order of the structure function, shown together with the corresponding variation of the compensated exponent for varying order of the structure function, corresponding to Fig. 13.The time interval is 250.0-254.1 s on the down-leg part of the flight.

Table 1 .
Input data for the numerical simulations.