A&A 417, 381-389 (2004)
DOI: 10.1051/0004-6361:20031765

Evolution of X-ray spectra in down-Comptonization. A comparison of the extended Kompaneets equation with Monte Carlo simulation and the Ross-McCray equation

D.-B. Liu1 - L. Chen1 - J. J. Ling1 - J. H. You1 - X. M. Hua2


1 - Institute for Space and Astrophysics, Department of Physics, Shanghai Jiao Tong University, Shanghai 200030, PR China
2 - Laboratory for High Energy Astrophysics, NASA/GSFC Code 661, Greenbelt, Maryland 20771, USA

Received 13 May 2003 / Accepted 18 November 2003

Abstract
The original Kompaneets equation fails to describe down-Comptonization, which is the most important radiative transfer process in hard X-ray and $\gamma$-ray astronomy, compared with up-Comptonization. In this paper, we improve our previous derivation of the extended Kompaneets equation and present it more clearly. The new equation can be used to describe a more general Comptonization process, including up- and down-Comptonization, suitable for any case, $h\overline{\nu}\ll kT_{\rm e}$, $h\overline{\nu}\gg kT_{\rm e}$ and $h\overline{\nu}\sim kT_{\rm e}$. The condition of the original Kompaneets equation $h\overline{\nu}\ll kT_{\rm e}$ is no longer necessary. Using the extended equation, we give some typical solutions in X-ray astronomy, and compare them with those of the prevailing Monte Carlo simulations and the Ross-McCray equation. The excellent consistency between the extended Kompaneets equation and Monte Carlo simulation or Ross-McCray equation confirms the correctness of our extended Kompaneets equation. The numerical solution of the extended Kompaneets equation is less expensive in terms of computational time than the Monte Carlo simulation. Another advantage of the equation method is the simplicity and the clarity in physics. The potential applications in X-ray and $\gamma$-ray astronomy are also emphasized.

Key words: scattering - radiative transfer - methods: analytical - X-rays, $\gamma$-rays: ISM

1 Introduction

The study of radiative transfer for photons passing through a plasma has been an important topic in both astrophysics and radiation physics. The transfer process significantly affects the properties of emergent radiation from plasmas, e.g. the spectral shape of the continuum, the profile and centroid energy of an emission or absorption line, the intensity ratio of different lines and the polarization, etc. (Fabian et al. 1995; Misra & Kembhavi 1998; Titarchuk 1994). The detailed consideration of the problem of energy transfer between the radiation and the medium is as important as that of the basic radiation mechanism, particularly in the X-ray band. Owing to the rapid development of X-ray astronomy in recent years, it seems necessary to reexamine and improve the available methods and theories for the study of transfer properties in radiation propagation (Kompaneets et al. 1957; Illarionov et al. 1979; Fabian et al. 1995; Sunyaev & Titarchuk 1985; Titarchuk 1994; Hua et al. 1997; Kazanas et al. 1997; Misra & Kembhavi 1998; Zhang et al. 2000). The specialty of the radiative transfer in X-ray band is that, for an almost fully ionized plasma with high temperature, which often occurs in high energy astrophysics, the dominant mechanism of energy exchange between radiation and plasma is non-elastic photon-electron scattering. This special transfer mechanism in X-ray astronomy is known as Comptonization. Although the change of the photon wavelength (or energy) in each individual collision is very small, the integrated variation of photon energy in multiple Compton scatterings in actual astronomical circumstances may be significant, particularly for the hard X-ray and $\gamma$-ray photons. The fractional change of the photon wavelength or frequency in each collision is given by the following well known formula (if the electron is approximately motionless before collision)
 
$\displaystyle \frac{\Delta \lambda}{\lambda_{0}}=\frac{\lambda^{\prime}-\lambda... ...ambda_{0}}\sin^{2}{\frac{\theta}{2}}\approx \frac{\lambda_{\rm C}}{\lambda_{0}}$     (1)

where $\lambda_{\rm C}\equiv h /m_{\rm e}c=0.024~ {\rm\AA}$ is the Compton wavelength and $\theta$ is the scattering angle. It is seen from Eq. (1) that the fractional change $\frac{\Delta \lambda}{\lambda_{0}}$ depends on the initial wavelenth $\lambda_{0}$. The shorter the initial wavelength $\lambda_{0}$ is, the larger the change $\frac{\Delta \lambda}{\lambda_{0}}$will be. For example, $\frac{\Delta \lambda}{\lambda_{0}}\approx 10^{-6}$ for $\lambda_{0}=5000~ {\rm\AA}$ in the optical band, but $\frac{\Delta \lambda}{\lambda_{0}}\approx 10^{-2}$ if $\lambda_{0}\simeq 0.5~ {\rm\AA}$ $(h\nu \sim 20~ \rm {keV}$, in the X-ray band). For a very hard $\gamma$-ray photon with energy $\sim$100 keV, $\frac{\Delta \lambda} {\lambda_{0}}\approx 10^{-1}$. This is why Comptonization is particularly important in hard X-ray and $\gamma$-ray astronomy.

There are two kinds of Comptonization processes in astrophysics. If the average energy of photons in radiation field $h\overline{\nu}$ is markedly larger than the average thermal energy of electrons $kT_{\rm e}$, $h\overline{\nu}\gg kT_{\rm e}$, which often occurs in hard X-ray and $\gamma$-ray astronomy, then the average energy of scattered photons decreases. This is known as down-Comptonization or Comptonization-softening. When $h\overline{\nu}\ll kT_{\rm e}$, which often occurs in radio and infra-red astronomy, the average energy of scattered photons increases, which we call up-Comptonization, or Comptonization-hardening. For the hardening process, if the condition $h\overline{\nu} \ll kT_{\rm e}\ll m_{\rm e}c^{2}$ is satisfied, the equation of radiative transfer, which is in the form of a diffusion equation, has been given by Kompaneets as follows,

 
$\displaystyle \frac{\partial n}{\partial t}~=~\frac{kT_{\rm e}}{m_{\rm e}c^{2}}... ... x}\left \{ x^{4}\left [ \frac{\partial n}{\partial x}+n(n+1)\right ] \right \}$     (2)

