Rotational evolution of the Crab pulsar in the wind braking model

Abstract

1 INTRODUCTION

The Crab pulsar (PSR B0531+21) is a young radio pulsar with a spin frequency ν = 30.2 Hz and frequency derivative |$\dot{\nu }=-3.86 \times 10^{-10} \, \rm Hz\,s^{-1}$| (Lyne, Pritchard & Smith 1993). Its characteristic magnetic field is about 7.6 × 10¹² G at the magnetic poles.¹ The Crab pulsar has been monitored continuously for decades by various telescopes since it was discovered in 1968 (Lyne et al. 1993, 2015; Wang et al. 2012). Its braking index is n = 2.51 ± 0.01 (Lyne et al. 1993). Different braking index values 2.45, 2.57 and 2.3 are also reported between glitches (Wang et al. 2012; Lyne et al. 2015). Observational details are given in Table 1.

Many mechanisms have been proposed in order to explain the braking index observations, e.g. the pulsar wind model (Xu & Qiao 2001; Wu, Xu & Gil 2003; Contopoulos & Spitkovsky 2006; Yue, Xu & Zhu 2007), a changing magnetic field strength (Lin & Zhang 2004; Chen & Li 2006; Espinoza et al. 2011), a changing inclination angle (Lyne et al. 2013), and additional torques due to accretion (Liu et al. 2014). However, these models are more or less not consistent with observations or cannot simulate the long-term evolution of pulsars. And these models cannot explain the different braking indices detected between glitches (Wang et al. 2012; Lyne et al. 2015). Furthermore, the effect of pulsar death (death line, Ruderman & Sutherland 1975; or death valley, Chen & Ruderman 1993; Zhang, Harding & Muslimov 2000) should be considered in modelling the long-term rotational evolution of the pulsar (Contopoulos & Spitkovsky 2006).

The pulsar wind model considers both the pulsar spin-down and the particle acceleration in the magnetosphere (Xu & Qiao 2001; Wu et al. 2003). In this paper, an updated pulsar wind model is built based on previous researches (Xu & Qiao 2001; Contopoulos & Spitkovsky 2006; Yue, Xu & Zhu 2007; Li, Spitkovsky & Tchekhovskoy 2012). It includes: (1) both magnetic dipole radiation and particle outflow, and their dependence on the inclination angle; (2) the particle outflow depends on the specific acceleration model; (3) the primary particle density may be much larger than the Goldreich–Julian charge density; (4) the effect of pulsar death is considered in modelling the rotational evolution of the pulsar. This model can calculate both the short-term and long-term evolution of pulsars. It is applied to the Crab pulsar which has the most detailed timing observations. Possible applications to other sources are also illustrated.

The pulsar wind model and model calculations of the Crab pulsar are presented in Section 2. Discussions and conclusions are given in Sections 3 and 4, respectively.

2 SPIN-DOWN OF THE CRAB PULSAR IN THE PULSAR WIND MODEL

2.1 Description of the pulsar wind model

2.2 Understanding the braking index of the Crab pulsar

Parameters of the Crab pulsar are calculated in this section. The VG(CR) model is taken as an example to show the calculation process. For pulsars with the observed ν, |$\dot{\nu }$|⁠, |$\ddot{\nu }$| and |$\stackrel{...}{\nu }$|⁠, their observational braking index n and second braking index m can be obtained by equations (2) and (3). Pulsar parameters such as magnetic field, inclination angle, and particle density can be calculated by these three equations (9), (10), and (11). The primary particle density may be 10³–10⁴ times of Goldreich–Julian charge density (Yue et al. 2007). According to the observational range of inclination angle and the characteristic magnetic field (which can be viewed as order of magnitude estimation of the true magnetic field), the range of particle density κ can be limited by equations (9) and (10). In our calculation, κ = 10³ is adopted,⁴ the inclination angle α and magnetic field B are about 55° and 8.1 × 10¹² G which are consistent with observational constraints. The second braking index calculated by equation (11) is 10.96, which is roughly consistent with the observed value 10.15 (Lyne et al. 1993), considering the observational uncertainties (Lyne et al. 1993, 2015).

