Classical and general relativistic post-Keplerian effects in binary pulsars hosting fast rotating main sequence stars

TL;DR

This paper investigates post-Keplerian effects in binary pulsar systems with fast rotating, distorted main sequence stars, using numerical simulations to assess the magnitude of relativistic and classical perturbations on pulsar timing over ten years.

Contribution

It introduces a numerical approach to quantify Newtonian and post-Newtonian perturbations in pulsar timing caused by distorted stellar gravitational fields in such binaries.

Findings

Classical oblateness effects cause timing shifts up to 150 seconds.

Post-Newtonian gravitoelectric shifts are up to 20 seconds.

Most pulsars in these systems have residuals less than 1 millisecond.

Abstract

We consider a binary system composed of a pulsar and a massive, fast rotating, highly distorted main sequence star as a potential scenario to dynamically put to the test certain post-Keplerian effects of both Newtonian and post-Newtonian nature. We numerically produce time series of the perturbations $\Delta\left(\delta\tau\right)$ of the R{\o}mer-like, orbital component of the pulsar's time delay $\delta\tau$ induced over 10 years by the pN gravitoelectric mass monopole, quadrupole, gravitomagnetic spin dipole and octupole accelerations along with the Newtonian quadrupolar one. We do not deal with the various propagation time delays due to the travelling electromagnetic waves. It turns out that, for a Be-type star with $M = 15\ \textrm{M}_\odot$, $R_\textrm{e} = 5.96\ \textrm{R}_\odot$, $\nu = 0.203$, $S = 3.41\times 10^{45}\ \textrm{J}\ \textrm{s}$, $J_2 = 1.92\times 10^{-3}$ orbited by a pulsar with an orbital period $P_\textrm{b}\simeq 40-70\ \textrm{d}$, the classical oblateness-driven effects are at the $\lesssim 4-150\ \textrm{s}$ level, while the pN shifts are of the order of $\lesssim 1.5-20\ \textrm{s}\ \left(GMc^{-2}\right)$, $\lesssim 10-40\ \textrm{ms}\ \left(GMR^2_\textrm{e} J_2 c^{-2}\right)$, $\lesssim 0.5 - 6\ \textrm{ms}\ \left(GSc^{-2}\right)$, $\lesssim 5 - 20\ \mu\textrm{s}\ \left(GSR^2_\textrm{e} \varepsilon^2 c^{-2}\right)$, depending on their orbital configuration. The root-mean-square (rms) timing residuals $\sigma_{\tau}$ of almost all the existing non-recycled, non-millisecond pulsars orbiting massive, fast rotating main sequence stars are $\lesssim\textrm{ms}$. Thus, such kind of binaries have the potential to become interesting laboratories to measure, or, at least, constrain, some Newtonian and post-Newtonian key features of the distorted gravitational fields of the fast rotating stars hosted by them [Abridged].

Tables

Table 1. Table 1: Binary pulsars with a massive, fast rotating main sequence star companion discovered so far. For each of them, we list the companion, the distance (when available), the mass M 𝑀 M , the equatorial and polar radii R e , R p subscript 𝑅 e subscript 𝑅 p R_{\textrm{e}},\,R_{\textrm{p}} (when available), the orbital period P b subscript 𝑃 b P_{\textrm{b}} , the orbital eccentricity e 𝑒 e , the spin period P 𝑃 P , the rms timing residuals s τ subscript s 𝜏 {\mathchoice{\mbox{s}}{\mbox{s}}{\mbox{s}}{\mbox{s}}}_{\tau} , the data analysis time span Δ ​ T Δ 𝑇 \Delta T . For the values of the listed parameters, see the references cited in the text and the online pulsar catalog at http://www.atnf.csiro.au/research/pulsar/psrcat/. For PSR B1256-63 we give the equatorial R e subscript 𝑅 e R_{\textrm{e}} and polar R p subscript 𝑅 p R_{\textrm{p}} radii as estimated by Negueruela et al . ( 2011 ) .

Pulsar	Companion	Distance (kpc)	$M\left({\text{M}}_{\odot }\right)$	$R_{\text{e}}\left({\text{R}}_{\odot }\right)$	$R_{\text{p}}\left({\text{R}}_{\odot }\right)$	$P_{\text{b}}$ (d)	$e$	$P$ (ms)	${\text{s}}_{\tau }$ (ms)	$\Delta T$ (yr)
PSR J0045-7319	B1V star	In SMC	$8.8\pm 1.8$	$6.4\pm 0.7$	$-$	$51.17$	$0.808$	$930$	$7.4$	$2$
PSR J1740-3052	Main sequence star	$-$	$>11$	$-$	$-$	$231$	$0.578$	$570$	$0.8$	$2.29$
PSR B1259-63	B2e star LS 2883	$2.75$	$\simeq 30$	$\simeq 9.7$	$\simeq 8.1$	$1237$	$0.870$	$48$	$0.46$	$13$
PSR J1638-4725	Be star	$-$	$>4$	$-$	$-$	$1940$	$0.95$	$764$	$5.3$	$4.35$
PSR J2032+4127	Be star MT91 213	$1.7$	$\simeq 15$	$-$	$-$	$8578$	$0.93$	$143$	$0.5-1$	$6$

Table 2. Table 2: ESTER models for Zero-Age Main Sequence (ZAMS) stars with an equatorial angular velocity at 70% of the critical angular velocity. R e subscript 𝑅 e R_{\textrm{e}} and V e subscript 𝑉 e V_{\textrm{e}} are the equatorial radius and velocity, respectively, S 𝑆 S is the total spin angular momentum, ν 𝜈 \nu is the flattening, and J 2 subscript 𝐽 2 J_{2} is the dimensionless mass quadrupole moment. The metallicity is Z=0.02 and the hydrogen mass fraction is X = 0.7 𝑋 0.7 X=0.7 . In our simulations, we use the parameters of the star with 15 ​ M ⊙ 15 subscript M direct-product 15\,\textrm{M}_{\odot} listed here.

