11institutetext: Lijing Shao 22institutetext: Kavli Institute for Astronomy
and Astrophysics, Peking University, Beijing 100871, China
National Astronomical Observatories, Chinese Academy of Sciences, Beijing
100012, China
22email<EMAIL_ADDRESS>
# Radio Pulsars as a Laboratory for Strong-field Gravity Tests

Lijing Shao

###### Abstract

General relativity offers a classical description to gravitation and
spacetime, and is a cornerstone for modern physics. It has passed a number of
empirical tests with flying colours, mostly in the weak-gravity regimes, but
nowadays also in the strong-gravity regimes. Radio pulsars provide one of the
earliest extrasolar laboratories for gravity tests. They, in possession of
strongly self-gravitating bodies, i.e. neutron stars, are playing a unique
role in the studies of strong-field gravity. Radio timing of binary pulsars
enables very precise measurements of system parameters, and the pulsar timing
technology is extremely sensitive to various types of changes in the orbital
dynamics. If an alternative gravity theory causes modifications to binary
orbital evolution with respect to general relativity, the theory prediction
can be confronted with timing results. In this chapter, we review the basic
concepts in using radio pulsars for strong-field gravity tests, with the aid
of some recent examples in this regard, including tests of gravitational
dipolar radiation, massive gravity theories, and the strong equivalence
principle. With more sensitive radio telescopes coming online, pulsars are to
provide even more dedicated tests of strong gravity in the near future.

##  1 Introduction

Pulsars are rotating magnetized neutron stars. On the one hand, due to their
large moment of inertia ($I\sim 10^{38}\,{\rm kg\,m}^{2}$) and usually small
external torque, their rotation is extremely stable. If a pulsar sweeps a
radiating beam in the direction of the Earth, a radio pulse could be recorded
using large-area telescopes for each rotation. As fundamentally known in
physics, such a periodic signal can be viewed as a clock. Therefore, pulsars
are famously recognized as astrophysical clocks in astronomy. Even better,
thanks to a sophisticated technique called pulsar timing Taylor:1992kea ,
pulsar astronomers can accurately record a number of periodic pulse signals.
These pulses’ times of arrival are compared with atomic clocks at the
telescope sites. Some of these observations can be carried out and last for
decades. From a large number of times of arrival of these pulse signals, the
physical properties of pulsar systems are inferred to a great precision
Lorimer:2005misc . For example, a recent study with sixteen years of timing
data of the Double Pulsar,111Currently, PSR J0737$-$3039A/B is the only
discovered double neutron star system whose two neutron stars were both
detected as pulsars Burgay:2003jj ; Lyne:2004cj ; Shao:2021iax , known as
Pulsar A and Pulsar B. PSR J0737$-$3039A/B, gives the rotational frequency of
pulsar A in the binary system Kramer:2021jcw ,

$\displaystyle\nu=44.05406864196281(17)\,{\rm Hz}\,.$ (1)

It has sixteen significant digits, and the numbers in the parenthesis give the
uncertainty of the last-two digits. Such a precision rivals the precision of
atomic clocks on the Earth Hobbs:2012apa , and also it possibly calls for an
extension of the usual use of floating numbers in computer numerics for future
precision pulsar timing experiments. Pulsars are truly precision clocks.

On the other hand, neutron stars are the densest objects known that are made
of standard-model materials. For such a compact object, gravity plays a vital
role in shaping its internal structure and affecting its external dynamics. As
explicitly demonstrated by Damour and Esposito-Farèse Damour:1993hw , if
gravity is described by an alternative theory to the general relativity—in
their case, a class of scalar-tensor gravity theories—nonperturbative phase-
transition-like behaviours might happen for neutron stars, resulting in large
deviations from general relativity in the strong field of neutron stars
Damour:1996ke ; Esposito-Farese:2004azw ; Sennett:2017lcx . These large
deviations will manifest in the timing data of pulsars in some way (cf. Sec.
2), and they could provide smoking-gun signals for gravity theories regarding
the strong-field properties. Combining the strong-field nature of neutron
stars and the precision measurements of times of arrival, radio pulsars are
truly ideal to test alternative theories of gravity Wex:2014nva ; Shao:2016ezh
; Wex:2020ald , augmenting what have been done in the weak field of the Solar
System Will:2014kxa , and complementing what are recently being performed with
gravitational waves LIGOScientific:2019fpa ; LIGOScientific:2020tif ;
LIGOScientific:2021sio  and black hole shadows EventHorizonTelescope:2019dse ;
EventHorizonTelescope:2022xnr ; EventHorizonTelescope:2022xqj .

Currently, more than three thousands of radio pulsars are
discovered222https://www.atnf.csiro.au/people/pulsar/psrcat/ Manchester:2004bp
. The most useful subset of pulsars in testing alternative gravity theories
are millisecond pulsars in clean binaries.333In one case, a pulsar in a triple
system, PSR J0337+1715, provides the best limit on the strong equivalence
principle Archibald:2018oxs ; Voisin:2020lqi ; Shao:2016ubu . Their times of
arrival at telescopes are imprinted with information from the following
sources:

  1. (i)

the Solar system dynamics which affect the motion of radio telescopes;

  2. (ii)

the binary dynamics which are resulted from the mutual gravitational
interaction between the two binary components; and

  3. (iii)

the interstellar medium which affects the propagation of radio waves in a
frequency-dependent way, in terms of dispersion, scattering, and so on.

A formalism, which includes the above effects and connects the proper time of
the pulse signals in the pulsar frame to the observed coordinate time at the
telescopes, is called a pulsar timing model. One of the widely used timing
models for binary pulsars is the Damour-Deruelle timing model Damour:1986ads .
It is a phenomenological model that applies to a large set of alternative
gravity theories which are possibly being the underlying theory for the
binary’s orbital motion.

In the Damour-Deruelle timing model, a handful of parameterized post-Keplerian
(PPK) parameters are introduced for generic Lorentz-invariant extensions of
gravity theories Damour:1991rd . The values of PPK parameters differ in
different gravity theories. Therefore, measurements of these PPK parameters
can be converted into constraints on parameters in the alternative gravity
theories. The most frequently used PPK parameters include $\dot{\omega}$,
$\dot{P}_{\rm b}$, $\gamma$, $r$, and $s$. The PPK parameter $\dot{\omega}$
describes the periastron advance of the binary orbit, the PPK parameter
$\dot{P}_{\rm b}$ describes the orbital period decay caused by the radiation
of gravitational waves, the PPK parameter $\gamma$ describes combined effects
from the Doppler time delay and gravitational time delay, and the PPK
parameters $(r,s)$ describe the Shapiro time delay imprinted by the spacetime
curvature of the companion star. The values of these five PPK parameters in
the general relativity are given by Damour:1986ads ; Lorimer:2005misc ,