The magnetic field and inclination angle are assumed to be constant in the pulsar wind model. The particle density decreases in long-term evolution of pulsars but remains unchanged in short duration at early time. Fig. 1 shows the braking index of the Crab pulsar as a function of spin period. As the pulsar period increases, the braking index decreases from 3 to 1. It illustrates that the Crab pulsar is initially braked by magnetic dipole (n → 3) and then by the particle wind (n → 1). We take n = 2 as the transition point (see Fig. 5, at n = 2, the slop of the curve is 0) and the corresponding spin period is about 55 ms. Since Table 1 has shown the timing observations of the Crab pulsar in 1969 (its age is 915 yr old by then), its initial period is P₀ ≈ 18.8 ms in AD 1054 by integrating equation (9). The calculations of all acceleration models are shown in Table 3. As shown in this table, the particle densities are different in different models, but are all more than 100 times of ρ_GJ. The inclination angles is within the observational range. The initial periods are all about 19 ms.

Figure 1.

The braking index of the Crab pulsar as function of spin period in the VG(CR) model. The dashed line is observational braking index of 2.51. The dotted line is braking index in the transition point (n = 2).

Open in new tabDownload slide

2.3 An increasing particle density results in a lower braking index during glitches

Fig. 2 shows the braking index as function of τ_κ for the Crab pulsar in the VG(CR) model. In order to better understand the observations of n < 2.51 (Lyne et al. 2015), we consider the τ_κ > 0 only. As we can see in this figure, the braking index is insensitive to τ_κ when it is larger than 10⁴ yr. But, when it is comparable with the characteristic age (about 10³ yr), the braking index decreases sharply. When braking index is 2.3, the changing rate of the particle density is |$\dot{\kappa }\approx 1.2\times 10^{-8} \, \rm s^{-1}$|⁠. In this interval (from MJD 51000 to MJD 53000, about 5 yr), the particle density has changed by 2 ρ_GJ. The change of particle density will result in changes of radiation energy loss rate and period derivative. Their changes, respectively, are |$\delta \dot{E}\approx 2\times 10^{35} \,\rm erg\,s^{-1}$| and |$\delta \dot{P}\approx 1.88\times 10^{-16} \, \rm s\,s^{-1}$|⁠.

Figure 2.

The braking index of the Crab pulsar as function of τ_κ in the VG(CR) model. The dashed line is the underlying braking index 2.51. The dotted line is the smaller braking index 2.3 measured during glitch activities (Lyne et al. 2015, fig. 7 there).

Open in new tabDownload slide

The pulsar wind model can be tested with observations (Jodrell Bank monthly ephemeris of the Crab pulsar⁶). The effect of glitches is not taken into consideration in the pulsar wind model. The changes of period and period derivative caused by glitches should be added (Lyne et al. 2015). Fig. 3 is the model calculation compared with observations (without modelling the transient variation caused by glitches). When the amplitude of glitches is relatively small, the change of particle density is not obvious, then the model calculations (the blue points) fit well with observations (the black points). However, when the amplitude of glitches is large, the effect of glitches must be added (the red points). Moreover, the effect of particle density change should be taken into consideration. A factor |$\delta \dot{P}=1.88\times 10^{-16}\, \rm s\,s^{-1}$| (which is calculated above) is added to the red points (from MJD 51804.75 afterwards), shown as the green points. The green points can match the general trend of the observations.

Figure 3.

Evolution of the Crab pulsar in |$P-\dot{P}$| diagram in VG(CR) model compared with observations. The black points are timing results from Jodrell Bank monthly ephemeris. The blue points are the rotational evolution of Crab pulsar in the VG(CR) model. The red ones are obtained by adding the effect of glitches (Lyne et al. 2015, table 3 there). And the green ones are the summation of red points and a factor |$\delta \dot{P}\approx 1.88\times 10^{-16} \, \rm s\,s^{-1}$| (from MJD 51804.75 afterwards), which may be caused by the change of outflow particle density. At first, the black, blue, red, and green points are coincide with each other. (Notes: the fitting does not model the transient process caused by glitches.)

