---

Non-Maxwellian electron distributions in clusters of galaxies

J. S. Kaastra^1,2 - A. M. Bykov³ - N. Werner⁴

1 - SRON Netherlands Institute for Space Research, Sorbonnelaan 2, 3584 CA Utrecht, The Netherlands
2 - Sterrenkundig Instituut, Universiteit Utrecht, PO Box 80000, 3508 TA Utrecht, The Netherlands
3 - A.F. Ioffe Institute of Physics and Technology, St. Petersburg 194021, Russia
4 - Kavli Institute for Particle Astrophysics and Cosmology, Stanford University, 452 Lomita Mall/mc 4085, Stanford, CA 94305, USA

Received 14 May 2009 / Accepted 28 May 2009

Abstract
Context. Thermal X-ray spectra of clusters of galaxies and other sources are commonly calculated assuming Maxwellian electron distributions. There are situations where this approximation is not valid, for instance near interfaces of hot and cold gas and near shocks.
Aims. The presence of non-thermal electrons affects the X-ray spectrum. To study the role of these electrons in clusters and other environments, an efficient algorithm is needed to calculate the X-ray spectra.
Methods. We approximate an arbitrary electron distribution by the sum of Maxwellian components. The decomposition uses either a genetic algorithm or an analytical approximation. The X-ray spectrum is then evaluated with a linear combination of those Maxwellian components.
Results. Our method is fast and leads to an accurate evaluation of the spectrum. The Maxwellian components allow the use of the standard collisional rates that are available in plasma codes such as SPEX. We give an example of a spectrum for the supra-thermal electron distribution behind a shock in a cluster of galaxies. The relative intensities of the satellite lines in such a spectrum are sensitive to supra-thermal electrons. These lines can only be investigated with high spectral resolution. We show that the instruments on future missions like Astro-H and IXO will be able to demonstrate the presence or absence of these supra-thermal electrons.

Key words: acceleration of particles - radiation mechanisms: general - X-rays: galaxies: clusters - X-rays: general

1 Introduction

Most of the visible baryonic matter in clusters of galaxies is in the form of a hot, tenuous gas. Usually it is assumed that this gas is in collisional ionisation equilibrium. In the central parts, the density is high enough and the timescales long enough to fulfil the conditions for collisional ionisation equilibrium. In the outer parts this may not always be the case. Freshly accreted gas may still be ionising. For a review of equilibration processes in such tenuous plasmas, see Bykov et al. (2008). Nevertheless, usually a Maxwellian electron distribution is still assumed to be valid; however, there are situations where this is not an obvious constraint. For instance, in cooling cores of clusters, electrons from hot gas may penetrate the colder parts. Also, when shocks are present, deviations from a Maxwellian distribution may occur associated with the temperature gradients or caused by particles accelerated by the shock. Shocks can stem from merger activity, AGN jets, and accretion shocks may occur in cluster outskirts. As the Coulomb thermal relaxation times increases with electron energy E as E^3/2, supra-thermal electrons are relatively long-lived and may yield a pressure contribution that is potentially interesting for mass profiles. The presence of such supra-thermal electrons can be revealed by excess emission in satellite lines.

In this paper we describe a way to model the emerging X-ray emission spectra for the case of a non-Maxwellian electron distribution. Although our main focus here is on clusters of galaxies, the same procedure can be applied to different circumstances, for instance in stellar coronae, supernova remnants, and the Galactic ridge; for the last case see also Masai et al. (2002) and Tatischeff (2003).

2 X-ray spectra for supra-thermal electron distributions

To calculate an X-ray spectrum, essentially two steps are needed: 1) determine the ionisation balance and; 2) calculate the corresponding spectrum. For more details see e.g. Kaastra et al. (2008).

We first consider the ionisation balance. To determine that, it is necessary to evaluate the collisional ionisation and recombination rates for all ions. These rates are derived by integrating the relevant cross sections over the electron distribution. For instance, the collisional ionisation rate is given by

where n_e and n_i are the densities of electrons and the ion i considered, is the energy-dependent collisional ionisation cross-section, is the kinetic energy of an electron colliding with the ion and is the electron velocity distribution.

In the most simple approximation, the cross section for collisional ionisation can be written as (Lotz 1967):

where E is the kinetic energy of the free electron, n_s is the number of electrons in a given atomic shell and the normalisation  m² keV². The cross section (2) can be inserted into (1), and the integration can be done analytically for simple electron distribution functions like Maxwellians or power laws, resulting in an expression proportional to E₁(I/kT)/IT^0.5, where I is the ionisation potential of the relevant shell, and is the exponential integral. Equation (2) is obviously too simple in realistic cases. But more sophisticated approximations to the collisional ionisation cross section can be made, e.g. Younger (1981), that still allow analytical integration for Maxwellian or power law electron distributions:

with and A, B, C and D adjustable parameters. Analytical integration has the obvious advantage of enhanced computational speed, as during the evaluation of a spectrum many atomic shells of hundreds of ions need to be taken into account.

A similar treatment can be made for other relevant rates, like the collisional excitation rates that are important for the line emission. However, some rates cannot be written in a form that allows analytical integration over a Maxwellian electron distribution. A good example are the radiative recombination rates. We recall that the cross section for radiative recombination can be expressed in terms of the photoionisation cross section through Milne's relation as

where is the photoionisation cross section at the photon energy . Although there exist analytical approximations to the photoionisation cross sections (for instance Verner & Yakovlev 1995), they are too complicated to allow for analytical integration of (4) over a Maxwellian electron distribution, in particular as they are more complicated than a simple power law near the ionisation edges (Cooper minima). Even for the hydrogenic case, where the cross section can be calculated analytically, its shape is too complex for analytical convolution with the electron distribution. Therefore, the integration must be done numerically. In practice, often the exact cross sections (including also resonances near the edge) are convolved numerically with a Maxwell distribution, and the resulting rates are then approximated by analytical functions.

It is clear that we need to follow an alternative approach for non-Maxwellian electron distributions. There have been attempts to make a generalisation of the Maxwell distribution. For instance, Porquet et al. (2001) combine a Maxwellian at low energies with a power-law at high energies, but they only consider the ionisation balance, not the resulting spectrum. Prokhorov et al. (2009) adopt a similar approximation for the electron distribution, but only consider continuum emission and the Fe XXV and Fe XXVI 1s-2p blends. Owocki & Scudder (1983) use the so-called kappa-distribution (Olbert et al. 1967):

This distribution tends to a Maxwellian for and it has a power law tail proportional to for high energies. However, the above approach has also clear disadvantages. Not every nonthermal electron distribution can be cast into the shape of (5), and not all rates can be convolved analytically with (5), in particular again the recombination rates. Owocki & Scudder (1983) had to assume very simple analytical approximations to the cross sections, like (2), in order to be able to evaluate analytically the plasma rates, which are given in terms of hypergeometric functions. Therefore we propose a different approach.

3 The multi-Maxwellian approach

Plasma codes generally may contain (ten)thousands of spectral lines, and at high spectral resolution, spectra contain many bins. Therefore, numerical integration for each rate over the electron distribution is a very time consuming task, and not of practical use, in particular when spectral fitting is performed and the model needs to be evaluated many times.

But all relevant rates in the case of the collisionally ionised plasmas that we consider here, are proportional to the product n_en_i, as for each of these rates the interaction (collision) of an electron with an ion is the driving process. Therefore, all rates are linearly proportional to the electron density n_e. If we decompose the electron distribution into a linear combination of elementary components, the total rate for any process is simply the sum of the rates for these individual elementary components. Now, for Maxwellian electron distributions all rates are well known and fast to evaluate analytically or with analytical approximations. Thus, if we decompose an arbitrary electron distribution f(E) into a linear combination of Maxwellians:

with A_i the normalisation constant for component i and g_M(E,T_i) the Maxwell distribution for temperature T_i, then the calculation of any rate (recombination, ionisation, excitation) is simple and straightforward. The problem is then reduced to obtaining the normalisations A_i and temperatures T_i for an arbitrary electron distribution f(E).

In principle, one might solve formally the integral equation

Using

(7) can be cast into the shape of a Laplace transform as

with

The calculation of the inverse Laplace transform is not always trivial. For instance, when we consider a mono-energetic beam, f(E) is essentially a delta-function, and the formal solution is a sine-wave with infinite frequency, which is not very practical. Fortunately, in most practical cases f(E) is a smooth function of energy, spanning a broad range of energies.