$\displaystyle\dot{\omega}$
$\displaystyle=3\left(\frac{P_{\mathrm{b}}}{2\pi}\right)^{-5/3}\left(T_{\odot}M\right)^{2/3}\left(1-e^{2}\right)^{-1}\,,$
(2) $\displaystyle\dot{P}_{\mathrm{b}}$
$\displaystyle=-\frac{192\pi}{5}\left(\frac{P_{\mathrm{b}}}{2\pi}\right)^{-5/3}\left(1+\frac{73}{24}e^{2}+\frac{37}{96}e^{4}\right)\left(1-e^{2}\right)^{-7/2}T_{\odot}^{5/3}m_{\rm
A}m_{\rm B}M^{-1/3}\,,$ (3) $\displaystyle\gamma$
$\displaystyle=e\left(\frac{P_{\mathrm{b}}}{2\pi}\right)^{1/3}T_{\odot}^{2/3}M^{-4/3}m_{\rm
B}(m_{\rm A}+2m_{\rm B})\,,$ (4) $\displaystyle r$
$\displaystyle=T_{\odot}m_{\rm B}\,,$ (5) $\displaystyle s$
$\displaystyle=x\left(\frac{P_{\mathrm{b}}}{2\pi}\right)^{-2/3}T_{\odot}^{-1/3}M^{2/3}m_{\rm
B}^{-1}\,,$ (6)

where $P_{\rm b}$ and $e$ are respectively the orbital period and orbital
eccentricity, $m_{\rm A}$ and $m_{\rm B}$ are the masses of the pulsar and its
companion in unit of the Solar mass ($M_{\odot}$), the total mass $M\equiv
m_{\rm A}+m_{\rm B}$, and $T_{\odot}\equiv GM_{\odot}/c^{3}=4.925490947\mu$s.
Equations (2)–(6) take different forms in alternative gravity theories, often
with dependence on the extra charges of the binary components in the theory,
e.g., these PPK parameters depend on scalar charges of the pulsar and its
companion in the scalar-tensor theory Damour:1996ke . In pulsar-timing
observation, each PPK parameter is independently measured. Eventually, for a
gravity theory to pass the tests from pulsar timing, it should give consistent
predictions to all the measured values of PPK parameters with a unique set of
physical parameters of the binary system. These consistency checks are often
illustrated in the mass-mass diagram. For an example, in Figure 1 the
measurements of three PPK parameters, $\dot{\omega}$, $\gamma$, and
$\dot{P}_{\rm b}$, from the Hulse-Taylor pulsar PSR B1913+16, give consistent
component masses when the general relativistic equations (2)–(4) are used
Weisberg:2016jye . Therefore, general relativity passes the tests posed by the
Hulse-Taylor pulsar Weisberg:2016jye .

Figure 1: Consistency of general relativity in describing three measured PPK
parameters ($\dot{\omega}$, $\gamma$, and $\dot{P}_{\rm b}$) from PSR B1913+16
in the mass-mass diagram Weisberg:2016jye .

In this following, we will give a few more concrete and recent examples where
binary pulsars play a key role in limiting alternative gravity theories,
including the gravitational dipolar radiation in the scalar-tensor gravity
(Sec. 2), two classes of massive gravity theories (Sec. 3), and the strong
equivalence principle (Sec. 4). These examples are by no means complete, and
certainly reflect the somehow biased topics that the author is interested in.
A short perspective discussion is given in Sec. 5. For more extensive reviews
on using radio pulsars for gravity tests, readers are referred to Refs.
Stairs:2003eg ; Wex:2014nva ; Manchester:2015mda ; Kramer:2016kwa ;
Shao:2016ezh ; Wex:2020ald .

##  2 Strong-field effects and gravitational dipolar radiation

Scalar-tensor gravity theories represent a well posed, healthy extension of
Einstein’s general relativity by including a nonminimally coupled scalar field
in the Lagrangian of gravity Brans:1961sx ; Berti:2015itd ; Will:2014kxa .
Shortly after the first discovery of the Hulse-Taylor binary pulsar, Eardley
Eardley:1975ads  pointed out that a gravitational dipolar radiation could be
used as a discriminant for such a class of gravity theories. An extra dipolar
radiation term can be tested with the PPK parameter $\dot{P}_{b}$.
Investigation along this line was boosted by the theoretical discovery that in
a slightly extended version of the original scalar-tensor gravity,
nonperturbative effects develop for certain neutron stars Damour:1993hw ;
Damour:1996ke . The so-called spontaneous scalarization introduces a much
enhanced gravitational dipolar radiation for a scalarized neutron star in a
binary. The dipolar radiation in principle can even dominate over the
quadrupolar radiation predicted by the general relativity in binary pulsar
observations [cf. Eq. (3)], but still keeping all weak-field gravity tests
satisfied. This enters the regime of strong-field gravity tests, where weak-
field tests have a rather limited power.

A general class of scalar-tensor gravity theories have the following action in
the Einstein frame,

$\displaystyle S=\frac{c^{4}}{16\pi
G_{*}}\int\frac{\mathrm{d}^{4}x}{c}\sqrt{-g_{*}}\left[R_{*}-2g_{*}^{\mu\nu}\partial_{\mu}\varphi\partial_{\nu}\varphi-V(\varphi)\right]+S_{\rm
m}\left[\psi_{\rm m};A^{2}(\varphi)g_{\mu\nu}^{*}\right]\,,$ (7)

where $g^{\mu\nu}_{*}$ and $R_{*}$ are the metric tensor and Ricci scalar
respectively, $\psi_{\rm m}$ collectively denotes standard-model matter
fields, $\varphi$ is an extra scalar field, and quantities with stars are in
the Einstein frame. The novel aspect lies in the fact that it is a conformal
metric $A^{2}(\varphi)g_{\mu\nu}^{*}$ instead of $g_{\mu\nu}^{*}$ itself that
couples to matter fields. Such a nonminimal coupling is important for the
discussions below.

The class of scalar-tensor gravity theories carefully examined by Damour and
Esposito-Farèse Damour:1993hw ; Damour:1996ke  has

$\displaystyle V(\varphi)$ $\displaystyle=0\,,$ (8) $\displaystyle A(\varphi)$
$\displaystyle=\exp\left(\beta_{0}\varphi^{2}/2\right)\,,$ (9)
$\displaystyle\alpha_{0}$ $\displaystyle=\beta_{0}\varphi_{0}\,,$ (10)