Open in new tabDownload slide

2.4 Long-term evolution of the Crab pulsar

The pulsar period and its fist derivative evolution with time in the VG(CR) model are shown in the first and second panel of Fig. 4. The curves evolve slowly in its early age and then rise sharply, indicating that the Crab pulsar is mainly braked by magnetic dipole radiation first and then mainly by particle wind. The transition period is 55 ms and age of the Crab pulsar will be 2771 calculated in the VG(CR) model (see Table 3). The braking index evolution with time in the VG(CR) model is shown in the third panel of Fig. 4. Braking index decreases from 3 to 6/7 – different minimum braking indices for different acceleration models. The bottom panel of Fig. 4 shows the second braking index evolution with time. The curve gradually decreases from 15 to 30/49. Fig. 5 shows the evolution of the Crab pulsar in |$P-\dot{P}$| diagram from birth to its 30 000 yr old. These points are, respectively, taken when the Crab pulsar is: 1, 915, 2771, 9000 and 30 000 yr old. Period derivative evolves with period by a function of |$\dot{P} \propto P^{2-n}$| (Espinoza et al. 2011). In the log–log plot, the evolution curve evolves with a slop (2 − n). As the braking index falling from 3 to 1, the pulsar evolve from magnetic dipole radiation dominated case (the left points) to the particle wind-dominated case (the right points) and the bottom is the transition point (n = 2). Clearly, the curve evolves more sharply after the transition point. The asymptotic behaviours are discussed in the Appendix B.

Figure 4.

Rotational evolution of the Crab pulsar from birth to 30 000 yr in the VG(CR) model. The first and second panels are, respectively, the period and period derivative evolution with time. The points are observations of period and its derivative (Lyne et al. 1993). The third and fourth panels are, respectively, the braking index and second braking index evolution with time. The dashed lines are observed values and the solid lines are, respectively, the minimum braking index and second braking index in the VG(CR) model.

Open in new tabDownload slide

Figure 5.

Evolution of the Crab pulsar in the P−|$\dot{P}$| diagram. Five points are chosen when the pulsar is: 1, 915, 2771, 9000 and 30 000 yr old.

Open in new tabDownload slide

2.5 The effect of pulsar death

As shown is the Fig. 5, the line evolves up right with a slope about 1 after the transition point just like PSR J1734−3333 whose braking index is 0.9 ± 0.2 (Espinoza et al. 2011). If the pulsar continually evolve in this trend, it would not be dead and might end itself as a magnetar. It will be hard to explain the observations for: (1) the number of magnetars are smaller than radio pulsars and most of the pulsars occupy the central region in the P−|$\dot{P}$| diagram (see Fig. 6); (2) only a small portion of magnetars are radio emitters (Olausen & Kaspi 2014); (3) the properties of the Crab pulsar, PSR J1734−3333 and other high magnetic field pulsars are significantly different from magnetars (Ng & Kaspi 2011). As pulsars slow-down, the density of outflow particles may decrease and the effect of particle wind will recede [for the SCLF cases, the number of particles (or κ) may close to constant, but the particle density still decreases because a reducing Goldreich–Julian charge density]. Therefore, the effect of pulsar death must be included in modelling the long-term rotational evolution of pulsars.

Figure 6.

Rotational evolution of the Crab pulsar in the VG(CR) model. The dashed curve is the evolution without considering pulsar death. The solid line is the one considering pulsar death. The |$P-\dot{P}$| diagram has listed all known radio pulsars, magnetars, and millisecond pulsars (updated from fig. 1 in Tong & Wang 2014). The dot–dashed line represents the fiducial death line (Ruderman & Sutherland 1975). Pulsars with meaningful braking indices (Espinoza et al. 2011; Lyne et al. 2015) are defined by large square. The arrows represent their evolution directions and each arrow indicate its motion in the next 10 000 yr (For the Vela pulsar, 20 000 yr is used; 5000 yr for PSR J1846−0258 and PSR J1119−6127).

