© ESO 2019

1. Introduction

Vela X is a middle-aged pulsar wind nebula (PWN) that is powered by the Vela pulsar (PSR B0833-45; Weiler & Panagia 1980), and located in the centre of the Vela supernova remnant (SNR). As one of the closest pulsars to the Earth, the Vela pulsar has a period of P = 89.3 ms (Large et al. 1968) and a period derivative of Ṗ = 1.25 × 10⁻¹³ s s⁻¹ (Lyne et al. 1996), which implies that the current spin-down luminosity of this pulsar has a value of L₀ = 6.9 × 10³⁶ erg s⁻¹ and a characteristic age of τ_c = 11.3 kyr (Taylor et al. 1993; Lyne et al. 1996). The value of distance to the Vela SNR is estimated to be 290 pc according to the radio and optical parallax measurements (Caraveo et al. 2001; Dodson et al. 2003), and the age is usually considered to be the characteristic age of the Vela pulsar (Reichley et al. 1970; Manzali et al. 2007). Lyne et al. (1996) estimated that the braking index for the Vela pulsar has a value of n = 1.4. However, based on the new long-term braking index measurements, Espinoza et al. (2017) suggested that the value of the braking index of the Vela pulsar is n = 1.7, which is adopted in this paper.

The photon emission from Vela X has been detected at radio, X-ray, and γ-ray bands, and showed that the nebula is a morphologically complex source. At radio and GeV energy ranges, Vela X shows a filamentary structure usually referred to as the halo (e.g. Abdo et al. 2010; Hinton et al. 2011; Grondin et al. 2013). In the radio band, Alvarez et al. (2001) revealed that the PWN has a power-law spectral index of 0.39 ± 0.03 in the 30 − 8400 MHz energy range. In the GeV γ-rays, according to the Fermi Large Area Telescope (LAT) observations, Abdo et al. (2010) revealed that the non-thermal spectrum of Vela X can be well fitted by a power law with a photon index of 2.41 ± 0.09 ± 0.15 in the 0.2 − 20 GeV energy range, and Grondin et al. (2013) pointed out that the photon spectrum above 0.2 GeV can be well described by a smoothly broken power law with spectral indices of Γ₁ = 1.83 ± 0.07 ± 0.27, Γ₂ = 2.88 ± 0.21 ± 0.06 and a break energy of E_b = 2.1 ± 0.5 ± 0.5 GeV. In the X-ray band, the ROAST observations revealed a much smaller emission region usually called the cocoon. (Markwardt & Ögelman 1995), and the Chandra X-ray images showed a double torus and jet-like structure at the centre of the PWN (Helfand et al. 2001). The ROAST and ASCA data found that the power-law photon index in the X-ray band has a value of 2.1 (Markwardt & Ögelman 1995, 1997), and the XMM-Newton data showed that the photon index is 2.3 (LaMassa et al. 2008). In the inner regions of the nebula, based on the XMM-Newton/MOS and BeppoSAX/MECS data, Mangano et al. (2005) pointed out that the power-law photon index and flux density in the 3 − 10 keV energy range have obvious changes as the radial distance increases. In the TeV γ-rays, the photon emission region, which has been detected by High Energy Stereoscopic System (H.E.S.S.), is coincident with the radio halo and X-ray cocoon. (Aharonian et al. 2006; Abramowski et al. 2012). Aharonian et al. (2006) have shown that the energy spectrum of Vela X in the 0.55 − 65 TeV energy range can be well fitted by an exponentially cut-off power-law function, with a photon index of 1.45 ± 0.09_stat ± 0.2_sys and an exponential cut-off at an energy of 13.8 ± 2.3_stat ± 4.1_sys TeV. While according to the deeper observation of the entire Vela X nebula region with H.E.S.S., Abramowski et al. (2012) revealed that the total differential energy spectrum in the energy of 0.75 − 70 TeV can be well described with a hard photon index of 1.32 ± 0.06_stat ± 0.12_sys and a cut-off energy of 14.0 ± 1.6_stat ± 2.6_sys TeV. Moreover, the H.E.S.S. Collaboration reported that the energy averaged surface brightness of Vela X decreases with the increase of radial distance in the TeV γ-ray band (Abramowski et al. 2012).

To interpret the non-thermal radiative properties of Vela X, various models have been proposed (e.g. Horns et al. 2006; de Jager et al. 2008; Zhang & Yang 2009; Abdo et al. 2010; Hinton et al. 2011; Mattana et al. 2011). Under the assumption that the particles in the nebula are composed of electrons and hadrons, Horns et al. (2006) proposed a leptonic-hadronic model, in which the X-ray emission is produced by the synchrotron radiation of relativistic electrons and the TeV photons are mainly from the decay of neutral pions produced in proton-proton collisions inside the nebula. In parallel, de Jager et al. (2008) suggested that there are two distinct electron populations existing in Vela X, one is responsible for the observed spectrum in radio and GeV γ-ray bands, and the other for the X-ray and TeV γ-ray emission. Both the leptonic-hadronic and leptonic models can well reproduce the non-thermal photon spectrum of Vela X. However, these two models cannot explain the spatial variations of the energy spectrum of the PWN observed at X-ray and TeV γ-ray bands, because the distribution of particles in the nebula is assumed to be uniform. In recent years, a spatially dependent particle transport model was proposed to explain the spatial variation properties of PWNe (e.g. Vorster & Moraal 2013; Porth et al. 2016; Lu et al. 2017a; Ishizaki et al. 2017; van Rensburg et al. 2018), in which the particle number density is spatially dependent and the processes of convection and diffusion are taken into account. Assuming that the PWN is a spherical symmetric system, Porth et al. (2016) studied the radial variations of the photon index and flux density of Vela X in the X-ray band with a steady-state particle transport model, and found that the spatial variations of photon index and flux density of the inner regions of the nebula can be well reproduced in the framework of their model. However, Porth et al. (2016) did not investigate the multi-band photon emission and spatially radiative properties of γ-rays of the PWN, because their main aims are to investigate the particle transport and the influence of the diffusion effect on the particle spectrum.