where $\varphi_{0}$ is the asymptotic value of $\varphi$ at infinity, and
$\alpha_{0}$ and $\beta_{0}$ are two theory parameters. This is the class of
scalar-tensor theories, sometimes denoted as $T_{1}(\alpha_{0},\beta_{0})$ and
called the Damour-Esposito-Farèse theory, that are most widely confronted with
pulsar observations Freire:2012mg ; Wex:2014nva ; Shao:2017gwu ; Zhao:2022vig
.

By integrating the modified Tolman-Oppenheimer-Volkoff equations derived from
theory (7), one gets a boost in a neutron star’s scalar charge when its mass
reaches a critical point. This phenomenon is understood from the viewpoint of
Landau’s phase transition theory when a tachyonic instability kicks in and a
new branch of neutron star solutions with scalar charges are energetically
favored Esposito-Farese:2004azw ; Sennett:2017lcx ; Khalil:2022sii . We define
the effective scalar charge of a neutron star Damour:1993hw ,

$\alpha_{\rm A}\equiv\frac{\partial\ln m_{\rm A}}{\partial\varphi_{0}}\,,$
(11)

which is a representative quantity characterizing the strength of deviation
from general relativity. In Figure 2, example curves for the effective scalar
charge as a function of neutron star mass are given in blue lines from top to
bottom for $\beta_{0}=-4.8,-4.6,-4.4,-4.2$, assuming the AP4 equation of state
and $|\alpha_{0}|=10^{-5}$. As we can easily seen, indeed that for certain
mass range of neutron stars, $|\alpha_{\rm A}|$ can be very large while
keeping its value very small in weak-gravity fields.

Figure 2: Blue curves show the effective scalar charge in the Damour-Esposito-
Farèse scalar-tensor gravity theory with $|\alpha_{0}|=10^{-5}$ and, from top
to bottom, $\beta_{0}=-4.8,-4.6,-4.4,-4.2$. The AP4 equation of state is
assumed in the calculation. Triangles show the observational bounds from
binary pulsars Shao:2017gwu ; Zhao:2022vig  and gravitational waves
TheLIGOScientific:2017qsa ; LIGOScientific:2018dkp  at the 90% confidence
level. The mass uncertainty for these neutron stars is indicated at the 68%
confidence level.

The emission of gravitational dipolar radiation in a binary pulsar is
proportional to the difference in the effective scalar couplings of the two
binary components A and B, and to the leading order, it contributes to an
additional decay rate of orbital period via Damour:1996ke ,

$\dot{P}_{\rm b}^{\text{dipole }}=-\frac{2\pi
G_{*}}{c^{3}}\left(1+\frac{e^{2}}{2}\right)\left(1-e^{2}\right)^{-5/2}\left(\frac{2\pi}{P_{\rm
b}}\right)\frac{m_{\rm A}m_{\rm B}}{M}\left(\alpha_{\rm A}-\alpha_{\rm
B}\right)^{2}\,.$ (12)

While neutron stars have significant scalar charges, white dwarfs, being weak-
field objects, are hardly different from their counterparts in general
relativity with a vanishingly small scalar charge $\alpha_{\rm
B}\simeq\alpha_{0}\to 0$, where $\alpha_{0}$ is well constrained by Solar
System weak-field tests Will:2014kxa . Therefore, neutron-star white-dwarf
binaries turn out to be the most sensitive probe in this regard Freire:2012mg
; Shao:2017gwu . Recently, a new study Zhao:2022vig  shows explicitly that
neutron-star neutron-star binaries with a significant difference in the masses
of binary components are also excellent laboratories. Therefore, to test the
gravitational dipolar radiation in scalar-tensor gravity, asymmetric binary
pulsars are needed444Unfortunately, we have not detected yet suitable neutron-
star black-hole binaries for this test, which are also potentially very good
testbeds Liu:2014uka . Wex:2014nva .

Some illustration for a specific equation of state, AP4, is given in Figure 2,
along with constraints on the gravitational dipolar radiation from seven
binary pulsars Zhao:2022vig : five neutron-star white-dwarf binaries (PSRs
J0348+0432, J1012+5307, J1738+0333, J1909$-$3744, and J2222$-$0137) and two
asymmetric neutron-star neutron-star binaries (PSRs J0737$-$3039A and
J1913+1102). For comparison, we also show a constraint from the first binary
neutron star merger observed via gravitational waves LIGOScientific:2018dkp .
In principle, the uncertainty in the superanuclear neutron-star matter is
entangled with strong-field gravity tests Shao:2019gjj . Nevertheless,
nowadays we have enough well-measured binary pulsar systems to populate the
whole mass range for neutron stars, and a combined study Shao:2017gwu ;
Zhao:2022vig  has verified that for each reasonable equation of state, the
possibility for spontaneous scalarization in the Damour-Esposito-Farèse
scalar-tensor gravity theory is very low. Following the method developed by
Shao et al. Shao:2017gwu , a dedicated Bayesian parameter-estimation study
combining the above-mentioned seven pulsar systems has basically closed the
possibility of developing spontaneous scalarization for an effective scalar
coupling larger than $10^{-2}$ for the theory given by Eqs. (7)–(10), no
matter of the underlying yet-uncertain equation of state for supranuclear
neutron-star matters.

It is worth to mention that, when performing Markov-chain Monte Carlo Bayesian
parameter estimation, the integration of the modified Tolman-Oppenheimer-
Volkoff equations needs to be carried out by more than millions of times on
the fly thus computationally expensive. Recently, reduced-order surrogate
models, which extract dominating features to represent accurate enough
integration results, were bulit to aid the speedup of the calculation
Zhao:2019suc ; Guo:2021leu . The codes of these reduced-order surrogate models
are publicly available at https://github.com/BenjaminDbb/pySTGROM and
https://github.com/mh-guo/pySTGROMX for community use.

Although the original Damour-Esposito-Farèse scalar-tensor gravity theory is
disfavored by binary pulsar timing results, in further extended, generic
scalar-tensor gravity theories, neutron stars can still be scalarized. This is
particularly true for a massive scalar-tensor theory with $V(\varphi)\sim
m^{2}\varphi^{2}$ when the Compton wavelength of the scalar field is smaller
than the orbital separation of the binary Ramazanoglu:2016kul ;
Yazadjiev:2016pcb ; Xu:2020vbs . Basically the modification with respect to
the general relativity in the orbital dynamics is suppressed exponentially in
a Yukawa fashion. Fortuitously, without giving much details, such kind of
massive scalar-tensor theories can be efficiently probed via the tidal
deformability measurement in gravitational waves Hu:2020ubl ; Hu:2021tyw ;
TheLIGOScientific:2017qsa . In this sense, a combination of pulsar timing data
and gravitational wave data is called for to probe a larger parameter space
for scalar-tensor gravity theories Shao:2017gwu .