$M\left({\text{M}}_{\odot }\right)$	$R_{\text{e}}\left({\text{R}}_{\odot }\right)$	$V_{\text{e}}\left(\text{km}{\text{s}}^{-1}\right)$	$S(\times {10}^{44}\text{J}\text{s})$	$\nu$	$J_2(\times {10}^{-3})$
10	4.74	444	15.3	0.201	1.63
15	5.96	485	34.1	0.203	1.92
30	8.89	562	125.	0.210	2.17

Table 3. Table 3: Same as in table  2 but for ESTER models of stars at mid-main-sequence, namely when the hydrogen mass fraction in the convective core is half of the initial one.

$M\left({\text{M}}_{\odot }\right)$	$R_{\text{e}}\left({\text{R}}_{\odot }\right)$	$V_{\text{e}}\left(\text{km}{\text{s}}^{-1}\right)$	$S(\times {10}^{44}\text{J}\text{s})$	$\nu$	$J_2(\times {10}^{-3})$
10	6.93	367	13.3	0.203	0.816
15	9.07	393	28.0	0.207	0.788
30	14.8	435	870.	0.227	0.495

Equations

$$U\left(\mathbf{r}\right)=-\frac{\mu }{r}\left[1-\sum _{\ell }{J}_{\ell }{\left(\frac{R_{\text{e}}}{r}\right)}^{\ell }{\mathcal{P}}_{\ell }\left(\xi \right)\right],$$

$$J_{2p}=-\frac{1}{M{R}_{\text{e}}^{2p}}{\int }_{\mathcal{V}}{r}^{2p}{\mathcal{P}}_{2p}\left(\xi \right)\rho \left(\mathbf{r}\right)d^3\mathbf{r}$$

$$J_2=-\frac{1}{M{R}_{\text{e}}^2}{\int }_{\mathcal{V}}{r}^2{\mathcal{P}}_2\left(\xi \right)\rho \left(\mathbf{r}\right)d^3\mathbf{r}\mathrm{and}S={\int }_{\mathcal{V}}{r}^2(1-{\xi }^2)\Xi \left(\mathbf{r}\right)\rho \left(\mathbf{r}\right)d^3\mathbf{r}$$

$$U\left(\mathbf{r}\right)=U_0+\Delta U_2=-\frac{\mu }{r}\left[1-{\left(\frac{R_{\text{e}}}{r}\right)}^2{J}_2{\mathcal{P}}_2\left(\xi \right)\right],$$

$${\mathbf{A}}^{\text{N}{J}_2}=-\nabla \Delta U_{J_2}=\frac{3\mu J_2{R}_{\text{e}}^2}{2r^4}\left[\left(5{\xi }^2-1\right)\hat{\mathbf{r}}-2\xi \hat{\mathbf{S}}\right].$$

$${\mathbf{A}}^{\textrm{1pN}M}=\frac{\mu}{c^{2}r^{2}}\left[\left(\frac{4\mu}{r}-{\mathrm{v}}^{2}\right)\mathbf{\hat{r}}+4{\mathrm{v}}_{r}{\color[rgb]{0,0,0}\mathbf{v}}\right],$$

$${\mathbf{A}}^{\textrm{1pN}S}=\frac{2GS}{c^{2}r^{3}}\left[3\xi\mathbf{\hat{r}}\mathbf{\times}{\color[rgb]{0,0,0}\mathbf{v}}+{\color[rgb]{0,0,0}\mathbf{v}}\mathbf{\times}\mathbf{\hat{S}}\right].$$

$${\mathbf{A}}^{\textrm{1pN}MJ_{2}}=\frac{\mu J_{2}R_{\textrm{e}}^{2}}{c^{2}r^{4}}\left\{\frac{3}{2}\left[\left(5\xi^{2}-1\right)\mathbf{\hat{r}}-2\xi\mathbf{\hat{S}}\right]\left({\mathrm{v}}^{2}-\frac{4\mu}{r}\right)-6\left[\left(5\xi^{2}-1\right){\mathrm{v}}_{r}-2\xi{\mathrm{v}}_{S}\right]{\color[rgb]{0,0,0}\mathbf{v}}-\frac{2\mu}{r}\left(3\xi^{2}-1\right)\mathbf{\hat{r}}\right\},$$

$$\mathbf{A}^{\textrm{1pN}SJ_{2}}=\frac{3GSR_{\textrm{e}}^{2}\varepsilon^{2}}{7c^{2}r^{5}}{\color[rgb]{0,0,0}\mathbf{v}}\mathbf{\times}\left\{5\xi\left[7\xi^{2}-3\right]\mathbf{\hat{r}}+3\left[1-5\xi^{2}\right]\mathbf{\hat{S}}\right\}.$$

$$\displaystyle{\mathbf{A}}^{\textrm{1pN}M}=\frac{\mu_{\textrm{tot}}}{c^{2}r^{2}}\left\{\left[\left(4+2\zeta\right)\frac{\mu_{\textrm{tot}}}{r}-\left(1+3\zeta\right){\mathrm{v}}^{2}+\frac{3}{2}\zeta{\mathrm{v}}_{r}^{2}\right]\mathbf{\hat{r}}+\left(4-2\zeta\right){\mathrm{v}}_{r}{\color[rgb]{0,0,0}\mathbf{v}}\right\},$$

$$\mathbf{S}\to \left(1+\frac{3}{4}\frac{M_{\text{p}}}{M}\right)\mathbf{S}+\left(1+\frac{3}{4}\frac{M}{M_{\text{p}}}\right){\mathbf{S}}_{\text{p}},$$

$${\mathcal{Q}}_2$$

$${\mathcal{Q}}_2^{\text{p}}$$

$$Q_2^{\text{p}}=-q\frac{M_{\text{p}}^3{G}^2}{c^4},$$

$$Q_2$$

$$Q_2^{\text{p}}$$

$${\dot{\Omega }}_{S_{\text{p}}}=\frac{3\mu {\sigma }_{\text{p}}\csc I\left({\hat{\sigma }}_{\text{p}}\cdot \hat{\mathbf{m}}\right)}{2c^2{a}^3{\left(1-e^2\right)}^{3/2}}.$$