Open in new tabDownload slide

The essential condition of radio emission for a pulsar is its pair production (Sturrock 1971; Ruderman & Sutherland 1975). The so called pulsar ‘death’ means the stopping of pulsar emission (Ruderman & Sutherland 1975; Chen & Ruderman 1993; Zhang et al. 2000; Contopoulos & Spitkovsky 2006). The death defined by Ruderman & Sutherland (1975) is that a pulsar dead when the maximum potential drop ΔΦ available from the pulsar is smaller than the acceleration potential needed to accelerate particles. Because this death line is model-dependent, Chen & Ruderman (1993) proposed death ‘valley’ by modifying the boundary conditions. And this death ‘valley’ was updated by considering different particle acceleration and photon radiation processes (Zhang et al. 2000). However these works are just to separate pulsars with radio emission from those without. Contopoulos & Spitkovsky (2006) included the effect of pulsar death when modelling the rotational evolution of pulsars. However, the braking indices there are always larger than 3. Considering particle acceleration and pulsar death simultaneously may give a comprehensive interpretation for both short- and long-term rotational evolution of pulsars.

The visual explanation for pulsar death in the pulsar wind model is that the particles in the magnetosphere are exhausted and there are no primary particles generated. For the physical understanding, we can refer to the equation (8), expressions ΔΦ ∝ Ω² and Δϕ (i.e. |$\Delta {\phi }_{\rm VGCR}=2.8\times 10^{13}R_{\rm 6}^{2/7}B_{12}^{-1/7}\Omega ^{-1/7} \rm \,V$|⁠; Xu & Qiao 2001), with the pulsar spinning down, the maximum potential of the pulsar (ΔΦ) drops, when the maximum potential cannot meet the need to accelerate particles (Δϕ > ΔΦ), the pulsar death. Besides, the acceleration potential is weakly dependent on the Ω, that is why we can make calculations under the constant acceleration potential (the CAP model; Yue et al. 2007). Fig. 6 shows the long-term evolutions of the Crab pulsar in the VG(CR) model. The dashed curve is the line shown in Fig. 5 and the solid is the one with pulsar death (according to equation 12). These two lines evolve similar at the early age, and divide obviously after the transition point. The second line falls down afterwards and moves down-right with a slop about −1 finally. The initial period and transition period can be calculated according the second line, which are same with the values given by the first line within precision. Indicating that the effect of pulsar death can be ignored in the early time. As mentioned before, the particle density decreases as pulsar spin-down. The gap between these two lines is caused by the reducing of particle density. When the spin period of Crab approaches to the death period 2.52 s, the radio emission tends to stop and the slop of the second line is very large. Crab evolution in various acceleration models are given in Fig. 7. The bottom line is the evolution curve in SCLF(ICS) model. There is no transition point, indicating that the particle wind effect of this model is very weak. The top line is the evolution of the Crab pulsar in the OG model. This model is for high-energy radiation, the death of this model may differ from radio emission models. The results here for the OG model are just crude approximations to show the effect of pulsar death. Evolution of the Crab pulsar in other models are similar with each other.

Figure 7.

Rotational evolution of the Crab pulsar in all acceleration models. Because the VG(CR) and SCLF(I) models are very similar, their lines are coincident. Only the line in the VG(CR) model is plotted. See also Fig. 6 for explanations.

Open in new tabDownload slide

3 DISCUSSIONS