where $x\equiv h\nu /kT_{\rm e}$ is the dimensionless photon energy in units of $kT_{\rm e}$, $\sigma_{T}$is the Thomson cross-section, $N_{\rm e}$ is the number density of the scattering electron gas, $n(x,t)\equiv n(\nu,t)$ is the frequency distribution function of the photon gas, which represents the "photon number'' in each photon state per unit volume at frequency $\nu$. Therefore the real number density of photons in the interval $\nu$- $\nu+{\rm d}\nu$ will be $n(\nu,t)\left (8\pi \nu^{2}/c^{3}\right ){\rm d}\nu$. For a tenuous radiation field, the "photon number'' $n(\nu, t)\ll 1$.

However, we emphasize that the Kompaneets Eq. (2) is restricted to a limited application area, only useful for the calculation of the up-Comptonization process owing to the severe condition $h\overline{\nu} \ll kT_{\rm e}\ll m_{\rm e}c^{2}$. For down-Comptonization, which could be more important in the hard X-ray or $\gamma$-ray astronomy where the condition $h\overline{\nu}\gg kT_{\rm e}$ is often satisfied, the adoption of Eq. (2) will inevitably bring a large error in the resultant calculation. Furthermore, even when the condition $h\overline{\nu}\ll kT_{\rm e}$ is satisfied in the initial stage of scattering process, the error in the calculation of the up-Comptonization process may also be significant if the scattering depth is large, $\tau_{\rm s}\gg 1$. In this case, the average energy $h\overline{\nu}$ of the scattered photons could be markedly increased to the level $h\overline{\nu}\sim kT_{\rm e}$, which destroys the necessary condition of Eq. (2), i.e. $h\overline{\nu}\ll kT_{\rm e}$.

Therefore, the original Kompaneets Eq. (2) should be extended and improved in order to describe a more general non-elastic scattering of photon-electron, particularly for the down-Comptonization process. Ross & McCray obtained another diffusion equation for the case $kT_{\rm e}\ll h\overline{\nu}\ll m_{\rm e}c^{2}$ to describe the down-Comptonization process.

 
$\displaystyle \left(\frac{\partial n}{\partial t}\right)_{\rm c}=\frac{kT_{\rm ... ...e}}{m_{\rm e} c^2}x^2\right) \frac{\partial n} {\partial x}\right]\right\}\cdot$     (3)

Based on the same diffusion approximation as given by Kompaneets (1957), Chen & You presented an extended Kompaneets equation under the much looser conditions $h\overline{\nu}\ll m_{\rm e}c^{2}$ and $kT_{\rm e} \ll m_{\rm e}c^{2}$ (Chen & You 1992, 1993; Chen et al. 1994),which is suitable to describe both the up- and down-Comptonizations ( $h\overline{\nu}\ll kT_{\rm e}$ or $h\overline{\nu}\gg kT_{\rm e}$, even $h\overline{\nu}\sim kT_{\rm e}$). In Sect. 2 of this paper, we present our improved primary derivation of this equation in a clearer way. In Sect. 3, we use the extended Kompaneets equation to calculate the Compton-evolution of some typical X-ray spectra, where we particularly emphasize the properties of the emergent X-ray spectrum in down-Comptonization. We compared the numerical solutions of our equation with the results of the prevailing Monte Carlo simulation, which confirmed the correctness and effectiveness of our extended Kompaneets equation. The calculated results using the Ross-McCray Eq. (3) (hereafter, RMC equation) are also presented in this paper for comparison with ours. Sect. 4 gives the conclusions and discussions.

2 The extended Kompaneets equation

Following Kompaneets (1956, 1957), we consider a "mixed gas'' consisting of an electron-gas (the fully ionized plasma) and a photon-gas (the radiation field) under a condition $kT_{\rm e} \ll m_{\rm e}c^{2}$ and $h\overline{\nu}\ll m_{\rm e}c^{2}$(Chen et al. 1994). We make two assumptions: (i) the common thermal equilibrium between the two gases has not yet been reached, but the electron-gas itself is already in thermal equilibrium as the interaction between electrons is the Coulomb long-range force. Therefore the Maxwellian distribution $f(p) =f_0\exp\left[-\frac{1}{kT_{\rm e}}\frac{p^2}{2m_{\rm e}}\right]$ can be used to describe the electron gas; (ii) owing to the frequent Compton scattering, the frequency distribution function of the photon-gas $n(\nu,t)$ is assumed to be isotropic, independent of the direction of the photon-wavevector $\vec{k}$, and $n(\nu,t)$will change with time t until the photon-gas reaches a thermal equilibrium with the electron-gas. In the following, the change of distribution function $n(\nu,t)$is treated as a diffusion of photon-gas in the "frequency space'' along the frequency axis $x\equiv h\nu /kT_{\rm e}$. The "diffusion rate'' of the "photon number'' $\partial n/\partial t$ in the "frequency space'' is determined by the diffusion Eq. (5) given as follows.

Consider an individual collision between an electron with momentum $\vec{p}$and a photon with frequency $\nu$. The energy and momentum conservations in the non-relativistic limit are