In the past few years, other variants of scalar-tensor gravity theories
triggered great enthusiasm. Some of them not only give scalarized neutron
stars, but also scalarized black holes, in contrast to the no-hair theorem. A
particularly interesting class of such theory includes a topological Gauss-
Bonnet term,

$\mathcal{G}=R_{\alpha\beta\gamma\delta}R^{\alpha\beta\gamma\delta}-4R_{\alpha\beta}R^{\alpha\beta}+R^{2}\,,$
(13)

coupling to the scalar field Doneva:2017bvd ; Silva:2017uqg ; Xu:2021kfh . In
Eq. (13), $R_{\alpha\beta\gamma\delta}$ and $R_{\alpha\beta}$ are the Riemann
tensor and Ricci tensor respectively. Preliminary constraints on the scalar-
Gauss-Bonnet gravity from binary pulsars are presented by Danchev et al.
Danchev:2021tew . This is a new field where observations of compact objects
including neutron stars and black holes are crucial to reveal the strong-field
information of gravitation.

##  3 Radiative effects in massive gravity theories

Radiative tests from binary pulsars are powerful, as the related PPK
parameter, $\dot{P}_{b}$, can be very well measured from a long-term timing
project on suitable pulsars Damour:1991rd . This parameter improves with
observational time span $T_{\rm obs}$ quite fast, as $T^{-5/2}_{\rm obs}$. The
orbital decay rate $\dot{P}_{b}$ is not only useful for constraining the
dipolar gravitational wave emission, but also in other radiative aspects of
gravitation, for example, in constraining the extra radiation caused by the
breaking down of the Lorentz symmetry Yagi:2013ava , or a nonzero mass of
gravitons Finn:2001qi ; Miao:2019nhf ; deRham:2012fw ; Shao:2020fka . Here we
give a brief introduction to the latter.

In general relativity, the hypothetical quantum particle for gravity,
graviton, is a massless spin-2 particle. However, massive gravity theories are
found to provide interesting phenomena related to the evolution of the
Universe, e.g. the accelerated expansion and dark energy deRham:2012fw .
Therefore, probing the upper bounds of the graviton mass is fundamentally
important to field theories and cosmology studies, and it is one of the
central topics in gravitational physics.

One of the early study of using binary pulsars to test the graviton mass was
performed by Finn and Sutton in 2002 Finn:2001qi . They investigated a
linearized gravity with a massive graviton with the action,

$\displaystyle S=\frac{1}{64\pi}\int\mathrm{d}^{4}x\bigg{[}$
$\displaystyle\partial_{\lambda}h_{\mu\nu}\partial^{\lambda}h^{\mu\nu}-2\partial^{\nu}h_{\mu\nu}\partial_{\lambda}h^{\mu\lambda}+2\partial^{\nu}h_{\mu\nu}\partial^{\mu}h$
$\displaystyle-\partial^{\mu}h\partial_{\mu}h-32\pi
h_{\mu\nu}T^{\mu\nu}+m_{\rm
g}^{2}\left(h_{\mu\nu}h^{\mu\nu}-\frac{1}{2}h^{2}\right)\bigg{]}\,,$ (14)

where the last term gives a unique graviton mass under certain
conditions555The conditions are that (i) the wave equation takes a standard
form  $\left(\square-m_{\rm g}^{2}\right)\bar{h}_{\mu\nu}+16\pi
T_{\mu\nu}=0\,,$ (15) and the theory recovers the general relativity in the
limit when $m_{\rm g}\to 0$, namely, there is no van Dam-Veltman-Zakharov
discontinuity Finn:2001qi . Finn:2001qi  while the others are just linearized
expansions from the Einstein-Hilbert action with $h_{\mu\nu}\equiv
g_{\mu\nu}-\eta_{\mu\nu}$ and $h\equiv h^{\mu}_{~{}~{}\mu}$. It was shown that
extra gravitational wave radiation exists in theory (3), which results in a
fractional change in the orbital decay rate, by Finn:2001qi

$\displaystyle\frac{\dot{P}_{\rm b}-\dot{P}_{\rm
b}^{\mathrm{GR}}}{\dot{P}_{\rm
b}^{\mathrm{GR}}}=\frac{5}{24}\frac{\left(1-e^{2}\right)^{3}}{1+\frac{73}{24}e^{2}+\frac{37}{96}e^{4}}\left(\frac{P_{\rm
b}}{2\pi\hbar}\right)^{2}m_{\rm g}^{2}\,.$ (16)

Here $\dot{P}_{\rm b}^{\mathrm{GR}}$ is the value predicted by the general
relativity in Eq. (3). Notice that the fractional change is proportional to
$\propto P_{\rm b}^{2}m_{\rm g}^{2}$. Therefore, if the precision of
$\dot{P}_{\rm b}$ is given, binary pulsars with larger orbits have a larger
figure of merit for the test. However, usually, the precision of $\dot{P}_{\rm
b}$ crucially depends on the orbital size, and it turns out that, still,
binary pulsars with smaller orbits have a larger figure of merit.

The most recent constraint in this Finn-Sutton framework was provided by a
combination of multiple best-timed binary pulsars with a Bayesian statistical
treatment. A collection of nine best-timed binary pulsars (PSRs J0348+0432,
J0737$-$3039, J1012+5307, B1534+12, J1713+0747, J1738+0333, J1909$-$3744,
B1913+16, and J2222$-$0137) provide a tight bound on the graviton mass,

$\displaystyle m_{g}<5.2\times 10^{-21}\,\mathrm{eV}/c^{2},\quad(90\%\text{
C.L.})\,,$ (17)

using a uniform prior in $\ln m_{\rm g}$ Miao:2019nhf . This limit is not the
strongest limit on the graviton mass deRham:2016nuf . However, from a
theoretical point of view, it is a bound from binary orbital dynamics,
complementary to, e.g. the kinematic dispersion-relation tests from the
LIGO/Virgo/KAGRA observation of gravitational waves LIGOScientific:2021sio .
It is worth mentioning that the theory (3) has some drawbacks including ghosts
and instability Finn:2001qi ; deRham:2016nuf , and here it is only used as a
strawman target for illustration.

It is interesting to note, that in different massive gravity theories, the
dependence of the extra radiation on the graviton mass is in general
different. It depends on the specifics of the illustrated gravity theory. This
is due to the deep fundamental principles in the designs of a number of
variants of massive gravity theories. For example, in a cosmologically
motivated massive gravity theory, known as the cubic Galileon theory with the
action deRham:2012fw ,