In the magnetic dipole radiation model, a pulsar is assumed as a magnetic dipole rotating in vacuum. It is just a low-order approximation of the true magnetosphere, and there must be particles acceleration and generation process so that the radio emission can be received. The actual magnetosphere surrounding the pulsar must be filled with plasma and corotate with pulsar because of magnetic freeze. According to the special relativity, the corotation magnetosphere cannot infinitely extend but only exist within the light cylinder. The magnetosphere is divided into two parts by the light cylinder: the closed magnetic field lines and the open magnetic field lines (including polar gap and outer gap). In the open field line region, particles are accelerated, generate radiation and flow out (thus taken away some of the rotational energy of the central neutron star). In the pulsar wind model, the pulsar is braked down by magnetic dipole radiation and particle wind. But in the electric current loss model, the electric current flow generates a spin-down torque (Harding, Contopoulos & Kazanas 1999; Beskin & Nokhrina 2007). Based on the same magnetospheric physics, the pulsar wind model and the current loss case are just different point of views, will produce the same results (e.g. section 3.2 in Tong et al. 2013). And they can be used to check each other. It is meaningful to take advantage of the existing information and explain more observations. The particle wind worked on pulsar spin-down was first proven by observations of intermittent pulsar PSR B1931+24 whose radio emission switches between ‘on’ and ‘of’ state (Kramer et al. 2006). Corresponding calculations for intermittent pulsars in pulsar wind model were made by Li et al. (2014).

As one of the best observed pulsar, observations of the Crab pulsar are relatively comprehensive (see Table 1 & timing results from Jodrell Bank website). It will be the best target to test radiation models. The pulsar parameters can be calculated by equations (9), (10), and (11) when model-dependent acceleration potential (or η) is given (see Table 2). The Goldreich–Julian charge density is taken as the primary particle density in these models of Xu & Qiao (2001) and Wu et al. (2003). However, Yue et al. (2007) constrained the primary particle density by observations of several pulsars and predicted that it should be 10³–10⁴ times of Goldreich–Julian charge density. In our calculations, the primary particle density is at least 100 times of Goldreich–Julian charge density. In the VG(CR) model, a particle density 10³ |ρ_GJ| is chosen. When applying the pulsar wind model to intermittent pulsars, a Goldreich–Julian charge density is assumed (Li et al. 2014). As the pulsar spin-down, its outflow particle density decreases. In the |$P-\dot{P}$| diagram, these intermittent pulsars are on the right of those young pulsars (see Fig. 6). This may explain the different particle density in the Crab pulsar and intermittent pulsars.

Other studies also found the need for a high particle density in the pulsar magnetosphere, e.g. when modelling the eclipse of the double pulsar system, giant pulses of the Crab pulsar, and magnetar X-ray spectra (Lyutikov 2008 and references therein), optical/ultraviolet excess of X-ray dim isolated neutron stars (Tong et al. 2010; Tong, Xu & Song 2011). Magnetohydrodynamical simulations of the Crab nebula also require a much higher electron flow from the magnetosphere (Buhler & Blandford 2014 and references therein). Besides, a super Goldreich–Julian current has been proposed in a series theoretical works about the pairs creations in the acceleration gap (Hirotani 2006; Zanotti, Morozova & Ahmedov 2012; Timokhin & Arons 2013). Combined with the results in this paper, it is possible that there are high-density plasmas (e.g. 10³–10⁴ times the Goldreich–Julian density) in both open and closed field line regions in the pulsar magnetosphere.

3.1 Applications to other sources

Here, the updated pulsar wind model is employed to compare with the observations of the Crab pulsar. It can also be applied to other sources. For eight pulsars with meaningful braking index measured (Espinoza et al. 2011; Lyne et al. 2015), their magnetospheric parameters can also be calculated, including magnetic field, inclination angle, particle density etc. But the higher second frequency derivatives of some pulsars (Hobbs, Lyne & Kramer 2010) may be attributed to the fluctuation of magnetosphere (Contopoulos 2007). During this process, more information will be helpful to constrain the model parameters, e.g. inclination angle, second braking index. Their evolution on the |$P-\dot{P}$| diagram can be calculated, which will be similar to Fig. 6. The pulsar PSR J1734−3333 has braking index n = 0.9 ± 0.2 (Espinoza et al. 2011). Its rotational energy loss rate will be dominated by the particle wind. Also from Fig. 6, the braking indices of young pulsars will naturally be divided into two groups. Braking indices of pulsars in the first group are close to 3, their evolution directions are down to the right in the |$P-\dot{P}$| diagram (⁠|$\dot{P}\propto P^{2-n}\sim P^{-1}$|⁠). Braking indices in the second group are close to one and their evolution directions are up to the right (⁠|$\dot{P}\propto P^{2-n}\sim P^{1}$|⁠). Generally speaking, sources in the second group will be older than the ones in the first group. At present, there are some hints for the evolution of pulsar braking index (Espinoza 2013). Future more observations will deepen our understanding on this aspect.