\begin{displaymath}\left(\frac{h\nu}{c}\right)\vec{n}+\vec{p} = \left(\frac{h\nu^{\prime}}{c}\right) \vec{n}^{\prime}+\vec{p}^{\prime} \end{displaymath}


 \begin{displaymath} h\nu+\frac{p^2}{2m_{\rm e}} = h\nu^{\prime}+\frac{p^{\prime 2}}{2m_{\rm e}} \end{displaymath} (4)

where $\vec{p}^{\prime}$ and $\nu^{\prime}$ represent respectively the momentum of the electron and the frequency of the photon after collision, $\vec{n}$ and $\vec{n}^{\prime}$ are the directions of the photon before and after collision, respectively. Such a collision process $(\vec{p},\nu,\vec{n})\longrightarrow(\vec{p}^{\prime}, \nu ^{\prime},\vec{n}^{\prime})$ leads to a decrease of the photon number $n(\nu,t)$.

In real astrophysical circumstances, the plasma is always tenuous enough. Therefore, as a good approximation, the electron-gas can be regarded as a classical system rather than a Fermi-gas, and the number density of electrons $N_{\rm e} f(\vec{p})~ {\rm d}^3\vec{p}$ in the interval $\vec{p}$- $\vec{p}+{\rm d}\vec{p}$can be discribed by the classical formula, i.e. $f(\vec{p})$ is the Boltzman distribution function. Denoting the transition probability of the elementary collision $(\vec{p},\nu,\vec{n})\longrightarrow(\vec{p}^{\prime}, \nu ^{\prime},\vec{n}^{\prime})$ by ${\rm d} W$, the total transition number in unit volume is $(1+n^{\prime})n N_{\rm e}f(\vec{p})~{\rm d}^3\vec{p}~{\rm d} W$, where $n\equiv n(\nu, t)$, $n^{\prime}\equiv n^{\prime}(\nu^{\prime}, t)$ are the photon numbers before and after collision, respectively. The occurrence of the factor $(1+n^{\prime})$ is due to the fact that the photon-gas is a boson system for which the transition number depends on both n and $n^{\prime}$. Similarly, the inverse-collision process $(\vec{p}^{\prime},\nu ^{\prime}, \vec{n}^{\prime})\longrightarrow(\vec{p},\nu,\vec{n})$ leads to an increase of $n(\nu,t)$. The corresponding transition number is $(1+n)n^{\prime}N_{\rm e} f(\vec{p}^{\prime})~{\rm d}^3\vec{p}^{\prime}~{\rm d} W$, where the transition probability ${\rm d} W$ is the same as that in the collision process $(\vec{p},\nu,\vec{n})\longrightarrow(\vec{p}^{\prime}, \nu ^{\prime},\vec{n}^{\prime})$, because the Klein-Nishina cross-section ${\rm d} \sigma_{\rm T}/{\rm d}\theta$ has the same value both for the scattering angle $\theta$ and $\pi -\theta$ (see Eq. (11) below). Therefore $\partial n/\partial t$ can be written as:

 
$\displaystyle \left(\frac{\partial n}{\partial t}\right)_{\rm c} =-N_{\rm e}\in... ...n(1+n^{\prime})f(\vec{p})-n^{\prime}(1+n) f(\vec{p}^{\prime})\right ]~{\rm d}W.$     (5)

The change of frequency is very small in each collision if $kT_{\rm e} \ll m_{\rm e}c^{2}$, $\vert\Delta\vert\equiv \vert\nu^{\prime}-\nu\vert\ll \nu$. Expanding $n^{\prime}=n(\nu^{\prime}, t)$and $f(\vec{p}^{\prime})$ in terms of the small quantity $\Delta$to second order under the condition $h\vert\Delta\vert\ll kT_{\rm e}$, and replacing the frequency $\nu$ by a convenient dimensionless frequency x, $x\equiv h\nu /kT_{\rm e}$, Eq. (5) becomes:
 
                          $\displaystyle \left(\frac{\partial n}{\partial t}\right)_{\rm c}$ = $\displaystyle \left[\frac{\partial n}{\partial x} +n\left(n+1\right)\right] \frac{N_{\rm e}h}{kT_{\rm e}}\int {\rm d}^3\vec{p}\int {\rm d} W f(\vec{p})\Delta$  
    $\displaystyle +\left[\frac{\partial^2n}{\partial x^2}+2(n+1)\frac{\partial n}{\partial x} +n(n+1)\right]$  
    $\displaystyle \times \frac{N_{\rm e}}{2}\left(\frac{h}{kT_{\rm e}}\right)^2\int{\rm d}^3 \vec{p}\int{\rm d} W f(\vec{p})\Delta^2.$ (6)

In the following derivations, we first calculate the second integral on the right-hand side of Eq. (6). Let:
 
$\displaystyle I \equiv h^2\int {\rm d}^3\vec{p}\int {\rm d} W f(\vec{p})\Delta^2.$     (7)

The expression $\Delta\equiv \nu^{\prime}-\nu$ in the integral can be obtained from Eq. (4). Retaining only the first order of the small quantity $\Delta$, we obtain:
 
$\displaystyle h\Delta=-\frac{hc\nu\vec{p}\cdot(\vec{n}-\vec{n}^{\prime})+(h\nu)... ...\cdot\vec{n}^{\prime}) -(\vec{p}\cdot\vec{n}^{\prime}/m_{\rm e} c)\right]}\cdot$     (8)

Equation (8) can be further simplified because in most cases in the hard X-ray and $\gamma$-ray astronomy, we have $h\nu\ll m_{\rm e}c^2$ and $kT_{\rm e} \ll m_{\rm e}c^{2}$. Therefore the value in the square bracket in the denominator of Eq. (8) approximately equals 1. Thus we have:
 
