Publications of the Astronomical Society of the Paciﬁc, 129:124201 (11pp), 2017 December https://doi.org/10.1088/1538-3873/aa8dee © 2017. The Astronomical Society of the Paciﬁc. All rights reserved. Printed in the U.S.A. On the Geometry of the X-Ray Emission from Pulsars. I. Model Formulation and Tests 1 2 Rigel Cappallo1,2, Silas G. T. Laycock1,2, and Dimitris M. Christodoulou1,3 Lowell Center for Space Science and Technology, 600 Suffolk Street, Lowell MA, 01854, USA Department of Physics and Applied Physics, University of Massachusetts Lowell, Lowell MA, 01854, USA; rigelcappallo@gmail.com, silas_laycock@uml.edu 3 Department of Mathematical Sciences, University of Massachusetts Lowell, Lowell MA, 01854, USA; dimitris_christodoulou@uml.edu Received 2017 April 10; accepted 2017 September 20; published 2017 October 24 Abstract X-ray pulsars are complex magnetized astronomical objects in which many different attributes shape the pulse proﬁles of the emitted radiation. For each pulsar, the orientation of the spin axis relative to our viewing angle, the inclination of the magnetic dipole axis relative to the spin axis, and the geometries of the emission regions all play key roles in producing its unique pulse proﬁle. In this paper, we describe in detail a new geometric computer model for X-ray emitting pulsars and the tests that we carried out in order to ensure its proper operation. This model allows for simultaneous tuning of multiple parameters for each pulsar and, by ﬁtting observed proﬁles, it has the potential to determine the underlying geometries of many pulsars whose pulse proﬁles have been cataloged and made public in modern X-ray databases. Key words: accretion, accretion disks – methods: numerical – pulsars: general – stars: magnetic ﬁeld – stars: neutron – X-rays: binaries Online material: color ﬁgures 1. Introduction This phenomenology makes XRPs a rich environment for the study of high-energy astrophysics. The pulse proﬁles of XRPs vary greatly between objects in the class and they are both energy and geometry dependent (Basko & Sunyaev 1976; Nagel 1981; Pavlov et al. 1994; Ferrigno et al. 2011). They display a variety of morphologies when viewed in different energy bands; and they often exhibit contrasting patterns due to variations of the emission regions near the surfaces of the NSs or projection effects (Karino 2007; Klus 2015). This phenomenology suggests strongly that various physical processes become more/less important at distinct energy thresholds (Hong et al. 2017). Many factors modifying pulse proﬁles have been investigated with a variety of empirical models built by many researchers in the past (Wang & Welter 1981; Mészáros & Nagel 1985a, 1985b; Mészáros & Riffert 1988; Riffert & Meszáros 1988; Parmar et al. 1989; Riffert et al. 1993; Kraus et al. 1995; Leahy & Li 1995; Beloborodov 2002, to name a few key investigations), but these models limit their parameter spaces in different ways by introducing various assumptions about the emitted beams, or relativistic effects, or the geometry of the emitting regions and their magnetic ﬁelds. In hopes of further understanding the structures of the emission regions and the various physical processes that dictate differing beam patterns at various energies in XRPs, we have developed a computer code that simulates a variety of orientations and emission geometries through the use of multiple free Pulsars are magnetized neutron stars (NSs) that rotate at various spin periods (PS ∼ 10−3–104 s) while emitting beams of radiation from regions believed to be coincident with their magnetic poles (e.g., Mereghetti 2001; Frank et al. 2002; Lyne & Graham-Smith 2012). If these beams happen to cross our line of sight (LOS) and the NS spin and magnetic axes are not aligned, then our telescopes detect periodic pulsations from these objects (Lyne & Graham-Smith 2012). Pulsars have been observed emitting in all areas of the electromagnetic spectrum, from radio waves to gamma-rays (Thompson 2000). These types of NSs all have differing mechanisms causing their particular emission energies (Lewin et al. 1995; Mereghetti & Stella 1995; DeDeo et al. 2001; Eggleton 2001; Mereghetti 2001; Cordes et al. 2004; Becker & Wolff 2005; Kylaﬁs et al. 2012). The focus of this paper is on modeling the pulse proﬁles of NSs that emit strongly in the X-ray regime, known as X-ray pulsars (XRPs), and especially on the subset of Be/X-ray pulsars (Coe et al. 2010; Reig 2011; Klus et al. 2014; Coe & Kirk 2015; Klus 2015; Haberl & Sturm 2016). These objects are found in binary systems where matter accretion from their stellar companions results in highly energetic phenomena. Their spin periods vary in time t, exhibiting both spin-up (dPS/dt < 0) and spin-down (dPS/dt > 0) states, as well as transitions between these states (Christodoulou et al. 2017a). 1 Publications of the Astronomical Society of the Paciﬁc, 129:124201 (11pp), 2017 December Cappallo, Laycock, & Christodoulou parameters to be determined by direct ﬁts to observed pulse proﬁles. The proﬁles created from this model will subsequently be ﬁtted to existing large catalogs of observed pulse proﬁles so that the underlying physical mechanisms responsible for the X-ray emission patterns from each source may become clearer. Similar models to our own have been developed previously by Wang & Welter (1981), Parmar et al. (1989), and Beloborodov (2002); yet, these models include various additional constraints that limit the number of free parameters and their applicability to the entire class of XRPs. For example, in Beloborodov (2002) the emission regions were assumed to be antipodal and isotropic sources of radiation, and in Parmar et al. (1989) the magnetic-dipole axis did not pass through the center of the NS, resulting in some new complex features in the modeled proﬁle. Our model is capable of reproducing results from these models by ﬁxing the same geometrical parameters that these authors held constant. But the new model also takes the next step by providing the ability to ﬁne-tune a larger number of free parameters. As such, it will enable us to ﬁt large databases of widely varying pulse shapes and it will ultimately give us a better understanding of the orientations and underlying geometries of the emission regions in many pulsars. In what follows, we describe our theoretical framework for modeling XRP emission (Section 2), the geometry and the computer code of our model (Section 3), and some preliminary tests and comparisons with previous landmark results from the literature (Sections 4 and 5). We conclude in Section 6 with a discussion of the applicability of our new model and future applications. In paper II of this series, we intend to present our detailed model ﬁts to the pulse proﬁles of the Be/X-ray Magellanic pulsars for which we have created a new library of all previous X-ray observations in the past 15 years (Yang et al. 2017). Figure 1. Orientations of pencil and fan beams in an accreting NS. AC represents the accretion columns and z is the rotation axis. The magnetic dipole axis is along the ACs. (A color version of this ﬁgure is available in the online journal.) Christodoulou et al. 2016, 2017b), exotic physics may possibly occur in this realm (Bachetti et al. 2014; Klus et al. 2014; Tendulkar et al. 2014). A better understanding of how these emission regions are structured and behave will eventually lead to a clearer grasp of the physics in very high-energy regimes. 2.1. Hot Spot Geometry It is believed that there are two underlying, fundamental components to the HS geometry. At lower luminosities (LX < 1037 ergs−1), and for certain magnetic ﬁeld strengths and orientations, there exists a “pencil-beam” component which consists of radiation directed along the local magnetic magnetic ﬁeld lines (Zavlin et al. 1995). At higher luminosities, a “fan-beam” component begins to emerge which exhibits maximum radiation perpendicular to the magnetic dipole axis. This orthogonal component explains the double-peaked structures that are found only in proﬁles from high luminosity (LX > 1037 ergs−1) sources (Ferrigno et al. 2011; Klus et al. 2014; Klus 2015). These double peaks appear when the local energy ﬂux exceeds the Eddington ﬂux (deﬁned as the ﬁrst-order moment of the radiation ﬁeld) by a factor of ∼100, forming a shock front above the NS surface with photons “leaking out” of the sides of the column of slowly falling plasma contained below the front (see Becker 1998, and our Figure 1). 2. Theoretical Framework X-ray binary (XRB) pulsars consist of a NS and a stellar companion supplying material (by wind outﬂow or Roche lobe overﬂow, or when the NS is immersed in circumstellar material) onto an accretion disk that forms, in many cases, around the NS. The plasma in the inner accretion disk redistributes its angular momentum and some of the matter with reduced angular momentum descends toward the surface of the NS following magnetic-ﬁeld lines. As this matter is ﬁnally decelerated on impact, it releases large amounts of gravitational potential energy in the form of high-energy radiation. For this reason, the electromagnetic emission regions of pulsars are thought to coincide with the positions of their magnetic poles. The geometry of such emission regions, known as hot spots (HSs), is believed to be complicated. Due to the high energies of the particles and the large magnitudes of the NS magnetic ﬁelds (∼1012–13 G; Stella et al. 1986; 2 Publications of the Astronomical Society of the Paciﬁc, 129:124201 (11pp), 2017 December Cappallo, Laycock, & Christodoulou Figure 2. Geometry in the light-bending approximation of Beloborodov (2002). (A color version of this ﬁgure is available in the online journal.) The alignment of the dipolar magnetic ﬁeld lines in the HSs is also an issue that must be considered. The ﬁeld lines could be normal to the NS surface or, alternatively, they could be aligned with the magnetic dipole axis, which may not be passing through the NS center (e.g., Burnard et al. 1988; Parmar et al. 1989). Only for an antipodal arrangement does the dipole axis pass through the NS center and then the above two ﬁeld line geometries coincide. This simpliﬁed arrangement has been used in many models, including Beloborodov (2002) and Wang & Welter (1981), but not in Parmar et al. (1989) who produced very different proﬁles resulting from such an offcentered magnetic ﬁeld (among other considerations as well). Some results from these models have been reproduced with our new code and they are discussed in Sections 3.4 and 5 below. Figure 3. Far side of a NS (outer circle). The shaded regions represent the unobserved sections of the surface and correspond to different ratios of rg/RNS. The circle colored in dark blue and black is for a ratio of 1/3 and the inner black circle represents the canonical NS with a ratio of 0.4135. (A color version of this ﬁgure is available in the online journal.) normal vector to the point of emission; and α is the angle between the initial direction of photon emission and the same normal vector (Figure 2). Equation (1) is used to determine α for a choice of ψ and the ratio rg/RNS. For the canonical NS, rg=4.135 km and the term in the ﬁrst parentheses in Equation (1) is reduced to a constant, viz. 2.2. Gravitational Considerations The gravitational ﬁeld of a typical NS is quite strong near its surface (g ∼ 1012 m s−2 for the canonical NS with mass 1.4 Me and radius 10 km) and the exterior space deviates markedly from Euclidean space. For this reason, gravitational lightbending effects onto the emitted radiation must be taken into account. For a relatively slowly rotating NS, frame dragging can be neglected and then the Schwarzschild metric is a sufﬁcient approximation for the surrounding space (Pechenick et al. 1983). Although the general relativistic equations governing electromagnetic wave propagation in the Schwarzschild metric cannot be solved analytically, a simple yet powerful approximation for light bending has been developed by Beloborodov (2002). This approximation is represented by the equation ⎛ rg ⎞ 1 - cos a = ⎜1 ⎟ (1 - cos y) , ⎝ RNS ⎠ 1 - cos a = 0.5865 (1 - cos y). (2 ) All of the proﬁles created by our model in this paper use Equation (2), but the model generally retains the ratio rg/RNS as a free parameter. This ratio is effectively related to the gravitational redshift zg of the emitted radiation, that is ⎛ rg ⎞- 2 z g = ⎜1 ⎟ - 1. ⎝ RNS ⎠ 1 (3 ) For comparison purposes, zg=0.3058 for the canonical NS with rg/RNS=0.4135. Equation (2) reproduces the results of the results of Pechenick et al. (1983) to within an accuracy of 3% for α90° and the error decreases to 1% for α75° (Beloborodov 2002). For our model, errors of 3% or less are acceptable when ﬁtting pulse proﬁles since larger errors are inherent in the observations. But we retain four signiﬁcant ﬁgures in the above equations to avoid introducing systematic errors of order of 1%. (1 ) where rg = 2GM c 2 is the Schwarzschild radius for a body of mass M, gravitational constant G, and light speed c; RNS is the radius of the NS; ψ is the angle between the LOS and the 3 Publications of the Astronomical Society of the Paciﬁc, 129:124201 (11pp), 2017 December Cappallo, Laycock, & Christodoulou Figure 5. Screen shot of the visual representation of our model. This speciﬁc geometry uses asymmetric HSs and produces the proﬁles shown in Figures 16 and 17 below. The two HS vectors are green and the rotation axis is red. (A color version of this ﬁgure is available in the online journal.) 3. Model and Computer Code The computer code, named Polestar,4 is written in Python 2.7 and VPython 3.0 and simulates radiation emission from a single NS with an arbitrary number of HSs radiating a combination of pencil and fan beams, along with the capability of incorporating more complex beam structures. The inclusion of multiple HSs facilitates the implementation of extended regions of emission on the surface of the NS rather than just two point sources (Ho 2007). The model NS rotates at a given angular velocity w , the magnetic ﬁeld B can be chosen to be be misaligned in a number of different ways, and the resulting pulse proﬁle of the radiation that reaches the observer is then computed as a function of phase. Finally, VPython 3.0 is used to create an animated, three-dimensional, visual representation of the geometry and the orientation of the source relative to our LOS, such as that shown in Figure 5. The mathematical framework that forms the basis of the Polestar code is described below. Figure 4. Pulse proﬁle generated by two antipodal isotropic HSs (a) without and (b) with gravitational light bending included. The differences are described in the text. As a consequence of light bending, more than half of the surface of a NS is observed at any time implying that there may be points in the pulse proﬁle where both poles are visible. The terminator that separates the region of the surface that cannot be observed is a circle on the far side of the NS (Figure 3) whose radius depends on RNS and the NS mass (Beloborodov 2002). The inclusion of gravitational light bending in the model is necessary because the effect alters pulse proﬁles in a dramatic fashion. Figure 4 shows an example of two antipodal isotropically emitting HSs with and without light bending. Whenever both poles are visible, the observed ﬂux remains constant (panel (b)) producing the plateaus that ﬁll the deep troughs which dominate in the nonrelativistic case (panel (a)). It should be noted however that this behavior emerges only for antipodal HS orientations. 3.1. Mathematical Framework Conceptually, the model is vector-based and consists of a line-of-sight unit vector L, an arbitrary number of HS normal unit vectors mi (i = 1, 2, ...), and a geometric beam (intensity) pattern for each HS (see Figure 6). The rest frame is taken to be the NS frame, so L is rotated in incremental steps Δξ through one full revolution about the rotation axis w of the NS (Figure 6). This axis is deﬁned as the z-axis in this system. At each Δξ step, the intensities Ii of the HSs are propagated in space obeying Equations (1) or (2), and then they are used to populate the resulting pulse proﬁle. 4 An interactive web interface for the Polestar code is publicly available at the address http://www.polestar.live. 4 Publications of the Astronomical Society of the Paciﬁc, 129:124201 (11pp), 2017 December Cappallo, Laycock, & Christodoulou Table 1 Input and Output Parameters in the Polestar Code Symbol Name Description Input Parameters z L ξ Δξ f i0 mi θi γi ψ rg/RNS α Ii Cp Cf Figure 6. Vectors of the model. The z axis is the NS rotation axis. Two additional angles are not marked for clarity: angle ψ in Equations (1) and(2) is the angle between the two vectors, and the inclination angle i0 is the colatitude of f. (A color version of this ﬁgure is available in the online journal.) Fi å Fi ξ i=1, 2, ... Free parameter (Figure 6) Free parameter (Figure 6) Figure 2 Equations (1) and(3) Equations (1) and(2) i=1, 2, ... Equations (11) and(12) Equations (11) and(12) Flux from each HS to observer Sum of ﬂuxes from all HSs to observer Phase of the proﬁle i=1, 2, ... For each Δξ step Modulo 2π L (x + Dx ) = R (Dx ) L (x ). (4 ) (7 ) Then a ﬂux value F is calculated by summing the two contributions from the HSs along the new LOS vector: F = (L · m1) I1(a1) + (L · m 2) I2 (a2) , (8 ) where I1 and I2 are the beam intensities of the two HSs. The intensities are free parameters of the model and they may also be functions of the angles αi between the LOS vector and each HS vector mi (Nagel 1981). The dot products must be positive. If any of them is negative, then the corresponding HS is not visible and the term is reset to zero. The above process is iterated for a total of 2π/Δξ steps in order to produce a complete pulse proﬁle, and then F is normalized by the maximum computed ﬂux value. A detailed description of the steps in the code is given in the following section and a complete listing of the input and the output parameters of Polestar is shown in Table 1. (5 ) The rotation matrix R (Dx ) is used to rotate L through an angle Δξ, where ⎛ cos (Dx ) -sin (Dx ) 0 ⎞ R (Dx ) = ⎜⎜ sin (Dx ) cos (Dx ) 0 ⎟⎟. ⎝ 0 0 1⎠ Modulo 2π Free parameter (Figure 6) Colatitude of f, Equation (13) In the typical case with two HSs and only pencil-beam emission, each iteration by Δξ involves the following sequence of steps: First, L(ξ) is rotated about the z-axis by an angle Δξ to L (x + Dx ): and the HS vectors as ⎛ cos qi cos gi ⎞ ⎜ ⎟ mi = ⎜ sin qi cos gi ⎟. ⎝ sin gi ⎠ This axis is ﬁxed Initially set on xz-plane Initially set to zero (Figure 6) Output Parameters The orientation of the LOS vector L is speciﬁed by the latitude f and the azimuthal angle x = å Dx (modulo 2π); this vector always begins on xz-plane (ξ = 0) of the rest frame. This is done for simplicity and manifests itself as a deﬁnition of phase which can, in turn, be shifted to any other time once the proﬁle is complete. Furthermore, the orientation of each HS vector mi is speciﬁed by a latitude γi and an azimuthal angle θi. Based on Figure 6, the LOS vector can be expressed as ⎛ cos x cos f ⎞ ⎜ ⎟ L = ⎜ sin x cos f ⎟ , ⎜ ⎟ ⎝ sin f ⎠ Rotation axis of NS LOS vector, ∣L∣ = 1 Azimuth of L and phase of the proﬁle Advance of ξ to new phase Latitude of LOS vector Inclination of LOS to NS spin axis z HS vectors, ∣ mi ∣ = 1 Azimuth of mi vector Latitude of mi vector Angle between LOS and mi Compactness parameter Angle between initial photon path and mi HS intensities Pencil-beam coefﬁcient of Ii Fan-beam coefﬁcient of Ii (6 ) The size of the increment Δξ is an input parameter and it allows for the model to have higher resolution than the observed pulse proﬁles. 5 Publications of the Astronomical Society of the Paciﬁc, 129:124201 (11pp), 2017 December Cappallo, Laycock, & Christodoulou 3.2. Computer Code The input parameters to Polestar (see Table 1) are held ﬁxed if a speciﬁc geometry is desired, or they are varied if the code is in the process of ﬁtting an observed proﬁle. There is also an option to ﬁx any number of the input parameters if the user wishes to explore the inﬂuence of the remaining parameters to the overall proﬁle. On initialization, the code calculates the initial angles between the HS unit vectors mi and the LOS unit vector L. After initialization, Polestar runs an iterative loop over angle ξ (modulo 2π) in which the calculation of the pulse proﬁle takes place. The iterative loop carries out the following sequence of computations: (1) It updates the angles ξ and ψ (Table 1). (2) It calculates angle α (Equation (1) or Equation (2)). (3) It calculates the dot and cross products between L and mi , that is L · m i = cos ai for each pencil beam and L ´ m i = sin ai for each fan beam, respectively. If any of these values is negative, then the corresponding beam is not visible to the observer and the product is reset to zero. (4) The pencil and fan intensity components from each HS are calculated and added together in order to obtain the total contribution from each HS (see Equation (11) below). (5) The total ﬂux at this phase is ﬁnally determined by summing the intensities of all the HSs. Figure 7. Emission pattern from a single HS with ﬂux proportional to either a cosine (simple pencil in red) or a sine (simple fan in blue) function. Here, I1=1 in Equations (9) and(10). The dashed line represents the magnetic axis through the two HSs. (A color version of this ﬁgure is available in the online journal.) of fan to pencil ﬂuxes shown in Figure 7, viz. F1,total (a1) = (Cp cos a1 + Cf sin a1) I1(a1) , where Cp and Cf are tuning coefﬁcients that obey the relation Cp + Cf = 1. In the next step, L is rotated through an angle Δξ and the computation runs again for the new phase. The simplest and most widely used approximation is that of emission from a single point on the NS surface (e.g., Beloborodov 2002; Karino 2007). This mimics a single pencil beam and can be achieved by taking the inner product between the LOS vector L and the HS vector m1; then the ﬂux is (9 ) where ∣ L ∣ = ∣ m1 ∣ = 1, angle α1 is determined from Equation (1) orEquation (2) for a given LOS, and cos a1 must be positive (otherwise it is reset to zero). On the other hand, the geometry of the fan-beam component can be represented by a sine function (Nagel 1981; Wang & Welter 1981; Karino 2007), in which case we can write its ﬂux as F1,fan (a1) = I1(a1) sin a1, (12) More complex beam architectures can also be accommodated through manipulation of the beam-generating function within the code. This function effectively controls the intensities Ii of any number of HSs and allows for more complicated radiation patterns than the simple sin ai and cos ai geometric dependencies of the ﬂuxes discussed above (see also Section 3.4). For example, narrowing of the emission peaks can be accomplished by using higher powers of these sine and cosine functions in Ii(αi), and the higher the power the narrower the peak. This is illustrated in Figure 8, where the solid lines show ﬂuxes F1 = cos10 a1 (Cp = 1) and F1 = sin10 a1 (Cf = 1) from a single HS; in comparison to the ﬂuxes in Figure 7 (plotted again in Figure 8 as dotted lines). 3.3. Hot Spot Emission F1,pencil (a1) = I1(a1) cos a1, (11) 3.4. Direction of Emission and Beaming Functions The model is ﬂexible in accommodating emission from an offset magnetic dipole (see Figure 9): by subtracting the two HS vectors m1 and m2 , we obtain a new vector along the chord joining the two poles. The direction of this vector deﬁnes the orientation of the offset magnetic-dipole axis (vertical blue line in Figure 9). The contribution to the observed ﬂux from each HS can then be calculated by unitizing this vector and by (10) where sin a1 must be positive (otherwise it is reset to zero). An overall emission pattern can then be realized by tuning the ratio 6 Publications of the Astronomical Society of the Paciﬁc, 129:124201 (11pp), 2017 December Cappallo, Laycock, & Christodoulou Figure 10. Antipodal conﬁguration with the magnetic (green vectors) and spin (red lines) axes offset by 90°. (A color version of this ﬁgure is available in the online journal.) Figure 8. Comparison of the simple fan and the simple pencil beams of Figure 7 (dotted blue and dotted red line, respectively) to the beams with ﬂuxes proportional to a sine and cosine raised to the tenth power (solid lines, respectively). Here, I1 = cos9 a1 and I1 = sin9 a1 in Equations (9) and(10), respectively. Note the ﬂattening of the trough near the magnetic axis (dashed line) in the fan-beam sine function (solid blue line) and the steepening of both beams around angles 0 and p 2 . (A color version of this ﬁgure is available in the online journal.) Figure 11. Flux from antipodal HSs with identical pencil beams. The dashed lines correspond to the individual contributions, proportional to cos ai , from each HS. (A color version of this ﬁgure is available in the online journal.) particular intensities used in the original work. So this model serves as a test of the geometric features of our code, but it does not test the beam functions that we commonly use in our modeling. In the tests described in Sections 2.2 and 3.3 above, we used some arbitrary geometric functions of αi for demonstration purposes. But for detailed modeling of observed pulse proﬁles, we rely on results of anisotropic radiative transfer calculations and beaming calculations (Nagel 1981; Kaminker et al. 1982; DeDeo et al. 2001). These results support intensity functions of the forms (Equations(2) and(3) in Nagel 1981): Ii µ cos2 ai (pencil); Ii µ sin2 ai (fan); and in some cases, even higher powers, depending on choices for the mean-free-path functions of the emitted photons. Figure 9. Diagram of a NS (circle) with an offset dipolar magnetic ﬁeld. The red line denotes the rotation axis through the NS center. The green vectors point to the HSs, and their difference deﬁnes the direction of the magnetic axis (vertical blue line). (A color version of this ﬁgure is available in the online journal.) projecting it and its opposite vector along the LOS vector (as in the models shown in Figures 16–18 below). This conﬁguration, along with more complicated beam functions, was used successfully by Parmar et al. (1989) in the detailed modeling of the pulse proﬁle of the 42 s X-ray pulsar EXO2030+375. We have reproduced this result with our code and an offset magnetic ﬁeld, but we had to also incorporate the 7 Publications of the Astronomical Society of the Paciﬁc, 129:124201 (11pp), 2017 December Cappallo, Laycock, & Christodoulou Figure 12. Flux from antipodal HSs with identical fan beams. The dashed lines correspond to the individual contributions, proportional to sin ai , from each HS. (A color version of this ﬁgure is available in the online journal.) Figure 14. As in Figure 12, but the ﬂuxes from the two fan beams are proportional to (sin ai )10 . (A color version of this ﬁgure is available in the online journal.) Figure 13. As in Figure 11, but the ﬂuxes from the two pencil beams are proportional to (cos ai )10 . (A color version of this ﬁgure is available in the online journal.) Figure 15. Combination of (sin ai )10 and (cos ai )10 ﬂuxes from two identical antipodal HSs with Cp=0.7 and Cf=0.3 giving rise to four peaks in the pulse proﬁle. The dashed lines correspond to the individual contributions from each HS. (A color version of this ﬁgure is available in the online journal.) 4. Variation of Some Free Parameters and Proﬁle Features (5)The compactness parameter rg/RNS (Figure 3) can also generally be a free parameter (although in this paper we hold it constant to its canonical value of 0.4135). The output of the code is the total ﬂux that reaches the observer as a function of phase at an arbitrary resolution (currently set to 1000 points per period), as well as a graphical pulse proﬁle and an animation of the geometry of the pulsar. In the following subsections we describe various proﬁle features that can be produced by variations of some of these free parameters. For simplicity, we assume an antipodal For the sake of clarity, we begin by listing the free parameters and the output products of our model (Table 1): (1)The location of each HS is designated by a latitudinal angle γi and a longitudinal angle θi (Figure 6). (2)Similarly, the location of the LOS is speciﬁed by a latitude f and a longitude ξ (Figure 6), where ξ represents the phase of the overall proﬁle as the NS spin axis remains ﬁxed in time. (3)The intensities of the HSs Ii(αi) are functions of the emission angles αi (Figure 2). (4)The coefﬁcients Cp and Cf control the pencil and beam contribution from each HS. 8 Publications of the Astronomical Society of the Paciﬁc, 129:124201 (11pp), 2017 December Cappallo, Laycock, & Christodoulou Figure 16. Pulse proﬁle generated by two HSs emitting pencil-beam radiation (Fi = I0 cos ai ). The individual HSs are shown in red and green, and their sum is shown in black (see also Figure 5 for visualization of the geometry). The HSs are asymmetric (f = 20°; γ1 = 70° and γ2 = −60°; and θ1 = 0 and θ2 = 90°). (A color version of this ﬁgure is available in the online journal.) Figure 18. Our model of the pulse proﬁle of Cen X-3 that was originally produced by Wang & Welter (1981) (Figure 7(j) in their paper). The model uses pure fan emission with Fi = I0 sin ai from two identical HSs and our LOS is inclined by i=31° to the pulsar’s spin axis, allowing for the secondary HS (green dashed line) to become visible for a small fraction of the phase (see also Figure 19 for visualization of the geometry). Model parameters are given in the text. (A color version of this ﬁgure is available in the online journal.) following are ﬁve examples of the inherent differences between pencil and fan behavior from the two identical HSs and the effects of using powers of the trigonometric functions for the intensities of the beams. 4.1. Pencil Versus Fan In Figure 11, the pencil-beam peaks are formed as each pole moves across the LOS. Since the emission comes from a pure cosine function, the total ﬂux is constant whenever both poles are in view, leading to the plateaus between the two peaks (see also Beloborodov 2002). In Figure 12, the fan-beam peaks are formed when the angle between the magnetic-ﬁeld axis and the LOS is 90°, with the eclipse-like features occurring when a given pole travels across the terminator on the far side of the NS. Note that the fan-beam peaks are shifted in Figure 12 by a quarter of the phase in comparison to the pencil proﬁle of Figure 11. This shift can account for quadruple-peaked structures in cases where the proﬁle is a combination of the two beams. Figure 17. As in Figure 16, but for fan-beam emission (Fi = I0 sin ai ). Note the eclipse-like behavior as the southern pole is brieﬂy obscured when traveling across the terminator (see also Figure 5 for visualization of the geometry, as well as Wang & Welter 1981). (A color version of this ﬁgure is available in the online journal.) arrangement with both the LOS and the magnetic axis offset initially by 90° from the spin axis (Figure 10). The model parameters for this conﬁguration are: f=0, implying that the inclination of the LOS to the NS spin axis is i 0 = 90, where i 0 º 90 - f ; 4.2. Narrowing of the Beam In Figure 13, the geometry is the same as that in Figure 11, but the cos ai function for each HS has been raised to the tenth power. Note that the peaks are now signiﬁcantly narrower and that the plateaus between the peaks have vanished—when both poles are visible, constant ﬂux is seen only in the strict Fi µ cos ai case. In this fashion, a geometric aspect of energy (13) γ1=γ2=0, θ1=0, and θ2=180°. This simple conﬁguration is chosen because it can be visualized easily. The 9 Publications of the Astronomical Society of the Paciﬁc, 129:124201 (11pp), 2017 December Cappallo, Laycock, & Christodoulou dependence can be accommodated in the code because the opening angle of the pencil beam is proportional to the square root of the energy (e.g., Kaminker et al. 1982; Pavlov et al. 1994; Ho & Lai 2001). Yet, other complicated processes (e.g., self-absorption, scattering, and reprocessing of the X-rays) that are energy-dependent and modify observed pulse proﬁles cannot be accounted for in this model. Similarly, in Figure 14 the sin ai function for the fan beam of each HS has been raised to the tenth power. Again the proﬁle peaks are narrower, and, in addition, the vertical lines seen in Figure 12 created by the terminator are here not as steep. This makes sense as the positions of these lines in the pulse proﬁle are dependent only on the geometry, and with the narrowing of the peaks there is less total ﬂux near the maxima. Figure 19. Visual representation of the geometry of Cen X-3 that produces the pulse proﬁle shown in Figure 18. (A color version of this ﬁgure is available in the online journal.) 5.2. Wang & Wekter (1981) Proﬁles Wang & Welter (1981) were the ﬁrst to describe the eclipselike cutoffs that one sees in proﬁles dominated by fan-beam emission (Fi = I0 sin ai ) from two identical HSs. Our model reproduces this behavior using pure fan-beam radiation: compare the sharp eclipse in Figure 17 (pure fan beam) to the gradual dip in Figure 16 (pure pencil beam). The eclipse and the dip are both due to the secondary HS traveling brieﬂy across the terminator. But the eclipse in Figure 17 occurs when the total ﬂux is near maximum, as the HS disappears from view abruptly. We note that, without gravitational light-bending, this eclipse would occur precisely at maximum ﬂux. Another test case is shown in Figure 18 in which we have reproduced a model ﬁt to the pulse proﬁle of the source Cen X-3 originally published by Wang & Welter (1981) (their Figure7(j)). The model parameters for this proﬁle are: θ1=44° and γ1=60° (primary HS); θ2=220° and γ2=−74° (secondary HS); f=59°, thus inclination of the LOS to the NS rotation axis i 0 = 31; Cp=0 and Cf=1; and Fi = I0 sin ai . The visual representation of this conﬁguration is shown in Figure 19. In the phase interval of approximately 0.3–0.5, the secondary HS is peaking out at the bottom of the NS and this causes the pronounced plateau seen in the pulse proﬁle near maximum ﬂux. 4.3. Four Peaks in the Proﬁle Figure 15 is an example showing a quadruple-peaked proﬁle that occurs when a combination of pencil and fan emission is used from two identical antipodal HSs. In this model with Cp=0.7 and Cf=0.3, the contribution from the pencil beams is dominant and not only due to the particular choice of the coefﬁcients. For the orthogonal orientation used in this model, the absolute maxima at phases 0 and 0.5 are produced by the pencil beam of each HS as it crosses the LOS and the ﬂux tapers off gradually as the HS travels into the far side of the NS; whereas the fan beam shuts off (comes up) abruptly as each HS moves into (out of) the terminator of the NS. On the other hand, both fan beams become fully visible at phases 0.25 and 0.75, where the contributions of the pencil beams are minimal and the total ﬂux rises to local maxima of »2Cf Cp = 0.857. 5. Tests and Comparisons with Previous Results Our code can reproduce several landmark results from the literature when we adopt the corresponding input parameters. Two speciﬁc test cases are described in this section. 5.1. Beloborodov (2002) Proﬁles 6. Discussion Beloborodov (2002) deﬁnes four different classes of pulse proﬁles, depending on the location of the LOS and the number of poles that are visible in each case (the latter depends on the inclination of the magnetic-ﬁeld axis to the rotation axis). His model assumes two antipodal, point-source HSs of equal intensity I0 emitting beamed radiation. In our model, this is duplicated by two antipodal HSs and pure pencil-beam radiation (Fi = I0 cos ai ). The resulting pulse proﬁles are identical to those published by Beloborodov (2002). For example, compare our Figure 4(b) to Figure 4 (class III) in that work. In this paper, we have described a new code and its tests for simulations of pulse proﬁles of pulsars. The code uses a wide variety of geometries for the relativistic emission regions and the dipolar magnetic ﬁeld, as well as different orientations for the observer. The parameter space is very large with typical simulations ﬁtting simultaneously 10 or more free parameters to an observed proﬁle. This ﬂexibility allows us to ﬁt the observed pulse proﬁles of high-mass X-ray binary (HMXB) pulsars and other XRBs with relative ease, but the results may also be subject to degeneracies (different regions of the parameter space effectively producing the same result; see, 10 Publications of the Astronomical Society of the Paciﬁc, 129:124201 (11pp), 2017 December Cappallo, Laycock, & Christodoulou e.g., the degenerate class IV ﬂat proﬁles in Beloborodov 2002). Our tests of the code have reproduced well-known results previously published by other researchers in the ﬁeld (Wang & Welter 1981; Parmar et al. 1989; Beloborodov 2002). The X-ray emission model of the new code is purely geometrical in nature; i.e., there is no radiative transfer (cf., Nagel 1981; Kaminker et al. 1982; Pavlov et al. 1994; Ho & Lai 2001) or input from physical properties of the emitting regions (Ho 2007); nevertheless, it is capable of reproducing many of the features observed in XRB pulse proﬁles; also the proﬁles of the same object viewed in different energy bands can be modeled in order to differentiate between proﬁles that have strong energy and X-ray luminosity dependence. We plan to use this code in conjunction with an automated ﬁtting algorithm in order to model the pulse proﬁles in a large database of the 80 known HMXB pulsars in the Magellanic Clouds (Coe et al. 2010; Coe & Kirk 2015; Klus 2015; Haberl & Sturm 2016; Yang et al. 2017). Even if some degeneracies are encountered in this type of analysis, such a large sample of pulsars will certainly be broken down to various subclasses each characterized by similar pulsar parameters (mass, radius, spin, magnetic ﬁeld) or proﬁle shapes. In conjunction with additional data such as X-ray spectra, the companions’ optical spectra, and typical magnetic-ﬁeld magnitudes such as those proposed by Christodoulou et al. (2017b), these subclasses can then provide information about the physical processes and the pulsar geometries that are responsible for the various types of pulsations. Christodoulou, D. M., Laycock, S. G. T., & Kazanas, D. 2017a, MNRAS, 470, L21 Christodoulou, D. M., Laycock, S. G. T., Yang, J., & Fingerman, S. 2016, ApJ, 829, 30 Christodoulou, D. M., Laycock, S. G. T., Yang, J., & Fingerman, S. 2017b, RAA, 17, 59 Coe, M. J., & Kirk, J. 2015, MNRAS, 452, 969 Coe, M. J., McBride, V. A., & Corbet, R. H. D. 2010, MNRAS, 401, 252 Cordes, J. M., Kramer, M., Lazio, T. J. W., et al. 2004, NewA, 48, 1413 DeDeo, S., Psaltis, D., & Narayan, R. 2001, ApJ, 559, 346 Eggleton, P. 2001, Evolutionary Processes in Binary and Multiple Stars (Cambridge: Cambridge Univ. Press) Ferrigno, C., Falanga, M., Bozzo, E., et al. 2011, A&A, 532, A76 Frank, J., King, A., & Raine, D. 2002, Accretion Power in Astrophysics (3rd ed.; Cambridge: Cambridge Univ. Press) Haberl, F., & Sturm, R. 2016, A&A, 586, A81 Ho, W. C. G. 2007, MNRAS, 380, 71 Ho, W. C. G., & Lai, D. 2001, MNRAS, 327, 1081 Hong, J., Antoniou, V., Zezas, A., et al. 2017, ApJ, arXiv:1705.09733 Kaminker, A. D., Pavlov, G. G., & Shibanov, Yu. A. 1982, Ap&SS, 86, 249 Karino, S. 2007, PASJ, 59, 961 Klus, H. 2015, PhD thesis, Univ. Southampton Klus, H., Ho, W. C. G., Coe, M. J., Corbet, R. H. D., & Townsend, L. J. 2014, MNRAS, 437, 3863 Kraus, U., Nollert, H.-P., Ruder, H., & Riffert, H. 1995, ApJ, 450, 763 Kylaﬁs, N. D., Contopoulos, I., Kazanas, D., & Christodoulou, D. M. 2012, A&A, 538, A5 Leahy, D. A., & Li, L. 1995, MNRAS, 277, 1177 Lewin, W. H. G., Paradijs, J., & van der Heuvel, E. P. J. (ed.) 1995, X-ray Binaries (Cambridge: Cambridge Univ. Press), 27 Lyne, A., & Graham-Smith, F. 2012, Pulsar Astronomy (4th ed.; Cambridge: Cambridge Univ. Press) Mereghetti, S. 2001, in Frontier Objects in Astrophysics and Particle Physics, ed. F. Giovannelli & G. Mannocchi (Italian Physical Society, Editrice Compositori, Bologna, Italy), 73, 239 Mereghetti, S., & Stella, L. 1995, ApJL, 442, L17 Mészáros, P., & Nagel, W. 1985a, ApJ, 298, 147 Mészáros, P., & Nagel, W. 1985b, ApJ, 299, 138 Mészáros, P., & Riffert, H. 1988, ApJ, 327, 712 Nagel, W. 1981, ApJ, 251, 278 Parmar, A. N., White, N. E., & Stella, L. 1989, ApJ, 338, 373 Pavlov, G. G., Shibanov, Yu. A., Ventura, J., & Zavlin, V. E. 1994, A&A, 289, 837 Pechenick, K. R., Ftaclas, C., & Cohen, J. M. 1983, ApJ, 274, 846 Reig, P. 2011, Ap&SS, 332, 1 Riffert, H., & Meszáros, P. 1988, ApJ, 325, 207 Riffert, H., Nollert, H.-P., Kraus, U., & Ruder, H. 1993, ApJ, 406, 185 Stella, L., White, N. E., & Rosner, R. 1986, ApJ, 308, 669 Tendulkar, S. P., Fürst, F., Pottschmidt, K., et al. 2014, ApJ, 795, 154 Thompson, D. J. 2000, AdSpR, 25, 659 Wang, Y. M., & Welter, G. L. 1981, A&A, 102, 97 Yang, J., Laycock, S. G. T., Christodoulou, D. M., et al. 2017, ApJ, 839, 119 Zavlin, V. E., Pavlov, G. G., Shibanov, Y. A., & Ventura, J. 1995, A&A, 297, 441 We are obliged to Malcolm Coe, Wynn Ho, and Helen Klus for their continuing assistance and advice during this project. We also thank the referee for a thorough review that has greatly improved and clariﬁed this paper. This work was supported in part by NASA grant NNX14-AF77G. References Bachetti, M., Harrison, F. A., Waltron, D. J., et al. 2014, Natur, 514, 202 Basko, M. M., & Sunyaev, R. A. 1976, MNRAS, 175, 395 Becker, P., & Wolff, M. 2005, ApJ, 630, 465 Becker, P. A. 1998, ApJ, 498, 790 Beloborodov, A. M. 2002, ApJL, 566, L85 Burnard, D. J., Klein, R. I., & Arons, J. 1988, ApJ, 324, 1001 11

1/--страниц