DIFFUSIVE COSMIC-RAY ACCELERATION IN SAGITTARIUS A*

The first 25 months of observations of the Galactic center region by Fermi-LAT have revealed the presence of a source, 1FGL J1745.6-2900, at energies above 10 GeV, coincident with the TeV central object HESS J1745-290 (Chernyakova et al. 2011). The combined Fermi-HESS spectrum appears to be inflected, with a relatively steep region at intermediate energies, connected to flatter components at both low and high energies. More specifically, although the spectrum in the 100 MeV–300 GeV energy range can be fit by a power law with slope α = 2.212  ±  0.005 and a flux normalization F(100 MeV) = (1.39  ±  0.02) × 10⁻⁸ cm⁻² s⁻¹ MeV⁻¹, a significantly better fit is obtained with a power-law index α = 2.196  ±  0.001 between 300 MeV and 5 GeV, and a separate power law with α = 2.681  ±  0.003 in the range 5 GeV–100 GeV. The latter is also significantly steeper than the spectrum reported by the HESS collaboration, where α ≈ 2.1, with a total flux of (1.87  ±  0.30) × 10⁻⁸ m⁻² s⁻¹ above 1 TeV (Aharonian et al. 2009).

The Fermi-HESS source appears to be the high-energy counterpart to Sgr A* (see Melia 2007 for a review). We note, however, that although the HESS source is coincident within ∼30'' of Sgr A*, its centroid is displaced roughly 7'' (∼0.4 pc) to the east of the Galactic center, an issue that we shall revisit below.

In earlier work (Liu et al. 2006; Ballantyne et al. 2007), we considered the interesting possibility that the TeV γ-rays are produced via π⁰ decays generated when relativistic protons, accelerated near the black hole, collide with hadrons farther out. Such a scenario is intriguing because our basic understanding of Sgr A* and its nearby environment precludes any possibility of this object producing a significant flux of TeV photons directly (Liu & Melia 2001). However, protons can be energized to TeV energies by stochastic acceleration in a magnetically dominated funnel within 20–30 Schwarzschild radii of the event horizon. Then, as these cosmic rays diffuse out through the surrounding medium, they may scatter with hydrogen nuclei in a shocked stellar-wind region and with the circumnuclear molecular torus surrounding Sgr A*.

But the latest ∼10 GeV Fermi-LAT observations, and more recent theoretical work with diffusive particle acceleration (see, e.g., Fatuzzo et al. 2010; Melia & Fatuzzo 2011; Fatuzzo & Melia 2012), render this basic picture inadequate for several reasons. First, the combined Fermi-HESS measurements show that a simple extrapolation of the TeV spectrum into the GeV range does not provide an acceptable explanation for the data. Second, it now appears unlikely that the energy of a cosmic ray remains constant as it diffuses out to larger radii. Thus, our earlier conclusion—that the incipient proton power-law index has to be ≈0.75 near the black hole to produce the observed photon index ≈2.25 farther out—is not valid when the cosmic rays undergo continued acceleration on their way out into the circumnuclear environment.

Attempting to explain the surprising Fermi-HESS spectrum of the source 1FGL J1745.6-2900, Chernyakova et al. (2011) modified the Ballantyne et al. (2007) treatment in several ways. First, instead of adopting a two-phase medium surrounding Sgr A*, in which the relatively cold inner molecular torus is surrounded by a hotter, more tenuous wind-shocked region, they assumed a single plasma with a 1/r² density profile. Second, they invoked a time-dependent cosmic-ray injection process which, by virtue of diffusion-induced and energy-dependent time delays, can produce the kind of inflected spectrum seen by Fermi-LAT and HESS. They were able to produce a reasonable fit to the data, but only under the assumption that Sgr A* is highly variable on timescales of several hundred to several thousand years.

However, Chernyakova et al.'s (2011) analysis appears to be incomplete and inconclusive because several important aspects of the particle diffusion used by them are unrealistic. First and foremost, the assumption of a constant particle energy following ejection is difficult to justify, given what we know about the physical conditions surrounding Sgr A*. One of the principal goals of this Letter is in fact to incorporate this important effect into the calculation of the spectrum.

Second, the supposition that protons are ejected episodically by Sgr A* is motivated by the detection of reflected X-ray emission from the Sgr B2 cloud complex some 100 pc away (and other similar clouds at various distances from the center). This iron-line fluorescence is often interpreted as the light echo of a flare produced by the supermassive black hole hundreds (or thousands) of years ago (see, e.g., Sunyaev et al. 1993; Fromerth et al. 2001; Revnivtsev et al. 2004; Terrier et al. 2010). But this kind of phenomenon requires a variation in the overall power of Sgr A* by some six orders of magnitude during a very short time compared to its age. Though one cannot decidedly argue against such an unprecedented change, it is far more likely that the source of light responsible for these echoes was the impact onto the 50 km s⁻¹ molecular cloud behind Sgr A* by the supernova shell of the explosion that produced Sgr A East (see Fryer et al. 2006).