$\displaystyle h\Delta \simeq -\frac{h\nu c}{m_{\rm e}c^2}\vec{p}\cdot(\vec{n}-\vec{n}^{\prime}) -\frac{(h\nu)^2}{m_{\rm e} c^2}(1-\vec{n}\cdot\vec{n}^{\prime}).$     (9)

Different from Kompaneets (1957), now we retain the second term in Eq. (9). The reason is that the order of the momentum of thermal electrons is $p\sim m_{\rm e}v\sim m_{\rm e}\sqrt{kT_{\rm e}/m_{\rm e}}$, therefore the ratio of the orders of the second to the first term is $\sim$ $ \sqrt{(\frac{h\nu}{m_{\rm e}c^{2}})\cdot (\frac{h\nu}{kT_{\rm e}})}$. If $h\nu$ is large enough and comparable with $kT_{\rm e}$, even $h\nu \gg kT_{\rm e}$, which often occurs in hard X-ray or $\gamma$-ray astronomy, the second term cannot be neglected anymore. Inserting Eq. (9) into Eq. (7), the integral can be reduced to:
 
                                    I $\textstyle \equiv$ $\displaystyle h^2\int {\rm d}^3\vec{p}\int {\rm d} W f(\vec{p})\Delta^{2}$  
  = $\displaystyle \int {\rm d}^3\vec{p}\int {\rm d} W f(\vec{p})\left[-\frac{hc\nu}{m_{\rm e} c^2} \vec{p}\cdot(\vec{n}-\vec{n}^{\prime})\right .$  
    $\displaystyle ~~~~~~~~~~~~~~~~~~~~~~~~~~~\left .-\frac{(h\nu)^2}{m_{\rm e} c^2} (1-\vec{n}\cdot\vec{n}^{\prime})\right]^2$  
  = I1+I2+I3 (10)

with

\begin{eqnarray*}I_1&=&\left(\frac{h\nu}{m_{\rm e}c}\right)^2\int {\rm d}^3\vec{... ...\cdot(\vec{n}-\vec{n}^{\prime})(1-\vec{n}\cdot\vec{n}^{\prime}). \end{eqnarray*}


Fixing $\vec{n}-\vec{n}^{\prime}$ as the z-axis and denoting the angle between $\vec{p}$ and $\vec{n}-\vec{n}^{\prime}$ as $\Theta$, we get:

\begin{eqnarray*}I_1&=&\left(\frac{h\nu}{m_{\rm e}c}\right)^2\int {\rm d}^3\vec{... ...)^2 \int {\rm d} W \vert\vec{n}-\vec{n}^{\prime}\vert^2\right]. \end{eqnarray*}


Inserting $f(\vec{p})=f_0 \exp(-p^2/2m_{\rm e}kT_{\rm e})$ into the expression I1, we have $\int p^2f(p)4\pi p^2{\rm d}p=3m_{\rm e}kT_{\rm e}$, therefore,

\begin{eqnarray*}I_1&=&\left(\frac{h\nu}{m_{\rm e} c}\right)^2m_{\rm e}kT_{\rm e... ...\rm e}kT_{\rm e} \int(1-\vec{n}\cdot\vec{n}^{\prime}){\rm d} W. \end{eqnarray*}


Similarly, we have:

\begin{eqnarray*}I_2&=&\left(\frac{h^2\nu^2}{m_{\rm e}c^2}\right)^2\int {\rm d}^... ...Theta\sin\Theta {\rm d}p {\rm d}\Theta {\rm d}\varphi\\ &=& 0. \end{eqnarray*}


Therefore, we have to calculate the integrals $\int(1-\vec{n}\cdot\vec{n}^{\prime}){\rm d}W$and $\int(1-\vec{n}\cdot\vec{n}^{\prime})^2{\rm d} W$ for all scattering directions $\theta$. In the non-relativistic limit, the Compton differential scattering cross section can be approximately expressed by the Thomson cross section, so the transition probability is:
 
$\displaystyle {\rm d} W =c{\rm d}\sigma_{\rm T}=c\frac{r^2_0}{2}(1+\cos^2\theta) 2\pi\sin\theta {\rm d}\theta$     (11)

where $r_0=e^2/m_{\rm e}c^2$ is the classical radius of the electron, $\theta$ is the scattering angle, $\cos\theta=\vec{n}\cdot\vec{n}^{\prime}$. From Eq. (11) we see that ${\rm d} W$ is symmetric for angles $\theta$ and $\pi -\theta$. Therefore,

\begin{eqnarray*}\int\vec{n}\cdot\vec{n}^{\prime}~{\rm d} W =\int\cos\theta {\rm d} W = 0 \end{eqnarray*}


and

\begin{eqnarray*}\int(1-\vec{n}\cdot\vec{n}^{\prime})~{\rm d} W =\int {\rm d} W = c\int {\rm d}\sigma_{\rm T}=\sigma_{\rm T} c \end{eqnarray*}



\begin{eqnarray*}\int(1-\vec{n}\cdot\vec{n}^{\prime})^2{\rm d} W =\pi c r_0^2\in... ...2\theta)^2 \sin\theta {\rm d}\theta=\frac{7}{5}\sigma_{\rm T} c. \end{eqnarray*}


Hence we can get:
 
$\displaystyle I=I_1+I_2=\frac{2kT_{\rm e}}{m_{\rm e}c^2}(h\nu)^2\sigma_{\rm T} c+ \frac{(h\nu)^4}{m^2_{\rm e} c^4}\frac{7}{5}\sigma_{\rm T} c.$     (12)

Now we calculate the first integral in Eq. (6). Let:
 
$\displaystyle H\equiv h\int {\rm d}^3\vec{p}\int {\rm d} W f(\vec{p})\Delta.$     (13)