$${\dot{\Omega }}_{\text{LT}}=\frac{2GS\csc I\left(\hat{\mathbf{S}}\cdot \hat{\mathbf{m}}\right)}{c^2{a}^3{\left(1-e^2\right)}^{3/2}}.$$

Full text

keywords gravitation $-$ binaries: general $-$ stars: rotation $-$ pulsars: general $-$ celestial mechanics

1 Introduction

In its weak-field and slow-motion approximation, general relativity predicts that, in addition to the time-honored post-Newtonian (pN) gravitoelectric and gravitomagnetic precessions induced by the mass $M$ (Schwarzschild) and the spin angular momentum $\mathbf{S}$ (Lense-Thirring) of the central body acting as source of the gravitational field, other pN gravitoelectric and gravitomagnetic orbital effects related to its oblateness arise as well (Soffel, 1989; Soffel et al., 1987; Heimberger, Soffel & Ruder, 1990; Brumberg, 1991; Will, 2014; Panhans & Soffel, 2014; Iorio, 2015; Meichsner & Soffel, 2015; Soffel & Frutos, 2016; Frutos-Alfaro & Soffel, 2018; Schanner & Soffel, 2018). So far, they have never been put to the test in any astronomical and astrophysical scenarios, despite some recent preliminary investigations pertaining the planet Jupiter in our solar system (Iorio, 2013, 2019); for some embryonic thoughts about an Earth-spacecraft scenario, see Iorio (2013, 2015). To the pN level, the oblateness of astronomical bodies modifies also the propagation of the electromagnetic waves in their deformed spacetime. About the perspectives of measuring the resulting deflection due to Jupiter with astrometric techniques, see, e.g., Crosta & Mignard (2006), Kopeikin & Makarov (2007), Le Poncin-Lafitte & Teyssandier (2008), Abbas, Bucciarelli & Lattanzi (2019), and references therein.

An analysis of the analytical expressions of the pN gravitoelectric and gravitomagnetic orbital precessions due to the asphericity of the primary (Iorio, 2015, 2019) shows that the key ingredients needed to enhance their magnitude are a strongly distorted central body, and a highly eccentric and close orbit of the moving particle.

Binaries composed by a pulsar and a main-sequence star (Wex, 1998) may offer, in principle, interesting natural laboratories to try to investigate such little known pN effects. Indeed, they are systems composed of a neutron star regularly emitting electromagnetic radio pulses orbiting an usually more massive main sequence star, which, in most cases, is highly distorted due to its fast rotation, along a generally elliptical orbit. A very accurate observable quantity in binary pulsars is represented by the measurement of the times of arrival (TOAs) $\tau$ of the pulses emitted by the neutron star which, in case it has a gravitationally bound companion, exhibit a regular variation $\delta\tau$ due to, among other things, the Keplerian motion about the common centre of mass: it is the Rømer-like time delay. The full variation of the pulses’ times of arrival is due to several other effects connected, e.g., with the propagation of the electromagnetic waves through the deformed spacetime of the system (Wex, 1998).

The first binary pulsar hosting a main sequence star to be discovered was PSR B1259-63 (Johnston et al., 1992; Shannon, Johnston & Manchester, 2014); it is characterized by a highly eccentric orbit ($e=0.870$) with an orbital period of $P_{\rm b}=1237\leavevmode\nobreak\ \textrm{d}=3.38\leavevmode\nobreak\ \textrm{yr}$. The pulsar’s non-degenerate companion is the fast spinning Be star LS 2883, whose equatorial velocity $V_{\textrm{e}}$ is about $280\leavevmode\nobreak\ \textrm{km\leavevmode\nobreak\ s}^{-1}$ corresponding to $\sim 70\%$ of its break-up velocity (Porter, 1996). Its mass $M$ and equatorial radius $R_{\textrm{e}}$ amounts to about 30 Solar masses $(\textrm{M}_{\odot})$ and $9.7$ Solar radii $(\textrm{R}_{\odot})$ (Negueruela et al., 2011). Later, the eccentric ($e=0.808$) binary PSR J0045-7319 was discovered (Kaspi et al., 1994). To date, it is the fastest orbiting system since it is $P_{\rm b}=51.17\leavevmode\nobreak\ \textrm{d}$. Its primary is a main sequence B-star spinning close to its break-up velocity (Lai, Bildsten & Kaspi, 1995; Kaspi et al., 1996). PSR J1638-4725, having an orbital period of $P_{\rm b}=1940\leavevmode\nobreak\ \textrm{d}=5.3\leavevmode\nobreak\ \textrm{yr}$ and $e=0.95$, was found by Lorimer et al. (2006). Its stellar companion should be a rapidly rotating Be star. PSR J1740-3052 (Stairs et al., 2001), with $P_{\rm b}=231\leavevmode\nobreak\ \textrm{d}$ and $e=0.578$, hosts most likely a B-type main-sequence star (Tam et al., 2010; Madsen et al., 2012). The most recently discovered main-sequence-star binary pulsar is the highly eccentric ($e=0.93$) PSR J2032+4127 (Lyne et al., 2015) characterized by $P_{\rm b}=8578\leavevmode\nobreak\ \textrm{d}=23.5\leavevmode\nobreak\ \textrm{yr}$. The companion of the neutron star is the massive Be star MT91 213 with $M\simeq 15\,\textrm{M}_{\odot}$.