However, Chernyakova et al. (2011) also considered a steady-state proton ejection and were able to account for the Fermi data, albeit under unusual circumstances. In order for the proton diffusion to produce the required stratified proton energy distribution, with the highest-energy protons moving more or less rectilinearly toward larger radii, they found that the diffusion coefficient must be strongly dependent on the energy (i.e., D(E)∝E^0.65). A second principal goal of this Letter is to calculate D(E) from first principles, using the physical conditions prevalent around Sgr A*, to see if such an energy dependence is realistic. As we shall see, our results do not support such a strongly variable diffusion coefficient.

What this means, of course, is that it may be too simplistic to do away with the two-phase gaseous environment surrounding Sgr A*. The physical conditions in the circumnuclear disk are drastically different from those in the rest of the tenuous, hot medium filling the inner 3 pc. Adopting an average particle density, and a concomitantly simplified magnetic field structure in lieu of the actual variation of these quantities between the molecular and ionized phases, significantly alters the particle diffusion. For these three principal reasons, it is therefore necessary to examine whether a steady-state cosmic-ray ejection by Sgr A* can provide a reasonable alternative explanation for its 100 MeV–100 TeV spectrum, but without having to invoke a diffusion coefficient that depends strongly on the particle energy.

In order to investigate how cosmic rays are energized as they propagate through the environment immediately surrounding Sgr A*, we adopt a numerical formalism based on the pioneering work of Giacalone & Jokipii (1994) that has been used to study the transport of cosmic rays in static turbulent magnetic fields (Casse et al. 2002; De Marco et al. 2007; Wommer et al. 2008; Fatuzzo et al. 2010), and has recently been extended to study the energy diffusion of cosmic rays in time-dependent turbulent magnetic fields in highly resistive environments. For this latter scenario, the added presence of a turbulent electric field that is everywhere perpendicular to the total magnetic field drives the stochastic acceleration of cosmic-ray particles (O'Sullivan et al. 2009; Fatuzzo & Melia 2012).

Following these recent works, we model the magnetic field within the inner few parsecs of the galaxy as a uniform field ${\bm B}_0 = B_0 \hat{z}$ superimposed by a turbulent component (likely produced by the stellar winds in this region) that is expressed as a sum of N randomly directed Alfvén waves with wavelengths λ_n = 2π/k_n:

$$\begin{equation} \delta {\bm B} = \sum _{n=1}^N {\bm A}_n \, e^{i({\bm k}_n \cdot {\bm r}-\omega _n t+\beta _n)}\,, \end{equation} \tag{ 1 }$$

where k₁ = k_min = 2π/λ_max and k_N = k_max = 2π/λ_min are, respectively, the wavenumbers corresponding to the maximum and minimum wavelengths associated with the turbulent field, and the phase β_n of each term is randomly selected. The dispersion relations ω_n = v_Ak_n|cos θ_n| and magnetic orientations $\hat{A}_n$ follow from linear MHD theory, where v_A is the Alfvén speed and θ_n is the angle between k_n and B₀. The total electric field δE that must accompany a time-varying magnetic field is found by invoking the MHD condition, which then guarantees that δE and B = B₀ + δB are everywhere perpendicular. Finally, the desired spectrum of the turbulent magnetic field is set through the appropriate choice of index Γ in the scaling

$$\begin{equation} A_n^2 = A_1 ^2\left[{k_n \over k_1} \right]^{-\Gamma } {\Delta k_n\over \Delta k_1} = A_1^2\left[{k_n \over k_1} \right]^{-\Gamma +1}, \end{equation} \tag{ 2 }$$

where A₁ is set by the ratio η of the turbulent field energy density to the underlying static field energy density (see Fatuzzo & Melia 2012 for a more complete discussion).

In theory, the dynamics of cosmic-ray protons diffusing through the central few parsecs of our Galaxy can be solved by numerically integrating the governing equations of motion (see, e.g., Fatuzzo & Melia 2012). In practice, however, this approach is too computationally taxing given the large difference between the radius of gyration (∼10⁻⁶ pc for TeV protons in the ∼ mG fields found at the galactic center) and the size of the central region. This point is further compounded by the large number of particles that must be tracked in order to obtain meaningful statistics at the highest energies.

In order to develop a formalism that is computationally viable, we perform a suite of test runs designed to inform a simple, random-walk model that adequately captures the essential features of particle diffusion and acceleration. Specifically, we obtain displacement values Δx, Δy, Δz, and Δγ for two successive time intervals Δt = λ_max/c by numerically solving the full set of equations of motion for a particle with initial Lorentz factor γ₀ moving through a specified turbulent field environment (as defined by B₀, v_A, Γ, λ_max, and η). The process is then repeated 10³ times for a new realization of the same field (i.e., same field parameters but different values of the random variables k_n and β_n) in order to produce a scatter plot of successive displacement pairs (e.g., Δγ_{n + 1} versus Δγ_n). This entire process is then repeated at several different particle energies, and the results are used to produce correlated distributions of displacement values that can be randomly sampled to produce a robust random-walk model for the spatial and energy diffusion of cosmic rays throughout the prescribed turbulent field.

While this process will be detailed in a later work, we note that diffusion coefficients obtained by integrating the full equations of motion for protons with 10⁴ < γ < 10⁸ and a turbulent field with Γ = 5/3, λ_max = 0.1 pc, B₀ = 1.4 mG, η = 0.55, and v_A = 1.1 × 10⁻³ c indicate that the often adopted scaling D_γ∝γ^{(2 − Γ)} is not valid at lower energies. Interestingly, these diffusion coefficients match well with those generated by our random-walk model for γ ∼ 10⁴–10⁵, with the latter scheme then yielding diffusion coefficients for γ < 10⁵ that scale as D_x = D_y = 9.8 × 10²⁴ γ^0.11 cm² s⁻¹, D_z = 1.1 × 10²⁶γ^0.13 cm² s⁻¹, and D_γ = 1.0 × 10⁻¹⁰γ^1.5 s⁻¹.

Our analysis indicates that good fits to the Fermi-HESS data can be obtained for reasonable model parameters. For simplicity, we treat the inner parsecs of the galaxy as a uniform "wind zone" of radius R_esc that encompasses a high-density "torus" with an inner radius of 1.2 pc and a thickness of 1 pc. The density distribution in the wind zone is due almost entirely to the interactions of stellar winds from the surrounding young stars. According to numerical simulations (Rockefeller et al. 2004; also guided by Ruffert & Melia 1994 and Falcke & Melia 1997), the average density in this region is 〈n^sw_H〉 = 121 cm⁻³ (outside of the clouds, defined by the condition n^sw_H < 3 × 10³ cm⁻³). The torus is comprised primarily of molecular gas and has an average density 〈n^cnd_H〉 = 233, 222 cm⁻³. A static magnetic field B₀ is assumed to be oriented parallel to the plane of the torus, with a magnitude B₀ = 1.4 mG, not unlike the ∼3 mG value obtained by invoking equipartition with the hot (kT = 1.3 keV) stellar-wind gas (Rockefeller et al. 2004). The turbulent magnetic field is set by adopting the values Γ = 5/3, λ_max = 0.1 pc, and η = 0.55.

Our Monte Carlo scheme is used to track the evolution of N_p protons injected at mildly relativistic energies near the black hole. We postulate that these protons are drawn from the extreme high-energy tail of a thermal distribution near the black hole, having sufficiently large energies for a similar stochastic mechanism to efficiently accelerate them out of their thermal state. The remaining particles near the black hole are then efficiently re-thermalized. We further postulate that none of the ambient particles in the wind zone are sufficiently energetic to be accelerated efficiently from their thermal state.

The initially mildly relativistic particles random walk through the surrounding medium until they either collide with an ambient proton in the stellar-wind region (where τ_pp = 1.5 × 10¹³ s), enter the torus, or escape once they diffuse beyond a radius R_esc. The probability of scattering in the wind zone for a given random-walk step Δt = λ_max/c is given by

$$\begin{eqnarray} &&p = 6.7 \times 10^{-7} \left({\lambda _{{\rm max}}\over 0.1\,{\rm pc}} \right) \left({\big\langle n_{\rm H}^{\rm sw}\big\rangle \over 121\,{\rm cm}^3}\right) \left({ \kappa \over 0.45}\right)\left({ \sigma _{pp}\over 40 \,{\rm mb}}\right),\nonumber \\ \end{eqnarray} \tag{ 3 }$$

but it should be noted that particles only diffuse a distance ∼0.02 pc during that time. In contrast, the probability of scattering in the torus for a similar step size is ∼1.3 × 10⁻³ given the much larger density. Since a particle would need ∼(0.5 pc/0.02 pc)² = 625 steps to traverse the torus (assuming similar field strengths and turbulence profiles as in the wind zone), particles entering the torus have a high probability of scattering. For simplicity, we therefore assume that any proton that enters the torus undergoes pp scattering within that region.

To estimate the escape radius R_esc of the wind zone, we note that ∼10% of a relativistic proton's energy goes into each of the two photons produced via the π₀ decay following an inelastic scattering event. The turnover in the HESS data therefore indicates that the maximum Lorentz factor attained by the underlying cosmic-ray population is γ_max ∼ 10⁵. We can therefore estimate R_esc by equating the acceleration time τ_acc ≡ γ²_max/D_γ(γ_max)  ≈  3 × 10¹² s required for particles to reach the maximum Lorentz factor to the escape time τ_esc ≡ R²_esc/D_z(γ_max) required for particles to diffuse out of the wind-zone region (since particles diffuse more readily in the direction parallel to the underlying field). This procedure yields an escape radius R_esc ≈ 10 pc. While somewhat larger than the HESS source point-spread function, this value is reasonable given the simplicity of our model.

A further complication for our Monte Carlo scheme arises from the fact that the Fermi-HESS data range over five orders of magnitude in energy, with the luminosity emitted at lower energies exceeding that at higher energies by ∼2 orders of magnitude. As a result, one can estimate that there are ∼10⁷ times more lowest-energy particles than highest-energy particles producing the observed emission. Obtaining good statistics for the entire energy range would then require that N_p ∼ 10⁹ particles be tracked. Since this would exceed our computational resources, we instead track N_p = 2 × 10⁶ particles, and extrapolate our results up to Lorentz factors γ ∼ 10⁵.

The resulting distributions in energy with which particles pp scatter are shown in Figure 1, where the dashed histogram denotes particles that scatter in the wind zone and the dotted histogram denotes particles that scatter in the torus. The solid lines show the linear fits to each distribution that are then used to extrapolate our results up to γ_max = 10⁵. Figure 1 shows that the majority of injected particles hit the torus. This result is not surprising given that particles diffuse a greater distance along the direction of the underlying field than across it, and are therefore naturally "funneled" into the inner edge of the torus for our adopted field geometry. However, a small fraction of particles are able to diffuse sufficiently far across the field to initially miss the torus. These particles thus diffuse throughout the surrounding wind-zone region until they either propagate back into the torus, pp scatter in the wind zone, or escape. This scenario naturally produces two distinct power-law distributions, and a total particle distribution that hardens from an index of ∼2.6 to ∼1.8 above ∼10¹³ eV.

Zoom In Zoom Out Reset image size

Upon pp scattering, these particles produce pions that either decay into photons (π⁰ → 2γ) or (in the case of π^±) via the muon channels into electrons and positrons. Each particle's contribution to the ensuing gamma-ray emissivity (via the π⁰ channel) can be calculated (see, e.g., Fatuzzo & Melia 2003), and because we are assuming a steady-state ejection of particles from Sgr A*, the total (unnormalized) emissivity is then found by summing over the full ensemble of protons that scatter. The (normalized) results shown in Figure 2 clearly indicate that steady-state stochastic acceleration within the central parsecs provides a viable explanation for the Fermi-HESS data (the lowest-energy data point in Figure 2 can be likely accounted for by the emission of the secondary leptons produced by the charged pion decays—see for example Figure 1 in Fatuzzo & Melia 2005).

Zoom In Zoom Out Reset image size

A full analysis of how the model parameters affect our results is beyond the scope of this work, and will be presented elsewhere. We note, however, that in general, the strength of the electric field increases as B₀, v_A, and η increase. As such, the particle distributions presented in Figure 1 would become harder as these parameters are increased. In contrast, particle acceleration becomes less efficient as λ_max increases, and is not sensitive to the value of Γ when η ≳ 1 (Fatuzzo & Melia 2012). In addition, particles diffuse more readily across the underlying magnetic field as η is increased, which has the effect of increasing the number of protons that "miss" the torus and produce the HESS spectrum.

We have lent support to the idea that—given the right environment—cosmic rays can be accelerated to TeV energies via stochastic processes in a turbulent magnetic field. Although this has been suspected for many years, we can only now demonstrate this quantitatively by comparing detailed numerical simulations with high-quality HESS and Fermi-LAT data. Previously, we showed that the cosmic-ray population permeating the inner several hundred parsecs of the Galaxy could be accounted for with this process, but only with the physical conditions encountered in that region (Fatuzzo & Melia 2012). As such, the fact that such a relativistic particle distribution is found only near the Galactic center can be understood as a product of the unique conditions found there.

In this Letter, we have focused on the particle acceleration occurring near Sgr A* itself. The combined data exhibited in Figure 1 present quite a challenge to any effort at understanding how these cosmic rays are produced. Related to this question is the issue of whether or not these observations provide evidence for strong variability in Sgr A*, as suggested by the fluorescent X-ray echoes detected in nearby molecular clouds.

But here we have shown that the use of a molecular torus and an interstellar medium filled with shocked stellar winds can provide just the right framework for steady-state stochastic acceleration to produce the observed two-component spectrum. We therefore suggest that Sgr A*'s high-energy spectrum does not necessarily provide evidence for the kind of rapid variability required to produce the fluorescent X-ray echoes.

We thank the referee for a careful review and useful comments that improved the quality of the writing M.F. is supported at Xavier University through the Hauck Foundation.