It is difficult to calculate Eq. (13) directly, but it can be deduced simply from the integral value I given by Eq. (12) from the consideration of conservation of total number of photons in the scattering process. The conservation equation in the "frequency space'' can be written as
 
$\displaystyle \left(\frac{\partial n}{\partial t}\right)_{\rm c}=-\nabla \cdot \vec{j} =-x^{-2}\frac{\partial(x^2j)}{\partial x}$     (14)

where $\vec{j}$ is the flux of photons defined in the "frequency space'' where the Cartesian coordinates are x1, x2 and x3 respectively, therefore $\sqrt{x^2_1+x^2_2+x^2_3}=x\equiv \frac{h\nu}{kT}$. Using spherical coordinates $(x,\theta,\varphi)$ to replace (x1,x2,x3), we have $\nabla\cdot\vec{j}=\frac{1}{x^2} \frac{\partial(x^2j)}{\partial x}$ and Eq. (14) is obtained. Equation (14) can be rewritten as:
 
$\displaystyle \left(\frac{\partial n}{\partial t}\right)_{\rm c}=-\frac{2}{x}j-\frac{\partial j}{\partial x}\cdot$     (15)

Comparing Eq. (15) with Eq. (6), we see that the flux of photons j has to be taken in the form[*]:
 
$\displaystyle j(x)=g(x)\left[\frac{\partial n}{\partial x}+n(1+n)\right].$     (16)

Inserting Eq. (16) into Eq. (15), we get:
 
$\displaystyle \left(\frac{\partial n}{\partial t}\right)_{\rm c}=-g(x)\left[\fr... ...rtial x}+\frac{2g}{x}\right]\left[\frac{\partial n} {\partial x}+n(1+n)\right].$     (17)

On the other hand, Eq. (6) can be written as:
 
$\displaystyle \left(\frac{\partial n}{\partial t}\right)_{\rm c}=\left[\frac{\p... ...\partial^2 n} {\partial x^2}+(2n+1)\frac{\partial n}{\partial x}+n(1+n)\right].$     (18)

Comparing Eq. (17) with Eq. (18) and noting that the coefficient of $\partial^{2}n/\partial x^{2}$ should be the same, g(x) is obtained as:
 
$\displaystyle g(x)=-\frac{N_{\rm e}}{2}(kT_{\rm e})^{-2}I=-Ax^2\left(1+Bx^2\right)$     (19)

where $A \equiv (kT_{\rm e}/m_{\rm e}c^2)N_{\rm e}\sigma_{\rm T} c$, ${B} \equiv \frac{7}{10}(kT_{\rm e}/m_{\rm e}c^2)$. Inserting Eq. (19) into Eq. (17) and comparing Eq. (17) with (18) again, we obtain
 
$\displaystyle H=\frac{kT_{\rm e}}{N_{\rm e}}{A}\left[4x-x^2+6{B}x^3-{B}x^4\right].$     (20)

Therefore, the extended Kompaneets equation under the condition $kT_{\rm e} \ll m_{\rm e}c^{2}$ and $h\nu\ll m_{\rm e}c^2$is obtained as follows:
 
$\displaystyle \left(\frac{\partial n}{\partial t}\right)_{\rm c}=\frac{kT_{\rm ... ...rm e} c^2}x^2\right)\left[\frac{\partial n} {\partial x}+n(1+n)\right]\right\}.$     (21)

Equation (21) is a significant improvement of the original Eq. (2) because the conditions $h\nu\ll m_{\rm e}c^2$ and $kT_{\rm e} \ll m_{\rm e}c^{2}$are much looser than those of Eq. (2). The comparison between $h\nu$ and $kT_{\rm e}$ is no longer necessary, i.e. Eq. (21) can be applied for various cases, $h\nu \ll kT_{\rm e}$, $h\nu \simeq kT_{\rm e}$, particularly for $h\nu \gg kT_{\rm e}$. It is obvious that Eq. (21) returns to the original Kompaneets Eq. (2) in the low frequency limit, i.e. when $x\equiv h\nu/kT_{\rm e}\ll 1$.

Using the diffusion Eq. (21), we studied the effect of Comptonization on the X-ray spectrum emergent from the scattering medium. The effect of down-Comptonization, for which Eq. (2) is difficult to use, is specially emphasized. Given the initial condition, i.e., the initial spectrum, and the diffusion time scale, we can solve Eq. (21) numerically. For simplicity, we assume that the X-ray source is stable, with a time-independent spectrum f(x). So the initial condition of Eq. (21) at t=0 can be expressed as

 
n(x,0)=f(x).     (22)

The diffusion time scale T, which represents the average time before a trapped photon escapes from the surface of the scattering region of the electron gas, is determined by both the geometric size L and the scattering depth $\tau_{\rm s}$ of the scattering region, where L denotes either the radius of a spherical scattering cloud surrounding the central X-ray source, or the thickness of a plane-parallel scattering slab, dependent on the geometry of the problem. According to the random walk theory (Ribicky & Lightman 1979), the average scattering number of a scattered photon before escaping is $N\approx \max(\tau_{\rm s},\tau_{\rm s}^{2})\approx \tau_{\rm s}^{2}$ for $\tau_{\rm s}\ga 1$, where $\tau_{\rm s}=N_{\rm e}\sigma_{\rm T}L$ is the scattering depth. The average trapping time $T=N\cdot \frac{l}{c}$, where $l=\frac{1}{N_{\rm e}\sigma_{\rm T}}$ is the mean free path of Compton scattering. Therefore we get:
 
$\displaystyle T\simeq \tau_{\rm s}^{2}\frac{1}{N_{\rm e}\sigma_{\rm T}c} =\frac{\tau_{\rm s} L}{c}, \qquad {\rm for~case }~~\tau_{\rm s} > 1.$     (23)