For the sake of completeness, we mention also a few other binary pulsars hosting a non-degenerate star, although they are not relevant for our purposes in view of the nature of their non massive and fast-rotating partners. They are PSR J1903+0327 (Arzoumanian et al., 2018), whose companion is a F5V-GOV $\sim 1\leavevmode\nobreak\ \textrm{M}_{\odot}$ star moving in $P_{\rm b}=95\leavevmode\nobreak\ \textrm{d}$ along a rather eccentric orbit with $e=0.44$, the transitional millisecond pulsar PSR J1023+0038 (Archibald et al., 2009) orbiting a low-mass ($0.2\leavevmode\nobreak\ \textrm{M}_{\odot}$) companion star in a circular path with $P_{\rm b}=4.75\leavevmode\nobreak\ \textrm{hr}$. The rms timing residuals of the aforementioned binary pulsars are all of the order of $\lesssim\textrm{ms}$, apart from PSR J1903+0327 which is at the \simeq{\color[rgb]{0,0,0}1}\,{\mathchoice{\mbox{m}}{\mbox{m}}{\mbox{m}}{\mbox{m}}}\textrm{s} level; more specifically, they are $\simeq 0.46\leavevmode\nobreak\ \textrm{ms}$ over $13\leavevmode\nobreak\ \textrm{yr}$ for PSR B1259-63 (Wang, Johnston & Manchester, 2004), $7.4\leavevmode\nobreak\ \textrm{ms}$ over $2\leavevmode\nobreak\ \textrm{yr}$ for PSR J0045-7319 (Kaspi et al., 1994), $\simeq 5.3\leavevmode\nobreak\ \textrm{ms}$ over $4.35\leavevmode\nobreak\ \textrm{yr}$ for PSR J1638-4725 (Lorimer et al., 2006), $\simeq 0.8\leavevmode\nobreak\ \textrm{ms}$ over $2.29\leavevmode\nobreak\ \textrm{yr}$ for PSR J1740-3052 (Stairs et al., 2001), $\simeq 0.5-1\leavevmode\nobreak\ \textrm{ms}$ over about $6\leavevmode\nobreak\ \textrm{yr}$ for PSR J2032+4127 (Lyne et al., 2015), $\simeq 1\leavevmode\nobreak\ \textrm{\mbox{m}s}$ over about $3\leavevmode\nobreak\ \textrm{yr}$ for PSR J1903+0327 (Arzoumanian et al., 2018), $0.1\leavevmode\nobreak\ \textrm{ms}$ over about $4\leavevmode\nobreak\ \textrm{yr}$ for PSR J1023+0038 (Archibald et al., 2013). Table 1 summarizes the key data for the binary pulsars hosting massive, fast rotating main sequence stars.

[FIGURE:]

About the achievable accuracy level in timing residuals, in the case of the pulsars orbiting a main sequence star, their timing seems doomed to stay at the $\simeq\textrm{ms}$ level. It is so because, for evolutionary reasons, they are not fully recycled (Srinivasan, 2010). Thus, their spinning periods $P$ are not at the millisecond level, and their TOAs are not measured with a precision of the order of $\simeq{\mathchoice{\mbox{m}}{\mbox{m}}{\mbox{m}}{\mbox{m}}}\textrm{s}$, and their timing is often contaminated by timing noise (Hobbs, Lyne & Kramer, 2010). The timing of the non-recycled pulsars is almost always less accurate than for the fully recycled pulsars. Indeed, PSR J1903+0327 rotates with a period $P=2.5\,\textrm{ms}$, and its rms timing residuals are as little as $\simeq 1\,{\mathchoice{\mbox{m}}{\mbox{m}}{\mbox{m}}{\mbox{m}}}\textrm{s}$. The spin period of PSR J1023+0038 is $P=1.7\,\textrm{ms}$, and its rms timing residuals are $0.1\,\textrm{ms}$. Incidentally, we mention the fact that, according to Table 2 of Arzoumanian et al. (2018), the rms timing accuracy of some fully recycled pulsars with $P=1.5-10\,\textrm{ms}$, isolated or with a white dwarf as companion, is of the order of $0.1-0.2\,{\mathchoice{\mbox{m}}{\mbox{m}}{\mbox{m}}{\mbox{m}}}\textrm{s}$. It is expected that future instrumental improvements may push the rms timing accuracy of some of them to the111A. Possenti, personal communication to L.I., April 2019. $\simeq 10\,\textrm{ns}$ level over time spans some yr long. For binary pulsars hosting another neutron star, the rms timing accuracy is of the order of $1-100\,{\mathchoice{\mbox{m}}{\mbox{m}}{\mbox{m}}{\mbox{m}}}\textrm{s}$.

Here, we will preliminarily investigate the size and the temporal patterns of the perturbations $\Delta\left(\delta\tau\right)$ induced on the Rømer-like orbital part of the pulsar’s time delay $\delta\tau$ by both the standard (Schwarzschild and Lense-Thirring) and the oblateness-driven pN accelerations felt by a fictitious neutron star orbiting a highly distorted, fast rotating B-type main sequence star in view of a possible detection in new binaries that may eventually be discovered in the future. However, caution is in order before inferring too optimistic conclusions from a straightforward comparison of our simulated time series with the rms timing residuals listed in Table 1. Even if the size of some pK signatures were to be larger than the $\simeq\textrm{ms}$ level, it does not necessarily mean that such effects will be measurable in the actual processing of the real observations. Indeed, careful, dedicated simulated data reductions and covariance analyses should be performed by explicitly modeling the signals of interest, estimating its characteristic parameters and inspecting the resulting correlations with the other parameters usually estimated. It should be kept in mind that, in principle, an unmodeled effect may be removed from the post-fit residuals, at least to a certain extent, being “absorbed” in the estimated values of the other parameters determined in the data reduction. Thus, our investigation should be deemed just as a sensitivity analysis able to preliminarily explore the potential of the scenarios considered. We will assume the validity of general relativity throughout the paper, which is organized as follows.

In Section 2, we discuss the magnitude of the angular momentum $S$ and the dimensionless quadrupole mass moment $J_{2}$ of typical fast rotating massive B-type stars. Section 3 is devoted to the numerical calculation of the perturbations $\Delta\left(\delta\tau\right)$ induced on the pulsar’s Rømer time delay by some post-Keplerian (pK) classical and pN accelerations. We do not calculate the propagation delays accounting, e.g., for the effects on the pulsar’s travelling electromagnetic waves through the deformed spacetime of the B-type star. Section 4 summarizes our findings and contains our conclusions. Once again, we stress the preliminary nature of our sensitivity investigation; we do not perform a full covariance analysis implying, e.g., the simulation of the pulsar’s TOAs and their reduction along with parameter estimation.

2 Quadrupole mass moment, flattening, and angular momentum of fast rotating main sequence B-type massive stars

