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

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].

An analysis of the analytical expressions of the pN gravitoelectric and gravitomagnetic orbital precessions due to the asphericity of the primary (Iorio 2015(Iorio , 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) τ of the pulses emitted by the neutron star which, in case it has a gravitationally bound companion, exhibit a regular variation δτ 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 b = 1237 d = 3.38 yr. The pulsar's non-degenerate companion is the fast spinning Be star LS 2883, whose equatorial velocity V e is about 280 km s −1 corresponding to ∼ 70% of its break-up velocity (Porter 1996). Its mass M and equatorial radius R e amounts to about 30 Solar masses (M ) and 9.7 Solar radii (R ) (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 b = 51.17 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 b = 1940 d = 5.3 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 b = 231 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 b = 8578 d = 23.5 yr. The companion of the neutron star is the massive Be star MT91 213 with M 15 M .
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 ∼ 1 M star moving in P b = 95 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 M ) companion star in a circular path with P b = 4.75 hr. The rms timing residuals of the aforementioned binary pulsars are all of the order of ms, apart from PSR J1903+0327 which is at the 1 µs level; more specifically, they are 0.46 ms over 13 yr for PSR B1259-63 (Wang, Johnston & Manchester 2004), 7.4 ms over 2 yr for PSR J0045-7319 (Kaspi et al. 1994), 5.3 ms over 4.35 yr for PSR J1638-4725 (Lorimer et al. 2006), 0.8 ms over 2.29 yr for PSR J1740-3052 (Stairs et al. 2001), 0.5 − 1 ms over about 6 yr for PSR J2032+4127 (Lyne et al. 2015), 1 µs over about 3 yr for PSR J1903+0327 (Arzoumanian et al. 2018), 0.1 ms over about 4 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.
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 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 µ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 ms, and its rms timing residuals are as little as 1 µs. The spin period of PSR J1023+0038 is P = 1.7 ms, and its rms timing residuals are 0.1 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 ms, isolated or with a white dwarf as companion, is of the order of 0.1 − 0.2 µs. It is expected that future instrumental improvements may push the rms timing accuracy of some of them to the 1 10 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 µs.
Here, we will preliminarily investigate the size and the temporal patterns of the perturbations ∆ (δτ) induced on the Rømer-like orbital part of the pulsar's time delay δτ 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 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 ∆ (δτ) 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 ν (R e − R p )/R e , expressed in terms of the equatorial and polar radii R e , R 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 M 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 are defined by the multipole expansion of the gravitational potential of a mass M, namely where µ GM is the star's gravitational parameter, G is the Newtonian constant of gravitation, ξ Ŝ ·r is the cosine of the angle θ between the directions of the body's spin axis and of an external point at r, P (ξ) is the Legendre polynomial of degree . We consider the mass distribution of the stellar models to be symmetric with respect to equator thus making the J of odd order all vanish. The remaining J 2p can be computed with the integral expression where the integration is over the volume 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 which are directly computed from the ESTER models; Ξ (r) 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 ∼ 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 2.1 × 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 {x, y} 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 v in the external field of an oblate body of mass M, equatorial and polar radii R e , R p , ellipticity ε 1 − R 2 p /R 2 e , angular momentum 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 r is, from Equation (1), where P 2 (ξ) = 3ξ 2 − 1 /2 is the Legendre polynomial of degree 2. The classical acceleration due to J 2 is 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) where c is the speed of light in vacuum, and v r v·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) 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;Iorio 2015) where v S v·Ŝ 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 moment 2 of a uniformly rotating homogenous oblate spheroid (Panhans & Soffel 2014;Meichsner & Soffel 2015;Iorio 2019) can be cast into the compact form The pK accelerations of Equations (5) to (9) perturb the otherwise Keplerian motion of the binary causing a change ∆ (δτ) of the regular variation δτ 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, ∆ (δτ) can be calculated by looking at the perturbations ∆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 p = 1.4 M having as a companion a Be-type star with M = 15 M , R e = 5.96 R , ν = 0.203, S = 3.41 × 10 45 J s, J 2 = 1.92 × 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 Ω, the argument of periastron ω, 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 φ of the projection of the stellar spin onto the plane of the sky. For a chosen pK acceleration A pK , in order to compute its perturbation ∆ (δτ) = ∆z (t) /c over, say, 10 yrs, we perform two runs sharing the same initial conditions with and without A pK , calculate the resulting time series of z (t) /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Ŝ 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 s δτ 50 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 s to −150 s. Instead, for the given value f 0 = 228 deg, the sensitivity to the other orbital parameters is rather modest, amounting to about 10 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 s to 25 s. The other orbital parameters have less impact since the resulting variation of the signals is of the order of about 2 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 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 ms. The ranges of variation induced by the other parameters are, instead, of the order of 1 − 5 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 µs, is due to the inclination I, the node Ω, and the longitude φ of the spin's projection onto the plane of the sky. We also checked that, for a pulsar orbiting in 20 − 30 d a star as massive as those in the last lines of Tables 2 to 3 with periastron distances of r min 1.2 − 1.1 R eq , the magnitude of the timing signatures would reach the 1 − 10 ms level.

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 A , M B is (Soffel 1989 where µ tot = GM tot = G (M A + M B ), and ζ M A M B /M 2 tot . By taking the standard value M p = 1.4 M for the mass of the pulsar, it is ζ = 0.078 for the Be-star assumed in Section 3.1. It turns out that the introduction of ζ in our numerical code changes the size of the time series of Figure 2 by 0.1 − 0.2 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 S of the Be-star in Equation (7) should be replaced by the sum (Barker & O'Connell 1975) where S p is the angular momentum of the pulsar. By assuming, e.g., S p = 3 × 10 40 J 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 7 × 10 −5 of that of the first one. Moreover, it is (3/4) M p /M = 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)  The effect of the pulsar's quadrupole mass moment J p 2 can be accounted for by the replacement M → M tot in Equation (5) and by writing in its right-hand-side another term for J p 2 analogous to that for the stellar oblateness J 2 . As a result, by introducing the dimensional quadrupole mass moment Q 2 −J 2 MR 2 e , the two terms in the right-hand-side of the resulting modified version of Equation (5) are weigthed by (Barker & O'Connell 1975) For a neutron star, it is (Laarakkers & Poisson 1999;Berti & Stergioulas 2004;Pappas & Apostolatos 2012;Bauböck et al. 2013) with 0.07 q 3.507 for a variety of Equations of State (EOSs). Thus, we have Q 2 = −9.8 × 10 47 kg m 2 , Q p 2 = −q 1.2 × 10 37 kg m 2 , so that Q p 2 10 −10 Q 2 . On the other hand, M 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 Q 2 in the numerical integration changes the size of the time series in Figure 1 by 1 s which may not be neglected, given the current level in σ τ .
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Ω S p of, say, the node Ω, averaged over one orbital revolution of the spinning pulsar in its motion around the massive Be-type star, assumed nonrotating, is (Iorio 2012) In this expression, σ p = S p /M p is the pulsar's spin angular momentum per unit mass, whilê m = {− cos I sin Ω, cos I cos Ω, sin I} 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) 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 Ω S p /Ω LT 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 10 −7 of, say, the ms-level Lense-Thirring signal.

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 S, equatorial and polar radii R e , R p , flattening ν, ellipticity ε and dimensionless quadrupole moment J 2 . Indeed, in addition to the usual 1pN Schwarzschild and Lense-Thirring effects due to the mass monopole ∝ GMc −2 and spin dipole ∝ GS c −2 moments, respectively, of the distorted stellar field, there are also other 1pN orbital effects, induced by the mass quadrupole ∝ GMR 2 e J 2 c −2 and spin octupole ∝ GS R 2 e ε 2 c −2 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 ∆ (δτ) induced over 10 yr by all of such pK perturbations on the otherwise Keplerian Rømer-type delay δτ in the pulsar's TOAs for a Be-type main sequence star characterized by M = 15 M , R e = 5.96 R , ν = 0.203, S = 3.41 × 10 45 J s, J 2 = 1.92 × 10 −3 orbited by a pulsar with an orbital geometry compatible (P b 40 − 70 d) 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 4 − 150 s, while the pN gravitoelectric quadrupolar one is 10 − 40 ms. The pN gravitoelectric (Schwarzschild-like) signal due to the stellar mass monopole can be as large as 1.5 − 20 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 ε is crucial, are of the order of 0.5 − 6 ms, and 5 − 20 µ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 ms over 2 − 13 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 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 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 1 − 10 ms level only for a very tight, eccentric binary (P b 20 − 30 d, r min 1.1 − 1.2 R e ) hosting a Be-star with 30 M , S = 125 − 870 × 10 44 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.   Negueruela et al. (2011). 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 and V e are the equatorial radius and velocity, respectively, S is the total spin angular momentum, ν is the flattening, and J 2 is the dimensionless mass quadrupole moment. The metallicity is Z=0.02 and the hydrogen mass fraction is X = 0.7. In our simulations, we use the parameters of the star with 15 M listed here.   Be star-pulsar J 2 timing shift   Be star-pulsar GMc -2 timing shift   Be star-pulsar GMJ 2 c -2 timing shift Be star-pulsar GMJ 2 c -2 timing shift   Be star-pulsar GSc -2 timing shift