In this paper, we will use Eq. (23) to estimate the diffusion time scale T to get the final emergent spectrum. In Eq. (23) we assumed a large scattering depth $\tau_{\rm s} >1$, optically thick, as otherwise the Comptonization process becomes unimportant. Combining Eqs. (21)-(23), we obtain the emergent Comptonization-spectrum n(x,T). Finally, the emergent spectral intensity Iis $I(\nu,T)\equiv I(x,T)\propto x^{3}n(x,T)$[*].

3 Some typical numerical solutions and comparison with Monte Carlo simulation and the Ross-McCray equation

In this section, we give some numerical solutions of the diffusion Eq. (21) under given conditions (22) and (23) to show the time-evolution behavior of a spectrum in the Comptonization, particularly in down-Comptonization process. We focus our attention on the high-frequency end of the primary incident hard X-ray or $\gamma$-ray spectrum, where $h\overline{\nu}\gg kT_{\rm e}$ is satisfied, for which the original Kompaneets Eq. (2) fails to give a correct results. We will compare our solutions with those of a Monte Carlo simulation to confirm the correctness and the effectiveness of our Eq. (21). Under the same conditions, we also calculated the evolution of the X-ray spectrum by use of the RMC Eq. (3) to compare with the results of our Eq. (21). For simplicity and to be specific, we assume that the scattering medium is homogeneous in space with constant electron density $N_{\rm e}$ and temperature $T_{\rm e}$. We fix $N_{\rm e}\approx 2\times 10^{16} ~\rm {cm^{-3}}$ and $kT_{\rm e}=1~\rm {keV}$ respectively, which are the typical values for the plasma in the accretion disk around a compact star, e.g., a neutron star, a singular star or a black hole with stellar mass (Treves et al. 1989; Sunyaev et al. 1973). Different from Chen et al. (1994) and Deng et al. (1998), in this paper Eq. (21) is solved numerically by use of the three-point finite difference method. For generality, we abandon the assumption of a "tenuous radiation field'', $n\ll 1$, adopted by Chen (1994) and Deng (1998), and retain the n2 term in Eq. (21). This means that our results can be used to discuss very luminous high-energy objects with a strong radiation field. For comparison, we present the Monte Carlo results using Hua's program of Monte Carlo simulation of Comptonization (Hua 1997).

3.1 Emission line with Gaussian profile

The expression of the spectral intensity for an emission line with normal Gaussian profile is given as:

 
$\displaystyle I(\nu) ~\sim ~ \exp\left [-\frac{4\ln 2}{(\Delta \nu)^{2}}\cdot (\nu-\nu_0)^2\right ].$     (24)

For definiteness, in the following calculation, we take the centroid energy of the line as $h\nu_{0}=10~{\rm keV}$, the FWHM of the line as $h\Delta \nu=1 ~{\rm {keV}}$. Using the relation $I_{\nu}=\frac{2h\nu^{3}}{c^{2}}n(\nu)$, the initial line profile of a stable source is given as
 
$\displaystyle n(x,0)=f(x)~\sim ~ x^{-3}\exp \left [-\frac{4\ln 2}{(\Delta x)^{2}}\cdot (x-x_0)^2\right ].$     (25)

The frequency range $2~\rm {keV}\leq h\nu \leq 14~\rm {keV}$ is taken in our calculation. The numerical solutions of our extended Eq. (21), the Monte Carlo simulation, and the RMC Eq. (3) are shown in Figs. 1a, b, c, respectively. The solid line represents the original line-profile before down-Comptonization; the dashed, the dash-dot and the dotted curve respectively represent the emergent line-profile with different diffusion time scale $T=1.0\times 10^{-2}~{\rm s}$, $5.0\times 10^{-2}~{\rm s}$ and $1.0\times 10^{-1}~{\rm s}$ under a given temperature $kT_{\rm e}=1.0 ~\rm {keV}$. The evolution behavior of the line profile with time T is obvious. And the difference of the resultant line-profiles between Figs. 1a and 1b, is very small, which confirms the correctness of Eq. (21). The line-profiles in Compton-evolution obtained by our extended equation are also consistent with those of the RMC equation (see Figs. 1a and 1c). From Figs. 1a, b or c we see, the main changes of the line profile in down-Comptonization are as follows: (i) the peak position of the line is shifted to lower energies. We call this "Compton redshift'' of the line. The longer the diffusion time T, the larger the line redshift will be; (ii) the line-width increases with T, which really shows a diffusion behavior in the "frequency space''. (iii) the line-profile becomes asymmetrical, steeper on the low-energy side and flatter on the high-energy side.
  \begin{figure} \par\includegraphics[width=8.3cm]{Ekg.EPS}\vspace*{6mm} \includeg... ...h=8.3cm]{Mcg.EPS}\vspace*{6mm} \includegraphics[width=8.3cm]{Rg.EPS}\end{figure} Figure 1: The down-Comptonization evolution behavior of a Gaussian line; a), b) and c) represent the calculated results by the extended Kompaneets equation, by the Monte Carlo simulation and by the RMC equation respectively. The solid line is the original spectrum, Dashed line, dash-dot-line and dotted line represent the profile of the emergent spectrum with diffusion time scale $T=1.0\times 10^{-2}~{\rm s}$, $5.0\times 10^{-2}~{\rm s}$ and $1.0\times 10^{-1}~{\rm s}$ ( $T_{\rm e}=1.0 ~\rm {keV}$), respectively.
Open with DEXTER

3.2 Black body spectrum

The Planckian formula of a black body spectrum gives

 
$\displaystyle I(\nu)=B_{\nu}(T_{\rm b})=\frac{2h\nu^{3}}{c^{2}}\cdot \frac{1}{\exp\left [h\nu/kT_{\rm b}\right ]-1}$     (26)