The binary pulsars on which we shall focus presumably own a main sequence B-type massive star. Such stars are mostly fast rotators Levato & Grosso (2013) and are hence distorted by the centrifugal acceleration. Such a distortion takes the mass distribution away from spherical symmetry and endows these stars with a quadrupolar and higher gravitational moments. In the case of Be stars (namely B stars with emission lines), rotation is believed to be almost critical, namely the rotational velocity is taken close ($>70\%$) to the Keplerian velocity at equator (e.g. Rivinius, Carciofi & Martayan, 2013). When critical rotation is reached, the surface distortion is maximum and the flattening $\nu\doteq(R_{\textrm{e}}-R_{\textrm{p}})/R_{\textrm{e}}$, expressed in terms of the equatorial and polar radii $R_{\mathrm{e}},\,R_{\mathrm{p}}$, respectively, is close to one third. For such stars the computation of the gravitational quadrupolar moment cannot be done perturbatively as it is the case for the Sun, which is a slow rotator (e.g. Rozelot, Godier & Lefebvre, 2001). Modelling these stars requires two-dimensional models. Fortunately, self-consistent 2D-models have recently been achieved with the ESTER code Espinosa Lara & Rieutord (2013); Rieutord, Espinosa Lara & Putigny (2016). Compared to previous 2D-models, ESTER models include self-consistently the baroclinicity of the stellar envelopes and can thus predict the associated differential rotation. They therefore provide unambiguously the total angular momentum of the star given, for instance, its equatorial velocity.

With the ESTER code, we computed the parameters of three stellar models of 10, 15, and 30 $\textrm{M}_{\odot}$ as they can represent the companions of PSR J0045-7319, PSR J2032+4127 and PSR B1256-63 respectively. Since the evolutionary status of the stars is unknown but presumably on or close to the main sequence, we computed their steady state at ZAMS (Zero-Age Main Sequence) and at half-main sequence to get an idea of the effects of evolution. We use standard galactic metallicity $Z=0.02$ with the solar mixture (which may be approximate for PSR J0045-7319 which is in the SMC, known to be less metallic than the Galaxy). Results of ZAMS and evolved models are displayed in Tables 2 to 3 respectively. There we give the total spin angular momentum $S$ and the dimensionless quadrupole mass moment $J_{2}$ along with other bulk parameters of the models. We recall that multipole gravitational moments $J_{\ell}$ are defined by the multipole expansion of the gravitational potential of a mass $M$, namely

[TABLE]

where $\mu\doteq GM$ is the star’s gravitational parameter, $G$ is the Newtonian constant of gravitation, $\xi\doteq\mathbf{\hat{S}}\mathbf{\cdot}\mathbf{\hat{r}}$ is the cosine of the angle $\theta$ between the directions of the body’s spin axis and of an external point at $\mathbf{r}$, $\mathcal{P}_{\ell}\left(\xi\right)$ is the Legendre polynomial of degree $\ell$. We consider the mass distribution of the stellar models to be symmetric with respect to equator thus making the $J_{\ell}$ of odd order all vanish. The remaining $J_{2p}$ can be computed with the integral expression

[TABLE]

where the integration is over the volume $\mathcal{V}$ of the star. The same expression is given in, for instance, Helled et al. (2011). Here we are especially interested in $J_{2}$ and $S$, namely in

[TABLE]

which are directly computed from the ESTER models; $\Xi\left(\mathbf{r}\right)$ is the local angular speed.

We choose an angular rotation rate of 70% of the actual critical (Keplerian) angular velocity of the star. Such a rotation rate is typical of the nearby fast rotating stars that have been measured by interferometry. Their flattening is typically $\sim 0.2$ (e.g. Domiciano de Souza et al., 2014), as our models.

From Tables 2 to 3, we clearly see that as evolution proceeds, namely as the hydrogen content of the core decreases, $J_{2}$ decreases as expected from the resulting contraction of the convective core. From the work of James (1964) we can compute $J_{2}$ for a polytrope of index $n=3$ with a similar flattening as the ESTER models. We find that $J_{2}\simeq 2.1\times 10^{-3}$ which is quite similar to the ZAMS models. For evolved models one should use polytropes with a higher polytropic index, as they are more centrally condensed, and typically $n=3.43$ matches the ESTER evolved models as far as $J_{2}$ is concerned. This result may be useful for simulating the orbital evolution of binary pulsars since polytropic models are much easier to compute.

3 The perturbations of the Rømer-type pulsar’s time delay due to some pK Newtonian and pN accelerations

Here, we will assume a coordinate system centered in the binary’s barycenter whose reference $z$-axis is directed along the line of sight from the binary to the observer, while the reference $\left\{x,\leavevmode\nobreak\ y\right\}$ plane spans the plane of the sky.

3.1 The pulsar as a structureless, pointlike particle

We will, first, consider the following pK accelerations experienced by a test particle moving with velocity {\color[rgb]{0,0,0}\mathbf{v}} in the external field of an oblate body of mass $M$, equatorial and polar radii $R_{\textrm{e}},\leavevmode\nobreak\ R_{\textrm{p}}$, ellipticity $\varepsilon\doteq\sqrt{1-R_{\textrm{p}}^{2}/R_{\textrm{e}}^{2}}$, angular momentum $\mathbf{S}$ and dimensionless quadrupole mass moment $J_{2}$. In Section 3.2, we will discuss the limits of validity of the point-particle approximation.

To the Newtonian level, the external potential of the distorted star at the position $\mathbf{r}$ is, from Equation (1),

[TABLE]

where $\mathcal{P}_{2}\left(\xi\right)=\left(3\xi^{2}-1\right)/2$ is the Legendre polynomial of degree 2. The classical acceleration due to $J_{2}$ is

[TABLE]

The 1pN gravitoelectric Schwarzschild-like acceleration affecting the motion of a test particle in the static field of a nonrotating, spherically symmetric body is (Petit, Luzum & et al., 2010)

[TABLE]

where $c$ is the speed of light in vacuum, and {\mathrm{v}}_{r}\doteq{\color[rgb]{0,0,0}\mathbf{v}}\mathbf{\cdot}\mathbf{\hat{r}} is the radial velocity of the test particle. Eq. (6) is responsible for the formerly anomalous perihelion precession of Mercury whose explanation was the first empirical confirmation of general relativity (Einstein, 1915).