$\displaystyle
S=\int\mathrm{d}^{4}x\left[-\frac{1}{4}h^{\mu\nu}(\mathcal{E}h)_{\mu\nu}+\frac{h^{\mu\nu}T_{\mu\nu}}{2M_{\mathrm{Pl}}}-\frac{3}{4}\left(\partial\varphi\right)^{2}\left(1+\frac{1}{3m_{\rm
g}^{2}M_{\mathrm{Pl}}}\square\varphi\right)+\frac{\varphi
T}{2M_{\mathrm{Pl}}}\right]\,,$ (18)

the specific way of the addition of the scalar field $\varphi$ introduces the
so-called screening mechanism, thus avoids the stringent constraints from the
Solar System, yet provides important changes to the cosmological evolution. In
the action (18), $\varphi$ is the Galileon scalar field, $T_{\mu\nu}$ is the
matter energy-momentum tensor, $T\equiv T^{\mu}_{~{}~{}\mu}$, $M_{\rm Pl}$ is
the Planck mass, and

$(\mathcal{E}h)_{\mu\nu}\equiv-\frac{1}{2}\square h_{\mu\nu}+\cdots$ (19)

is the Lichnerowicz operator. For a central massive body with mass $M$, the
screening radius is $r_{\star}=\big{(}M/16m_{\rm g}^{2}M_{\rm
Pl}^{2}\big{)}^{1/3}$, within which, the theory exhibits strong couplings and
it reduces to the canonical gravity.

Figure 3: Cumulative probability for the graviton mass with two different
priors in the cubic Galileon theory Shao:2020fka . Shaded regions show the
excluded graviton mass values at the 95% confidence level.

According to de Rham et al. deRham:2012fw , though with a screening mechanism
to suppress modification at the high density region within $r_{\star}$, this
cubic Galileon theory predicts a different scaling behaviour for the
gravitational radiation. For a system with a typical length scale $L$, the
fifth-force suppression factor is $\sim\big{(}L/r_{\star}\big{)}^{3/2}$, and
the suppression factor for the gravitational radiation is $\sim(P_{\rm
b}/r_{\star})^{3/2}$. As for a binary system, one has $L\sim vP_{\rm b}$ where
$v$ is a characteristic velocity. Therefore, the gravitational radiation is,
compared with the fifth force, less suppressed by a factor of $v^{3/2}$, and
it provides a valuable window to look for evidence of this theory via
radiative channels, for example, in binary pulsar systems.

Analytic radiative powers were worked out by de Rham et al. deRham:2012fw ,
and the extra radiative channels include monopolar radiation, dipolar
radiation, and quadrupolar radiation. For binary pulsar systems with different
orbital periods and orbital eccentricities, the dominate radiation channel can
be different Shao:2020fka . For the current set of binary pulsars, the
quadrupole radiation is the dominating factor among the extra channels
deRham:2012fw ; Shao:2020fka .

The most up-to-date constraint from binary pulsars is

$\displaystyle m_{g}<2\times 10^{-28}\,\mathrm{eV}/c^{2},\quad(95\%\text{
C.L.})\,,$ (20)

for the cubic Galileon theory, and the cumulative probability distributions of
the graviton mass are given in Figure 3 for two different priors Shao:2020fka
. Such a tight constraint was obtained from the combination of fourteen best-
timed binary pulsar systems, including PSRs J0348+0432, J0437$-$4715,
J0613$-$0200, J0737$-$3039, J1012+5307, J1022+1001, J1141$-$6545, B1534+12,
J1713+0747, J1738+0333, J1756$-$2251, J1909$-$3744, B1913+16, and
J2222$-$0137\. One should keep in mind that, the limit (20) is theory
specific, and in this situation, only applies to the cubic Galileon theory
given in Eq. (18). Nonetheless, it provides an interesting example that for a
gravity theory designed for cosmological purposes at corresponding
lengthscales, binary pulsar systems with astronomical lengthscales still
provide intriguing and useful bounds. It is an illustration of using binary
pulsars in the studies of cosmology by examining the modification to binary
orbits brought by a cosmologically-motivated modified gravity.

##  4 Strong equivalence principle and dark matters

Binary pulsars are not only useful for the radiative tests introduced in the
above sections, they also provide superb limiting power in the conservative
aspects of gravitational dynamics for orbital evolutions. Below we introduce
an example of examining the strong equivalence principle via the conservative
dynamics of binary pulsars Damour:1991rq ; Zhu:2018etc , and its extension to
test certain interesting properties of dark matters Wagner:2012ui ;
Shao:2018klg .

As discovered by Damour and Schäfer Damour:1991rq , a perturbed binary orbit
with an equivalence-principle-violating abnormal acceleration has a
characteristic evolution in its orbital elements. The notable change is the
appearance of a vectorized superposition of two eccentricity vectors for the
real orbital eccentricity. It provides a graphical understanding of the
underlying dynamics for a binary in presence of equivalence principle
violations. The real orbital eccentricity vector, $\bm{e}(t)$, is an addition
of a rotating normal eccentricity vector, $\bm{e}_{\rm R}(t)$, in its post-
Newtonian fashion, and an extra abnormal eccentricity vector,
$\bm{e}_{\Delta}$, which is time independent and whose length is proportional
to the Eötvös parameter, $\Delta$, describing the violation of the equivalence
principle. If $\Delta=0$, the abnormal eccentricity vector $\bm{e}_{\Delta}=0$
and it returns to the precessing case in the general relavitity. A graphical
illustration is given in Figure 4. As we discussed in Sec. 1, the pulsar
timing technique is very sensitive to tiny changes in the orbit, and such a
change can be captured in pulsar timing data Damour:1991rq .

Figure 4: Graphical illustration of the time-varying orbital eccentricity
vector, $\bm{e}(t)$, for a binary pulsar, in the presence of strong
equivalence principle violation Damour:1991rq . The orbital eccentricity
vector evolves according to $\bm{e}(t)=\bm{e}_{\Delta}+\bm{e}_{\rm R}(t)$,
where $\bm{e}_{\rm R}(t)$ is the usual precessing eccentricity vector in the
general relativity, and the constant abnormal eccentricity is in the direction
of $\bm{a}_{\perp}$, which is the projection of external Galactic acceleration
in the orbital plane.

At the beginning, such a scenario was applied to a few binary pulsars in a
statistical sense by marginalizing over some unknown angles to obtain
constraints on the violation of the equivalence principle Damour:1991rq .
Later it was implemented to a handful of binary pulsars with an improved
statistical methodology to better account for the movements of binary pulsars
in the Milky Way Wex:2014nva . Then, with better data and more information
about binary pulsar systems, a direct method was developed Freire:2012nb . The
direct method not only can constrain the equivalence principle violation, but
in principle can detect it if it exists.