Glitch is a phenomenon of the pulsar suddenly spinning up and then slowing down by a larger rate. In the pulsar wind model, the lower braking index 2.3 during glitches (of the Crab pulsar) is due to the larger outflow particle density. The increasing particle density is a possible performance associated with glitches. Glitches induced magnetospheric activities are very common in magnetars (Kaspi, Gavriil & Woods 2003; Dib & Kaspi 2014). For the high magnetic field pulsar PSR J1846−0258, it also shows a smaller braking index after glitch (Livingstone et al. 2011). The cause of a smaller braking index may be similar to the Crab pulsar case. Due to a larger particle outflow, the pulsar should experience some net spin-down after glitch. This has already been observed (Livingstone, Kaspi & Gavriil 2010). Some primary calculations for PSR J1846−0258 has been done previously (Tong 2014, based on a wind braking model designed for magnetars). The smaller braking index of PSR J1846−0258 after glitch is consistent with the results here.

From Figs 6 and 7, the pulsar will go up-right in the |$P-\dot{P}$| diagram in the wind braking dominated case. Therefore, for pulsars in the up right corner of the |$P-\dot{P}$| diagram (e.g. high magnetic field pulsars and magnetars), their true dipole magnetic field may be much lower than the characteristic magnetic field. For magnetar SGR 1806−20, its characteristic magnetic field is about 5 × 10¹⁵ G (Tong et al. 2013). This dipole magnetic field is too high to understand comfortably (Vigano et al. 2013). Using the physical braking mechanism in this paper, the dipole magnetic field of SGR 1806−20 will be much lower. In the late stage of a pulsar, it will bent and go downward in the |$P-\dot{P}$| diagram. The Crab pulsar and other pulsars will not evolve into the cluster of magnetars. Furthermore, for magnetars, they also will go downwards in the |$P-\dot{P}$| diagram when they are aged enough. This may corresponds to the case of low magnetic field magnetars (Tong & Xu 2012).

4 CONCLUSIONS

An updated pulsar wind model which considered the effect of particle density and pulsar death is developed. It is employed to calculate parameters and simulate the evolution of the Crab pulsar. The braking index n < 3 is the combined effect of magnetic dipole radiation and particle wind. The lower braking indices measured between glitches are caused by the increasing particle density. By adding the glitch parameters and a change of period derivative caused by the change of particle density, the theoretical evolution curve is well fitted with observations (Fig. 3). This may be viewed as glitch induced magnetospheric activities in normal pulsars. Giving the model-dependent acceleration potential drop, the magnetic field, inclination angle and particle density of the Crab pulsar can be calculated. The evolution of braking index and the Crab pulsar in |$P-\dot{P}$| diagram can be obtained. In the |$P-\dot{P}$| diagram, the Crab pulsar will evolve towards the death valley, not to the cluster of magnetars (Fig. 6). Different acceleration models are also considered. The possible application of the present model to other sources (pulsars with braking index measured, and the magnetar population) is also mentioned.

ACKNOWLEDGEMENTS

The authors would like to thank the Referee Barsukov D.P. for helpful comments and R.X.Xu, J.B.Wang for discussions. HT is supported Xinjiang Bairen project, West Light Foundation of CAS (LHXZ201201), Qing Cu Hui of CAS, and 973 Program (2015CB857100).

REFERENCES

APPENDIX A: CONSIDERATION OF PULSAR ‘DEATH’

APPENDIX B: ASYMPTOTIC BEHAVIORS IN FIGS 4 AND 5