The imprint of pulsar parameters on the morphology of Pulsar Wind Nebulae
The morphology of young Pulsar Wind Nebulae (PWN) is largely determined by the properties of the wind injected by the pulsar. We have used a recent parametrization of the wind obtained from Force Free Electrodynamics simulations of pulsar magnetospheres to simulate nebulae for different sets of pulsar parameters. We performed axisymmetric Relativistic Magnetohydrodynamics simulations to test the morphology dependence of the nebula on the obliquity of the pulsar and on the magnetization of the pulsar wind. We compare these simulations to the morphology of the Vela and Crab PWN. We find that the morphology of Vela can be reproduced qualitatively if the pulsar obliquity angle is and the magnetization of the wind is high (). A morphology similar to the one of the Crab Nebula is only obtained for low magnetization simulations with . Interestingly, we find that Kelvin-Helmholtz instabilities produce small scale turbulences downstream of the reverse shock of the pulsar wind.
keywords:instabilities – MHD – shock waves – pulsars: general – pulsars: individual: Crab, Vela – ISM: supernova remnants
Most of the rotational energy lost by a pulsar is transferred to a relativistic particle wind. These particles are by numbers predominantly electrons and positrons (referred to together as electrons in the following). The wind is thought to be cold, meaning that its thermal energy is much less than its bulk kinetic and magnetic energies. When the wind interacts with ambient material, the particles become isotropised and radiate (Arons, 2012). This is what is then seen as a Pulsar Wind Nebula (PWN) in the sky. To date, around 100 of these systems have been found (Kargaltsev et al., 2015). In the X-rays, several of them show a torus morphology, with a jet emerging perpendicular to it (Kargaltsev & Pavlov, 2008).
Young PWN, with an age smaller than 10000 yrs, have not yet been distorted by the reverse shock of the stellar explosion (Gaensler & Slane, 2006; Kargaltsev et al., 2015). Their morphology is still closely related to the properties of the pulsar wind. Particularly in two cases, the Crab and Vela PWN, the plasma outflow can be resolved observationally in great detail, down to spatial scales of ly. These systems therefore provide a test bed to study the behaviour of relativistic plasma, which is also of relevance for other non-thermal sources as Active Galactic Nuclei (Netzer, 2014; Massaro et al., 2015) or Gamma-ray Bursts (Gehrels & Mészáros, 2012; Berger, 2014).
The properties of the pulsar wind – as its magnetic field, particle and velocity distributions – are not known with certainty today. However, over the past years there has been great progress in this respect. Several groups have performed Force Free Electrodynamics (FFE) (Spitkovsky, 2006; Kalapotharakos et al., 2012; Tchekhovskoy et al., 2013) and Particle in Cell (PIC) simulations of pulsar magnetospheres (Philippov & Spitkovsky, 2014; Cerutti et al., 2016). These simulations allow one to trace the properties of the wind out to several light cylinder radii . Due to the relativistic speed of the wind, it is expected that the wind does not have time to re-arrange itself on larger scales afterwards. In particular, the latitude dependence of its energy flux is expected to remain unchanged as the wind moves out into the nebula (Tchekhovskoy et al., 2016).
Recently, the first analytic parametrization of the latitude dependent luminosity of the pulsar wind has been derived from FFE simulations (Tchekhovskoy et al., 2016). The main parameter determining the wind properties is the obliquity angle between the pulsar spin axis and its magnetic moment . The latter can not be measured directly. Constraints from pulsar light curve modelling usually differ vastly between pulsar emission models (Pierbattista et al., 2015). The other important unknown parameter is the magnetization of the wind. It is thought that most of the wind energy is in its magnetic fields () at . When and where this energy is transferred to kinetic particle energy is not known. FFE simulations of pulsar magnetospheres do not include non-thermal particle acceleration, this question cannot be addressed by them.
In this paper we will study the dependence of the nebula morphology on and . Both of these parameters strongly affect the forces acting on the wind plasma. They are therefore expected to shape the morphology of the resulting nebula. We performed Relativistic Magnetohydrodynamic (RMHD) simulations to scan the parameter space of different obliquity angles and a high and low magnetization of the wind. RMHD simulations of PWN performed in the past have primarily focused on the Crab nebula (Hester, 2008; Bühler & Blandford, 2014). Qualitatively, the toroidal structure and the jet are well reproduced in 2D axisymmetric simulations (Komissarov & Lyubarsky, 2004; Del Zanna et al., 2004; Del Zanna et al., 2006; Volpi et al., 2009; Bucciantini et al., 2011). Dynamically, the motion of thin filaments, so called “wisps” is also reproduced (Camus et al., 2009). Recently, the first 3D simulations of the Crab nebula showed that axisymmetric simulations overpredict the strength of the jet (Porth et al., 2013, 2014). In addition, compared to 2D simulations, significant turbulence emerges far downstream of the wind termination shock in 3D. This potentially enhances the magnetic dissipation, allowing for larger magnetizations of the wind to reproduce the Crab’s morphology.
Unfortunately, performing several 3D simulations to scan the phase space of pulsar wind parameters is still computationally too expensive. We therefore performed 2D axisymmetric simulations. In contrast to most previous studies, we simulate both hemispheres. As was shown by Porth et al. (2014), this enhances the magnetic dissipation also in the axisymmetric case in the equatorial regions. Nevertheless, we will keep the axisymmetric limitation of our simulations in mind and will come back to it in the discussion of the simulation results in section 3.
We chose the length scales and spin-down power of the pulsar to values appropriate for the Vela PWN (Pavlov et al., 2003; Durant et al., 2013). To our knowledge this system has not been simulated in RMHD to date. We expect that apart from scaling factors, the PWN morphology does not depend strongly on this choice. Qualitatively, we expect the simulated morphologies to be similar also in other young PWN. We will use cgs units throughout, except for length scales, which for convenience will be given in light years.
2 Pulsar Wind Nebula simulations
In our simulations, a pulsar wind is injected into a homogeneous ambient medium at rest. The later has only the role to confine the pulsar wind. The regions of interest for the high-energy emission are in the inner nebula, close to the reverse shock of the pulsar wind. Once the pulsar wind has blown a sufficiently large bubble inside of its surrounding medium, the properties of the inner region do not depend strongly on the ambient environment. Instead, the inner nebula morphology is predominantly determined by the wind properties. As mentioned before, the latter are determined by the obliquity angle of the pulsars and by the wind magnetization. We performed six simulations with different combinations of these parameters shown in table 1. The simulation setup will be described in more detail in the following.
2.1 Simulation setup
The simulations were run with the RMHD module of the PLUTO111http://plutocode.ph.unito.it code version 4.2 (Mignone et al., 2007). The RMHD equations were evolved in time using an HLLC solver. We assumed a polytropic equation of state with an adiabatic index of . We applied spherical coordinates, with a linear binning in the polar angle and a logarithmic binning in radius . Adaptive Mesh Refinement (AMR) was used with four refinement levels, with a factor two in cell sizes between the levels (Mignone et al., 2012). The first grid was divided in 88 radial bins, from ly to ly, and 32 polar angle bins. At a typical distance of the reverse shock of ly the resolution is ly at the highest AMR level. To speed up the simulations, AMR is not activated in the unshocked wind region for ly. Reflective boundary conditions were applied along the symmetry axis and continuous outflow boundary conditions were applied at . More details on the simulation parameters can be found in the PLUTO configuration files shown in Appendix A.
The pulsar wind was injected into an ambient medium at . In order to speed up the simulations, we used a relatively thin ambient medium with a density of g cm. The work required to blow the pulsar wind bubble into the surrounding medium is thereby reduced. In order to further accelerate the simulations, the region within ly was filled with the unperturbed pulsar wind parameters at the start of the simulations.
We setup the pulsar wind parameters following the prescriptions given in Porth et al. (2014). The total energy flux density of the wind, , as a function of the polar angle and the obliquity angle of the pulsar is given by:
The angular dependence of the wind is given by the function , and is added for numerical reasons to avoid vanishing energy flux at the poles. is a scaling factor that normalizes the total energy flux to the spin-down luminosity of the Vela pulsar ergs s (Manchester et al., 2005). The angular dependence of the wind is obtained by averaging the parametric solutions given by Tchekhovskoy et al. (2016) over the azimuth angle :
The energy flux density is divided between a magnetic term and a kinetic term :
where is the wind magnetization. To assure that the Poynting flux vanishes at the poles, the magnetization saturates at:
, where is one of the simulation parameters reported in table 1. We assume that dissipation of magnetic fields in the striped wind zone changes the magnetization to
with . This assumes a perfect dissipation of the magnetic energy in opposite field lines in the stripped wind region. It is dynamically not important if the magnetic dissipation happens before the wind termination or right after it. Efficient dissipation after the wind termination has been shown in PIC simulation (Sironi & Spitkovsky, 2011). The full annihilation of opposite magnetic field lines in the stripped region implies that the average wind magnetization depends strongly on . In particular, for large values of the average magnetization will be low, independently of , as can be seen in table 1. The resulting energy flux distributions for the different simulations are shown in figure 1.
The magnetic field in the simulations is assumed to be purely toroidal, and inverts its direction in the equatorial plane:
Asymptotically one expects the poloidal and radial components of the magnetic field to decrease rapidly with distance with , while the toroidal component decreases with . For distances one therefore expects the toroidal component to dominate. Observationally, the toroidal magnetic field structure is supported by the fact that polarization measurements of PWN show an aligned magnetic structure, which is perpendicular to the symmetry axis of the nebula (Moran et al., 2014).
One complication in RMHD simulations of PWN is that the wind is expected to have a Lorentz factor . However, the code becomes numerically unstable for . Fortunately, in PWN the global plasma dynamics are primarily determined by the total kinetic energy (Porth et al., 2014). We therefore followed the approach taken by all previous RMHD simulation of PWN and compensated for the lower by increasing the plasma density:
The simulations were performed over a time period of yrs. After this time, the PWN has reached the self similar expansion phase, where spatial scales increase approximately linearly (Komissarov & Lyubarsky, 2004). This can be seen in fig. 2 for the position of the reverse shock: some initial oscillations are followed by a rapid expansion. Around yrs, the reverse shock position begins to expand more slowly and approximately linear in time for all simulations. The average shape of the shock surface as a function of latitude is shown in figure 3.
2.2 Radiative model
Synchrotron emission of relativistic electrons energetically dominates the radiative output of PWN. Over a wide range in frequency , non-thermal photon spectra are observed. The photon flux is typically well described by a power-law function . A typical value for the photon index in the X-rays is , which we will assume throughout.
In order to properly describe the synchrotron emission of PWN, the energy distribution of electrons as a function of position and energy needs to be known. In RMHD emission models of the Crab nebula, particles are typically injected at the wind termination shock. However, it is not known what the latitude dependence of the injected spectrum is (Olmi et al., 2014). Furthermore, it is well possible that non-thermal particle acceleration also happens in the body of the nebula, far downstream of the wind termination shock (Gallant & Arons, 1994; Amato & Arons, 2006; Bucciantini et al., 2011). Given these uncertainties, we proceed using a simpler approach. Following Komissarov & Lyubarsky (2004) and Del Zanna et al. (2004), we assume that the photon emissivity depends primarily on the Doppler factor and magnetic field perpendicular to the line of sight n in the plasma rest frame:
and is the velocity vector normalized to the speed of light. Previous studies have shown that and are the most important parameters in determining the synchrotron emissivity (Del Zanna et al., 2006; Camus et al., 2009; Bucciantini et al., 2011; Porth et al., 2014). We will therefore use the prescription above to gain a qualitative picture of the PWN morphology and to identify the most likely radiation regions. After having described the simulation procedure, we proceed to discuss their results.
In the following, we will focus our discussion on the snapshot at the end of the simulations at yrs. The nebula is already close to a self similar expansion at this time.
3.1 The plasma flow
The global plasma flow patterns of the different simulations are shown in figure 4. Common patterns are observed for all simulations: the pulsar wind is decelerated abruptly at an oblate reverse shock. Plasma at higher latitudes is transported to the polar regions due to the magnetic hoop stress and a jet emerges. In the equatorial region, plasma continues to move radially outward, creating a torus in the equatorial plane. The torus region is highly turbulent and the radial flow pattern is lost. The polar flow becomes stronger compared to the equatorial one with increasing average wind magnetization. As expected, the most relevant parameter is and not . This can be seen clearest for the case, where the flow patterns are very similar for and .
A zoom in on the flow pattern in the region close to the reverse shock is shown in figure 5. The region before the reverse shock, referred to as “wind region” in the following, becomes smaller with increasing . This can also be seen in the shock position at different latitude angles shown in figure 3. The wind region also becomes more oblate with increasing . This is expected, as the increased hoop stress results in increasing pressure in the polar regions (Lyutikov et al., 2016b). Despite these dependencies, the flow patterns are similar. The exception is the simulation of high magnetization and a small obliquity angle ( and ). The size of the wind region is greatly reduced and the polar flow greatly enhanced compared to the other simulations.
It is apparent from figure 5 that the degree of plasma turbulence increases with decreasing . The flow is very regular for the and simulation. There is also a regular plasma flow in the downstream of the reverse shock for and . These are the two simulations with the largest . In contrast, all other simulations show a high degree of turbulence, which emerges almost directly downstream of the termination shock.
Interestingly, loop like patterns from magnetic Kelvin-Helmholtz (KH) instabilities are observed downstream of the shear flow regions. An example is shown in figure 6. As the shear flow is reduced for high magnetizations, KH loops are predominantly observed in the low-sigma simulations. KH instabilities in PWN have been studied by Bucciantini & Del Zanna (2006) for local features within nebulae. To our knowledge, KH instabilities have not been reported in previous publications of global PWN simulations. However, its emergence had already been seen in the simulations discussed in the PhD thesis of Camus (2009). The reason that most previous simulations did not reveal this instability is likely that their spatial resolution was a factor larger compared to the ones presented here (Del Zanna et al., 2006; Camus et al., 2009; Bucciantini et al., 2011). As will be discussed in more detail in section 3.4, the presence of the KH instability could lead to interesting radiative signatures.
3.2 Synchrotron emission maps
In this section, we will look at the morphology expected from the synchrotron emission of the plasma flow described in the previous section. The emission of the inner regions of PWN is very anisotropic due to the strong Doppler boosting. The viewing angle therefore plays an important role in determining the observed morphology. For our discussion, we will assume a viewing angle of . This is approximately the angle under which we observe the Vela and Crab nebulae from Earth (Weisskopf et al., 2000; Helfand et al., 2001; Ng & Romani, 2004; Weisskopf et al., 2012). These sources will be discussed in more detail in sections 3.3 and 3.4. Maps for viewing angles of and are shown in appendix B.
The synchrotron maps for simulations with are shown in figure 7. All of them show a torus, which comes from the equatorial region of the PWN. The synchrotron maps for and are qualitatively similar, with wider rings than the ones seen in simulation with . The reason for this is that in the latter case more of the emission comes from higher latitudes. This is a result of the higher average magnetization and higher wind energy flux for (see figure 1).
Figure 8 shows the synchrotron maps calculated for the simulations done with . In this case, the difference between different obliquities is stronger. The main reason for this is that for and the average magnetization is close to unity. This results in elongated synchrotron nebulae. In the case of , the torus is barely visible anymore as the emission is dominated by the higher latitudes. The jet is also clearly visible in this case.
Before we continue to confront these simulations with observations, we would like to recall two caveats: (1) Our simulations are axisymmetric. 3D simulations have shown that axisymmetry is a good approximation close to the wind reverse shock, in particular in the equatorial regions (Porth et al., 2014). However, at high latitudes and further away from the reverse shock, 3D simulations show significantly different emission patterns. It follows from this that emission associated to the equatorial region in the inner nebula can be expected to be more robust. Polar regions on the other hand, and in particular the jet, are less trustworthy and will therefore not be a focus of discussion. (2) The emission model we apply here ignores spatial differences in the particle distribution function. The lack of strong spectral variation in the emission from the inner regions of PWN indicates that this is likely a good approximation (Mori et al., 2004). Nevertheless, this is certainly oversimplified; e.g. individual structures in the Crab nebula are known to be visible only at particular wavebands (Hester et al., 2002). Keeping these caveats in mind, we will proceed to compare our simulation results to the Vela and Crab nebulae.
3.3 Vela Pulsar Wind Nebula
The Vela PWN is embedded in the Vela Supernova Remnant, also known as G263.9-3.3. Several regions of non-thermal emission are known within this remnant. The brightest one in radio is labelled Vela-X. The Vela pulsar is located in Vela-X (Horns et al., 2006; de Jager et al., 2008; Abramowski et al., 2012; Grondin et al., 2013). It is one of the closest and brightest gamma-ray pulsar known to date. Its distance has been determined via parallax to ly (Caraveo et al., 2001; Dodson et al., 2003). Due to its proximity the Vela pulsar might contribute significantly to the local cosmic-ray electron flux (Hinton et al., 2011). In the region surrounding the pulsar, the Chandra X-ray Observatory revealed a double ring structure, which is likely related to the reverse shock of the pulsar wind (Helfand et al., 2001; Pavlov et al., 2003; Durant et al., 2013). It is this innermost emission which we aim to reproduce in the simulations here. Unfortunately, the observational data is restricted to the X-ray band, as the ring structure has not been detectable so far at other wavebands (Moran et al., 2014; Marubini et al., 2015).
The three dimensional structure of the Vela rings has been interpreted as two equal torii, which are ontop of each other (Helfand et al. (2001), for an high contrast image of the Vela rings see figure 2 in Pavlov et al. (2003)). This results in a morphology which is more elongated along the symmetry axis compared to the Crab Nebula. The rings were also found to be closer to the pulsar as the innermost ring in the Crab nebula. It was suggested that this might be due to a larger wind magnetization of order unity (Helfand et al., 2001). Indeed, our simulations confirm that solutions with higher magnetization resemble the ring structure better. Solutions with a lower – as all simulations with and both simulations with – are to wide and the torii are too narrow. Interestingly, several bright small-scale structures are found along the symmetry axis of the nebula. These features are related to the highly beamed right downstream of the reverse shock. Such features are also observed in the X-ray data (Levenfish et al., 2013).
The best morphology agreement is found for the and simulation. In this simulation rings of similar size originate downstream of the reverse shock. One lower ring from the equatorial regions right behind the wind termination and an upper ring from the torus. In contrast, all other solutions result in rings of increasing size as one moves away from the pulsar. The exceptions are the simulations for . However, for a the nebula is far too elongated and for the case of no double ring structure is obtained. The value of agrees well with the value of inferred from polarization measurements of the pulsar profile (Johnston et al., 2005). Models of the pulsed gamma-ray emission typically give higher values – (Abdo et al., 2010).
3.4 Crab Pulsar Wind Nebula
The Crab is the most studied PWN. Its torus and jet structure can be observed in great detail from the radio to X-ray band (Hester, 2008; Bühler & Blandford, 2014). It has therefore been the primary target for RMHD studies of PWN. Several authors have qualitatively reproduced the morphology of the inner nebula in axisymmetric simulations assuming a low wind magnetization (Komissarov & Lyubarsky, 2004; Del Zanna et al., 2006; Volpi et al., 2009; Camus et al., 2009; Porth et al., 2014). We confirm this findings in our simulations. A good agreement is found for case of and and for both simulations with . Also in agreement with previous studies, a small bright feature is found just below the pulsar position (Lyutikov et al., 2016b; Yuan & Blandford, 2015). This “inner knot” was proposed to be the site of the recently discovered gamma-ray flares (Buehler et al., 2012; Mayer et al., 2013). However, no observational evidence for this has been found to date (Rudy et al., 2015).
It is puzzling, that the innermost ring of the Crab nebula observed in X-rays does not show a brightness profile as expected from Doppler beaming. The back side of the inner ring has a brightness which is comparable to its front side. In RMHD simulations the ring is found to be much fainter than in observations (Porth et al., 2014). In addition, the ring is composed of a series of knots, in contrast to the smooth profiles found in the simulations. Our simulations show a similar disagreement. However, the observed KH instabilities at the shear flow downstream of the reverse shock lead to Doppler boosted emission regions also on the back side of emission rings close to the reverse shock. This results in a brighter emission from the receding part of the flow than expected in radial flow models. Our simulations show that this effect is not strong enough to result in a ring of equal brightness. However, the increased turbulence in the KH loops might trigger increased particle acceleration via magnetic reconnection (Cerutti et al., 2013; Sironi & Spitkovsky, 2014) and magnetoluminescense (Blandford et al., 2014; East et al., 2015; Nalewajko et al., 2016; Yuan et al., 2016; Lyutikov et al., 2016a). Stochastic Fermi acceleration might also occur (Rieger et al., 2006). In either case, fresh injection of high energy particles in back-flowing plasma could lead to a decreased Doppler asymmetry only for the highest synchrotron frequencies. Indeed, the inner ring is only observed in X-rays and has not been detected at lower frequencies to date (Hester et al., 2002). As the development of the KH instability is reduced for higher plasma magnetization, the latter could be constrained if this interpretation is correct. Testing this idea quantitatively requires dedicated simulations which include particle acceleration, which is beyond the scope of this paper.
4 Summary and outlook
We have performed axisymmetric RMHD simulations of PWN to scan the parameter space for different pulsar wind properties. We have simulated the wind emerging for different pulsar obliquities, which has recently been derived from FFE simulations (Tchekhovskoy et al., 2016). In addition we have tested different wind magnetizations. In general, we find that the average wind magnetization is the most important parameter in determining the PWN morphology. The main effect of increasing obliquity angle is to increase the size of the striped wind region, where we have assumed a perfect dissipation of opposite magnetic field lines. We found that the wind region upstream of the reverse shock is smaller in size and becomes more oblate with increasing .
With the exception of the wind morphology for and , all simulations showed a torus in their emission maps, which emerges from the equatorial region downstream of the reverse shock. We have compared the morphologies of the different simulations to the Vela and Crab PWN. For Vela, we found that the simulation with the parameters and gives the best match to the observed morphology. For the Crab nebula, all simulations with a low match the observed morphology ( and ; and ; for and ).
We found that, particularly for low magnetizations, KH instabilities develops at the downstream at the shear flow of the reverse shock. The KH loops have the effect to increase the emission of the receding side of the nebula compared from what is expected from Doppler boosting of a radial flow. We suggest that this effect might help to explain that the innermost ring observed in X-rays in the Crab nebula has almost constant brightness.
We have pointed out the caveat that these conclusions rely on axisymmetric simulations. It would be desirable to confirm these findings with 3D simulations in the future. For a more quantitative comparison of observations and simulations, it will also be important to include a model for particle acceleration in the simulations. From the observational side, the detection of the rings observed in the Vela PWN outside of the X-ray band would be crucial to constrain the electron energy distribution. Taken together these steps provide the prospects in understanding the plasma flow quantitatively in PWN. This would be the first time this is achieved for relativistic plasmas flows and would likely have implications also for other sources as GRBs or AGN.
We thank the referee Oliver Porth for very helpful comments and discussions. We believe they greatly improved the quality of the article. We also want to thank Andrea Mignone and Claudio Zanni very much for their help with the PLUTO code.
- Abdo et al. (2010) Abdo A. A., et al., 2010, The Astrophysical Journal, 713, 154
- Abramowski et al. (2012) Abramowski A., et al., 2012, Astronomy & Astrophysics, 548, A38
- Amato & Arons (2006) Amato E., Arons J., 2006, The Astrophysical Journal, 653, 325
- Arons (2012) Arons J., 2012, Space Science Reviews, 173, 341
- Berger (2014) Berger E., 2014, Annual Review of Astronomy and Astrophysics, 52, 43
- Blandford et al. (2014) Blandford R., Simeon P., Yuan Y., 2014, Nuclear Physics B (Proceedings Supplements), Volume 256, p. 9-22., 256, 9
- Bucciantini & Del Zanna (2006) Bucciantini N., Del Zanna L., 2006, Astronomy and Astrophysics, 454, 393
- Bucciantini et al. (2011) Bucciantini N., Arons J., Amato E., 2011, Monthly Notices of the Royal Astronomical Society, 410, 381
- Buehler (2016) Buehler R., 2016, PWN Imprint webpage, http://www.rolfbuehler.net/notes/PWNimprint.html
- Buehler et al. (2012) Buehler R., et al., 2012, The Astrophysical Journal, 749, 26
- Bühler & Blandford (2014) Bühler R., Blandford R., 2014, Reports on Progress in Physics, 77, 066901
- Camus (2009) Camus N. F., 2009, Relativistic Magnetohydrodynamical Models of Pulsar Wind Nebulae, %****␣PWNsim.bbl␣Line␣75␣****https://www1.maths.leeds.ac.uk/applied/Research/agf/thesis/nicholas_camus.pdf
- Camus et al. (2009) Camus N. F., Komissarov S. S., Bucciantini N., Hughes P. a., 2009, Monthly Notices of the Royal Astronomical Society, 400, 1241
- Caraveo et al. (2001) Caraveo P. A., De Luca A., Mignani R. P., Bignami G. F., 2001, The Astrophysical Journal, 561, 930
- Cerutti et al. (2013) Cerutti B., Werner G. R., Uzdensky D. A., Begelman M. C., 2013, The Astrophysical Journal, Volume 782, Issue 2, article id. 104, 15 pp. (2014)., 782
- Cerutti et al. (2016) Cerutti B., Philippov A., Spitkovsky A., 2016, Monthly Notices of the Royal Astronomical Society, 457, 2401
- Del Zanna et al. (2004) Del Zanna L., Amato E., Bucciantini N., 2004, Astronomy and Astrophysics, 421, 1063
- Del Zanna et al. (2006) Del Zanna L., Volpi D., Amato E., Bucciantini N., 2006, Astronomy and Astrophysics, 453, 621
- Dodson et al. (2003) Dodson R., Legge D., Reynolds J. E., McCulloch P. M., 2003, The Astrophysical Journal, 596, 1137
- Durant et al. (2013) Durant M., Kargaltsev O., Pavlov G. G., Kropotina J., Levenfish K., 2013, The Astrophysical Journal, 763, 72
- East et al. (2015) East W. E., Zrake J., Yuan Y., Blandford R. D., 2015, Physical Review Letters, Volume 115, Issue 9, id.095002, 115
- Gaensler & Slane (2006) Gaensler B. M., Slane P. O., 2006, Annual Review of Astronomy and Astrophysics, 44, 17
- Gallant & Arons (1994) Gallant Y. A., Arons J., 1994, The Astrophysical Journal, 435, 230
- Gehrels & Mészáros (2012) Gehrels N., Mészáros P., 2012, Science (New York, N.Y.), 337, 932
- Grondin et al. (2013) Grondin M.-H., Romani R. W., Lemoine-Goumard M., Guillemot L., Harding A. K., Reposeur T., 2013, The Astrophysical Journal, 774, 110
- Helfand et al. (2001) Helfand D. J., Gotthelf E. V., Halpern J. P., 2001, The Astrophysical Journal, 556, 380
- Hester (2008) Hester J. J., 2008, Annual Review of Astronomy and Astrophysics, 46, 127
- Hester et al. (2002) Hester J. J., et al., 2002, The Astrophysical Journal, 577, L49
- Hinton et al. (2011) Hinton J. A., Funk S., Parsons R. D., Ohm S., 2011, The Astrophysical Journal, 743, L7
- Horns et al. (2006) Horns D., Aharonian F., Santangelo A., Hoffmann A. I. D., Masterson C., 2006, Astronomy and Astrophysics, 451, L51
- Johnston et al. (2005) Johnston S., Hobbs G., Vigeland S., Kramer M., Weisberg J. M., Lyne A. G., 2005, Monthly Notices of the Royal Astronomical Society, 364, 1397
- Kalapotharakos et al. (2012) Kalapotharakos C., Contopoulos I., Kazanas D., 2012, Monthly Notices of the Royal Astronomical Society, 420, 2793
- Kargaltsev & Pavlov (2008) Kargaltsev O., Pavlov G. G., 2008, in AIP Conference Proceedings. AIP, pp 171–185, doi:10.1063/1.2900138
- Kargaltsev et al. (2015) Kargaltsev O., Cerutti B., Lyubarsky Y., Striani E., 2015, Space Science Reviews, 191, 391
- Komissarov & Lyubarsky (2004) Komissarov S. S., Lyubarsky Y. E., 2004, Monthly Notices of the Royal Astronomical Society, 349, 779
- Levenfish et al. (2013) Levenfish K. P., Bykov A. M., Durant M., Kargaltsev O. Y., Kropotina Y. A., Pavlov G. G., Krassilchtchikov A. M., Uvarov Y. A., 2013, Memorie della Societa Astronomica Italiana, 84, 588
- Lyutikov et al. (2016a) Lyutikov M., Sironi L., Komissarov S., Porth O., 2016a, p. 140
- Lyutikov et al. (2016b) Lyutikov M., Komissarov S. S., Porth O., 2016b, Monthly Notices of the Royal Astronomical Society, 456, 286
- Manchester et al. (2005) Manchester R. N., Hobbs G. B., Teoh A., Hobbs M., 2005, The Astronomical Journal, 129, 1993
- Marubini et al. (2015) Marubini T. E., Sefako R. R., Venter C., de Jager O. C., 2015, Proceedings of SAIP2012, p. 240
- Massaro et al. (2015) Massaro F., Thompson D. J., Ferrara E. C., 2015, The Astronomy and Astrophysics Review, 24, 2
- Mayer et al. (2013) Mayer M., Buehler R., Hays E., Cheung C. C., Dutka M. S., Grove J. E., Kerr M., Ojha R., 2013, The Astrophysical Journal, 775, L37
- Mignone et al. (2007) Mignone A., Bodo G., Massaglia S., Matsakos T., Tesileanu O., Zanni C., Ferrari A., 2007, The Astrophysical Journal Supplement Series, 170, 228
- Mignone et al. (2012) Mignone A., Zanni C., Tzeferacos P., van Straalen B., Colella P., Bodo G., 2012, The Astrophysical Journal Supplement Series, 198, 7
- Moran et al. (2014) Moran P., Mignani R. P., Shearer A., 2014, Monthly Notices of the Royal Astronomical Society, 445, 835
- Mori et al. (2004) Mori K., Burrows D. N., Hester J. J., Pavlov G. G., Shibata S., Tsunemi H., 2004, The Astrophysical Journal, 609, 186
- Nalewajko et al. (2016) Nalewajko K., Zrake J., Yuan Y., East W. E., Blandford R. D., 2016, eprint arXiv:1603.04850
- Netzer (2014) Netzer H., 2014, Annual Review of Astronomy and Astrophysics, 53, 150619171245001
- Ng & Romani (2004) Ng C., Romani R. W., 2004, The Astrophysical Journal, 601, 479
- Olmi et al. (2014) Olmi B., Del Zanna L., Amato E., Bandiera R., Bucciantini N., 2014, Monthly Notices of the Royal Astronomical Society, 438, 1518
- Pavlov et al. (2003) Pavlov G. G., Teter M. A., Kargaltsev O., Sanwal D., 2003, The Astrophysical Journal, 591, 1157
- Philippov & Spitkovsky (2014) Philippov A. A., Spitkovsky A., 2014, The Astrophysical Journal, 785, L33
- Pierbattista et al. (2015) Pierbattista M., Harding A. K., Grenier I. A., Johnson T. J., Caraveo P. A., Kerr M., Gonthier P. L., 2015, Astronomy & Astrophysics, 575, A3
- Porth et al. (2013) Porth O., Komissarov S. S., Keppens R., 2013, Monthly Notices of the Royal Astronomical Society, 431, L48
- Porth et al. (2014) Porth O., Komissarov S. S., Keppens R., 2014, Monthly Notices of the Royal Astronomical Society, 438, 278
- Rieger et al. (2006) Rieger F. M., Bosch-Ramon V., Duffy P., 2006, Astrophysics and Space Science, Volume 309, Issue 1-4, pp. 119-125, 309, 119
- Rudy et al. (2015) Rudy A., et al., 2015, The Astrophysical Journal, 811, 24
- Sironi & Spitkovsky (2011) Sironi L., Spitkovsky A., 2011, The Astrophysical Journal, 741, 39
- Sironi & Spitkovsky (2014) Sironi L., Spitkovsky A., 2014, The Astrophysical Journal Letters, Volume 783, Issue 1, article id. L21, 6 pp. (2014)., 783
- Spitkovsky (2006) Spitkovsky A., 2006, The Astrophysical Journal, 648, L51
- Tchekhovskoy et al. (2013) Tchekhovskoy A., Spitkovsky A., Li J. G., 2013, Monthly Notices of the Royal Astronomical Society, 435, L1
- Tchekhovskoy et al. (2016) Tchekhovskoy A., Philippov A., Spitkovsky A., 2016, Monthly Notices of the Royal Astronomical Society, 457, 3384
- Volpi et al. (2009) Volpi D., Del Zanna L., Amato E., Bucciantini N., 2009, preprint, 0903.4120
- Weisskopf et al. (2000) Weisskopf M. C., et al., 2000, The Astrophysical Journal, 536, L81
- Weisskopf et al. (2012) Weisskopf M. C., Elsner R. F., Kolodziejczak J. J., O’Dell S. L., Tennant A. F., 2012, The Astrophysical Journal, 746, 41
- Yuan & Blandford (2015) Yuan Y., Blandford R. D., 2015, Monthly Notices of the Royal Astronomical Society, 454, 2754
- Yuan et al. (2016) Yuan Y., Nalewajko K., Zrake J., East W. E., Blandford R. D., 2016, p. 23
- de Jager et al. (2008) de Jager O. C., Slane P. O., LaMassa S., 2008, The Astrophysical Journal, 689, L125
Appendix A Pluto configuration file
The relevant parts of the pluto.ini file are listed below.
The relevant parts of the definitions.h file are listed below.