The most stringent limit using binary pulsars comes from a precisely timed
long-orbital-period binary pulsar, PSR J1713+0747 Zhu:2018etc , as larger
orbits have higher figures of merit in such a test Damour:1991rq . Using the
improved direct method, the limit on the Eötvös parameter from PSR J1713+0747
is Zhu:2018etc ,

$\displaystyle|\Delta|<2\times 10^{-3}\,,\quad(95\%\text{ C.L.})\,.$ (21)

Though it is much less limiting than the constraint earlier obtained from the
Solar System Wagner:2012ui ; Will:2014kxa , the limit (21) encodes strong-
field effects. For example, in the case of the aforementioned scalar-tensor
gravity, the strong-field version of Eötvös parameter will be very different
from its weak-field counterpart Freire:2012nb . Therefore, such a limit from
neutron stars is a standalone bound and applicable to the strong version of
equivalence principle Wex:2014nva ; Shao:2016ezh .

The limit (21) is not only interesting to gravitational physics, it also has
its value when we look at it from a different angle. As we now know, the
binary pulsar is actually immersed in the ocean of dark matters in the Milky
Way. As we have not really understood what the very nature of dark matter is,
the above method for testing the equivalence principle provides a non-
traditional probe to dark matter’s properties. Shao et al. Shao:2018klg
proposed a method where such a limit, with a proper handle, can be converted
to the interaction properties between dark matters and ordinary matters.

If there is a long-range fifth force between dark matter particles and
ordinary matter fields, as many field theories will suggest Wagner:2012ui , it
is likely to introduce an apparent violation of the strong equivalence
principle if we have not taken the fifth force into account in our standard
assumptions. The role of the Galactic acceleration in Figure 4, whose
projection on the orbital plane is $\bm{a}_{\perp}$, is replaced by the
attraction of dark matters to the binary system. The difference in the
acceleration to two binary components (a neutron star and a white dwarf in the
case of PSR J1713+0747), described by $\Delta$, is replaced by a quantity
related to the long-range fifth-force between dark matters and ordinary
standard-model matters Shao:2018klg .

Detailed analysis of PSR J1713+0747 Shao:2018klg  took into consideration of
the Galactic distribution of dark matters, and gave a very different bound in
nature that could be obtained from terrestrial experiments Wagner:2012ui . The
current observational data of PSR J1713+0747 already imply that, if there is
such a long-range fifth force between dark matters and ordinary matters, its
magnitude should be no more than 1% of the gravitational force between them.
Such a limit provides a useful complement to other types of dark-matter
experiments, which are usually looking for short-range forces between the
hypothesized dark-matter particles and the standard-model particles
Wagner:2012ui , including the searches in underground laboratories, particle
colliders, and X-ray/$\gamma$-ray observations via high-energy satellites.

##  5 Summary

In this chapter, we present some basic concepts of using binary pulsars as
fundamental clocks in a curved spacetime to probe various types of
modifications to the binary orbits. These modifications could have been caused
by a modified gravity theory or some other new physics like a long-range fifth
force between dark matters and ordinary matters. As pulsar timing provides us
with very accurate measurements, it puts constraints on tiny changes caused by
an alternative gravity theory other than the general relativity. Moreover,
neutron stars are intrinsically strong-gravity objects, and nonperturbative
aspects of the strong-field gravity can also be studied via radio pulsar
experiments. Actually, quite many strong-field limits are still best provided
by pulsar timing experiments, even nowadays in presence of new types of
observations like gravitational waves and black hole shadows. A careful study
shows that the limits from pulsar timing are actually complementary to those
from gravitational wave detections and black hole shadows Shao:2017gwu ;
EventHorizonTelescope:2022xqj . Proper combinations of these strong-gravity
experiments could provide a more complete landscape to gravitation in the
strong-field .

Solely focusing on the radio pulsar side, the timing experiments can be
carried out for decades, in particular for some interesting systems like the
Hulse-Taylor pulsar PSR B1913+16 Weisberg:2016jye  and the Double Pulsar PSR
J0737$-$3039A/B Kramer:2021jcw . Long-term observations improve the precision
of PPK parameters with the observational time span $T_{\rm obs}$. For
examples, the precision in the orbital decay parameter, $\dot{P}_{\rm b}$,
improves very fast, as $T_{\rm obs}^{-5/2}$, and the precision in the
periastron advance rate, $\dot{\omega}$, improves as $T_{\rm obs}^{-3/2}$.
Furthermore, the sensitivity of radio telescopes is also improving, notably
with the Five-hundred-meter Aperture Spherical Telescope in China
Jiang:2019rnj ; Lu:2019gsr  and the Square Kilometre Array in South Africa and
Australia Shao:2014wja ; Bull:2018lat . The former has already been operating
for a couple of years, while the latter has also entered the construction
phase recently. The improvement in the sensitivity of radio telescopes
directly converts to improvements in the timing precision. Therefore, the real
improvement for PPK parameters is faster than the theoretical power law
predictions. Last but not the least, radio telescopes are also continuously
discovering new pulsar systems, and some of these systems with suitable system
properties will contribute to strong-field gravity tests. We are even looking
forward to discovering yet-undetected binary pulsar systems like neutron-star
black-hole binaries with short orbital periods $P_{\rm b}\lesssim 1$ day or
pulsars around the Sgr A∗ black hole with orbital periods $P_{\rm b}\lesssim
10$ years Liu:2014uka ; Bower:2018mta ; Bower:2019ads , which will provide
completely new gravity tests in the strong-field regimes Kramer:2004hd .

###### Acknowledgements.

We are grateful to the 740th WE-Heraeus-Seminar “Experimental Tests and
Signatures of Modified and Quantum Gravity”, organized by Christian Pfeifer
and Claus Lämmerzahl. LS was supported by the National SKA Program of China
(2020SKA0120300), the National Natural Science Foundation of China (11975027,
11991053, 11721303), and the Max Planck Partner Group Program funded by the
Max Planck Society.