To be specific, in the following calculation, we take $kT_{{\rm b}}=2~{\rm {keV}}$and $kT_{\rm e}=1.0 ~\rm {keV}$. $kT_{\rm e}\leq kT_{{\rm b}}$ to ensure the possibility of down-Comptonization, particularly at the high energy side of the black body spectrum. The initial condition now becomes:
 
$\displaystyle n(x,0)=f(x)~\sim ~ \frac{1}{\exp\left [x\cdot (T_{\rm e}/T_{\rm b})\right ]-1}\cdot$     (27)

In the following model calculation, we take the frequency range $0.1 ~{\rm {keV}}\leq h\nu \leq 26 ~{\rm {keV}}$, and $kT_{\rm e}=1.0 ~\rm {keV}$. The diffusion time scales are $T=1.0\times 10^{-2}$, $5.0\times 10^{-2}$ and $1.0\times 10^{-1}~{\rm s}$respectively. Figures 2a, b, c represent the time-evolution of the initial black body spectrum, calculated by the extended Kompaneets Eq. (21), by the Monte Carlo simulation, and by the RMC Eq. (3), respectively, which further indicates the equivalence of these methods. As we expected, the main changes occur in the high frequency part of the initial black body spectrum where $h\nu \gg kT_{\rm e}$ and down-Comptonization is very efficient. The peak position of spectrum is shifted to lower frequencies, as expected in down-Comptonization. Finally, from Fig. 2 we see a tendency of evolution towards a new black body spectrum with new temperature $kT_{\rm e}\sim 1.0 ~\rm {keV}$. This is easy to understand becuase both the radiation field and the plasma will achieve a common thermal equilibrium with temperature $T_{\rm e}$ if the diffusion time $T\rightarrow \infty$.
  \begin{figure} \par\includegraphics[width=7.9cm]{Ekb.EPS}\vspace*{6mm} \includeg... ...=7.9cm]{Mcb.EPS}\vspace*{6mm} \includegraphics[width=7.9cm]{Rbb.EPS}\end{figure} Figure 2: The down-Comptonization evolution of a blackbody spectrum; a), b) and c) represent the calculated results using the extended Kompaneets equation, the Monte Carlo simulation, and the RMC equation respectively.
Open with DEXTER

3.3 Power-law spectrum

The power-law spectrum which often occurs in X-ray astronomy has the form:

 
$\displaystyle I(\nu)~\sim ~ \nu^{-\alpha} \qquad (\nu_{{\rm min}}\leq \nu \leq \nu_{{\rm\max}}).$     (28)

In the following calculations, we adopt a typical observation value of spectral index $\alpha =0.7$. For simplicity, we treat the power-law spectrum with sharp cutoffs. We take $h\nu_{{\rm {min}}}=0.1~\rm {keV}$ and $h\nu_{{\rm {max}}}=20~\rm {keV}$as the lower and the upper frequency cut, respectively. The initial distribution function is:
 
$\displaystyle n(x,0)=f(x)~\sim ~ \left \{\begin{array}{ll}x^{-3}x^{-\alpha}~~~~... ...{max}}} \\ 0~~~~~~~~~~~~~~x< x_{{\rm min}},~x>x_{\rm max}\end{array}\right . .$     (29)

Again, we adopt the time parameters as $T=1.0\times 10^{-2}$, $5.0\times 10^{-2}$, $1.0\times 10^{-1}~{\rm s}$ and $kT_{\rm e}=1.0 ~\rm {keV}$. Calculated results by Eq. (21), by the Monte Carlo method, and by Eq. (3) are shown in Figs. 3a, b, c, respectively. The main change in the evolution is the steepening of the power-law spectrum at high frequencies, where the condition $h\nu \gg kT_{\rm e}$is satisfied and the down-Comptonization is efficient. On the other hand, in the low-frequency band where $h\nu \ll kT_{\rm e}$ (because $h\nu_{\rm {min}}=0.1~\rm {keV} \ll kT_{\rm e}$), up-Comptonization occurs, therefore we expect that the low-frequency photons will inevitably be hardened, and shifted towards the higher frequency region, thus causing a drop at the low frequency side, as shown in Fig. 3.
  \begin{figure} \par\includegraphics[width=8.1cm]{Ekp.EPS}\vspace*{6mm} \includeg... ...h=8.1cm]{Mcp.EPS}\vspace*{6mm} \includegraphics[width=8.1cm]{Rp.EPS}\end{figure} Figure 3: The down-Comptonization evolution of a power-law spectrum; a), b) and c) represent the calculated results using the extended Kompaneets equation, the Monte Carlo simulation, and the RMC equation respectively.
Open with DEXTER

3.4 Thermal bremsstrahlung spectrum

In this section, we show the evolution of the thermal bremsstrahlung radiation originating from an optically thin plasma with temperature $kT_{\rm {ff}}= 5~{\rm {keV}}$. When the radiation passes through a "cold'' plasma with a lower temperature $kT_{\rm e}\approx 1.0 ~{\rm keV}$, down-Comptonization occurs because of $kT_{\rm e}< kT_{{\rm ff}}$, particularly for the high frequency portion of the incident bremsstrahlung spectrum where the condition $h\nu \gg kT_{\rm e}$ is satisfied. The bremsstrahlung spectrum is given as (Rybicki & Lightman 1979; You 1998)

 
$\displaystyle I(\nu)~\sim ~ \overline{g}_{{\rm {ff}}}(\nu, T_{{\rm {ff}}})\cdot \exp\left [-h\nu/kT_{{\rm {ff}}} \right ]$     (30)