The 1pN gravitomagnetic Lense-Thirring acceleration in the stationary field due to the rotating primary is (Petit, Luzum & et al., 2010)

[TABLE]

The gravitomagnetic field of the Earth was unambiguously measured for the first time by the Gravity Probe B (GP-B) mission (Everitt et al., 2011). Tests of the Lense-Thirring orbital precessions with some terrestrial geodetic satellites are ongoing; see, e.g. Renzetti (2013), and Lucchesi et al. (2015) and references therein for comprehensive reviews.

The 1pN gravitoelectric acceleration felt by a test particle in the field of an oblate body is (Soffel et al., 1987; Soffel, 1989; Brumberg, 1991; Will, 2014; Iorio, 2015)

[TABLE]

where {\mathrm{v}}_{S}\doteq{\color[rgb]{0,0,0}\mathbf{v}}\mathbf{\cdot}\mathbf{\hat{S}} is the component of the particle’s velocity along the direction of the primary’s spin. Note that the parameter $J_{2}$ in Equation (8) is the same entering Equation (5), as per Equations 1 to 2 of Soffel et al. (1987).

The 1pN gravitomagnetic acceleration imparted to a test particle by the spin octupole moment222It will be shown that its effects are small enough to justify order-of-magnitude calculations, without need of detailed stellar models. of a uniformly rotating homogenous oblate spheroid (Panhans & Soffel, 2014; Meichsner & Soffel, 2015; Iorio, 2019) can be cast into the compact form

[TABLE]

The pK accelerations of Equations (5) to (9) perturb the otherwise Keplerian motion of the binary causing a change $\Delta\left(\delta\tau\right)$ of the regular variation $\delta\tau$ of the TOAs due to the relative orbital motion of the pulsar and the massive companion. It can be modeled as the ratio of the projection of the barycentric orbit of the pulsar onto the line of sight to $c$ (Damour & Schaefer, 1991; Konacki, Maciejewski & Wolszczan, 2000). Thus, $\Delta\left(\delta\tau\right)$ can be calculated by looking at the perturbations $\Delta z$ induced by Equations (5) to (9) on the $z$-component of the pulsar’s barycentric orbital motion. We do that by numerically integrating the equations of motion of a fictitious pulsar with mass $M_{\mathrm{p}}=1.4\,\mathrm{M}_{\odot}$ having as a companion a Be-type star with $M=15\leavevmode\nobreak\ \textrm{M}_{\odot},\leavevmode\nobreak\ R_{\textrm{e}}=5.96\leavevmode\nobreak\ \textrm{R}_{\odot},\leavevmode\nobreak\ \nu=0.203,\leavevmode\nobreak\ S=3.41\times 10^{45}\leavevmode\nobreak\ \textrm{J}\leavevmode\nobreak\ \textrm{s},\,J_{2}=1.92\times 10^{-3}$, as per Table 2, for different values of its orbital configuration, determined by the initial values of the semimajor axis $a$, the eccentricity $e$, the orbital inclination $I$ to the plane of the sky, the longitude of the ascending node $\Omega$, the argument of periastron $\omega$, the true anomaly at epoch $f_{0}$, and of the stellar spin axis characterized by its inclination $i$ to the line of sight, and the longitude $\phi$ of the projection of the stellar spin onto the plane of the sky. For a chosen pK acceleration ${\mathbf{A}}^{\textrm{pK}}$, in order to compute its perturbation $\Delta\left(\delta\tau\right)=\Delta z\left(t\right)/c$ over, say, 10 yrs, we perform two runs sharing the same initial conditions with and without ${\mathbf{A}}^{\textrm{pK}}$, calculate the resulting time series of $z\left(t\right)/c$ and take their difference. Figures 1 to 5 depict our results for a given orbital configuration and the aforementioned Be-type main sequence star; in the panels of each Figure, we vary the Keplerian orbital elements and the orientation of $\mathbf{\hat{S}}$ in order to investigate the sensitivity to the parameter space of the adopted binary system. For the sake of a comparison, the Keplerian delay for the adopted reference orbital configuration lies in the range $-250\,\mathrm{s}\lesssim\delta\tau\lesssim 50\,\mathrm{s}$ over one orbital revolution.

Figure 1 displays the Newtonian signatures due to the star’s $J_{2}$. It can be noted that the resulting signal is strongly dependent on the initial value of the true anomaly, spanning the range from $50\,\textrm{s}$ to $-150\,\textrm{s}$. Instead, for the given value $f_{0}=228\,\textrm{deg}$, the sensitivity to the other orbital parameters is rather modest, amounting to about $10\,\textrm{s}$.

The 1pN gravitoelectric Schwarzschild-like signatures due to the stellar mass monopole are reproduced in Figure 2. Also in this case, the initial value of the true anomaly induces a marked variability in the decadal time series which ranges from to $-5\,\textrm{s}$ to $25\,\textrm{s}$. The other orbital parameters have less impact since the resulting variation of the signals is of the order of about $2\,\textrm{s}$.

The 1pN gravitoelectric time series induced by the quadrupole mass moment $J_{2}$ of the star are the subject of Figure 3. In this case, the ranges of variation due to all the orbital parameters are rather similar, amounting to about $20-50\,\textrm{ms}$.

Figure 4 shows the 1pN gravitomagnetic Lense-Thirring signatures due to the spin dipole moment of the star. They are mainly sensitive to the orbital inclination $I$ and to the spin’s inclination $i$ to the line of sight which induces a variability as large as $10-12\,\textrm{ms}$. The ranges of variation induced by the other parameters are, instead, of the order of $1-5\,\textrm{ms}$.