In this paper, we study the multi-band photon emission and high energy spatially radiative properties of Vela X in the framework of a spatially dependent particle transport model proposed by Lu et al. (2017a). In the model, the object is assumed to be a spherical symmetric system. The X-ray and very high energy (VHE) γ-ray observations have shown that the Vela X nebula is a highly asymmetric PWN (Markwardt & Ögelman 1995; Aharonian et al. 2006; Abdo et al. 2010). However, the X-ray observations of the inner region of Vela X show a bright X-ray nebula that is approximately spherically symmetric (Mangano et al. 2005; Mattana et al. 2011). Besides, the multi-wavelength maps of the analysis region of Vela X are usually treated as a circular region (e.g. Aharonian et al. 2006; Abdo et al. 2010; Abramowski et al. 2012; Tibaldo et al. 2018). Thus, as a first-order approximation, the spatially dependent and spherical symmetric model can help us to study the spatial variations of azimuthally integrated spectra in the Vela PWN. The organization of the paper is as follows. In Sect. 2, the spatially dependent particle transport model is simply described and the modelling results are presented. The main conclusions and discussions are given in Sect. 3.

2. Model and results

In this section, the spatially dependent particle transport model of Lu et al. (2017a) will be briefly reviewed, and the modelling results of Vela X, including the spectral energy distributions (SEDs), radial variations of photon index, and flux density of the X-ray emission in the inner regions of the nebula, as well as the spatial variations of photon indices and surface brightness in the GeV and TeV bands, will be presented.

2.1. Model description

In a spherically symmetric system, the evolution of electron distribution n_e = n_e(r, γ, t) inside the PWN can be described with the particle transport equation

$$\begin{array}{cc}\hfill \frac{\partial n_{\mathrm{e}}}{\partial t}=& \kappa \frac{{\partial }^2{n}_{\mathrm{e}}}{\partial r^2}+[\frac{1}{r^2}\frac{\partial }{\partial r}(r^2\kappa )-V]\frac{\partial n_{\mathrm{e}}}{\partial r}-\frac{1}{r^2}\frac{\partial }{\partial r}\left[r^{2}V\right]n_{\mathrm{e}}\hfill \\ \hfill & +\frac{\partial }{\partial \gamma }\left[\dot{\gamma }{n}_{\mathrm{e}}\right]+Q_{\mathrm{inj}},\hfill \end{array}$$(1)

where the processes of convection, diffusion, adiabatic loss, and radiative losses are taken into account. In the above equation, r represent the radial position, γ is the Lorentz factor of electrons, κ = κ(r, γ, t) is the particle diffusion coefficient, V = V(r) is the convection velocity, $\dot{\gamma }$ is the total energy loss rate of electrons, and Q_inj is the electrons injected at the termination shock. The number of particles per energy interval that flow through the termination shock should be equal to the number of particles injected from the central pulsar, which implies the inner boundary condition can be expressed as (e.g. Ng & Gleeson 1975; Vorster & Moraal 2013)

$$\begin{array}{c}\hfill V_0{n}_{\mathrm{e}}-\kappa (r_0(t),\gamma ,t)\frac{\partial n_{\mathrm{e}}}{\partial r}=\frac{Q_{\mathrm{inj}}(\gamma ,t)}{4\pi r_0^2(t)},\end{array}$$(2)

where r₀(t) is the radius of the termination shock at time t and V₀ represents the particle convection velocity at the shock.

According to de Jager et al. (2008), the particles in Vela X are assumed to be composed of two distinct electron populations, one at low energies responsible for the radio and GeV γ-ray emission and the second higher energy population for the X-ray and TeV emission. Thus, the distribution of injected electrons at the termination shock is given by