where $\overline{g}_{{\rm ff}}(\nu, T_{{\rm {ff}}})\approx 1$ is the average Gaunt factor at temperature $kT_{\rm {ff}}$. In the following calculation, we adopt the same parameter group as above, $T=1.0\times 10^{-2}$, $5.0\times 10^{-2}$and $1.0\times 10^{-1}~{\rm s}$, the frequency range in our calculation is $0.1 ~{\rm {keV}}\leq h\nu \leq 20 ~{\rm {keV}}$. The initial condition is now expressed as:

 
$\displaystyle n(x,0)=f(x)~\sim ~ x^{-3}\cdot \overline{g}_{\rm {ff}}\cdot \exp \left [-x\cdot(T_{\rm e}/T_{\rm {ff}})\right ].$     (31)

The evolution spectra calculated by Eq. (21), by the Monte Carlo method, and by Eq. (3) are shown in Figs. 4a, b, c, respectively. From Fig. 4 we see that the evolution characteristics of the bremsstrahlung spectrum are quite similar to those of the power-law spectrum, shown in Fig. 3.


  \begin{figure} \par\includegraphics[width=8cm]{Ekbr.EPS}\vspace*{6mm} \includegr... ...dth=8cm]{Mcbr.EPS}\vspace*{6mm} \includegraphics[width=8cm]{Rbr.EPS}\end{figure} Figure 4: The down-Comptonization evolution of a thermal bremsstrahlung spectrum; a), b) and c) represent the calculated results using the extended Kompaneets equation, the Monte Carlo simulation, and the RMC equation respectively.
Open with DEXTER

4 Conclusion and discussion

The original Kompaneets equation fails to describe down-Comptonization, which is the most important radiative transfer process in hard X-ray and $\gamma$-ray astronomy. Therefore in this paper, we present an extended Kompaneets Eq. (21), which can be used to describe a more general Comptonization process, including down- and up-Comptonization, suitable for any case, $h\overline{\nu}\ll kT_{\rm e}$, $h\overline{\nu}\gg kT_{\rm e}$ and $h\overline{\nu}\sim kT_{\rm e}$. In other words, the condition $h\overline{\nu}\ll kT_{\rm e}$for the original Kompaneets Eq. (2) is no longer necessary. Using Eq. (21), we give some typical solutions in X-ray astronomy, and compare them with Monte Carlo simulation results. The excellent consistency of the Compton-evolution spectra obtained by the use of the two different methods confirmed the correctness of our extended Kompaneets equation. Though the prevailing approach to a quantitative analysis of the Comptonization process is still the Monte Carlo simulation, it has a remarkable defect in practice, namely the great quantity of computations. In addition to the simplicity and the clarity in physics, a more important advantage of solving Eq. (21) is that it is less expensive in terms of computational time. The solutions of the extended Kompaneets equation are also consistent with those of the RMC equation particularly for down-Comptonization. Therefore Eq. (21) can be regarded as an important improvement of the original Kompaneets Eq. (2), particularly for hard X-ray and $\gamma$-ray astronomy.

Our calculations show that the change of the emergent spectrum in Comptonization, particularly in down-Comptonization, depends on both the difference between $h\overline{\nu}$ and $kT_{\rm e}$, and the scattering depth $\tau_{\rm s}$(or equivalently, the diffusion time scale T). The larger the difference $h\overline{\nu} - kT_{\rm e}$, and/or the larger the depth $\tau_{\rm s}$(equivalently, the time scale T), the larger the change of the emergent spectrum will be.

The structure of Eq. (21) has two marked characteristics which can be regarded as a criterion of the correctness of our Eq. (21). (i) It has the form of $\frac{1}{x^2}\frac{\partial}{\partial x} \left\{x^4 \left(1+\frac{7}{10}\frac{k... ...rm e} c^2}x^2\right) \left[\frac{\partial n} {\partial x}+n(1+n)\right]\right\}$, which ensures the conservation of number of photons (see the conservation Eq. (14)). The invariance of the total number of photons is a basic requirement for the electron-photon scattering process. (ii) Equation (21) contains a factor $\left [\frac{\partial n}{\partial x}+n(n+1)\right ]$ which ensures $\partial n/\partial t=0$ when the photon gas reaches the thermal equilibrium distribution $n=({\rm e}^{x}-1)^{-1}$. This is also a necessary requirement for the correctness of this equation because the thermal equilibrium will inevitably be reached through the scattering and $\partial n/\partial t \rightarrow 0$, which implies that the diffusion finally stops. Ross et al. (1978) also noticed the restriction of the original Kompaneets Eq. (2) and suggested an alternative equation to replace Eq. (2) to suit the case $h\overline{\nu}\gg kT_{\rm e}$. Their equation has a form $\frac{1}{x^2}\frac{\partial}{\partial x} \left\{x^4 \left[n+\left(1+\frac{7}{10... ...{\rm e}}{m_{\rm e} c^2}x^2\right) \frac{\partial n} {\partial x}\right]\right\}$, which is quite similar to Eq. (21) but deviates from the necessary form $\left [\frac{\partial n}{\partial x}+n(n+1)\right ]$. According to their equation, the diffution will never stop, i.e., $\partial n/\partial t\neq 0$even when the thermal equilibrium, $n=({\rm e}^{x}-1)^{-1}$, has been reached.

Acknowledgements
We would like to thank the anonymous referee for his helpful suggestions. We are also grateful to Prof. McCray for his helpful discussions and suggestions when Dr. You visited JILA, Boulder in 2001. The work of D.B.L. is supported by the National Natural Science Foundation of China (Grant No. 10373010), and the Research Fund for the Doctoral Program of Higher Education of China (Grant No. 20030248012).

References


Copyright ESO 2004