The 1pN gravitomagnetic time series caused by the spin octupole moment of the star are depicted in Figure 5. The widest range of variability, $30-50\,{\mathchoice{\mbox{m}}{\mbox{m}}{\mbox{m}}{\mbox{m}}}\textrm{s}$, is due to the inclination $I$, the node $\Omega$, and the longitude $\phi$ of the spin’s projection onto the plane of the sky. We also checked that, for a pulsar orbiting in $\simeq 20-30\,\textrm{d}$ a star as massive as those in the last lines of Tables 2 to 3 with periastron distances of $r_{\textrm{min}}\simeq 1.2-1.1\,R_{\textrm{eq}}$, the magnitude of the timing signatures would reach the $\simeq 1-10\,\textrm{ms}$ level.

3.2 The effects of the mass and spin multipoles of the pulsar

Until now, we have considered the neutron star as a structureless, point particle moving around its more massive companion. In fact, a pulsar is an extended body with its own mass and spin multipole moments which, at least in principle, may have an impact on its orbital motion in a full two-body framework.

The modification of Equation (6) for two finite bodies of masses $M_{\textrm{A}},\,M_{\textrm{B}}$ is (Soffel, 1989)

[TABLE]

where $\mu_{\textrm{tot}}=GM_{\textrm{tot}}=G\left(M_{\textrm{A}}+M_{\textrm{B}}\right)$, and $\zeta\doteq M_{\textrm{A}}M_{\textrm{B}}/M_{\textrm{tot}}^{2}$. By taking the standard value $M_{\textrm{p}}=1.4\,\textrm{M}_{\odot}$ for the mass of the pulsar, it is $\zeta=0.078$ for the Be-star assumed in Section 3.1. It turns out that the introduction of $\zeta$ in our numerical code changes the size of the time series of Figure 2 by $\simeq 0.1-0.2\,\textrm{s}$, while their temporal patterns remain essentially unchanged. Given the current level of accuracy in the timing residuals, such a discrepancy might be significative, and Equation (10) should be used instead of Equation (6).

In regard to the angular momentum, the spin $\mathbf{S}$ of the Be-star in Equation (7) should be replaced by the sum (Barker & O’Connell, 1975)

[TABLE]

where ${\mathbf{S}}_{\textrm{p}}$ is the angular momentum of the pulsar. By assuming, e.g., $S_{\textrm{p}}=3\times 10^{40}\,\textrm{J}\,\textrm{s}$ as for PSR J0737-3039A (Burgay et al., 2003; Lyne et al., 2004; Kramer, 2012; Kehl et al., 2017; Kramer, 2018), it turns out that the magnitude of the second term in Equation (11) amounts to about $\simeq 7\times 10^{-5}$ of that of the first one. Moreover, it is $\left(3/4\right)\left(M_{\textrm{p}}/M\right)=0.07$. Thus, as far as the 1pN gravitomagnetic Lense-Thirring effect is concerned, our scenario can be well approximated by a restricted two-body system with a spinning primary, and Equation (7) is substantially adequate. Indeed, it turns out that rescaling the star’s angular momentum in our numerical simulations as dictated by Equation (11) slightly modifies the size of the time series of Figure 4 by just $\simeq 1-1.5\,\textrm{ms}$.

The effect of the pulsar’s quadrupole mass moment $J_{2}^{\textrm{p}}$ can be accounted for by the replacement $M\rightarrow M_{\textrm{tot}}$ in Equation (5) and by writing in its right-hand-side another term for $J_{2}^{\textrm{p}}$ analogous to that for the stellar oblateness $J_{2}$. As a result, by introducing the dimensional quadrupole mass moment $Q_{2}\doteq-J_{2}MR_{\textrm{e}}^{2}$, the two terms in the right-hand-side of the resulting modified version of Equation (5) are weigthed by (Barker & O’Connell, 1975)

[TABLE]

For a neutron star, it is (Laarakkers & Poisson, 1999; Berti & Stergioulas, 2004; Pappas & Apostolatos, 2012; Bauböck et al., 2013)

[TABLE]

with $0.07\lesssim q\lesssim 3.507$ for a variety of Equations of State (EOSs). Thus, we have

[TABLE]

so that $\mathcal{Q}_{2}^{\textrm{p}}\simeq 10^{-10}\,\mathcal{Q}_{2}$. On the other hand, $M_{\textrm{p}}/M=0.09$ in Equation (12). Thus, for the Newtonian signature of the quadrupole mass moment, the restricted two-body scenario with an oblate primary is adequate in the present case, and the use of Equation (5) is justified provided that the stellar quadrupole moment $Q_{2}$ is replaced by Equation (12). Indeed, it turns out that the introduction of $\mathcal{Q}_{2}$ in the numerical integration changes the size of the time series in Figure 1 by $\simeq 1\,\textrm{s}$ which may not be neglected, given the current level in ${\mathchoice{\mbox{s}}{\mbox{s}}{\mbox{s}}{\mbox{s}}}_{\tau}$.

Despite Equation (8) and Equation (9) were derived so far only for the motion of a test particle around a spinning, oblate mass, there are no doubts that they are adequate to the scenario considered here.

Actually, general relativity predicts that, in general, there is also a self-force due to the spin angular momentum of an extended rotating body in motion in an external gravitational field which acts on it modifying its orbit through a spin-orbit coupling (Papapetrou, 1951; Dixon, 1974; Barker & O’Connell, 1979; Mashhoon & Singh, 2006; Mathisson, 2010; Iorio, 2012). In order to quickly evaluate the possible impact of such an effect in our case, let us proceed as follows. To the 1pN level, the precession $\dot{\Omega}_{S_{\textrm{p}}}$ of, say, the node $\Omega$, averaged over one orbital revolution of the spinning pulsar in its motion around the massive Be-type star, assumed nonrotating, is (Iorio, 2012)

[TABLE]

In this expression, ${\mathbf{\sigma}}_{\textrm{p}}={\mathbf{S}}_{\textrm{p}}/M_{\textrm{p}}$ is the pulsar’s spin angular momentum per unit mass, while $\mathbf{\hat{m}}=\left\{-\cos I\sin\Omega,\,\cos I\cos\Omega,\,\sin I\right\}$ is a unit vector in the orbital plane perpendicular to the line of the nodes. Let us compare Equation (17) to the analogous precession induced by some of the effects previously considered in which the pulsar was treated as a point particle. The Lense-Thirring node rate, induced by Equation (7) which is responsible for the $\simeq\textrm{ms}$ time series of Figure 4, is (Iorio, 2012)