## References

  * [1] B. P. Abbott et al.  GW170817: Observation of gravitational waves from a binary beutron star inspiral.  Phys. Rev. Lett., 119(16):161101, 2017.
  * [2] B. P. Abbott et al.  Tests of General Relativity with GW170817.  Phys. Rev. Lett., 123(1):011102, 2019.
  * [3] B. P. Abbott et al.  Tests of General Relativity with the Binary Black Hole Signals from the LIGO-Virgo Catalog GWTC-1.  Phys. Rev. D, 100(10):104036, 2019.
  * [4] R. Abbott et al.  Tests of general relativity with binary black holes from the second LIGO-Virgo gravitational-wave transient catalog.  Phys. Rev. D, 103(12):122002, 2021.
  * [5] R. Abbott et al.  Tests of general relativity with GWTC-3.  2021\.  arXiv:2112.06861.
  * [6] K. Akiyama et al.  First M87 Event Horizon Telescope Results. I. The Shadow of the Supermassive Black Hole.  Astrophys. J. Lett., 875:L1, 2019.
  * [7] K. Akiyama et al.  First Sagittarius A* Event Horizon Telescope Results. I. The Shadow of the Supermassive Black Hole in the Center of the Milky Way.  Astrophys. J. Lett., 930(2):L12, 2022.
  * [8] K. Akiyama et al.  First Sagittarius A* Event Horizon Telescope Results. VI. Testing the Black Hole Metric.  Astrophys. J. Lett., 930(2):L17, 2022.
  * [9] A. M. Archibald, N. V. Gusinskaia, J. W. T. Hessels, A. T. Deller, D. L. Kaplan, D. R. Lorimer, R. S. Lynch, S. M. Ransom, and I. H. Stairs.  Universality of free fall from the orbital motion of a pulsar in a stellar triple system.  Nature, 559(7712):73–76, 2018.
  * [10] E. Berti et al.  Testing general relativity with present and future astrophysical observations.  Class. Quant. Grav., 32:243001, 2015.
  * [11] G. C. Bower et al.  Galactic Center Pulsars with the ngVLA.  ASP Conf. Ser., 517:793, 2018.
  * [12] G. C. Bower et al.  Fundamental Physics with Galactic Center Pulsars.  Bulletin of the American Astronomical Society, 51(3):438, May 2019\.
  * [13] C. Brans and R. H. Dicke.  Mach’s principle and a relativistic theory of gravitation.  Phys. Rev., 124:925–935, 1961.
  * [14] M. Burgay et al.  An increased estimate of the merger rate of double neutron stars from observations of a highly relativistic system.  Nature, 426:531–533, 2003.
  * [15] T. Damour and N. Deruelle.  General relativistic celestial mechanics of binary systems. II. The post-Newtonian timing formula.  Ann. Inst. Henri Poincaré Phys. Théor., 44:263–292, 1986\.
  * [16] T. Damour and G. Esposito-Farèse.  Nonperturbative strong field effects in tensor-scalar theories of gravitation.  Phys. Rev. Lett., 70:2220–2223, 1993.
  * [17] T. Damour and G. Esposito-Farèse.  Tensor-scalar gravity and binary pulsar experiments.  Phys. Rev. D, 54:1474–1491, 1996.
  * [18] T. Damour and G. Schaefer.  New tests of the strong equivalence principle using binary pulsar data.  Phys. Rev. Lett., 66:2549–2552, 1991.
  * [19] T. Damour and J. H. Taylor.  Strong field tests of relativistic gravity and binary pulsars.  Phys. Rev. D, 45:1840–1868, 1992.
  * [20] V. I. Danchev, D. D. Doneva, and S. S. Yazadjiev.  Constraining scalarization in scalar-Gauss-Bonnet gravity through binary pulsars.  2021\.  arXiv:2112.03869.
  * [21] C. de Rham, J. T. Deskins, A. J. Tolley, and S.-Y. Zhou.  Graviton Mass Bounds.  Rev. Mod. Phys., 89(2):025004, 2017.
  * [22] C. de Rham, A. J. Tolley, and D. H. Wesley.  Vainshtein mechanism in binary pulsars.  Phys. Rev. D, 87(4):044025, 2013.
  * [23] D. D. Doneva and S. S. Yazadjiev.  New Gauss-Bonnet black holes with curvature-induced scalarization in extended scalar-tensor theories.  Phys. Rev. Lett., 120(13):131103, 2018.
  * [24] D. M. Eardley.  Observable effects of a scalar gravitational field in a binary pulsar.  The Astrophysical Journal, 196:L59–L62, Mar. 1975.
  * [25] G. Esposito-Farèse.  Tests of scalar-tensor gravity.  AIP Conf. Proc., 736(1):35–52, 2004.
  * [26] L. S. Finn and P. J. Sutton.  Bounding the mass of the graviton using binary pulsar observations.  Phys. Rev. D, 65:044022, 2002.
  * [27] P. C. C. Freire, M. Kramer, and N. Wex.  Tests of the universality of free fall for strongly self-gravitating bodies with radio pulsars.  Class. Quant. Grav., 29:184007, 2012.
  * [28] P. C. C. Freire, N. Wex, G. Esposito-Farèse, J. P. W. Verbiest, M. Bailes, B. A. Jacoby, M. Kramer, I. H. Stairs, J. Antoniadis, and G. H. Janssen.  The relativistic pulsar-white dwarf binary PSR J1738+0333 II. The most stringent test of scalar-tensor gravity.  Mon. Not. Roy. Astron. Soc., 423:3328, 2012.
  * [29] M. Guo, J. Zhao, and L. Shao.  Extended reduced-order surrogate models for scalar-tensor gravity in the strong field and applications to binary pulsars and gravitational waves.  Phys. Rev. D, 104(10):104065, 2021.
  * [30] G. Hobbs et al.  Development of a pulsar-based timescale.  Mon. Not. Roy. Astron. Soc., 427:2780–2787, 2012.
  * [31] H. Hu, M. Kramer, N. Wex, D. J. Champion, and M. S. Kehl.  Constraining the dense matter equation-of-state with radio pulsars.  Mon. Not. Roy. Astron. Soc., 497(3):3118–3130, 2020.
  * [32] Z. Hu, Y. Gao, R. Xu, and L. Shao.  Scalarized neutron stars in massive scalar-tensor gravity: X-ray pulsars and tidal deformability.  Phys. Rev. D, 104(10):104014, 2021.
  * [33] P. Jiang et al.  Commissioning progress of the FAST.  Sci. China Phys. Mech. Astron., 62(5):959502, 2019.
  * [34] M. Khalil, R. F. P. Mendes, N. Ortiz, and J. Steinhoff.  Effective-action model for dynamical scalarization beyond the adiabatic approximation.  2022\.  arXiv:2206.13233.
  * [35] M. Kramer.  Pulsars as probes of gravity and fundamental physics.  Int. J. Mod. Phys. D, 25(14):1630029, 2016.
  * [36] M. Kramer, D. C. Backer, J. M. Cordes, T. J. W. Lazio, B. W. Stappers, and S. Johnston.  Strong-field tests of gravity using pulsars and black holes.  New Astron. Rev., 48:993–1002, 2004.
  * [37] M. Kramer et al.  Strong-field gravity tests with the double pulsar.  Phys. Rev. X, 11(4):041050, 2021.
  * [38] K. Liu, R. P. Eatough, N. Wex, and M. Kramer.  Pulsar–black hole binaries: prospects for new gravity tests with future radio telescopes.  Mon. Not. Roy. Astron. Soc., 445(3):3115–3132, 2014.
  * [39] D. R. Lorimer and M. Kramer.  Handbook of Pulsar Astronomy.  Cambridge University Press, Cambridge, England, 2005.
  * [40] J. Lu, K. Lee, and R. Xu.  Advancing Pulsar Science with the FAST.  Sci. China Phys. Mech. Astron., 63(2):229531, 2020.
  * [41] A. G. Lyne et al.  A double-pulsar system: A rare laboratory for relativistic gravity and plasma physics.  Science, 303:1153–1157, 2004.
  * [42] R. N. Manchester.  Pulsars and gravity.  Int. J. Mod. Phys. D, 24(06):1530018, 2015.
  * [43] R. N. Manchester, G. B. Hobbs, A. Teoh, and M. Hobbs.  The Australia Telescope National Facility pulsar catalogue.  Astron. J., 129:1993, 2005.
  * [44] X. Miao, L. Shao, and B.-Q. Ma.  Bounding the mass of graviton in a dynamic regime with binary pulsars.  Phys. Rev. D, 99(12):123015, 2019.
  * [45] F. M. Ramazanoğlu and F. Pretorius.  Spontaneous scalarization with massive fields.  Phys. Rev. D, 93(6):064005, 2016.
  * [46] N. Sennett, L. Shao, and J. Steinhoff.  Effective action model of dynamically scalarizing binary neutron stars.  Phys. Rev. D, 96(8):084019, 2017.
  * [47] L. Shao.  Testing the strong equivalence principle with the triple pulsar PSR J0337+1715.  Phys. Rev. D, 93(8):084023, 2016.
  * [48] L. Shao.  Degeneracy in Studying the Supranuclear Equation of State and Modified Gravity with Neutron Stars.  AIP Conf. Proc., 2127(1):020016, 2019.
  * [49] L. Shao.  General relativity withstands double pulsar’s scrutiny.  APS Physics, 14:173, 2021.
  * [50] L. Shao et al.  Testing Gravity with Pulsars in the SKA Era.  In Advancing Astrophysics with the Square Kilometre Array, volume AASKA14, page 042. Proceedings of Science, 2015.
  * [51] L. Shao, N. Sennett, A. Buonanno, M. Kramer, and N. Wex.  Constraining nonperturbative strong-field effects in scalar-tensor gravity by combining pulsar timing and laser-interferometer gravitational-wave detectors.  Phys. Rev. X, 7(4):041025, 2017.
  * [52] L. Shao and N. Wex.  Tests of gravitational symmetries with radio pulsars.  Sci. China Phys. Mech. Astron., 59(9):699501, 2016.
  * [53] L. Shao, N. Wex, and M. Kramer.  Testing the universality of free fall towards dark matter with radio pulsars.  Phys. Rev. Lett., 120(24):241104, 2018.
  * [54] L. Shao, N. Wex, and S.-Y. Zhou.  New graviton mass bound from binary pulsars.  Phys. Rev. D, 102(2):024069, 2020.
  * [55] H. O. Silva, J. Sakstein, L. Gualtieri, T. P. Sotiriou, and E. Berti.  Spontaneous scalarization of black holes and compact stars from a Gauss-Bonnet coupling.  Phys. Rev. Lett., 120(13):131104, 2018.
  * [56] I. H. Stairs.  Testing general relativity with pulsar timing.  Living Rev. Rel., 6:5, 2003.
  * [57] J. H. Taylor.  Pulsar timing and relativistic gravity.  Phil. Trans. A. Math. Phys. Eng. Sci., 341(1660):117–134, 1992\.
  * [58] G. Voisin, I. Cognard, P. C. C. Freire, N. Wex, L. Guillemot, G. Desvignes, M. Kramer, and G. Theureau.  An improved test of the strong equivalence principle with the pulsar in a triple star system.  Astron. Astrophys., 638:A24, 2020.
  * [59] T. A. Wagner, S. Schlamminger, J. H. Gundlach, and E. G. Adelberger.  Torsion-balance tests of the weak equivalence principle.  Class. Quant. Grav., 29:184002, 2012.
  * [60] J. M. Weisberg and Y. Huang.  Relativistic measurements from timing the binary pulsar PSR B1913+16.  Astrophys. J., 829(1):55, 2016.
  * [61] A. Weltman et al.  Fundamental Physics with the Square Kilometre Array.  Publ. Astron. Soc. Austral., 37:e002, 2020.
  * [62] N. Wex.  Testing relativistic gravity with radio pulsars.  In S. M. Kopeikin, editor, Frontiers in Relativistic Celestial Mechanics: Applications and Experiments, volume 2, page 39. Walter de Gruyter GmbH, Berlin/Boston, 2014.
  * [63] N. Wex and M. Kramer.  Gravity tests with radio pulsars.  Universe, 6(9):156, 2020.
  * [64] C. M. Will.  The confrontation between general relativity and experiment.  Living Rev. Rel., 17:4, 2014.
  * [65] R. Xu, Y. Gao, and L. Shao.  Strong-field effects in massive scalar-tensor gravity for slowly spinning neutron stars and application to X-ray pulsar pulse profiles.  Phys. Rev. D, 102(6):064057, 2020.
  * [66] R. Xu, Y. Gao, and L. Shao.  Neutron stars in massive scalar-Gauss-Bonnet gravity: Spherical structure and time-independent perturbations.  Phys. Rev. D, 105(2):024003, 2022.
  * [67] K. Yagi, D. Blas, E. Barausse, and N. Yunes.  Constraints on Einstein-Æther theory and Hořava gravity from binary pulsar observations.  Phys. Rev. D, 89(8):084067, 2014.  [Erratum: Phys.Rev.D 90, 069902 (2014), Erratum: Phys.Rev.D 90, 069901 (2014)].
  * [68] S. S. Yazadjiev, D. D. Doneva, and D. Popchev.  Slowly rotating neutron stars in scalar-tensor theories with a massive scalar field.  Phys. Rev. D, 93(8):084038, 2016.
  * [69] J. Zhao, P. C. C. Freire, M. Kramer, L. Shao, and N. Wex.  Closing a spontaneous-scalarization window with binary pulsars.  Class. Quant. Grav., 39(11):11LT01, 2022.
  * [70] J. Zhao, L. Shao, Z. Cao, and B.-Q. Ma.  Reduced-order surrogate models for scalar-tensor gravity in the strong field regime and applications to binary pulsars and GW170817.  Phys. Rev. D, 100(6):064034, 2019.
  * [71] W. W. Zhu et al.  Tests of gravitational symmetries with pulsar binary J1713+0747.  Mon. Not. Roy. Astron. Soc., 482(3):3249–3260, 2019.