There are various ways to solve (9). There is software available that can do the inverse Laplace transform numerically (for example d'Amore et al. 1999), but in practice this is difficult, because these algorithms only work if the left-hand-side of (9) is an analytical function. However, in practice the electron distributions are given in tabular form, and the algorithms used to interpolate such tables do not produce a formal analytical function (the requirement that the function is infinitely differentiable is usually violated). There are two solutions to this problem that we elaborate in the next subsections: direct approximation of the electron distribution by an analytical function, or direct decomposition of the electron distribution as a sum of Maxwellians.

3.1 Approximation of the electron distribution by an analytical function

In practice, apart from the most pathological cases, electron distributions contain a core that is close to Maxwellian, plus a high energy tail. For an example see Sect. 4. After some experimentation, we found that the electron distribution of that example can be approximated by the following series:

with x the dimensionless momentum defined by

and where T₀ is a characteristic temperature corresponding to the electron distribution, and the parameters c_k and a_k can be adjusted to give the best fit.

Equation (11) is an analytical function, and the inverse Laplace transform in this case can be done analytically, yielding

where is the -function and Dirac's delta-function.

For our example, the normalisations c₀-c₄ are 0.755, 0.0609, 2.54, 13.3, and 17.58, respectively; the scale factors a₀-a₄ are 0.483, 152, 6843, 57.4 and 12.4. The fit is not perfect, but better than 0.5% for x<30, and better than 1% for x<700. The relative difference to the exact distribution are shown in Fig. 2.

3.2 Direct decomposition of the electron distributions as the sum of Maxwellians

We seek an approximation to f(E) that can be written as a linear combination of n Maxwellians, with n a given number. The solution can be represented as a set of pairs (T_i,A_i) with T_i the temperatures and A_i the normalisations. We will choose for convenience a set of increasing temperatures ( T_i+1>T_i), and we will demand that (physically allowed components). Most realistic electron distributions can be represented as a Maxwellian with high-energy tails. Therefore, the first temperature T₁ should be close to the temperature T₀ of this dominant Maxwellian component.

This then leads to the following approach. We define a set of nlogarithmically spaced temperature intervals (T_1,i,T_2,i) with T_1,1 equal to T₀ and for all i, we take T_2,i=sT_1,i and T_1,i=T_2,i-1. For the parameter s>1 we adopt a value of 1.5. The total temperature range spanned by these n contiguous intervals is therefore T₀ - sⁿ T₀. For each interval, we choose an arbitrary temperature T_i within this interval ( ) and an arbitrary normalisation . We scale the A_i in such a way that their sum corresponds to the total electron density. We then check how close the resulting electron distribution is to the electron distribution f(E) that we want to approximate. By varying all values T_i within their allowed ranges (the above mentioned intervals) and also the values for A_i, we try to get the best possible approximation.

We use the genetic algorithm PIKAIA developed by Charbonneau (1995) to find the best solution. This algorithm needs a formal function of the parameters (T_i, A_i) that needs to be maximised. For this function we choose the inverse of the sum of the squared relative deviations of the approximation to the true electron distribution f(E).

After some experimenting, we obtained satisfactory convergence with the following initial parameters: number of individuals 500, number of generations 500. All other parameters in PIKAIA were kept to their default values. The number of individuals represents here the number of sets (T_i,A_i) from which the iteration starts; the number of generations is the number of iterations during which the individual sets evolve towards the optimum solution.

We have tested our algorithm on various different electron distributions, but give here only one illustrative example.

4 Example of decomposition of an electron distribution into Maxwellians

We have used electron distributions based on the models of Bykov & Uvarov (1999) for collisionless MHD shocks taking into account particle acceleration and a shock precursor. We consider here a case for a pre-shock temperature of T₀=10⁸ K, electron density 10³ m^-3, magnetic field 10^-10 T, and shock size 100 kpc. In the simulation we assumed a pure Maxwellian distribution of electrons in the far upstream region of the flow and thus only the direct injection of electrons from the upstream thermal pool to Fermi type shock acceleration is responsible for the non-thermal electron population. The electron distribution in the immediate post-shock region is shown in Fig. 1 for three different Mach numbers. The model distribution was calculated taking account of both Fermi-type acceleration in a collisionless shock and Coulomb losses of the electrons in both the upstream and downstream regions. It is normalised to the momentum p₀ corresponding to the typical energy kT₀. Note that the supra-thermal electron distribution well away from the shock front differs from that shown in Fig. 1, because the efficient Coulomb losses wash-out the non-relativistic supra-thermal electrons.

There is a variant of the model where mildly relativistic electrons (with Lorentz factors 30) comprising a putative long-lived cosmic ray electron population in the ICM are re-accelerated by the MHD-shocks. The energy efficiency problem is alleviated in that model in comparison with the case of only direct particle injection from the thermal plasma. However, in this paper we concentrate mostly on the diagnostic of non-relativistic electrons and thus we consider only a direct injection model.

For our example, we take the case of M=2.2. In our decomposition with the genetic algorithm we have taken n=32 Maxwellians. The relative differences of the approximation compared to the exact electron distribution is shown in Fig. 2. The temperatures and relative normalisations of the solution are shown in Fig. 3. We have also used the analytical approximation (13) and binned it into similar temperature bins as the solutions from the genetic algorithm.

5 Spectrum for the supra-thermal electron distribution

We have adjusted all models in the spectral fitting package SPEX that involve emission or absorption from a hot plasma. All these models now include an option to account for supra-thermal electrons. They have an additional parameter, which is the name of a file containing the temperatures and relative emission measures of the Maxwellian components. As SPEX does not use pre-calculated tables but calculates spectra on the fly, all relevant rates (ionisation, recombination, excitation) are simply calculated by adding the contributions from the Maxwellian components. Obviously, this process is done in two steps: first the composite multi-Maxwellian electron distribution is applied to determine the ionisation balance, and with the resulting non-equilibrium ion concentrations, we calculate the X-ray spectrum for the non-equilibrium electron distribution. It is tacitly assumed here that we consider a time-independent, steady situation. For the example given in the previous section, only a few dozen Maxwellian components are needed. This allows fast and accurate evaluation of the spectrum, without the need to make simplifications to the atomic physics.

The most important effects of non-thermal electrons on the spectrum is then a shift of the ionisation balance towards higher ionisation, the production of a non-thermal Bremsstrahlung tail on the continuum spectrum, and enhanced satellite line emission. These satellite lines (for instance in the Fe-K band) can be detected easily with high-resolution spectrometers that will fly on future missions such as Astro-H and IXO. Their relevance as indicators for non-thermal electrons was already indicated by Gabriel & Phillips (1979).

To illustrate the effects of such a supra-thermal electron distribution on data, we simulated an XMM-Newton EPIC/pn spectrum extracted from a circular region with a radius of 1 centred on the core of a bright cluster with a 0.3-10 keV luminosity of  W within the extraction region, at an assumed redshift of z=0.055. In the simulation of the spectrum, we assumed the above mentioned post-shock downstream electron distribution for the Mach number of M=2.2 and pre-shock temperature of kT=8.62 keV (10⁸ K). We assume a deep 100 ks observation. The resulting spectrum has very high statistics, which should in principle allow to detect any non-isothermality of the plasma. The best fit temperature of the simulated spectrum is  keV (  K), consistent with the post-shock temperature of the plasma given by the Rankine-Hugoniot jump condition. An isothermal model fits the data extremely well (reduced ) and the non-thermal tail of the electron distribution cannot be detected in the spectrum.

We also simulated a spectrum with the same input parameters as observed during a deep 100 ks observation with the X-ray micro-calorimeter on the proposed International X-ray Observatory (IXO). We fitted the simulated spectrum with an isothermal model and obtained a temperature of  keV, about 1 keV lower than the expected post-shock temperature. In Fig. 4 we show the 6.9-7.0 keV part of the spectrum (rest-frame energies) which shows the Fe XXVI Ly lines and the Fe XXV j-satellite line. Enhanced equivalent widths of satellite lines are good indicators of non-thermal electrons. The satellite line in the simulated spectrum in Fig. 4 is clearly stronger than that predicted by the thermal model with a Maxwellian electron distribution. This exercise illustrates, that in order to observationally reveal non-Maxwellian tails in the electron distributions, we will need high-resolution spectra obtained by future satellites with a large effective area.

However, even before IXO, Astro-H (expected launch 2014) will be able to detect supra-thermal electrons. We simulated the same spectrum for the same extraction region as above for the main instruments of Astro-H, again for 100 ks exposure time. The two hard X-ray telescopes detect the source up to 75 keV, and the spectrum is well approximated ( for 223 degrees of freedom) by an isothermal model with measured temperature of  keV. Thus, the presence of such a small amount of supra-thermal electrons cannot be revealed as a hard tail, but it can be revealed in high-resolution spectra. Figure 4 shows the simulated spectrum for the Soft X-ray Spectrometer (SXS) of Astro-H. The excess flux at the Fe XXV j-satellite has a significance. Obviously, longer exposure times will enhance the significance.

Interestingly, in all our simulations above, the best-fit iron abundance for the isothermal model is about 30% higher than the actual abundance that we have put into our model spectrum with supra-thermal electrons. This holds also if we restrict our fit only to the Fe L-shell or Fe K-shell band. This is because the higher-energy electrons have less efficient line emissivity relative to the continuum, compared to lower-energy electrons. Without knowing the amount of supra-thermal electrons, which can be determined only from high-resolution spectra, this abundance bias cannot be resolved and will result into incorrect interpretations.

6 Concluding remarks

In practice most effort goes into finding a good decomposition of an electron distribution into Maxwellians. We have indicated and illustrated in this paper two different methods: fits to simple analytical models that allow analytical inversion of the Laplace transform, and a direct decomposition with a genetic algorithm.

Finally we note that all relevant plasma rates that are used in the SPEX code are calculated with non-relativistic approximations. For instance, ionisation cross sections are approximated by analytical functions like (3) that loose their validity for relativistic energies. Thus, for electron distributions containing a significant fraction of relativistic electrons, the results will be less accurate. Fortunately, in most situations this is not a problem. For instance, the electron distribution of Fig. 1 has a high-energy tail roughly proportional to p^-3. Inserting this for instance into (1) and using (2) for the high-energy limit, shows that the integrand scales, apart from a logarithmic term, proportional to E^-3, and therefore the highest energy electrons do not contribute much to the rates. Therefore, even though the approximations made to the decomposition of the electron distribution, in particular the genetic algorithm, are not always very accurate at high energies (see Fig. 2), this affects the final spectrum to a much lesser extent.

Only in astrophysical situations with a significant number of relativistic electrons our method will not apply. Bykov (2002) has given an example of this in the context of supernova remnants, considering only line fluorescence due to collisional ionisation. Good approximations to the relativistic collisional ionisation cross section of K-shell and L-shell electrons are available (see references in Bykov 2002), but for a full plasma model relativistic corrections to all rates would be needed, which is beyond the scope of the present paper.

Acknowledgements

The Netherlands Institute for Space Research is supported financially by NWO, the Netherlands Organisation for Scientific Research. A.M. Bykov was supported by RBRF grant 09-02-12080. Support for this work was provided by the National Aeronautics and Space Administration through Chandra Postdoctoral Fellowship Award Number PF8-90056 issued by the Chandra X-ray Observatory Center, which is operated by the Smithsonian Astrophysical Observatory for and on behalf of the National Aeronautics and Space Administration under contract NAS8-03060.

References

• Bykov, A. M. 2002, A&A, 390, 327 [NASA ADS] [CrossRef] [EDP Sciences] (In the text)
• Bykov, A. M., & Uvarov, Y. A. 1999, Sov. J. Exper. Theor. Phys., 88, 465 [NASA ADS] [CrossRef] (In the text)
• Bykov, A. M., Paerels, F. B. S., & Petrosian, V. 2008, Space Sci. Rev., 134, 141 [NASA ADS] [CrossRef] (In the text)
• Charbonneau, P. 1995, ApJS, 101, 309 [NASA ADS] [CrossRef] (In the text)
• d'Amore, L., Laccetti, G., & Murli, A. 1999, ACM Trans. Math. Softw., 25, 306 [CrossRef] (In the text)
• Gabriel, A. H., & Phillips, K. J. H. 1979, MNRAS, 189, 319 [NASA ADS] (In the text)
• Kaastra, J. S., Paerels, F. B. S., Durret, F., Schindler, S., & Richter, P. 2008, Space Sci. Rev., 134, 155 [NASA ADS] [CrossRef] (In the text)
• Lotz, W. 1967, Z. Phys., 206, 205 [NASA ADS] [CrossRef] (In the text)
• Masai, K., Dogiel, V. A., Inoue, H., Schönfelder, V., & Strong, A. W. 2002, ApJ, 581, 1071 [NASA ADS] [CrossRef] (In the text)
• Olbert, S., Egidi, A., Moreno, G., & Pai, L. G. 1967, Trans. AGU, 48, 177 (In the text)
• Owocki, S. P., & Scudder, J. D. 1983, ApJ, 270, 758 [NASA ADS] [CrossRef] (In the text)
• Porquet, D., Arnaud, M., & Decourchelle, A. 2001, A&A, 373, 1110 [NASA ADS] [CrossRef] [EDP Sciences] (In the text)
• Prokhorov, D. A., Durret, F., Dogiel, V., & Colafrancesco, S. 2009, A&A, 496, 25 [NASA ADS] [CrossRef] [EDP Sciences] (In the text)
• Tatischeff, V. 2003, in EAS Publications Series, ed. C. Motch & J.-M. Hameury, EAS SP, 7, 79 (In the text)
• Verner, D. A., & Yakovlev, D. G. 1995, A&A, 109, 125 [NASA ADS] (In the text)
• Younger, S. M. 1981, J. Quant. Spectr. Rad. Trans., 26, 329 [NASA ADS] [CrossRef] (In the text)

All Figures

---

Copyright ESO 2009