[TABLE]

It turns out that, in the case of a pulsar like, e.g., PSR J0737-3039A orbiting the Be-type star of the Figures 1 to 5, it is $\left|\dot{\Omega}_{S_{\textrm{p}}}/\dot{\Omega}_{\textrm{LT}}\right|\simeq 10^{-5}$. Thus, we conclude that the spin-orbit self-force experienced by the pulsar is completely negligible in our scenario. In fact, there is also a further self-force acting on an extended rotating body moving in an external gravitational field due to the spin-spin coupling between the angular momenta of the source and of the orbiter itself (Barker & O’Connell, 1979; Iorio, 2012). By following the same strategy for the spin-orbit coupling, the results in Iorio (2012) allow to conclude that, in our case, such an effect is even smaller than the previous one, being of the order of $\simeq 10^{-7}$ of, say, the $\simeq\textrm{ms}$-level Lense-Thirring signal.

4 Summary and conclusions

We preliminarily explored the possibility of putting to the test several pK features of motion of Newtonian and pN origin in binaries hosting a pulsar and a massive, fast rotating, highly distorted main sequence star characterized by mass $M$, angular momentum $\mathbf{S}$, equatorial and polar radii $R_{\textrm{e}},\leavevmode\nobreak\ R_{\textrm{p}}$, flattening $\nu$, ellipticity $\varepsilon$ and dimensionless quadrupole moment $J_{2}$. Indeed, in addition to the usual 1pN Schwarzschild and Lense-Thirring effects due to the mass monopole $\left(\propto GMc^{-2}\right)$ and spin dipole $\left(\propto GSc^{-2}\right)$ moments, respectively, of the distorted stellar field, there are also other 1pN orbital effects, induced by the mass quadrupole $\left(\propto GMR^{2}_{\textrm{e}}J_{2}c^{-2}\right)$ and spin octupole $\left(\propto GSR^{2}_{\textrm{e}}\varepsilon^{2}c^{-2}\right)$ moments, whose magnitudes may, perhaps, lie above the sensitivity threshold of the pulsar timing residuals in yet-to-be-discovered close binaries. However, the Newtonian perturbations due to $J_{2}$ are larger than the pN ones.

In order to perform a preliminary sensitivity analysis, we numerically investigated the orbital shifts $\Delta\left(\delta\tau\right)$ induced over 10 yr by all of such pK perturbations on the otherwise Keplerian Rømer-type delay $\delta\tau$ in the pulsar’s TOAs for a Be-type main sequence star characterized by $M=15\leavevmode\nobreak\ \textrm{M}_{\odot},\leavevmode\nobreak\ R_{\textrm{e}}=5.96\leavevmode\nobreak\ \textrm{R}_{\odot},\leavevmode\nobreak\ \nu=0.203,\leavevmode\nobreak\ S=3.41\times 10^{45}\leavevmode\nobreak\ \textrm{J}\leavevmode\nobreak\ \textrm{s},\,J_{2}=1.92\times 10^{-3}$ orbited by a pulsar with an orbital geometry compatible $\left(P_{\rm b}\simeq 40-70\leavevmode\nobreak\ \textrm{d}\right)$ with some of the tightest binaries of this kind out of those discovered so far. We also investigated the sensitivity of the pK timing shifts to the whole system’s parameter space by varying both the orientation of the stellar spin axis and the orbital elements of the pulsar’s orbit. It turns out that the magnitude of the Newtonian signature due to $J_{2}$ can be as large as $\lesssim 4-150\leavevmode\nobreak\ \textrm{s}$, while the pN gravitoelectric quadrupolar one is $\lesssim 10-40\leavevmode\nobreak\ \textrm{ms}$. The pN gravitoelectric (Schwarzschild-like) signal due to the stellar mass monopole can be as large as $\lesssim 1.5-20\leavevmode\nobreak\ \textrm{s}$ level. The pN gravitomagnetic shifts due to the spin dipole (Lense-Thirring) and octupole moments, for the evaluation of whose size the knowledge of the stellar spin angular momentum $S$ and ellipticity $\varepsilon$ is crucial, are of the order of $\lesssim 0.5-6\leavevmode\nobreak\ \textrm{ms}$, and $\lesssim 5-20\leavevmode\nobreak\ {\mathchoice{\mbox{m}}{\mbox{m}}{\mbox{m}}{\mbox{m}}}\textrm{s}$, respectively. The rms of the timing residuals of all the non-recycled, non-millisecond pulsars like those having a fast rotating main sequence star companion discovered so far are $\lesssim\textrm{ms}$ over $2-13\leavevmode\nobreak\ \textrm{yr}$. It seems difficult that they can be substantially improved in the future. This implies that, in principle, all the pK effects considered fall within the potential measurability domain, except the pN spin octupole which is $\simeq 1-2$ orders of magnitude weaker, at least for the orbital configurations and the star considered in this study. The pN Lense-Thirring signatures are just at the $\simeq\textrm{ms}$ level. The different temporal patterns characteristic of the signals investigated may be helpful in separating them. It turns out that the tiniest pN effect may reach the $\simeq 1-10\,\textrm{ms}$ level only for a very tight, eccentric binary $(P_{\textrm{b}}\simeq 20-30\,\textrm{d},\,r_{\textrm{min}}\simeq 1.1-1.2\,R_{\textrm{e}})$ hosting a Be-star with $30\,\textrm{M}_{\odot},\,S=125-870\times 10^{44}\,\textrm{J\,s}$.

Finally, we stress the preliminary nature of our sensitivity analysis. To this aim, we remark that we did not compute the other kinds of time delay connected, e.g., with the propagation of the electromagnetic waves in the deformed spacetime. Moreover, we did not perform a full covariance analysis implying a simulation of the pulsar’s TOAs, their reduction, and parameter estimation.

Acknowledgements

We are grateful to A. Possenti for useful information. MR and ADS also acknowledge the strong support of the French Agence Nationale de la Recherche (ANR), under grant ESRR (ANR-16-CE31-0007-01).