Tsallis statistics as a tool for studying interstellar turbulence
DD RAFT VERSION , O
CTOBER
24, 2018
Preprint typeset using L A TEX style emulateapj v. 08/22/09
TSALLIS STATISTICS AS A TOOL FOR STUDYING INTERSTELLAR TURBULENCE
A. E
SQUIVEL AND
A. L
AZARIAN Instituto de Ciencias Nucleares, Universidad Nacional Autónoma de México, Apartado Postal 70-543, 04510 México D.F., México Astronomy Department, University of Wisconsin-Madison, 475 N.Charter St., Madison, WI 53706-1582, USA
Draft version, October 24, 2018
ABSTRACTWe used magnetohydrodynamic (MHD) simulations of interstellar turbulence to study the probability distri-bution functions (PDFs) of increments of density, velocity, and magnetic field. We found that the PDFs are welldescribed by a Tsallis distribution, following the same general trends found in solar wind and Electron MHDstudies. We found that the PDFs of density are very different in subsonic and supersonic turbulence. In order toextend this work to ISM observations we studied maps of column density obtained from 3D MHD simulations.From the column density maps we found the parameters that fit to Tsallis distributions and demonstrated thatthese parameters vary with the sonic and Alfvén Mach numbers of turbulence. This opens avenues for usingTsallis distributions to study the dynamical and perhaps magnetic states of interstellar gas.
Subject headings:
ISM: general — MHD — turbulence INTRODUCTIONTurbulence is known to be a crucial ingredient in manyastrophysical phenomena that take place in the interstellarmedium (ISM): such as star formation, cosmic ray dispersion,formation and evolution of the (also very important) mag-netic field, and virtually every transport process (see for in-stance the reviews of Elmegreen & Scalo 2004; Mac Low &Klessen 2004; Ballesteros-Paredes et al. 2007; McKee & Os-triker 2007, and references therein). This has brought manyresearchers to the study of turbulence, from both observa-tional and theoretical perspectives.From the observational point of view there are several tech-niques aimed to study properties of ISM turbulence. Scintil-lation studies (see for instance Narayan & Goodman 1989;Spangler & Gwinn 1990) have been very successful to char-acterize ISM turbulence, however they are limited to den-sity fluctuations in only ionized media. For neutral mediathe density fluctuations can also be obtained from columndensity maps (see Padoan et al. 2003). Other avenue of re-search is provided by radio spectroscopic observations, usingfor instance the centroids of spectral lines (von Hoerner 1951;Münch 1958; Dickman & Kleiner 1985; Kleiner & Dickman1985; O’Dell & Castaneda 1987; Falgarone et al. 1994; Mi-esch & Bally 1994; Miesch & Scalo 1995; Lis et al. 1998; Mi-esch et al. 1999), or linewidth maps (Larson 1981, 1992; Scalo1984, 1987). Spectroscopic observations have the advantageof containing information about the underlying turbulent ve-locity field. However, they are also sensitive to fluctuations indensity, and the separation of the two contributions alone hasproved to be a difficult problem (see Lazarian 2006b).Much of the effort to relate observations with models ofturbulence has focused on obtaining the spectral index (log-log slope of the power-spectrum) of density and/or velocity.To this purpose several methods have been developed andtested through numerical simulations. Among them, we canmention the “Velocity Channel Analysis” (VCA, Lazarian &Pogosyan 2000; Esquivel et al. 2003; Lazarian & Pogosyan2004; Lazarian et al. 2001; Padoan et al. 2003; Chepurnov& Lazarian 2009), the “Spectral Correlation Function” (SCF
Electronic address: [email protected],[email protected]
Rosolowsky et al. 1999; Padoan et al. 2001) , “Velocity co-ordinate spectrum” (Lazarian & Pogosyan 2008, 2006), or“Modified Velocity Centroids” (Lazarian & Esquivel 2003;Esquivel & Lazarian 2005; Ossenkopf et al. 2006; Esquivelet al. 2007; Chepurnov et al. 2006; Chepurnov & Lazarian2009; Padoan et al. 2009) Different methods are suitable indifferent situations, for instance VCA is adequate to retrievevelocity information from supersonic turbulence, while cen-troids are only useful in subsonic (or at most, transonic) tur-bulence.More sophisticated techniques for the analysis of the data,including using wavelets instead of Fourier transforms (Hen-riksen 1994; Stutzki et al. 1998) have been attempted. Arather sophisticated example of such analysis is a use ofthe Principal Component Analysis (PCA) (Heyer & Schloerb1997; Brunt & Heyer 2002).The techniques above provide the spectrum of turbulent ve-locity and/or density. Spectra, however, do not provide a com-plete description of turbulence. For example, Chepurnov et al.(2008) showed that fields with a very different distribution ofdensity could have the same power spectrum. Intermittencyand topology are some of the turbulence features that requireadditional statistical tools. Measures of turbulence intermit-tency include probability distribution functions (PDFs, see forinstance Falgarone et al. 1994; Lis et al. 1996, 1998; Klessen The VCA and SCF uses the same sort of two point statistics. An insub-stantial difference is that VCA stresses the use of fluctuation spectra, whileSCF uses structure functions. A more fundamental difference stem from thefact that the VCA is based on the analytical theory relating the observed spec-trum of intensity fluctuations and the underlying spectra of velocity and den-sity, while the choice of structure function normalization adopted in the SCFprevents the use of the developed theoretical machinery for a quantitativeanalysis (Lazarian 2009). We should mention, however, that the use of integral transforms, otherthan the the Fourier transform, provide a new way of getting the same in-formation. For instance, using Fourier transforms one can also attempt to getthe empirical relation between the underlying velocity and the statistics of thefluctuations. However, the VCA-obtained relations are very difficult to matchempirically because of the strong non-linear nature of the transforms. There-fore it is not surprising that the empirically obtained relations for the PCAdo not reproduce important cases discovered by the VCA. The PCA analysisfailed to reveal the dependence of the Position-Position-Velocity (PPV) fluc-tuations on the density, which should be the case for the shallow spectrum offluctuations. See more on the testing the PPV dependencies in Chepurnov &Lazarian (2009). a r X i v : . [ a s t r o - ph . GA ] D ec Esquivel & Lazarian2000), higher-order She-Lévêque exponents (She & Lévêque1994; Padoan et al. 2003; Falgarone et al. 2005; Hily-Blantet al. 2008). In addition, genus (see also Lazarian et al. 2002;Kim & Park 2007; Chepurnov et al. 2008) was used to studythe topological structure of the diffuse gas distribution.Astrophysical turbulence is very complex, e.g. it can hap-pen in the multi-phase media with numerous energy sources.Therefore the synergy of different techniques provide an in-valuable insight into the properties of the ISM turbulenceand the self-organization of interstellar medium. For in-stance, Kowal et al. (2007) and Burkhart et al. (2009) haveapplied several statistical measures, including a measure ofbispectrum, to a wide set of MHD simulations and synthetic-observational data. More recent work applying the same setof techniques to the Small Magellanic Cloud (SMC) data re-vealed the advantages of such an approach (Burkhart et al.2010). It is the consistency of the results obtained by differ-ent techniques, which made the study of turbulence propertiesin SMC more reliable.Thus it is important to proceed with the quest for new tech-niques for astrophysical studies (see Lazarian 2009 and ref.therein). In this paper we add the Tsallis statistics to the ar-senal of tools available for ISM studies of turbulence. Thisis a description of the PDFs that have been successfully usedpreviously in the context of solar wind observations and itshares similarities to previous studies of ISM turbulence (e.g.Klessen 2000).The earlier solar wind relation between the Tsallis statisticsand turbulence was sought on the basis of 1D MHD simula-tions. We know that the properties of turbulence strongly de-pend on its dimensions. Thus we use 3D MHD simulations.We analyze a reduced set of MHD simulations, including col-umn density maps which are the most easily available datafrom observations.This paper is organized as follows. In section 2 we brieflypresent the Tsallis distribution, which we will use to char-acterize magnetohydrodynamic (MHD) simulations. The de-scription of the numerical models and the results (for 3D dataas well as column density) can be found in sections 3–5. Wediscuss the significance of our findings in 6. Finally, in section7 we provide our summary. TSALLIS STATISTICSBurlaga & Viñas (2004a,b, 2005a,b); Burlaga, Ness, &Acuña (2006a); Burlaga, Viñas, Ness, & Acuña (2006b);Burlaga & Viñas (2006); Burlaga, Viñas, & Wang (2007b);Burlaga, Ness, & Acuña (2007a, 2009), have studied thePDFs of the variation of the magnetic field strength in thesolar wind, measured by the
Voyager 1 and spacecrafts. Inparticular, they analyzed the fluctuations of magnetic field in-crements [ B ( t + τ m ) − B ( t )] as a function of the timescale τ m ,and showed that for a wide range of scales their PDFs couldbe described by a Tsallis (or q -Gaussian) distribution (Tsallis1988). The Tsallis PDF of increments ( ∆ f ) of an arbitraryfunction has the form: R q = A (cid:34) + ( q − ∆ f ( t , τ ) w (cid:35) − / ( q − , (1)where q is the so called “entropic index” or “non-extensivityparameter” and is related to the size of the tail of the dis-tribution, w provides a measure of the with of the PDF (re-lated to the dispersion of the distribution), while A measures the amplitude. Notice that the Tsallis distribution is symmet-ric . Varying the parameter q in the Tsallis distribution oneobtain distributions that range from Gaussian to “peaky” dis-tributions with large tails. The parameter q is closely relatedto the kurtosis (fourth order one-point moment) of the PDF,and similarly the parameter w is related to the variance of thePDF. However, they are not the same, while fitting a Tsallisdistribution to a PDF, the two parameters are calculated si-multaneously, therefore w will depend on the value of q andvice versa, while the dispersion of a distribution is indepen-dent of the kurtosis. In the next section we discuss how thisdifference might play in our favor. Furthermore, the Tsallisdistribution is often used in the context of non-extensive sta-tistical dynamics, it was originally derived (Tsallis 1988) froman entropy generalization to extend the traditional Boltzmann-Gibbs statistics to multifractal systems (such as the ISM). TheTsallis distribution reduces to the classical Boltzmann-Gibbs(Gaussian) distribution in the limit of q →
1. It provides aphysical foundation to other functions that have been usedto model empirically velocity PDFs in space plasmas (suchas the kappa function, see Maksimovic et al. 1997; Leubner2002). However, the physical interpretation of ISM simula-tions/observations PDFs in the light of such statistical dynam-ics is beyond the scope of this paper.Turbulence intermittency manifests itself in non-Gaussiandistributions, thus one can hope to establish a link betweenthe properties of turbulence and the Tsallis distributions pa-rameters obtained from simulations and/or observations. Thisprovides an alternative/complementary method to higher or-der statistics (Burkhart et al. 2009), or to other PDF analysis(such as those in Falgarone et al. 1994; Klessen 2000, wherethe kurtosis and/or skewness of the distributions is measured),all of which are sensitive to the non-Gaussian nature of realturbulence.Following the procedure presented by Burlaga and collabo-rators in what follows we make histograms of the distributionof fluctuations from our simulations. But, in contrast to theirstudies, where only the magnetic field variations were studiedat different times, we consider the spatial variations of density,velocity and magnetic field intensity at a given time, as a func-tion of a spatial scale (or lag) r . To make the fits and presentthe results we remove the mean value of the increments andscale them to units of the standard deviation of each PDF. Inother words, we have used ∆ f ( r ) = [ f ( r ) − (cid:104) f ( r ) (cid:105) x ] /σ f , with f ( x ) = ρ ( x + r ) − ρ ( x ), v x ( x + r ) − v x ( x ), v y ( x + r ) − v y ( x ), v z ( x + r ) − v z ( x ), B x ( x + r ) − B x ( x ), B y ( x + r ) − B y ( x ), and B z ( x + r ) − B z ( x ),where (cid:104) ... (cid:105) x stands for spatial average. We have neglectedthe anisotropy introduced by the magnetic field (see Cho &Lazarian 2009), the lag r is considered along each of the threecardinal directions and the resulting histograms include vari-ations in all of them. MHD SIMULATIONSWe have taken a reduced subset of the MHD simulationsof fully-developed (driven) turbulence of Kowal et al. (2007).We use the data output of the simulations to study the result-ing PDFs of density, magnetic field and velocity in differentturbulence scenarios.The code used (described in detail in Kowal et al. 2007) Measurements of velocity increments in the solar wind show asymmetricPDFs, which has lead Burlaga & Viñas (2004a) to add a cubic term to the ex-pression in eq. (1), in order to form a “generalized” (non-symmetric) Tsallisdistribution that can be used to describe skewed PDFs. sallis statistics 3
TABLE 1P
ARAMETERS OF THE SIMULATIONS
Model p B ext M s a M A b Description1 0.1 0.1 2.0 2.0 Supersonic & super-Alfvénic2 0.1 1.0 2.0 0.7 Supersonic & sub-Alfvénic3 1.0 0.1 0.7 2.0 Subsonic & super-Alfvénic4 1.0 1.0 0.7 0.7 Subsonic & sub-Alfvénic a Sonic Mach number, defined as M s = (cid:104) v / c s (cid:105) (averaged over the entirecomputational domain). b Alfvén Mach number, defined as M A = (cid:104) v / v A (cid:105) (averaged over the en-tire computational domain). solves the ideal MHD equations in conservative form: ∂ρ∂ t + ∇· ( ρ v ) = 0 , (2) ∂ρ v ∂ t + ∇· (cid:20) ρ vv + (cid:18) p + B π (cid:19) I − π BB (cid:21) = f , (3) ∂ B ∂ t − ∇× ( v × B ) = , (4)along with the additional constraint ∇ · B = 0 in a Cartesian,periodic domain. An isothermal equation of state p = c s ρ isused, where ρ is the mass density, v the velocity, B the mag-netic field, p the gas pressure, and c s the isothermal soundspeed. The term f in the equation of momentum (3) is alarge scale driving. The driving is purely solenoidal, ap-plied in Fourier space at a fixed wave-number k = 2 . . v is inunits of the rms velocity of the system, the Alfvén velocity v A = B / (4 πρ ) / is also in the same units. The magnetic fieldis of the form B = B ext + b , that is a uniform background field B ext plus a fluctuating part b (initially equal to zero). As it iscustomary, the units are normalized such that the Alfén veloc-ity v A = B ext / (4 πρ ) / = 1 and the mean density ρ = 1. Thisway, the controlling parameters are the sound, and the Alfvénspeeds (i.e. the gas pressure and the magnetic pressure, re-spectively), yielding any combination of subsonic or super-sonic, with sub-Alfvénic or super-Alfvénic regimes of turbu-lence. For this paper we explore the different combinations ofparameters summarized in Table 1.We should note that it has been recently found that a com-pressive component of the driving can have an important ef-fect on the statistical properties of the flow (Schmidt et al.2008; Federrath et al. 2008; Schmidt et al. 2009). This is aninteresting effect that require further study, but for the timebeing we will restrict ourselves to the usual divergence-freedriving. STATISTICS OF THREE-DIMENSIONAL DATABelow, we analyze the results of our simulations to testwhether the statistics of the PDFs of the fluctuations of thevariables can be described by the Tsallis formalism. To tothis we fit to each histogram a Tsallis distribution with theLevenberg-Marquardt algorithm (Press et al. 1992), the his-tograms and the fits are presented in Figure 1, the symbolsare the data from the simulations, and the lines are the corre-sponding fits. For visual purposes we have shifted verticallythe results by successive factors of 100 for each different lag.Similarly to the results from solar wind observations (Burlaga& Viñas 2004a,b, 2005a,b; Burlaga et al. 2006a,b; Burlaga & Viñas 2006; Burlaga et al. 2007a,b, 2009), the magnetic fieldcomponents, as well as the density and velocity components,show kurtotic PDFs that become smoother as we increase thescale-length (lag). One difference is the fact that the velocityincrements in our simulations are symmetric, as opposed tothe skewed distributions observed in the solar wind at smallscales. It is also interesting that (only) for supersonic turbu-lence density distributions are quite distinct from other quan-tities. Such behavior resembles the power spectrum, wherethe magnetic field and velocity for supersonic turbulence havesteep spectra regardless of sonic Mach number, while densitybecomes shallower as the Mach number increases (see for in-stance Esquivel & Lazarian 2005). Such shallower slopes ofthe density power-spectrum have been also reported when thedriving changes from solenoidal to compressive (potential)for a fixed ( >
1) sonic Mach number (Federrath et al. 2009).This means that shocks in supersonic turbulence, which createmore small-scale structure in density, yield not only shallowerspectra, but also a larger degree of intermittency. At the sametime, such shocks can be formed easier if the turbulence driv-ing mechanism has compressive modes.We can observe from Figure 1 that the Tsallis distributionsfit remarkably well the PDFs in the simulations. Only smalldepartures are evident at the tails of the PDFs, in particularfor density in models 1 and 2 (supersonic turbulence). Thevalues of the q and w parameters from the fits are plotted as afunction of the separation r in Figures 2, and 3, respectively.In Figures 2 and 3 we see the same general trends foundin solar wind studies (Burlaga & Viñas 2004a,b, 2005a,b;Burlaga et al. 2006a,b; Burlaga & Viñas 2006; Burlaga et al.2007a,b, 2009) and electron MHD simulations (Cho & Lazar-ian 2009): a q parameter that decreases, along with a w valuethat increases as we take larger separations. In other words,the PDFs become more kurtotic, and narrower, as we go tosmall scales, while at large scales the PDFs are more Gaus-sian and wider. The same behavior (i.e. a decreasing kurtosisfor an increasing lag) has been found also for velocity cen-troids in simulations and observations (see Falgarone et al.1994; Lis et al. 1998; Klessen 2000; Federrath et al. 2009). Infact, we have computed the kurtosis and variance of the PDFs(not shown here) and found a similar monotonic decrease andincrease, respectively. Both the kurtosis and q tend to theirGaussian value (3 and 1, respectively) as the separation tendsto the injection scale.It is noticeable that the density PDFs in the supersonicregime have a very distinct shape than those of the velocityand magnetic field components. In the models of supersonicturbulence they are much narrower, and more peaked than insubsonic turbulence. They also show a shallower slope in the q and w vs r plots. We can also notice slightly larger valuesof q (and lower w ) for the velocity components with respectto the magnetic field components, more pronounced for highAlfén Mach numbers. STATISTICS OF OBSERVABLE DATA: COLUMNDENSITYThe measures of non-Gaussianity presented in the previoussection are interesting in their own right, and might provideus with insights to the intermittency of MHD turbulence fromsimulations. However, observations of the ISM cannot pro-vide direct three-dimensional information of the density, ve-locity, and magnetic fields. For instance, from spectroscopicobservations we can access only the velocity distribution ofthe emitting material at a given position (or positions) in the Esquivel & Lazarian F IG . 1.— PDFs of fluctuations for different spatial lags r , for all the models. The symbols represent data from the numerical simulations and the lines are thefits with a q-Gaussian (Tsallis) distribution. The different MHD variables are color coded according to the legend at the upper left corner, and the value of the lagis indicated at the left of each set of curves. plane of the sky.With the simulations presented in the previous section wehave constructed (2D) maps of column density, assuming thatthe emitting material is optically thin, and with an emissiv-ity linearly proportional to the density (e.g. H I ). From thetwo-dimensional column density maps we constructed PDFsof column density increments, and following the same pro-cedure we described before, we made fits to Tsallis distribu-tions. The PDFs and the fits are presented in Figure 4.From the figure we can confirm the same general behaviorof the PDFs, which become more Gaussian as we take largerseparations. We can see that the PDFs are still consistent withTsallis distributions, although the dependence on the lag ismore subtle than with the 3D fields used in the previous sec-tion.In Figure 5 we present the parameters q and w as obtainedfrom the fitting procedure. Both panels of Figure 5 confirmthe decrease of q along with an increase of the width of thePDFs ( w ) with the lag. From the figure it is also evidentthat the fitting parameters can clearly distinguish betweenthe supersonic and subsonic models. However, this can bedone with other measures, for instance the variance of densityor column density (Padoan et al. 1997; Passot & Vázquez-Semadeni 2003; Federrath et al. 2008, 2009), or the skewnessand kurtosis of density and column density PDFs (Kowal et al.2007; Burkhart et al. 2009). A more difficult task is to extractinformation of the Alfvén Mach number. One can see fromFigure 5 that the Tsallis parameters depend strongly on thesonic Mach number, and less so on the Alfvén Mach num-ber. However, an encouraging clear dependence on the AlfvénMach number is evident in the supersonic models, in partic- ular for w (lower panel of Figure 5), where all the M A = 2points lie above the M A = 0 . M s = 0 .
7) it not clear at all, and the behavior withthe Alfén Mach number is the opposite, higher w values forsub-alfvénic models. DISCUSSIONWe found that the Tsallis expression [see Eq. (1)] describeswell the statistics of the PDFs of our 3D-MHD numerical sim-ulations. This is welcome news for both the ISM, where wepropose to use Eq. (1) as a new tool to study turbulence, andfor the solar wind studies, where the Tsallis expression wasused, but its numerical testing was limited.We have to add that the present work is intended only as aproof of concept. While the Tsallis PDFs do provide physicalgrounds to interpret the distributions obtained in MHD simu-lations and/or observations, it is not clear whether the formal-ism of non-extensive statistical dynamics is applicable to theturbulent ISM. This, we feel deserve a separate study. How-ever, the fact that the expression provides both the fit to nu-merical MHD simulations and experimental data, and more-over seems to distinguish the Alfvén Mach number in somecases, is encouraging.In regard to observational quantities, we should note thatPDFs of velocity centroids (as opposed to column density,which is what we ave used in the present paper) have beenemployed to characterize the intermittency both in numeri-cal models (Falgarone et al. 1994; Klessen 2000; Federrathet al. 2009) and molecular cloud observations (for instance ρ -Ophiuchi in Lis et al. 1998; Polaris and Taurus in Hily-Blantet al. 2008). There is always a trade off while choosing dif-sallis statistics 5 F IG . 2.— Values of the q parameter from the fits shown in Figure 1, the symbols (and colors) indicate the MHD variable used, according to the legend at theleft of the plots. ferent quantities to measure turbulence. On the one hand ve-locity is a dynamical quantity that is easy to relate with turbu-lence theories, however (at least in the case of power-spectra)is trickier to obtain, and in fact (for power-spectra) velocitycentroids fail to trace the velocity statistics in strongly super-sonic turbulence (see for instance Esquivel & Lazarian 2005).On the other hand density fluctuations are not only producedfrom turbulence, but are more robustly obtained from obser-vations (e.g. Lazarian & Pogosyan 2008). We must stress thatthe healthy approach is to use these tools, and others as well,as complementary.The most important finding of this paper is the dependenceof the q and w parameters entering Eq. (1) on both sonic andAlfvén Mach number of the turbulence. These numbers areessential for many astrophysical processes from star forma-tion (see McKee & Ostriker 2007), magnetic diffusion en-abled by reconnection (Lazarian 2005), propagation of heat(Narayan & Medvedev 2001; Lazarian 2006a), cosmic rays(Yan & Lazarian 2004) and dust grains acceleration (Yan et al.2004).By now several ways of estimating the Alfvén and sonicMach numbers have been proposed (see Kowal et al. 2007),all with their own limitations and uncertainties. Thus the newapproach here based on Eq. (1) is a welcome addition to theexisting tools.We should add however, that a number of turbulence pa-rameters, which are as important as these Mach numbers, arenot being explored in this work. The scale of turbulence, forinstance, has important consequences for star formation (e.g.Klessen et al. 2000; Mac Low & Klessen 2004) as well as for dynamo theory (e.g. Brandenburg & Subramanian 2005b,a;Hanasz et al. 2009) SUMMARYWe have used three-dimensional MHD simulations to studythe distribution functions of increments of density, velocity,and magnetic field. In addition, we studied synthetic columndensity maps, which can be directly compared with observa-tions. The results can be summarized as follows.1. The PDFs of each of all the quantities is well describedby a Tsallis distribution. Such distributions are sym-metric, and characterized basically by three parameters:an amplitude A , a width w , and the entropic index q thatrelates to the shape (large values yield peaked distri-butions, while with small values they approximate to aGaussian).2. The general trend of the PDFs for different incrementsare consistent with previous studies, of solar wind ob-servations. That is, the PDFs are close to a Gaussiandistribution for large separations, and become morekurtotic and narrower (i.e. a larger degree of intermit-tency) for small separations.3. While the PDFs of (3D) density are quite similar tothose of magnetic field and velocity in subsonic turbu-lence (regardless of Alfvén Mach number), for super-sonic turbulence (also independently of M A ) they arenoticeable narrower and more peaked. That is, the fitsof the PDFs to a Tsallis distribution have smaller w , andlarger q parameters. Esquivel & Lazarian F IG . 3.— Values of the w parameter from the fits shown in Figure 1, the symbols (and colors) indicate the MHD variable used, according to the legend at theleft of the plots.