$$\begin{array}{c}\hfill Q_{\mathrm{inj}}(\gamma ,t)=\{\begin{array}{cc}{Q}_{0,\mathrm{R}}(t){\left(\frac{\gamma }{{\gamma }_{\max ,\mathrm{R}}}\right)}^{-{\alpha }_1}\hfill & {\gamma }_{\min ,\mathrm{R}}<\gamma <{\gamma }_{\max ,\mathrm{R}},\hfill \\ Q_{0,\mathrm{X}}(t){\left(\frac{\gamma }{{\gamma }_{\max ,\mathrm{X}}}\right)}^{-{\alpha }_2}\hfill & {\gamma }_{\min ,\mathrm{X}}<\gamma <{\gamma }_{\max ,\mathrm{X}},\hfill \end{array}\end{array}$$(3)

where Q_0, i(t) is the normalization constant and can be calculated by

$$\begin{array}{c}\hfill {\displaystyle {\int }_0^{\infty }{Q}_{0,i}(t){(\gamma /{\gamma }_{\max ,i})}^{-{\alpha }_i}\gamma m_{\mathrm{e}}{c}^2\mathrm{d}\gamma ={\eta }_{i}L(t).}\end{array}$$(4)

Here, η_i (i = R, X) is the fraction of spin-down luminosity L(t) of a pulsar converted to relativistic particle energy in each respective component. Assuming that the pulsar is a pure dipole radiator with a braking index of n, the time evolution of the spin-down luminosity of a pulsar can be given by

$$\begin{array}{c}\hfill L(t)=L_0{(1+\frac{t}{{\tau }_0})}^{-\frac{n+1}{n-1}},\end{array}$$(5)

where L₀ is the initial luminosity and τ₀ = 2τ_c/(n − 1)−T_age is the initial spin-down timescale, and T_age is the age of the pulsar. In Eq. (3), γ_max, i is the respective maximum energy of the two distinct electron populations, γ_min, i is the minimum energy of electrons, which is assumed to be 1.0 in this paper according to de Jager et al. (2008), and α₁ and α₂ = 2.0 are the corresponding spectral indices.

Following Kennel & Coroniti (1984a,b), the magnetic field in the PWN is assumed to be toroidal, and then the time and spatial evolutions of magnetic field in the nebula can be described as

$$\begin{array}{c}\hfill B(r,t)=B_0(t)\left[\frac{r_0(t)}{r\upsilon (r)}\right],\end{array}$$(6)

where υ(r) is the radial profiles of the convection velocity, and B₀(t) is the magnetic field at the termination shock that is obtained by $W_{\text{B}}(t)={\int }_{r_0(t)}^{R_{\text{pwn}}(t)}{B}^2(r,t)r^2\text{d}r/2$, where R_pwn(t) is the radius of the PWN. The function W_B(t) is calculated with dW_B(t)/dt = η_BL(t)−[W_B(t)/R_pwn(t)][dR_pwn(t)/dt] and η_B is the fraction of pulsar spin-down luminosity converted to magnetic field energy (e.g. Pacini & Salvati 1973). According to Lerche & Schlickeiser (1981), the evolution of the diffusion coefficient is given by

$$\begin{array}{c}\hfill \kappa (r,\gamma ,t)={\kappa }_0(t)\left(\frac{\gamma m_{\mathrm{e}}{c}^2}{1\mathrm{TeV}}\right){\left[\frac{r_0(t)}{r\upsilon (r)}\right]}^{-1},\end{array}$$(7)

and κ₀(t)∝1/B(t) is the diffusion at the termination shock.

In the model, the evolution of non-thermal photons arising from the radiation processes of relativistic electrons of the PWN is described as a photon conservation equation

$$\begin{array}{c}\hfill \frac{\partial n_{\gamma }}{\partial t}+\frac{n_{\gamma }}{t_{\mathrm{esc}}}=Q_{\gamma }^{\mathrm{syn}}+Q_{\gamma }^{\mathrm{ics}}+L_{\gamma }^{\mathrm{ssa}}+L_{\gamma }^{\gamma \gamma }.\end{array}$$(8)

In the above equation, n_γ = n_γ(r, x, t) is the number density of photons with the dimensionless photon energy x = E/m_ec², t_esc represents the photon escape timescale, $Q_{\gamma }^{\mathrm{syn}}$ and $Q_{\gamma }^{\mathrm{ics}}$ are the emissivity of synchrotron radiation and inverse Compton process, respectively, $L_{\gamma }^{\mathrm{ssa}}$ is the loss rate for photons caused by the process of synchrotron self-absorption, and $L_{\gamma }^{\gamma \gamma }$ is the lose rate of photon–photon pair production inside the PWN. A more detailed description and numerical calculation of the spatially dependent model are shown in Lu et al. (2017a); Lu et al. (2017b), and the radial variations of surface brightness and photon index at different bands are calculated according to Holler et al. (2012).

The dynamical evolution of the PWN is simulated according to the model presented by Bucciantini et al. (2011). In their model, the object is assumed to be a spherically symmetric system and the proper motion of the pulsar as well as the radiation losses of particles are negligible, and then the dynamical evolution of a PWN is divided into three phases: free-expansion, reverberation, and Sedov–Taylor phases. In the initial free-expansion phase, the radius of the PWN is given by

$$\begin{array}{c}\hfill R_{\mathrm{pwn}}(t)={\left(\frac{E_{\mathrm{sn}}^3{L}_0^2}{M_{\mathrm{ej}}^5}\right)}^{\frac{1}{10}}{\left(\frac{t}{1+\frac{t}{{\tau }_0}}\right)}^{\frac{6}{5}}\frac{1}{1-s}\sum _{i=0}^{\infty }{c}_i{s}^i,\end{array}$$(9)

where E_sn is the energy of the supernova explosion, M_ej is the mass of the supernova ejecta, s = (t/τ₀)/(1 + t/τ₀), and the coefficients c_i are given in Bucciantini et al. (2004). In the reverberation phase, the dynamical evolution radius is calculated with

$$\begin{array}{c}\hfill M_{\mathrm{sw}}{\ddot{R}}_{\mathrm{pwn}}(t)=4\pi R_{\mathrm{pwn}}^2(t)[P_{\mathrm{in}}(t)-P_{\mathrm{out}}(t)].\end{array}$$(10)

Here, M_sw is the mass of material swept up by the PWN when it expands inside the SNR, P_in is the pressure of the PWN shell, and P_out is the pressure of the SNR ejecta at the position of the PWN shell. In the Sedov–Taylor phase, the radius is determined by

$$\begin{array}{c}\hfill R_{\mathrm{pwn}}^4(t_{\mathrm{Sedov}})P_{\mathrm{in}}(t_{\mathrm{Sedov}})=R_{\mathrm{pwn}}^4(t)P(t),\end{array}$$(11)

where t_Sedov is the time at which the Sedov–Taylor phase begins, $P(t)={\rho }_{\text{ism}}{v}_{\text{fs}}^2/({\gamma }_{\text{snr}}+1)$ is the pressure in the SNR forward shock, where v_fs is the speed of the forward shock, and ρ_ism is the number density of the interstellar medium (ISM). Finally, the dynamical evolution radius of the termination shock is estimated by $r_0(t)=\sqrt{L(t)/4\pi c{P}_{\mathrm{ts}}(t)}$, where P_ts(t) is the pressure of the gas at the shock (Bucciantini et al. 2004), and the evolution of SNR are calculated with the semi-analytic model provided by Truelove & McKee (1999).

It should be noted that the SNR photon fields may prevent high energy γ-rays from leaving the PWN via γ + γ → e⁺ + e⁻ (e.g. Murase et al. 2015). In this paper, the γ-ray attenuation must also be taken into account and the photon fields of SNR are assumed to be a thermal spectrum for simplicity. According to Murase et al. (2015), the optical depth to pair production is approximated as

$$\begin{array}{c}\hfill {\tau }_{\gamma \gamma -\mathrm{sn}}=\frac{3}{16}{\sigma }_{\mathrm{T}}{n}_{\gamma }^{\mathrm{sn}}{R}_{\mathrm{ej}}(t)\mathcal{G}(y),\end{array}$$(12)

where R_ej(t) is the radius of the SNR ejecta, and $n_{\gamma }^{\text{sn}}=2\zeta (3){(k_{\text{B}}{T}_{\text{sn}})}^3/({\pi }^2{\hslash }^3{c}^3)$ is the photon density of SNR with a temperature of T_sn. Here k_B is the Boltzmann constant and σ_T represents the Thomson cross section. The function 𝒢(y) = ℱ(y)/ζ(3) with y = m_ec²/xk_BT_sn, in which ℱ(y) is given by Dermer et al. (2012) and ζ(3) = 1.202 is the Riemann’s zeta function. In order to consider the γ-ray attenuation by SNR photon fields, the photon spectrum is modified by multiplying by a factor of e^{−τ_γγ−sn}.

2.2. Results

Vela X is a dynamical system that shows a complex morphology at different observation bands (e.g. Aschenbach et al. 1995; Helfand et al. 2001; Aharonian et al. 2006; Grondin et al. 2013). It is generally believed that the PWN flow has been crushed by the SNR reverse shock, which leads to Vela PWN having a very complex observation morphology (e.g. Abdo et al. 2010; Hinton et al. 2011). With the parameters listed in Table 1, we calculated the dynamical radii of the Vela SNR system, including the radius of the SNR (R_SNR), reverse shock (R_rs), PWN (R_pwn), and termination shock (R_ts) in the framework of our model. The results are shown in Fig. 1, where the stars represent the current radius and termination shock radius of Vela PWN and imply that the nebula is in the phase of reverberation. From the dynamical modelling results, we obtain that the radius of Vela X has a value of R_pwn(T_age)=1.25° and the termination shock radius is R_ts(T_age)=37″ at the age of T_age = 11.3 kyr, which is consistent with the values R_pwn(T_age)=1.2° shown in Abramowski et al. (2012) and R_ts(T_age)=25″ derived by Ng & Romani (2004).

Fig. 1. Dynamical evolution radius of SNR (Rsnr, purple line), reverse shock (Rrs, green line), PWN (Rpwn, black line) and the termination shock (Rts, blue line). The stars represent the current radius and terminate shock radius of Vela PWN at Tage = 11.3 kyr.

Since the photon SEDs of Vela X at different energy bands of different emission regions and the radial variations of photon index in the X-ray band as well as the surface brightness in the TeV γ-rays have been observed, all of the observations should be reproduced by our model to obtain the best-fit parameters. In the calculation, we first reproduced the photon SEDs of different emission regions roughly, and then adjusted the parameters of the initial diffusion coefficient and magnetic fraction η_B to reproduce the radial variations of the photon index and flux in the X-ray band as well as surface brightness in the TeV γ-rays, which due to the radial profiles are mainly determined by the competition of synchrotron radiation loss and particle diffusion (e.g. Vorster & Moraal 2013), and then adjusted the fitted parameters again to reproduce the photon SEDs. We repeated the above steps and we obtained the parameters of Vela X, which are listed in Table 1.

The photon SEDs of Vela X calculated with our model are shown in Fig. 2, where the left panel is the comparison of our model results of the PWN with the observed data from the radio to γ-ray bands, and the right panel is the comparison of model results with the data observed in the γ-ray band. The soft photon fields, which include three components in the calculation, are as follows: the cosmic microwave background (CMB, temperature 2.73 K, photon density 0.26 eV cm⁻³), the galactic far-infrared background (FIR, temperature 25 K, photon density 12.0 eV cm⁻³), and the near-infrared and optical photon field (NIR) due to the stars (temperature 6500 K, photon density 0.4 eV cm⁻³). The χ² value of the model results with the observed data in the GeV band is 2.2 and in the TeV for the inner emission region of r = 48 arcmin we obtain 2.5. The luminosity of Vela X at different bands of the different emission regions respectively derived from the model and the observations are shown in Table 2, revealing that the luminosities at different bands obtained by our model are consistent with the results of observations with the exception of the TeV γ-ray luminosity. The TeV γ-ray luminosity obtained by our model is less than the result estimated by Aharonian et al. (2006), but the integral flux above 1 TeV derived from our model of 1.9 × 10⁻¹¹ cm⁻² s⁻¹ is consistent with the value 2.1 × 10⁻¹¹ cm⁻² s⁻¹ obtained by Abramowski et al. (2012) according to the deeper observation of H.E.S.S.

Fig. 2. Comparison of modelling results of Vela X with the observed data at different emission regions of r = 6 arcmin (dark-yellow line), r = 15 arcmin (red line), and r = 48 arcmin (blue line), as well as the total flux (black line). Left panel: modelling SEDs from radio to γ-ray bands. Right panel: modelling SEDs in the γ-ray band. The dash-dotted lines are the unabsorbed spectra, and the solid lines consider γ-ray attenuation by the SNR photon fields. The observed data are taken from Abdo et al. (2010) for radio band (black, full circles), Mattana et al. (2011) and Mangano et al. (2005) for X-rays of the inner regions of r = 6 arcmin (dark-yellow, solid thick line) and r = 15 arcmin (red, solid thick line), Grondin et al. (2013), Ajello et al. (2017), and Tibaldo et al. (2018) for GeV γ-rays (black, open circles), and Abramowski et al. (2012) for TeV γ-rays emission within r = 48 arcmin (blue, full circles), respectively.

As seen in Fig. 2, the observed SEDs of Vela X, including X-ray photon spectra of the inner emission regions and γ-ray spectra in the GeV and TeV bands, can be well reproduced in the framework of our model with the parameters shown in Table 1, while the radio photon spectrum cannot. It should be pointed out that we did not find any set of appropriate parameters to reproduce the SED of Vela X in the radio band based on the spherically symmetric spatially dependent model. From the modelling results, we can also see that the process of γ-ray attenuation by the SNR photon field is important for modifying the VHE γ-ray spectrm of Vela X (see the right panel of Fig. 2), and a temperature of T_sn = 42 K of Vela SNR photon field is obtained. We note that such a photon field around the nebula has not been detected by the current observations of Vela X, but this value is interestingly consistent with the supernova (SN) formed dust temperature for the young PWNe of 50 − 60 K estimated by Temim et al. (2012).

The photon field produced by synchrotron emission of the PWN also may prevent high energy photons from leaving the nebula. In order to study the effect of photon–photon annihilation of the synchrotron photon inside the PWN on the photon SEDs, we calculated the optical depths to the pair production with the parameters listed in Table 1 and the results are given in Fig. 3, where the optical depth to the SNR photon field is calculated with Eq. (12) and the optical depth to the synchrotron photon field of the PWN are determined by (e.g. Aharonian et al. 2008)

$$\begin{array}{c}\hfill {\tau }_{\gamma \gamma -\mathrm{pwn}}={\displaystyle {\int }_{r_0}^{R_{\mathrm{pwn}}}{\int }_0^{\infty }{\sigma }_{\gamma \gamma }{n}_{\gamma }(r,x,t)\mathrm{d}x\mathrm{d}r,}\end{array}$$(13)

Fig. 3. Present-day optical depths to pair production of synchrotron photon field of the PWN (red line) and SNR photon field (black line), with the parameters listed in Table 1.

where n_γ(r, x, t) is the number density of non-thermal photons and σ_γγ represents the photon–photon annihilation cross section. It can be seen from Fig. 3 that the process of photon–photon annihilation of synchrotron photons inside the PWN can be neglected due to the optical depth τ_γγ−pwn  ≪  τ_γγ−sn, and that the VHE γ-ray spectrum of Vela X is mainly modified by the photon–photon annihilation process of the SNR photon field.

To analyse the effect of the processes of synchrotron self-absorption and pair production on the photon spectrum, we also calculated the timescales of photon escape, synchrotron self-absorption, and pair production inside the PWN. The results are shown in Fig. 4. In the calculation, the photon escape timescale is calculated with t_esc = R_pwn(t)/c, the timescale of photon–photon pair production t_γγ inside the PWN is calculated by Eq. (28) of Lu et al. (2017a), and the synchrotron self-absorption timescale is given by (e.g. Mastichiadis & Kirk 1995)

$$\begin{array}{c}\hfill t_{\mathrm{ssa}}={[-{\sigma }_{\mathrm{T}}c\frac{\pi }{6{\alpha }_f}{x}^{-1/2}{b}^{-1/2}{\left[\frac{\partial }{\partial \gamma }\left(\frac{n_{\mathrm{e}}}{{\gamma }^2}\right)\right]}_{\gamma ={(x/b)}^{1/2}}]}^{-1},\end{array}$$(14)

Fig. 4. Current timescales of photon escape (black lines), pair production (upper panel), and synchrotron self-absorption (bottom panel) inside the PWN at different radial distances of r = 3.9 arcmin (blue lines) and r = 66.8 arcmin (red lines).

where α_f is the fine-structure constant and b = B(r, t)/B_cr with the critical magnetic field B_cr. The results reveal that although the processes of synchrotron self-absorption and pair production inside the PWN are considered in the model, they can be neglected because of the timescales of t_γγ  ≪  t_esc and t_ssa  ≪  t_esc in all of the emission regions. Therefore, the processes of synchrotron self-absorption and pair production inside the PWN only have little effect on the photon spectrum of Vela X.

With the parameters listed in Table 1, the high energy radiative spectra of Vela X are well reproduced at the age of 11.3 kyr. At this age, the present-day spin-down luminosity of the Vela pulsar is estimated to be 6.9 × 10³⁶ erg s⁻¹, which is consistent with the value given by Taylor et al. (1993). According to the model results, a fraction of pulsar spin-down power converted to magnetic energy is derived to be η_B = 0.0087, which implies that the current magnetic field of Vela X at the termination shock has a value of B₀(T_age)=11.8 μG and a spatial average of $\overline{B}(T_{\mathrm{age}})=8.7\mu \mathrm{G}$. The averaged magnetic field strength in our model is larger than the values of 5 μG obtained by de Jager et al. (2008) with a spatially independent model and 5.8 μG given by Porth et al. (2016) with a steady-state particle transport model. Moreover, with an initial diffusion coefficient of 9.0 × 10¹⁸ cm² s⁻¹, the current diffusion coefficient of Vela X at the termination shock is estimated to be κ₀(T_age)=6.13 × 10²³ cm² s⁻¹ and the spatial averaged diffusion coefficient has a value of $\overline{\kappa }(T_{\mathrm{age}})=1.18\times {10}^{24}{\mathrm{cm}}^2{\mathrm{s}}^{-1}$ with an electron energy of 1 TeV. These results are smaller than the values ${\kappa }_0(T_{\mathrm{age}})=\overline{\kappa }(T_{\mathrm{age}})=1.4\times {10}^{26}{\mathrm{cm}}^2{\mathrm{s}}^{-1}$ obtained by Porth et al. (2016), because the particle diffusion is energy independent in their model while it is energy dependent in our model.

The X-ray observations revealed that the photon index and flux density of Vela X in the 3 − 10 keV energy range increase with the increase of radial distance (Mangano et al. 2005). Thus, we calculated the variations of photon index and flux density with radial distance for the nebula according to the parameters listed in Table 1 and the results are shown in Fig. 5, where the fit residual is calculated with residual=(observation-model)/model. In the calculation, the observed data in the energy range of 3 − 10 keV of the emission region within r = 15 arcmin of the PWN are obtained from Mangano et al. (2005). It can be seen from Fig. 5 that the radial variations of photon index and flux density in this energy range can be well reproduced in the framework of our model when the emission regions are farther away from the central pulsar. While in the vicinity of the pulsar, the observed data cannot be reproduced well with our model, especially for the radial variations of the flux density, which have significant deviations between the data and the model. From the right panel of Fig. 5, it is clear that the flux density calculated with our model is less than the observed flux density of the inner emission regions. On the one hand, the particle distribution in the nebula is assumed to be isotropic and a toroidal magnetic field is adopted in our model, while the particle and magnetic field in the PWN may have more complicated distributions (e.g. Del Zanna et al. 2004; Porth et al. 2014). On the other hand, it could be that the observed data of the emission regions in the vicinity of the pulsar include photons that emit from the central pulsar, which leads to higher values of the observed flux density. The modelling results also show that the photon index and flux density of Vela X in the X-ray band have more significantly radial variations in the vicinity of the Vela pulsar, because the radial variations of magnetic field and diffusion in the inner regions are more obvious (Kennel & Coroniti 1984a; Lu et al. 2017b), and the SED shape in the X-ray band is determined by the competition between the radiation losses and particle diffusion (Vorster & Moraal 2013).

Fig. 5. Comparisons of radial variations of photon index (left panel) and flux density (right panel) in the energy range of 3 − 10 keV of the Vela X emission within r = 15 arcmin. The observed data are taken from Mangano et al. (2005).

Based on deeper observations of the entire Vela X nebula region with H.E.S.S., Abramowski et al. (2012) have shown that the energy averaged surface brightness of the γ-ray emission of the nebula in the 0.75 − 70 TeV energy range has obvious changes as the radial distance increases. We calculated the differential γ-ray spectra of the inner region of r = 48 arcmin and ring extension of 48 <  r <  72 arcmin of Vela X, as well as the averaged surface brightness in the 0.75 − 70 TeV energy range with the same model parameters. The results are shown in Fig. 6 and the observed data are taken from Abramowski et al. (2012). From the left panel of Fig. 6, it is clear that the γ-ray spectra of the inner region and ring extension of Vela X can be well reproduced in our model. The comparison of radial variations of the energy averaged surface brightness of the PWN in the 0.75 − 70 TeV energy range are shown in the right panel of Fig. 6, indicating that our model can also reproduce the decreased surface brightness. In our model, the nebula is assumed to be a spherical symmetric system and the particle distribution inside the nebula is considered to be isotropic, while the TeV γ-ray observations of H.E.S.S. revealed that Vela X is an asymmetric PWN, which leads the radial variations of surface brightness in the TeV band to have some deviations between the observations and the model.

Fig. 6. Left panel: comparisons of differential γ-ray spectrum of the inner region of r = 48 arcmin (red line) and ring extension of 48 <  r <  72 arcmin (black line) of Vela X in the TeV band. Right panel: radial variations of surface brightness in the 0.75 − 70 TeV energy range. The observed data are taken from Abramowski et al. (2012).

The central position of VHE γ-ray emission of Vela X is not located at the location of the Vela pulsar because the view fields in the γ-rays revealed that the significant emissions are mainly centred on the south of the pulsar (Aharonian et al. 2006; Abramowski et al. 2012), which is different from that shown in the X-ray band. The observations have shown that the central position of X-ray emission is located at the location of the central pulsar (Mangano et al. 2005). We also note that the central position of GeV γ-ray emission is also not located at the location of the Vela pulsar but located at the south of the pulsar, and that it is not the same as the central position of the TeV γ-ray emission because the observed image in the GeV band is different from that in the TeV band (e.g. Abdo et al. 2010; Tibaldo et al. 2018). In our model, the particles are injected from the central pulsar, that is, the central position of photon emission is located at the location of the central pulsar. Therefore, the TeV γ-ray spectra of the inner region and ring extension of Vela X, as well as the radial variations of surface brightness in 0.75 − 70 TeV energy range, can be reproduced by our model, but they are somewhat different from the observations. In order to investigate the spatial radiative properties of Vela X more accurately at different bands, a non-spherically symmetric spatially dependent model that involves the proper motion of the pulsar should be presented, though it is outside of the scope of our paper.

Abramowski et al. (2012) have pointed out that the energy spectrum of Vela X in the 0.75 − 70 TeV energy range can be well described by an exponentially cut-off power-law function, dN/dE ∝ (E/1 TeV)^−Γe^−E/E_cut, where E is the energy of photons, Γ is the power-law photon index, and E_cut is the cut-off energy. In order to investigate the radial variations of the spectrum in the TeV band, we also calculated the photon index and cut-off energy of Vela X in the 0.75 − 70 TeV energy range with an exponentially cut-off power-law function, and the results are shown in Fig. 7. It can be seen from the left panel of Fig. 7 that the power-law photon index in the 0.75 − 70 TeV energy range remains almost unchanged as the radial distance increases. The radial variations of the cut-off energy are shown in the right panel of Fig. 7, showing that the cut-off energy decreases with the increase of radial distance in the inner regions, while there are only slight changes in the outer regions. In fact, the shape of the photon spectrum in the TeV band is determined by the competition between the radiation losses and particle diffusion (Vorster & Moraal 2013; Lu et al. 2016). In our model, the magnetic field decreases while the diffusion coefficient increases with the increasing time (see Eqs. (6) and (7)). For Vela X with an age of 11.3 kyr, the current competition between the radiation losses and particle diffusion has transitioned to a diffusion-dominated regime. Particle diffusion is important for reducing the amount of radiation losses (Vorster & Moraal 2013), which leads to the radial variations of photon index and cut-off energy in the TeV energy range as shown in Fig. 7.

Fig. 7. Radial variations of photon index (left panel) and cut-off energy (right panel) of Vela X in the 0.75 − 70 TeV energy range.

Finally, we calculated the radial variations of the power-law photon index and surface brightness of Vela X in the 0.2 − 20 GeV energy range, and the results are shown in Fig. 8. Abdo et al. (2010) have pointed out that the total GeV γ-ray spectrum in the energy range of 0.2 − 20 GeV can be well described by a power-law function with a photon index of 2.41 ± 0.09 ± 0.15. Thus, the radial photon index of our model in the GeV energy range is also calculated with a power-law function. From the left panel of Fig. 8, it can be seen that the photon index of the inner regions of Vela X in the GeV band has a significant increase as the radial distance increases, while there is only a slight changes in the outer regions. In the calculation, the injected particles of Vela X are composed of two distinct electron populations (see Eq. (3)) and GeV γ-rays are produced by the lower energy population, indicating that the power-law photon index in the 0.2 − 20 GeV energy range is also determined by the competition between the radiation losses and particle diffusion. In our model, the particle diffusion is energy dependent and the lower energy particles have a slow diffusion, which give rise to the radial variations of the photon index of Vela X in the GeV band being mainly determined by the radiation losses and leading to the variations of photon index in the 0.2 − 20 GeV energy range as shown in the left panel of Fig. 8. We would like to point out that the radial variations of the photon index of Vela X in the GeV band is different from the results of MSH 15-52, G21.5-0.9, and 3C58 obtained by Lu et al. (2017a,b). This is because the injected particles of Vela X are assumed to be composed of two distinct electron populations, while a broken power-law form is assumed for MSH 15-52, G21.5-0.9, and 3C58 in Lu et al. (2017a,b). The radial variations of surface brightness in the 0.2 − 20 GeV energy range are given in the right panel of Fig. 8. Similar to the radial variations of energy averaged surface brightness in the TeV band, the surface brightness of Vela X in the GeV energy range also decreases with the increase of radial distance. It should be pointed out that although the Fermi-LAT observations have revealed that the photon spectra of Vela X in the southern and northern wings have some differences (Grondin et al. 2013), the spatial variations of photon index and surface brightness at the GeV band have not been reported in the literature. Therefore, the spatial variations of Vela X in the GeV band predicted by our model should be tested in the further observations.

Fig. 8. Radial variations of power-law photon index (left panel) and surface brightness (right panel) of Vela X in the energy range of 0.2 − 20 GeV.

3. Summary and discussion

In this paper, non-thermal photon emission and high energy spatially radiative properties of Vela X have been investigated in the framework of a spatially dependent particle transport model. The main difference between the model presented here and that given by Lu et al. (2017a) is that the injected particles are assumed to consist of two distinct electron populations and the γ-ray attenuation by the SNR photon fields is taken into account.

With the parameters listed in Table 1, the photon SEDs of Vela X observed at X-ray, GeV, and TeV γ-ray bands can be reproduced well in this model (see Fig. 2), and the process of γ-ray attenuation by the SNR photon fields is important for correcting the VHE γ-ray spectrum of Vela X (see right panel of Fig. 2), and the temperature of Vela SNR photon field is estimated to be 42 K.

The spatial variations of photon SEDs of Vela X also have been investigated. From our calculations, on the one hand, the radial variations of photon index and flux density of the emission region within r = 15 arcmin of Vela X in the 3 − 10 keV energy range are well reproduced in the model, with the exception of the data in the regions of the vicinity of the central pulsar (see Fig. 5). The difference between model results and the observed data may be that the observed data of the vicinity of the pulsar include the photons whose emission from the central pulsar and the particle distribution in the nebula is assumed to be isotropic and a toroidal magnetic filed is adopted in our model, while the particle and magnetic field in the PWN probably have more complicated distributions. On the other hand, the differential γ-ray spectra of the inner region of r = 48′ and ring extension of 48′< r <  72′ as well as the radial variations of averaged surface brightness of Vela X in the 0.75 − 70 TeV energy range can be reproduced well in the model (see Fig. 6).

The radial variations of the γ-ray spectrum in the 0.75 − 70 TeV energy range have been studied with an exponentially cut-off power-law function, indicating that the photon index of Vela X in the TeV band only has a slight change as the radial distance increases, and the cut-off energy decreases with the increase of radial distance in the inner regions but only has a slight change in the outer regions (see Fig. 7). Furthermore, the radial variations of the power-law photon index and surface brightness of Vela X in the 0.2 − 20 GeV energy range are also predicted, and the results show that the power-law photon index increases while the surface brightness decreases with the increase of radial distance (see Fig. 8). We note that the increased photon index of Vela X in the GeV band is different from the results of MSH 15-52, G21.5-0.9, and 3C58 obtained by Lu et al. (2017a); Lu et al. (2017b), because the injected particles of Vela X are assumed to be composed of two distinct electron populations and GeV γ-rays are produced by the lower energy population.

It should be pointed out that the particles in the model are injected from the central pulsar, which leads to the central position of photon emission of our model being located at the location of the pulsar. The observations showed that the central positions of GeV and TeV γ-ray emissions of Vela X are located at the south of the Vela pulsar (Aharonian et al. 2006; Abdo et al. 2010; Abramowski et al. 2012; Grondin et al. 2013; Tibaldo et al. 2018), and the central positions shown in the GeV and TeV bands are different. Therefore, although our model can roughly account for the high energy spatially radiative properties of Vela X, a non-spherically symmetric spatially dependent model, which includes the proper motion of the pulsar, should be presented to investigate the radiative properties of the PWN more precisely.

Acknowledgments

We are very grateful to the anonymous referee for his or her very constructive comments. This work is partially supported by the National Natural Science Foundation of China (NSFC 11433004, U1738211, 11173020, 11803027), Top Talents Program of Yunnan Province, and Yunnan local colleges applied basic research projects (2017FH001-17).

References

All Tables

All Figures

	Fig. 1. Dynamical evolution radius of SNR (Rsnr, purple line), reverse shock (Rrs, green line), PWN (Rpwn, black line) and the termination shock (Rts, blue line). The stars represent the current radius and terminate shock radius of Vela PWN at Tage = 11.3 kyr.
In the text

Fig. 2. Comparison of modelling results of Vela X with the observed data at different emission regions of r = 6 arcmin (dark-yellow line), r = 15 arcmin (red line), and r = 48 arcmin (blue line), as well as the total flux (black line). Left panel: modelling SEDs from radio to γ-ray bands. Right panel: modelling SEDs in the γ-ray band. The dash-dotted lines are the unabsorbed spectra, and the solid lines consider γ-ray attenuation by the SNR photon fields. The observed data are taken from Abdo et al. (2010) for radio band (black, full circles), Mattana et al. (2011) and Mangano et al. (2005) for X-rays of the inner regions of r = 6 arcmin (dark-yellow, solid thick line) and r = 15 arcmin (red, solid thick line), Grondin et al. (2013), Ajello et al. (2017), and Tibaldo et al. (2018) for GeV γ-rays (black, open circles), and Abramowski et al. (2012) for TeV γ-rays emission within r = 48 arcmin (blue, full circles), respectively.

In the text

	Fig. 3. Present-day optical depths to pair production of synchrotron photon field of the PWN (red line) and SNR photon field (black line), with the parameters listed in Table 1.
In the text

	Fig. 4. Current timescales of photon escape (black lines), pair production (upper panel), and synchrotron self-absorption (bottom panel) inside the PWN at different radial distances of r = 3.9 arcmin (blue lines) and r = 66.8 arcmin (red lines).
In the text

	Fig. 5. Comparisons of radial variations of photon index (left panel) and flux density (right panel) in the energy range of 3 − 10 keV of the Vela X emission within r = 15 arcmin. The observed data are taken from Mangano et al. (2005).
In the text

	Fig. 6. Left panel: comparisons of differential γ-ray spectrum of the inner region of r = 48 arcmin (red line) and ring extension of 48 <  r <  72 arcmin (black line) of Vela X in the TeV band. Right panel: radial variations of surface brightness in the 0.75 − 70 TeV energy range. The observed data are taken from Abramowski et al. (2012).
In the text

	Fig. 7. Radial variations of photon index (left panel) and cut-off energy (right panel) of Vela X in the 0.75 − 70 TeV energy range.
In the text

	Fig. 8. Radial variations of power-law photon index (left panel) and surface brightness (right panel) of Vela X in the energy range of 0.2 − 20 GeV.
In the text