measurement from the Hulse-Taylor binary system. The second-largest contribution to the observed orbital period decay is related to the effective mass loss of pulsar A, which now also has to be taken into account. Improving the precision of the double pulsar tests even further will, from now on, be constrained by our ability to correct for kinematic effects.

We pointed out that a direct comparison of tests of PN inspiral phase coefficients with different compact objects (BHs vs NSs) as well as different gravity regimes (mildly relativistic vs highly relativistic strong field) comes with certain caveats. Stating this, it is nevertheless obvious that our high-precision timing tests superbly complement the LIGO/Virgo tests, which currently are less precise at low-GW orders but allow us to probe higher-order contributions to the GW emission. This can be seen more clearly in approaches based on different theories of gravitation. For instance, in DEF gravity (Sec. VII A), the constraints from the double pulsar are considerably tighter than those from the double-NS merger GW170817 (binary BH mergers do not provide any constraints for DEF gravity). Section VII describes two alternatives to GR, namely, DEF gravity and Bekenstein's TeVeS, in some detail. We show that the new double pulsar results constrain effects that one would typically expect from certain modifications of GR including dipolar radiation and a periodic change of the MOI due to a varying local gravitational constant along the orbit of pulsar A. The double pulsar observations presented here lead to further constraints in the two-parameter space of DEF gravity and result in an additional falsification of TeVeS, which is qualitatively different to that from GW170817 [44], by improving on the limits in Ref. [208]. We also show that the double pulsar can be used to improve constraints on a time-varying gravitational constant, in particular, for effects related to strong gravitational fields. However, a more complete analysis is still needed.

The double pulsar system is, so far, the only double-NS system where the orientations of both spins relative to the total angular momentum vector can be measured. Combined with a superb knowledge about the systemic velocity and the very high precision measurements of both NS masses, the system is a unique laboratory to study binary evolution, the mechanisms of core-collapse supernovae, and the formation and structure of NSs. Obviously, the double pulsar is a rich laboratory for a wide variety of physics and astrophysics but especially for testing GR and its alternatives. For the latter application, past and future observations of the system provide constraints that are highly complementary to other methods, such as observations with GW detectors, the study of (supermassive) BHs via VLBI imaging [227], or orbits of stars [228] and flares [229]. Figure 21 demonstrates one way of using a two-parameter space to illustrate the complementarity of different gravity experiments and put the experiments presented in this paper into context. The first parameter is the

<figure><img src="image_24.png" /><figcaption>FIG. 21. Parameter space for putting the experiments of this paper into context with other gravity experiments. $\Phi$ is the potential of the gravitational interaction, and $\xi_{max}$ is the maximum spacetime curvature in the system. Experiments involving BHs are indicated by black filled circles, whereas other experiments are reported in blue. For pulsar experiments we show the Double Pulsar ('J0737' and 'J0737 (Shapiro)') and the universality-of-free-fall test with the triple system pulsar ('J0337') [76,209]. Besides these, the following experiments are shown: solar system (cf. [41]), LIGO/Virgo mergers (cf. Section VI B 2), experiments with the supermassive BHs Sgr A* and M87* [227-229,233-235], and compact double white dwarf systems ('WD-WD'; e.g., [236]). Experiments that probe the propagation of photons in curved spacetime are highlighted by red circles. Details can be found in the text.</figcaption></figure>

potential of the gravitational interaction $\Phi$. This is typically the (external) potential probed by a photon or a mass (e.g., a pulsar) in the gravitational field of another mass. $\Phi/c^2$ represents a quantity that typically enters the PN and post-Minkowskian approximation schemes. As a second parameter, we have chosen the maximum spacetime curvature in the system, $\xi_{max}$ defined as the square root of the Kretschmann scalar. Among material bodies, $\xi_{max}$ allows distinguishing between weakly and strongly self-gravitating masses. For a given EOS, a NS shows a monotonic relation between $\xi_{max}$ and, e.g., its surface potential $\Phi_s$. For NSs, one finds $\Phi_s/c^2 \sim -0.2$, which underlines the strong-field aspect of pulsar tests. For BHs, the maximum spacetime curvature (causally connected to its environment) is the one at the horizon, which is a measure for the size and mass of the BH ($\xi_{max} \propto M_{BH}^{-2}$, for a nonrotating BH). This is, for instance, of relevance for BHs in alternatives to GR, where in certain cases the magnitude of the modification decreases with increasing mass and, therefore, with decreasing $\xi_{max}$ (see, e.g., Ref. [230]).
Figure 21 illustrates well how the double pulsar with its strongly self-gravitating components probes the mildly relativistic strong field regime. The double pulsar appears twice in Fig. 21, once for the experiments related to orbital dynamics, like GW damping (Sec. VI B 2) and periastron precession (Sec. VI B 3), and a second entry (with label “Shapiro”) for the test related to photon propagation (Sec. VI B 5). It is evident that, in terms of coupling between gravitational and electromagnetic fields, the Shapiro delay test in Sec. VI B 5 is the precision timing experiment which probes the strongest spacetime curvature. When comparing the different gravity experiments in the parameter space in Fig. 21, one has to keep in mind the qualitative difference between them, for instance, BHs vs material bodies (cf. the discussion at the end of Sec. VI B 2), stationary vs dynamical or radiative situations, etc. Needless to say, a two-parameter plot cannot capture all quantities relevant for characterizing gravity tests and, therefore, always gives an incomplete comparison (for alternative parameter spaces, see, e.g., Refs. [34,228,231,232]). In conclusion, it is clear that studies of the double pulsar will continue to be extremely useful, with new applications undoubtedly awaiting us. We believe that this work presents an important milestone in this endeavor.

## ACKNOWLEDGMENTS

We acknowledge the contribution of our dear colleague and collaborator Nichi D'Amico, who sadly passed away far too early. We are indebted to a number of colleagues, who have contributed to our discussions or helped with the observations. These include in alphabetical order: David Champion, Gilles Esposito-Farèse, Aditya Parthasarathy, Nataliya Porayko, Roman Rafikov, Gerhard Schäfer, and Vivek Venkatraman Krishnan. We are especially grateful to Jumei Yao for exploring the YMW16 electron density model for our new astrometric constraints, to Michael Keith for his contributions to optimization of Tempo2 for operation with large datasets, to Victoria Grandy for calibrations of much of the GBT data, and to Lijing Shao for carefully reading the manuscript and many valuable comments that helped to improve it. We thank the staff of all observatories for their continued help over the years. The Green Bank Observatory is a facility of the National Science Foundation (NSF) operated under cooperative agreement by Associated Universities, Inc. The National Radio Astronomy Observatory is a facility of the NSF operated under cooperative agreement by Associated Universities, Inc. The Parkes radio telescope is part of the Australia Telescope National Facility (grid.421683.a) which is funded by the Commonwealth of Australia for operation as a National Facility managed by CSIRO. The Nançay Radio Observatory is operated by the Paris Observatory, associated with the French Centre National de la Recherche Scientifique (CNRS). This publication is also based on observations with the 100-m telescope of the

Max-Planck-Institut für Radioastronomie at Effelsberg. The Westerbork Synthesis Radio Telescope is operated by the Netherlands Institute for Radio Astronomy (ASTRON) with support from The Netherlands Foundation for Scientific Research (NWO). Pulsar research and access to the Lovell telescope is supported by consolidated grants from the Science and Technology Facilities Council in the United Kingdom. We acknowledge financial support from the Action Fédératrice PhyFOG funded by Paris Observatory and from the “Programme National Gravitation, Références, Astronomie, Métrologie” (PNGRAM) funded by CNRS/INSU and CNES, France. Pulsar research at U. B. C. is supported by an NSERC Discovery Grant and by the Canadian Institute for Advanced Research. J. W. M. is supported by a CITA Postdoctoral Fellowship. This work was supported by the Natural Sciences and Engineering Research Council of Canada (NSERC), funding reference No. CITA 490888-16. A. T. D. is supported by an Australian Research Council Future Fellowship (FT150100415). D. R. L. acknowledges support from the NSF Grants No. AAG-1616042 and No. OIA-1458952. P. C. C. F. gratefully acknowledges financial support by the European Research Council, under the European Unions Seventh Framework Program (FP/2007-2013) Grant Agreement No. 279702 (BEACON), and continuing support from the Max Planck Society. This work is supported by the Max-Planck Society as part of the “LEGACY” collaboration on low-frequency gravitational wave astronomy. M. K., N. W., and G. D. acknowledge financial support by the European Research Council for the ERC Synergy Grant BlackHoleCam (ERC-2013-SyG, Grant Agreement No. 610058). Finally, we thank three anonymous referees for their kind comments and suggestions, which helped us to improve the manuscript.

## APPENDIX A: VLBI OBSERVATIONS, DATA REDUCTION, AND VERIFICATION

A total of 18 observations were made using the VLBA between 2016 October and 2018 May (VLBA project code BD193). Each observation is one hour in duration and records 256 MHz of dual-polarization data with a center frequency of 1.56 GHz. A single scan per observation on the source ICRF J074533.0+101112 is utilized for calibrating the instrumental bandpass. The source ICRF J073038.2-320820 is used as the primary phase reference calibrator, with NVSS J073709-302710 observed in beam and contemporaneously with the pulsar to refine the calibration solutions and to define a relative position reference point. Four other background radio sources are also in beam and observed simultaneously and are subsequently used as astrometric check sources. NVSS J073709-302710 is separated by 15' from PSR J0737-3039A/B on the sky and by 5' to 15' from the other background sources.

The calibration procedure we employ makes use of the pipeline described in detail in Ref. [215]. Because of the southerly location of PSR J0737-3039A/B, we relax
<table><caption>TABLE VI. VLBI astrometric results for PSR J0737–3039A/B and in-beam background sources. Listed are reference position in right ascension (R.A.) and declination (Dec.), proper motion in each coordinate, parallax, and reference epoch. The top half of the table shows results when proper motion is fixed to expected values, while the second half shows results when proper motion is freely fitted. See the text for details.</caption><thead><tr><th></th><th>J0737–3039A/B</th><th>J0736–3027</th><th>J0736–3037</th><th>J0737–3018</th><th>J0738–3025</th></tr></thead><tbody><tr><td>RA (J2000)</td><td>07<sup>h</sup>37<sup>m</sup>51<sup>s</sup>.247(1)</td><td>07<sup>h</sup>36<sup>m</sup>45<sup>s</sup>.378(1)</td><td>07<sup>h</sup>36<sup>m</sup>51<sup>s</sup>.234(1)</td><td>07<sup>h</sup>37<sup>m</sup>09<sup>s</sup>.378(1)</td><td>07<sup>h</sup>38<sup>m</sup>19<sup>s</sup>.793(1)</td></tr><tr><td>Dec (J2000)</td><td>-30<sup>°</sup>39'40".68(1)</td><td>-30<sup>°</sup>27'18".77(1)</td><td>-30<sup>°</sup>37'46".16(1)</td><td>-30<sup>°</sup>18'52".68(1)</td><td>-30<sup>°</sup>25'04".92(1)</td></tr><tr><td>Parallax, π<sub>v</sub> (mas)</td><td>1.30<sup>+0.13</sup><sub>-0.11</sub></td><td>0.02<sup>+0.06</sup><sub>-0.07</sub></td><td>-0.24<sup>+0.18</sup><sub>-0.21</sub></td><td>-0.03<sup>+0.11</sup><sub>-0.11</sub></td><td>-0.07<sup>+0.11</sup><sub>-0.11</sub></td></tr><tr><td>μ<sub>α</sub> (mas yr<sup>-1</sup>, fixed)</td><td>-2.567 ± 0.030</td><td>0.0</td><td>0.0</td><td>0.0</td><td>0.0</td></tr><tr><td>μ<sub>δ</sub> (mas yr<sup>-1</sup>, fixed)</td><td>2.082 ± 0.038</td><td>0.0</td><td>0.0</td><td>0.0</td><td>0.0</td></tr><tr><td>Position epoch (MJD)</td><td>58000</td><td>58000</td><td>58000</td><td>58000</td><td>58000</td></tr><tr><td>RA (J2000)</td><td>07<sup>h</sup>37<sup>m</sup>51<sup>s</sup>.247(1)</td><td>07<sup>h</sup>36<sup>m</sup>45<sup>s</sup>.378(1)</td><td>07<sup>h</sup>36<sup>m</sup>51<sup>s</sup>.234(1)</td><td>07<sup>h</sup>37<sup>m</sup>09<sup>s</sup>.378(1)</td><td>07<sup>h</sup>38<sup>m</sup>19<sup>s</sup>.793(1)</td></tr><tr><td>Dec (J2000)</td><td>-30<sup>°</sup>39'40".68(1)</td><td>-30<sup>°</sup>27'18".77(1)</td><td>-30<sup>°</sup>37'46".16(1)</td><td>-30<sup>°</sup>18'52".68(1)</td><td>-30<sup>°</sup>25'04".92(1)</td></tr><tr><td>Parallax, π<sub>v</sub> (mas)</td><td>1.43<sup>+0.14</sup><sub>-0.13</sub></td><td>-0.05<sup>+0.07</sup><sub>-0.08</sub></td><td>-0.35<sup>+0.27</sup><sub>-0.34</sub></td><td>-0.12<sup>+0.17</sup><sub>-0.17</sub></td><td>-0.10<sup>+0.13</sup><sub>-0.14</sub></td></tr><tr><td>μ<sub>α</sub> (mas yr<sup>-1</sup>, fitted)</td><td>-2.23<sup>+0.22</sup><sub>-0.21</sub></td><td>-0.20<sup>+0.15</sup><sub>-0.14</sub></td><td>-0.32<sup>+0.42</sup><sub>-0.48</sub></td><td>-0.25<sup>+0.26</sup><sub>-0.25</sub></td><td>-0.12<sup>+0.21</sup><sub>-0.22</sub></td></tr><tr><td>μ<sub>δ</sub> (mas yr<sup>-1</sup>, fitted)</td><td>2.79<sup>+0.45</sup><sub>-0.43</sub></td><td>0.01<sup>+0.11</sup><sub>-0.11</sub></td><td>0.06<sup>+0.61</sup><sub>-0.57</sub></td><td>-0.50<sup>+0.30</sup><sub>-0.32</sub></td><td>-0.03<sup>+0.48</sup><sub>-0.46</sub></td></tr><tr><td>Position epoch (MJD)</td><td>58000</td><td>58000</td><td>58000</td><td>58000</td><td>58000</td></tr></tbody></table>

the elevation-based flagging to times when a source is below an elevation of 15°, although we test cutoff elevations in the range 12°–18° and find minimal variation. After the application of flagging and calibration, the calibrated visibility dataset for each of the in-beam sources is divided by a source model (derived using the concatenated observations from all 18 epochs), and the source position and uncertainty at that epoch are estimated with an image-plane fit to the resultant dataset. The division by a source model removes any source structure from the image to be fitted, leaving a pointlike source whose centroid position is unaffected by changes in resolution due to the absence of different VLBA antennas in some epochs. In addition to the statistical uncertainty recorded from the image-plane fit, an estimate for the systematic position uncertainty is made using Eq. (1) of Ref. [215], and this uncertainty is added in quadrature to form a total estimated positional uncertainty for each source at every epoch. As discussed in Ref. [215], this empirical estimator of systematic position uncertainties encapsulates likely contributions such as the differential troposphere and ionosphere between the sight lines to the pulsar and the nearby calibrators.

Again following Ref. [215], we estimate the astrometric parameters for PSR J0737–3039A/B and the in-beam background sources using a bootstrap method. The orbital reflex motion of the double pulsar is negligibly small and is not included in the fit. We repeat this process twice—once utilizing the proper motion as a free parameter in the fit and once fixing the proper motion to the expected value (zero for the background sources, and the pulsar timing result reported in Table IV for the pulsar). In this latter case, a proper motion value for right ascension and declination is randomly drawn for each bootstrap trial using the corresponding mean and standard deviation from the pulsar timing fit; i.e., the uncertainty in the pulsar timing proper motion is accounted for.

The best-fit VLBI results for PSR J0737–3039A/B and the background radio sources used as astrometric checks are shown in Table VI. Given that the background radio sources are expected to be distant radio-emitting active galactic nuclei, the parallax and proper motion measured for these sources are expected to be consistent with zero. Uncertainties for proper motion and parallax are derived from the bootstrap procedure described above. The uncertainties on the reference position are dominated by the uncertainty in the absolute position of the in-beam reference source NVSS J073709–302710. This quantity is relatively poorly constrained, being derived purely from phase referencing from the primary calibrator ICRF J073038.2–320820, and we estimate a nominal uncertainty of 10 mas based on the approximately 2° angular separation between these two sources. Accordingly, we assign an uncertainty of 10 mas to the reference position for all sources in Table VI. While we focus on the proper motion and parallax (which are unaffected by errors in the absolute position that are constant in time) from this point forward, we do note that comparing the reference positions for PSR J0737–3039A/B obtained from timing and VLBI at a common epoch shows agreement within this nominal 10-mas uncertainty.

As our reference case, we utilize the fit in which the proper motion is constrained based on the timing value (top section of Table VI) but note that the best-fit value is consistent to within $1\sigma$ regardless of whether the proper motion is fixed or not. The results are similarly insensitive to other potential choices in the data reduction, such as a scaling factor applied to the estimated systematic uncertainty added to the position measurements or the elevation flagging limit. As an illustration of the results, we show the offset in the position in right ascension as a function of time for both PSR J0737–3039A/B (where the parallax signature is clearly seen) and the in-beam background source NVSS J073709–301853 (which shows no parallax signature, as
<figure><img src="image_25.png" /><figcaption>FIG. 22. Bootstrap fit plots of right ascension vs time, with the fixed proper motion subtracted. Each line shows the result of a single bootstrap trial. The left panel shows PSR J0737-3039A/B, and the right panel shows the in-beam background source NVSS J073709-301853, which shows a parallax consistent with zero as expected.</figcaption></figure>

expected) after the subtraction of proper motion for our reference case in Fig. 22. Right ascension is shown for three reasons: first, because the elongated synthesized VLBA beam leads to greater precision in the right ascension axis than in declination; second, the parallax signature is larger in right ascension than declination; and third, because of the first two reasons, observations are scheduled at times of maximum parallax displacement in right ascension (and, hence, minimal parallax displacement in declination).

As we are primarily concerned with the parallax (and the corresponding uncertainty) for PSR J0737-3039A/B, the consistency of the astrometric results for the contemporaneously observed background radio sources with expectations is of considerable interest. As shown in Table VI, the fitted parallax and proper motion are consistent with expectations at the $1\sigma$ level in 75% of cases and at the $2\sigma$ level in 100% of cases. This offers a high degree of reassurance that the parallax uncertainty for PSR J0737-3039A/B is well estimated. It is, however, necessary to note that the background sources (and PSR J07373-3039A/B) span a range of angular separations from 5' to 15' to the in-beam calibrator, and at 15' the angular separation between PSR J0737-3039A/B and the in-beam calibrator is the equal highest of all source pairs considered. While systematic errors are likely to scale with angular separation [215], the bootstrap technique employed both here and in the study made by Deller et al. [215] offers a high degree of robustness in incorporating these into the estimated uncertainties.

The parallax probability distributions resulting from the bootstrap fits described above are not perfectly Gaussian, exhibiting a slight skew with a tail toward larger values.

We now consider the potential effects of refractive image wander on the VLBI results. Refraction caused by large-scale gradients in the ionized ISM can cause the apparent

position of a radio source to shift with a characteristic timescale usually on the order of months to years, depending on the observing frequency and relative pulsar-ISM-Earth velocity [237]. For most radio sources, this astrometric offset would be uncorrectable but, as discussed in Sec. VIII B, for radio pulsars the time variability of the pulsar DM can be related back to gradients in the ionized ISM and, hence, used to estimate (a component of) the necessary correction. Figure 24 shows both the DM variations and the derived image wander (see Appendix B 2 for details) over the entire observing span.

We apply these corrections to the VLBA positions presented above and refit the dataset for parallax and reference position. When doing so, the best-fit parallax (from our reference-case fit with the proper motion constrained by the timing values) becomes $1.27_{-0.10}^{+0.11}$ mas—a negligible change. This result is unsurprising when considering Fig. 24, in which the derived offsets can be seen to be small at the times sampled by our VLBA observations (MJD 57670–58260.)

Since the refractive wander corrections that we can derive based on DM variations are necessarily incomplete, we also estimate a “worst-case” impact of refractive wander by sampling the corrections that would have been applied had the VLBA observations began in earlier years. We sample corrections from the VLBA observation dates shifted by an integer number of years in the past, generating a total of ten potential sets of corrections. Across these ten sets of potential corrections, the mean change to parallax is 0.03 mas, with a maximum change of 0.09 mas—smaller in all cases than the currently estimated uncertainty.

These refractive wander corrections can also be applied to historical VLBI datasets, in particular, to the Long Baseline Array work from Ref. [93]. As can be seen from Fig. 24, the impact during the time of these observations (MJD 53970–54500) is considerably greater than during
<figure><img src="image_26.png" /><figcaption>FIG. 23. Probability distribution for parallax resulting from VLBI with proper motion fixed to the timing value (black), from pulsar timing (red), and the weighted mean of the two (blue).</figcaption></figure>

our VLBA observations. Moreover, the LBA parallax fit is particularly sensitive to the observation at MJD 54308 (which is the sole point taken near the parallax maximum in declination), and the estimated refractive shift is particularly pronounced at this time. After correcting the estimated refractive offsets, the best-fit parallax from the LBA dataset changes appreciably, from 0.87 to 0.96 mas. We also note that the original fit to the LBA data [93] did not have the advantage of the now well-constrained proper motion and uses a least-squares fit to estimate the parallax uncertainty rather than a bootstrap fit as we do for the newer VLBA data. Repeating the analysis of the LBA data using the same approach taken here yields a 68% confidence interval for the parallax of 0.61–1.17 mas, which becomes consistent at the approximately $1\sigma$ level with the VLBA results (but still favors a larger distance, in contrast to the timing, which favors a smaller distance). We note that the increased number of observations (by a factor of almost 3), a more compact primary calibrator source, and the presence of in-beam background sources that can be used to check for unmodeled systematic astrometric offsets all act to add extra confidence in the accuracy of the uncertainty estimation for the VLBA dataset, compared to the LBA dataset.

Finally, we combine the VLBA parallax $\pi_v$ with the timing parallax $\pi_t$ into a single estimate $\pi_c$. As $\pi_v$ and $\pi_t$ are independent, a weighted mean, with weights equal to the inverse variance of each measurement, is the best linear unbiased estimator of $\pi_c$. The variances are well estimated in both cases. In a weighted mean, one would normally compute a $\chi^2$ value as a goodness of fit measure. A large $\chi^2$ would prompt a search for an error. In this case, we have only two estimates to combine, so we simply take the difference $D = \pi_v - \pi_t = 0.85$. The variance of $D$ is the sum of the variances of $\pi_v$ and $\pi_t$, which gives a standard deviation of 0.51. Thus, the measured $D$ represents a $1.7\sigma$ event—not improbable enough to cause a search for errors in $\pi_v$ and $\pi_t$. We can put it as $\chi^2(1) = (0.85/0.51)^2$ and use the $\chi^2(1)$ distribution to show that the probability of $D > 0.85 = 10\%$, confirming that the two independent estimates are adequately consistent.

The probability distribution for $\pi_c$ is derived using a Monte Carlo approach, in which we randomly draw a timing parallax and a VLBI parallax from their distributions and then take a weighted average to form a combined sample. This process is repeated 100 000 times to form the combined parallax probability distribution. The final result for the parallax is $\pi_c = 1.36_{-0.10}^{+0.12}$ mas (68% confidence levels), a shift of approximately $0.5\sigma$ from the VLBA-only estimate. For each iteration, we also compute the distance. Because of the compensating effects of the tail on the high side of the parallax distribution and the bias resulting from the parallax inversion, the distance probability distribution is close to symmetric, giving our best estimate for the distance of 735 ± 60 pc. We show the parallax probability distributions from our VLBA observations, from pulsar timing, and for the weighted mean in Fig. 23.

# APPENDIX B: ANALYSIS OF ISM ELECTRON-DENSITY FLUCTUATIONS

## 1. PSD and SF simulations

The measured fluctuations DM($t$) (Sec. IV B) are distorted by the fitting technique in two ways. First, the Tempo2 analysis fits a piecewise linear model of DM($t$) with a sampling of 100 d. This is equivalent to convolving the data with a triangle of baseline ±100 d. Second, the Monte Carlo scheme adjusts the phase of the DM($t$) sampling and then bins it in 24-d bins. This is equivalent to convolving the data by a 24-d rectangle. To estimate the spectral exponent of the original DM($t$) variations, we simulate various pure power-law processes, convolve them by the two measurement effects, and compare their power spectra with that of the observations.

We make the comparison in the spectral domain (rather than the SF), because the errors in the spectral harmonics are independent and we can adjust the level and the exponent to minimize the number of harmonics in the observed spectrum that exceed the 90% confidence limits. The samples of the SF, like those of a covariance function, are highly correlated, and this makes them appear smoother than the errors would suggest. The mean and 90% confidence limits of the best match are shown in Fig. 18 in red. We also show the mean of the expected Kolmogorov exponent in blue. Note that the spectra are relatively steep, so the PSD must be computed with prewhitening and postdarkening to avoid a bias due to spectral leakage. We prewhiten the time series with a first difference, compute the power spectrum, and correct the spectrum with the Fourier transform of the difference operator [89].

The match is good from $2 \times 10^{-4} \text{ d}^{-1}$ up to about $7 \times 10^{-3} \text{ d}^{-1}$, where white noise begins to dominate the spectrum. It is very clear that a Kolmogorov spectrum cannot fit the data. The PSD also shows a peak near the annual frequency (at a period of 335 d). This is probably
due to Earth's orbital motion through phase gradients in the ISM. As the ISM is also moving, these phase gradients can essentially Doppler shift the annual motion to a higher frequency. The same simulation is used to compare with the observed structure function in more detail. This is shown in Fig. 19.

## 2. Refractive image wander

The basic phenomenon causing intensity scintillations is multipath propagation, also referred to as scattering. This both increases the apparent angular diameter of the pulsar and broadens its pulse width. Interference between the scattered rays causes the intensity scintillation. The scattering contribution to the pulse width can be estimated directly from the bandwidth of the intensity scintillations [92] as $\tau_p = 1/(2\pi\nu_{0.5}) \sim 3 \mu s$ at 820 MHz for the double pulsar. Both diffractive and refractive intensity scintillations cause fluctuations in the TOAs of the pulses on short timescales, and phase gradients cause fluctuations of the TOAs on longer scales. However, all of these fluctuations are comparable with, but less than, the pulse width [238], so we do not correct for them in this work. Phase gradients do have an important effect, because they displace the apparent position of the pulsar by $\theta_r = \nabla\phi/k$, where $k = 2\pi/\lambda$ is the propagation constant. This displacement can be comparable with the uncertainty in parallax measurements, so we examine the effect of refractive displacement in detail.

The phase gradients for PSR J0737-3039A/B are estimated directly from intensity scintillations [92]. Unfortunately, few of our observations have sufficient frequency resolution to measure the dynamic spectra of intensity scintillations accurately. However, we can use a daily interpolation of the observed DM variations to estimate their temporal gradient, which can be scaled to phase by

$$ \phi = 2\pi \times 10^6 \text{DM}/2.41 \times 10^{-4} \nu_{\text{MHz}}. \quad (\text{B1}) $$

We then scale the temporal gradient to a spatial gradient in the direction of the velocity by dividing by the velocity $v(t)$. We must assume that there is a comparable gradient perpendicular to the velocity, since the turbulence is known to be roughly isotropic [92]; consequently, the rms gradient should be direction independent. In Fig. 24, we show the displacement $\theta_r$ at 1.56 GHz in the direction of $v$. The rms($\theta_r$) = 0.125 mas, so the total rms image wander should be $\sqrt{2}$ times greater, approximately 0.177 mas.

The dates of the VLBI observations used to estimate the parallax (and, thus, the pulsar distance) are marked as blue x symbols in Fig. 24. One can see that the later VLBI observations, which are used in this paper, are at

<figure><img src="image_27.png" /><figcaption>FIG. 24. Top: dispersion measure of the double pulsar as a function of time. The black line shows the interpolated DM variations, and the blue line shows these interpolated variations with the worst-case solar-wind contribution subtracted. Bottom: the image wander in the direction of the velocity $\theta_r$ (mas) at 1.56 GHz. The times of VLBI observations are marked with the blue x symbols, with the early observations being made with the LBA and the more recent observations with the VLBA and presented in this paper.</figcaption></figure>

a time of relatively low phase gradients (at least in the direction of the velocity). To estimate the effect on parallax, we resolve the vector displacement $\theta_r$ into components parallel to the R.A. and Dec. axes. The actual effects of these displacements are already discussed in Appendix A.

Refractive image wander also potentially has an effect on the astrometric parameters derived from pulsar timing. The timing signature of a position offset is an annual sinusoid, with the phase determined by the ratio of the perpendicular components of the offset (in R.A. and Dec., for example). For proper motion the signature is an annual sinusoid of uniformly increasing amplitude, and for parallax it is a six-month sinusoid. Since the refractive image wander (Fig. 24) has a strong annual component, attributed to the annual motion of the LOS over a gradient in the screen electron density (Appendix B 1), there is a possibility of a significant covariance between image wander and the astrometric parameters.

To assess the significance of this effect, we take the variations in refractive angle shown in Fig. 24, interpolated to daily values, and transform these to daily values of TOA offset using

$$ \Delta t_{\text{refr}} = 0.5 \theta_r^2 \frac{d s}{c(1-s)}, \quad (\text{B2}) $$

where $\theta_r$ is the refractive angle in radians, $d$ is the pulsar distance (735 pc), $c$ is the velocity of light, and $s$ is the
distance of the scattering screen from Earth as a fraction of the pulsar distance. These daily TOA offsets are linearly interpolated to the MJD of each TOA in the 4-min dataset and then scaled by the ratio $(f_{TOA}/f_{Ref})^{-2}$, where $f_{TOA}$ is the radio frequency for the TOA and $f_{Ref}$ is the reference frequency for the refractive angles shown in Fig. 24, 1.56 GHz. These time offsets are then added to the corresponding TOA. Comparison of the results of a standard Tempo2 analysis solving for the astrometric parameters with and without the refractive delay corrections shows that the uncertainties in the astrometric parameters are unaffected by the refractive image wander and that the changes in the values are all less than or about $0.2\sigma$.

In this discussion, we assume that DM($t$) is entirely due to the ISM. However, the observed DM($t$) includes a solar-wind component. The worst-case solar-wind component would be a spherical solar wind with $n_e \propto R_{Sun}^{-2}$ and $n_e = 10 \text{ cm}^{-3}$ at Earth. We subtract this from the observations and plot the result as a blue line in Fig. 24. Evidently, the worst-case solar-wind contribution is negligible, and our assumption is justified.

## APPENDIX C: PK PARAMETERS IN SCALAR-TENSOR GRAVITY

In this section, we provide expressions for the PK parameters valid for a wide class of scalar-tensor gravity (STG), including DEF gravity [40,42] used in Sec. VII A of this paper. More details can be found in reviews like Refs. [34,35,41,63,203], and the corresponding references given therein. Here, we closely follow the representation used in Ref. [42]. We restrict the discussion to PK parameters that we actually use for the tests in Sec. VII, i.e., $k$, $\gamma_E$, $s$, and $\dot{P}_b$.

### 1. Advance of periastron: $k$

For the advance of periastron in the presence of two strongly self-gravitating bodies, one finds in fully conservative, boost-invariant gravity theories

$$k = (2 + 2\gamma_{AB} - X_B\beta_{BB}^A - X_A\beta_{AA}^B)\frac{\hat{\beta}_O^2}{1 - e^2} + \mathcal{O}(\hat{\beta}_O^4), \quad (C1)$$

where $\hat{\beta}_O \equiv (G_{AB}Mn_b)^{1/3}/c$. $G_{AB}$ is the effective gravitational constant, entering the orbital dynamics of the binary system already at the Newtonian level. In GR, $\hat{\beta}_O$ is equal to $\beta_O$ of Eq. (11). In STG, $G_{AB} = G_*(1 + \alpha_A\alpha_B)$, where $G_*$ is the bare gravitational constant and $\alpha_a$ ($a$ = A, B) denotes the effective scalar coupling of body a. $\alpha_a$ is a measure for the change of the inertial mass of body a with respect to the (asymptotic) scalar field $\phi_0$, i.e., $\alpha_a = \partial \ln m_a / \partial \phi_0$, where the number of baryons of the NS is kept fixed when taking the partial derivative. The parameters $\gamma_{AB}$, $\beta_{BB}^A$, and $\beta_{AA}^B$ represent

three body-dependent strong-field generalizations of Eddington's two weak-field PPN parameters $\gamma^{PPN}$ and $\beta^{PPN}$ [204]. They enter the modified Einstein-Infeld-Hoffmann equations of motion for a binary system consisting of strongly self-gravitating bodies at the 1PN level [63,137]. In GR, as a result of the effacement principle [239], $G_{AB} = G$ and $\gamma_{AB} = \beta_{BB}^A = \beta_{AA}^B = 1$. In STG, they are given by [63]

$$\gamma_{AB} = 1 - \frac{2\alpha_A\alpha_B}{1 + \alpha_A\alpha_B}, \quad (C2)$$

$$\beta_{BB}^A = 1 + \frac{\beta_A \alpha_B^2}{2(1 + \alpha_A \alpha_B)^2}, \quad (C3)$$

$$\beta_{AA}^B = 1 + \frac{\beta_B \alpha_A^2}{2(1 + \alpha_A \alpha_B)^2}. \quad (C4)$$

The quantity $\beta_a$ ($a$ = A, B) is calculated according to $\beta_a = \partial\alpha_a/\partial\phi_0$, where the number of baryons of the NS is kept fixed when taking the partial derivative with respect to the asymptotic scalar field $\phi_0$. Depending on the coupling parameters $\alpha_0$ and $\beta_0$ of DEF gravity, $\beta_a$ can assume rather large values. More generally, in a system consisting of two NSs, these generalized PPN parameters can be very different from the corresponding weak-field PPN parameters. Even if there are only very small (or no) deviations from GR in the weak field, one can have order unity deviations in a double-NS system (see, e.g., Ref. [40]).

As discussed in Sec. VII A, although the $\mathcal{O}(\hat{\beta}_O^4)$ corrections in Eq. (C1) do become relevant in the analysis of the double pulsar, presently it is still sufficient to use the corresponding GR expressions in the gravity tests presented here.

### 2. Amplitude of the Einstein delay: $\gamma_E$

In STG, the amplitude $\gamma_E$ of the time dilation (Einstein delay) gets modified as well. Apart from changes in the expressions for the time dilation, there is a modulation of the spin period of the pulsar caused by a periodic variation of the local gravitational constant, at the location of the pulsar as it moves around its companion. This effect has the same orbital dependence as the time dilation and, hence, can be absorbed in $\gamma_E$. The resulting (total) expression for A is [42]

$$\gamma_E = \frac{e}{n_b} X_B \left( X_B + \frac{1 + k_A \alpha_B}{1 + \alpha_A \alpha_B} \right) \hat{\beta}_O^2, \quad (C5)$$

where $k_A \equiv -\partial \ln I_A / \partial \phi_0$. Here again, the number of baryons of A is kept fixed when taking the partial derivative. Note that $\alpha_A = \alpha_B = k_A = 0$ in GR.
### 3. Shapiro shape parameter: s

The Shapiro shape parameter s can simply be identified with $sin i$. When expressed as a function of the masses of the binary system, one uses the mass function, which is linked to the theory-dependent third Kepler law. To leading order, one finds in STG [42]

$$s = \frac{xn_b}{X_B \hat{\beta}_O}. \qquad (C6)$$

NLO corrections to the mass function depend on $\gamma_{AB}$ [cf. Eq. (3.9) in Ref. [35]], where $\epsilon = 2\gamma_{AB} + 1$. Since such NLO corrections are at the one-sigma level of the measured s, it is sufficient to use the corresponding GR expression in Eq. (19) in the PK-parameter test.

### 4. GW damping: $\dot{P}_b$

Concerning GW damping, the (dynamical) scalar field leads to various modifications of the expression for $\dot{P}_b$. In terms of PN order, the leading contribution comes from scalar dipole radiation and enters already at the 1.5PN level, i.e., $\hat{\beta}_O^3$. In STG, the leading-order term for dipolar GW damping is given by

$$\dot{P}_b^{\text{dipole}} = -2\pi X_A X_B \hat{\beta}_O^3 \frac{1 + \frac{1}{2}e^2}{(1 - e^2)^{5/2}} \frac{(\alpha_A - \alpha_B)^2}{1 + \alpha_A \alpha_B}. \qquad (C7)$$

In our test, we also use those NLO corrections which are of the order of $(\alpha_A - \alpha_B) \times \hat{\beta}_O^5$, even though $(\alpha_A - \alpha_B)$ has to be small, given the tight agreement with GR as seen, e.g., in Eq. (48). The rather lengthy expression for these NLO corrections can be found in Eq. (6.52b) in Ref. [203]. There are also terms of the order of $(\alpha_A - \alpha_B)^2 \times \hat{\beta}_O^5$, which can be safely ignored.

At the 2.5PN order, one finds a modified quadrupole formula that combines quadrupolar contributions from the tensor and the scalar field. For STG, it is given by [203]

$$\dot{P}_b^{\text{quadrupole}} = -\frac{192\pi}{5} X_A X_B \hat{\beta}_O^5 \frac{1 + \frac{73}{24}e^2 + \frac{37}{96}e^4}{(1 - e^2)^{7/2}} \times \frac{1 + \frac{1}{6}(\alpha_A X_B + \alpha_B X_A)^2}{1 + \alpha_A \alpha_B}. \qquad (C8)$$

In addition, the leading contribution from the scalar monopole radiation also enters at the 2.5PN order and in STG is given by [203]

$$\dot{P}_b^{\text{monopole}} = -3\pi X_A X_B \hat{\beta}_O^5 \frac{e^2(1 + \frac{1}{4}e^2)}{(1 - e^2)^{7/2}} \times \frac{[\alpha_A(1 + \frac{2}{3}X_B) + \alpha_B(1 + \frac{2}{3}X_A) + \frac{\beta_A \alpha_B + \beta_B \alpha_A}{1 + \alpha_A \alpha_B}]^2}{1 + \alpha_A \alpha_B}. \qquad (C9)$$

Note that in the double pulsar this contribution is greatly suppressed, compared to the other $O(\hat{\beta}_O^5)$ contributions, as it is proportional to $e^2 \approx 0.008$.

## APPENDIX D: GLOSSARY OF FREQUENTLY USED SYMBOLS AND ABBREVIATIONS

For the convenience of the reader, we provide a list of frequently used symbols and their meaning. We also list abbreviations frequently used in the paper.

### 1. Symbols

* $\alpha$: Right ascension (RA) in equatorial coordinates; spectral index

* $\delta$: declination (DEC) in equatorial coordinates

* $l$: Galactic longitude

* $b$: Galactic latitude

* $\mu_\alpha$: proper motion in RA

* $\mu_\delta$: proper motion in DEC

* $\mu_l$: proper motion in $l$

* $\mu_b$: proper motion in $b$

* $v_T$: transverse velocity

* $\pi$: annual parallax

* $\nu$: pulsar spin frequency

* $\dot{\nu}$: first spin frequency derivative

* $\ddot{\nu}$: second spin frequency derivative

* $\dddot{\nu}$: third spin frequency derivative

* $\dddot{\ddot{\nu}}$: fourth spin frequency derivative

* $\phi$: rotational phase

* $N_0$: pulse number at a reference epoch $t_0$

* $c_i$: profile evolution “FD parameters”

* $D$: Doppler factor between SSB and pulsar system

* $P_b$: orbital period

* $n_b$: orbital angular frequency, $2\pi/P_b$

* $x$: projected semimajor axis

* $e_T$: eccentricity (Kepler equation)

* $e_\theta$, $e_r$: additional eccentricities [see Eqs. (4) and (6)]

* $T_0$: epoch of periastron

* $\omega$, $\omega_0$: longitude of periastron, at $T_0$

* $\Delta_R$: Rømer delay

* $\Delta_E$: Einstein delay

* $\Delta_S$: signal propagation delay

* $\Delta_A$: aberration delay

* $u$: relativistic eccentric anomaly

* $\theta$: true anomaly

* $\beta_O$: dimensionless orbital velocity $(GMn_b)^{1/3}/c$

* $k$: periastron advance parameter

* $\dot{\omega}$: secular periastron advance ($\dot{\omega} \equiv n_b k$)

* $\dot{P}_b$: change of orbital period

* $\gamma_E$: Einstein delay amplitude

* $r$: range parameter of the Shapiro delay

* $s$: shape parameter of the Shapiro delay; scintillation scale; fractional scattering screen distance

* $t$: time

* $z_s$: logarithmic Shapiro shape, $z_s \equiv -\ln(1-s)$
$\Lambda_u$: argument of the Shapiro delay [see Eq. (24)]

$\lambda_{LT}$: scaling parameter of Lense-Thirring precession

$q_{NLO}$: NLO factor for signal propagation

$\delta_0$: relativistic deformation of orbit

$A, B$: aberration coefficients [see Eq. (28)]

$\dot{x}$: change of projected semimajor axis

$\dot{e}_r$: change of eccentricity

$\Omega_B^{\text{spin}}$: rate of relativistic spin precession of B

$i$: orbital inclination

$M$: total system mass

$m_A$: mass of pulsar A

$m_B$: mass of pulsar B

$X_{A,B}$: $m_{A,B}/M$

$R$: mass ratio $m_A/m_B$

$I_A$: MOI of pulsar A

$\alpha_0, \beta_0$: coupling strength parameters in STG

$T_1(\alpha_0, \beta_0)$: Damour–Esposito-Farèse STG

$\alpha_a$: effective scalar coupling of body $a$ ($a = A, B$)

$\beta_a$: derivative of $\alpha_a$ with respect to scalar field

$k_A$: derivative of $-\ln I_A$ with respect to scalar field

$\gamma^{\text{PPN}}$: spatial-curvature parameter of PPN formalism

$\beta^{\text{PPN}}$: nonlinearity parameter of PPN formalism

$\gamma_{ab}$: body-dependent strong-field generalization of $\gamma^{\text{PPN}}(a, b \in \{A, B, 0\})$

$\beta_{ab}^c$: body-dependent strong-field generalization of $\beta^{\text{PPN}}(a, b, c \in \{A, B\})$

## 2. Abbreviations

BH: black hole

DM: dispersion measure

EOS: equation of state

GLSQ: generalized least-squares timing method

GR: general relativity

GW: gravitational wave

ISM: interstellar medium

LOS: line of sight

LT: Lense-Thirring

MJD: modified Julian date

MOI: moment of inertia

NLO: next-to-leading-order

NS: neutron star

PK: post-Keplerian

PN: post-Newtonian

PPN: parametrized post-Newtonian

PSD: power spectral density

SEP: strong equivalence principle

SSB: Solar System barycenter

SF: structure function

S/N: signal-to-noise ratio

STG: scalar-tensor gravity

TOA: time of arrival

VLBA: Very Long Baseline Array

VLBI: Very Long Baseline Interferometry

[1] B. P. Abbott, R. Abbott, T. D. Abbott, M. R. Abernathy, F. Acernese, K. Ackley, C. Adams, T. Adams, P. Addesso, R. X. Adhikari, V. B. Adya, C. Affeldt, M. Agathos, K. Agatsuma, N. Aggarwal *et al.* (LIGO Scientific Collaboration and Virgo Collaboration), *Observation of Gravitational Waves from a Binary Black Hole Merger*, Phys. Rev. Lett. **116**, 061102 (2016).

[2] J. H. Taylor, L. A. Fowler, and P. M. McCulloch, *Measurements of General Relativistic Effects in the Binary Pulsar PSR 1913+16*, Nature (London) **277**, 437 (1979).

[3] J. H. Taylor and J. M. Weisberg, *A New Test of General Relativity: Gravitational Radiation and the Binary Pulsar PSR 1913+16*, Astrophys. J. **253**, 908 (1982).

[4] D. R. Lorimer and M. Kramer, *Handbook of Pulsar Astronomy* (Cambridge University Press, Cambridge, England, 2004), Vol. 4.

[5] A. Hewish, S. J. Bell, J. D. H. Pilkington, P. F. Scott, and R. A. Collins, *Observation of a Rapidly Pulsating Radio Source*, Nature (London) **217**, 709 (1968).

[6] T. Gold, *Rotating Neutron Stars as the Origin of the Pulsating Radio Sources*, Nature (London) **218**, 731 (1968).

[7] M. I. Large, A. E. Vaughan, and B. Y. Mills, *A Pulsar Supernova Association*, Nature (London) **220**, 340 (1968).

[8] One parsec (pc) is approximately $3.086 \times 10^{16}$ m, One kpc = 1000 pc.

[9] J. W. T. Hessels, S. M. Ransom, I. H. Stairs, P. C. C. Freire, V. M. Kaspi, and F. Camilo, *A Radio Pulsar Spinning at 716 Hz*, Science **311**, 1901 (2006).

[10] A. Lyne and F. Graham-Smith, *Pulsar Astronomy* (Cambridge University Press, Cambridge, England, 2012).

[11] N. Andersson, B. Haskell, G. L. Comer, and L. Samuelsson, *The Dynamics of Neutron Star Crusts: Lagrangian Perturbation Theory for a Relativistic Superfluid-Elastic System*, Classical Quantum Gravity **36**, 105004 (2019).

[12] P. B. Demorest, T. Pennucci, S. M. Ransom, M. S. E. Roberts, and J. W. T. Hessels, *A Two-Solar-Mass Neutron Star Measured Using Shapiro Delay*, Nature (London) **467**, 1081 (2010).

[13] J. Antoniadis *et al.*, *A Massive Pulsar in a Compact Relativistic Binary*, Science **340**, 448 (2013).

[14] H. T. Cromartie *et al.*, *Relativistic Shapiro Delay Measurements of an Extremely Massive Millisecond Pulsar*, Nat. Astron. **4**, 72 (2020).

[15] E. Fonseca *et al.*, *Refined Mass and Geometric Measurements of the High-Mass PSR J0740+6620*, Astrophys. J. Lett. **915**, L12 (2021).

[16] S. K. Greif, K. Hebeler, J. M. Lattimer, C. J. Pethick, and A. Schwenk, *Equation of State Constraints from Nuclear Physics, Neutron Star Masses, and Future Moment of Inertia Measurements*, Astrophys. J. **901**, 155 (2020).

[17] R. A. Hulse and J. H. Taylor, *Discovery of a Pulsar in a Binary System*, Astrophys. J. **195**, L51 (1975).

[18] L. L. Smarr and R. Blandford, *The Binary Pulsar: Physical Processes, Possible Companions and Evolutionary Histories*, Astrophys. J. **207**, 574 (1976).

[19] T. Damour and R. Ruffini, *Certain New Verifications of General Relativity Made Possible by the Discovery of a*
*Pulsar Belonging to a Binary System*, Acad. Sci. Paris C. R. Ser. Sci. Math. **279**, 971 (1974).

[20] B. M. Barker and R. F. O'Connell, *Relativistic Effects in the Binary Pulsar PSR 1913+16*, Astrophys. J. **199**, L25 (1975).

[21] M. Burgay, N. D'Amico, A. Possenti, R. N. Manchester, A. G. Lyne, B. C. Joshi, M. A. McLaughlin, M. Kramer, J. M. Sarkissian, F. Camilo, V. Kalogera, C. Kim, and D. R. Lorimer, *An Increased Estimate of the Merger Rate of Double Neutron Stars from Observations of a Highly Relativistic System*, Nature (London) **426**, 531 (2003).

[22] A. G. Lyne, M. Burgay, M. Kramer, A. Possenti, R. N. Manchester, F. Camilo, M. A. McLaughlin, D. R. Lorimer, N. D'Amico, B. C. Joshi, J. Reynolds, and P. C. C. Freire, *A Double-Pulsar System: A Rare Laboratory for Relativistic Gravity and Plasma Physics*, Science **303**, 1153 (2004).

[23] J. H. Taylor, *Pulsar Timing and Relativistic Gravity*, Phil. Trans. R. Soc. A **341**, 117 (1992).

[24] D. Nice, P. Demorest, I. Stairs, R. Manchester, J. Taylor, W. Peters, J. Weisberg, A. Irwin, N. Wex, and Y. Huang, *Tempo: Pulsar Timing Data Analysis*, 2015, ascl:1509.002, https://ui.adsabs.harvard.edu/abs/2015ascl.soft09002N.

[25] G. B. Hobbs, R. T. Edwards, and R. N. Manchester, *Tempo2, a New Pulsar-Timing Package—I. An Overview*, Mon. Not. R. Astron. Soc. **369**, 655 (2006); http://ascl.net/1210.015.

[26] A. Einstein, *Erklärung der Perihelbewegung des Merkur aus der Allgemeinen Relativitätstheorie*, Sitzungsber. K. Preuß. Akad. Wiss. XLVII, 831 (1915).

[27] N. Wex, in *Frontiers in Relativistic Celestial Mechanics: Applications and Experiments*, Vol. 2, edited by S. M. Kopeikin (Walter de Gruyter GmbH, Berlin, 2014), Chap. 2, p. 39.

[28] T. Damour and N. Deruelle, *General Relativistic Celestial Mechanics of Binary Systems. I. The Post-Newtonian Motion*, Ann. Inst. Henri Poincaré, A **43**, 107 (1985).

[29] T. Damour and N. Deruelle, *General Relativistic Celestial Mechanics of Binary Systems. II. The Post-Newtonian Timing Formula*, Ann. Inst. Henri Poincaré, A **44**, 263 (1986).

[30] J. H. Taylor and J. M. Weisberg, *Further Experimental Tests of Relativistic Gravity Using the Binary Pulsar PSR 1913+16*, Astrophys. J. **345**, 434 (1989).

[31] M. Kramer, I. H. Stairs, R. N. Manchester, M. A. McLaughlin, A. G. Lyne, R. D. Ferdman, M. Burgay, D. R. Lorimer, A. Possenti, N. D'Amico, J. M. Sarkissian, G. B. Hobbs, J. E. Reynolds, P. C. C. Freire, and F. Camilo, *Tests of General Relativity from Timing the Double Pulsar*, Science **314**, 97 (2006).

[32] M. Kramer and I. H. Stairs, *The Double Pulsar*, Annu. Rev. Astron. Astrophys. **46**, 541 (2008).

[33] M. Kramer and N. Wex, *TOPICAL REVIEW: The Double Pulsar System: A Unique Laboratory for Gravity*, Classical Quantum Gravity **26**, 073001 (2009).

[34] N. Wex and M. Kramer, *Gravity Tests with Radio Pulsars*, Universe **6**, 156 (2020).

[35] T. Damour and J. H. Taylor, *Strong-Field Tests of Relativistic Gravity and Binary Pulsars*, Phys. Rev. D **45**, 1840 (1992).

[36] T. Damour and G. Schäfer, *Higher-Order Relativistic Periastron Advances and Binary Pulsars*, Nuovo Cimento Soc. Ital. Fis. **101B**, 127 (1988).

[37] B. M. Barker and R. F. O'Connell, *Gravitational Two-Body Problem with Arbitrary Masses, Spins, and Quadrupole Moments*, Phys. Rev. D **12**, 329 (1975).

[38] A. T. Deller, J. M. Weisberg, D. J. Nice, and S. Chatterjee, *A VLBI Distance and Transverse Velocity for PSR B1913+16*, Astrophys. J. **862**, 139 (2018).

[39] M. J. Reid and M. Honma, *Microarcsecond Radio Astronomy*, Annu. Rev. Astron. Astrophys. **52**, 339 (2014).

[40] T. Damour and G. Esposito-Farèse, *Nonperturbative Strong-Field Effects in Tensor-Scalar Theories of Gravitation*, Phys. Rev. Lett. **70**, 2220 (1993).

[41] C. M. Will, *Theory and Experiment in Gravitational Physics*, second ed. (Cambridge University Press, Cambridge, England, 2018).

[42] T. Damour and G. Esposito-Farèse, *Tensor-Scalar Gravity and Binary-Pulsar Experiments*, Phys. Rev. D **54**, 1474 (1996).

[43] J. D. Bekenstein, *Relativistic Gravitation Theory for the Modified Newtonian Dynamics Paradigm*, Phys. Rev. D **70**, 083509 (2004).

[44] S. Boran, S. Desai, E. O. Kahya, and R. P. Woodard, *GW170817 Falsifies Dark Matter Emulators*, Phys. Rev. D **97**, 041501(R) (2018).

[45] J. M. Cordes and T. J. W. Lazio, *NE2001.I. A New Model for the Galactic Distribution of Free Electrons and Its Fluctuations*, arXiv:astro-ph/0207156.

[46] J. M. Yao, R. N. Manchester, and N. Wang, *A New Electron-Density Model for Estimation of Pulsar and FRB Distances*, Astrophys. J. **835**, 29 (2017).

[47] B. Willems, J. Kaplan, T. Fragos, V. Kalogera, and K. Belczynski, *Formation and Progenitor of PSR J0737-3039: New Constraints on the Supernova Explosion Forming Pulsar B*, Phys. Rev. D **74**, 043003 (2006).

[48] I. H. Stairs, S. E. Thorsett, R. J. Dewey, M. Kramer, and C. A. McPhee, *The Formation of the Double Pulsar PSR J0737-3039A/B*, Mon. Not. R. Astron. Soc. **373**, L50 (2006).

[49] T. M. Tauris, M. Kramer, P. C. C. Freire, N. Wex, H. T. Janka, N. Langer, P. Podsiadlowski, E. Bozzo, S. Chaty, M. U. Kruckow, E. P. J. van den Heuvel, J. Antoniadis, R. P. Breton, and D. J. Champion, *Formation of Double Neutron Star Systems*, Astrophys. J. **846**, 170 (2017).

[50] R. N. Manchester, M. Kramer, A. Possenti, A. G. Lyne, M. Burgay, I. H. Stairs, A. W. Hotan, M. A. McLaughlin, D. R. Lorimer, G. B. Hobbs, J. M. Sarkissian, N. D'Amico, F. Camilo, B. C. Joshi, and P. C. C. Freire, *The Mean Pulse Profile of PSR J0737-3039A*, Astrophys. J. **621**, L49 (2005).

[51] R. D. Ferdman, I. H. Stairs, M. Kramer, R. P. Breton, M. A. McLaughlin, P. C. C. Freire, A. Possenti, B. W. Stappers, V. M. Kaspi, R. N. Manchester, and A. G. Lyne, *The Double Pulsar: Evidence for Neutron Star Formation without an Iron Core-Collapse Supernova*, Astrophys. J. **767**, 85 (2013).

[52] N. Pol, M. McLaughlin, M. Kramer, I. Stairs, B. B. P. Perera, and A. Possenti, *A Direct Measurement of Sense of Rotation of PSR J0737-3039A*, Astrophys. J. **853**, 73 (2018).
[53] M. A. McLaughlin, A. G. Lyne, D. R. Lorimer, A. Possenti, R. N. Manchester, F. Camilo, I. H. Stairs, M. Kramer, M. Burgay, N. D'Amico, P. C. C. Freire, B. C. Joshi, and N. D. R. Bhat, *The Double Pulsar System J0737-3039: Modulation of A by B at Eclipse*, Astrophys. J. 616, L131 (2004).

[54] R. P. Breton, V. M. Kaspi, M. Kramer, M. A. McLaughlin, M. Lyutikov, S. M. Ransom, I. H. Stairs, R. D. Ferdman, F. Camilo, and A. Possenti, *Relativistic Spin Precession in the Double Pulsar*, Science 321, 104 (2008).

[55] Note that the angle $\theta$ in Ref. [54] is (to a good approximation) $180^\circ - \delta_B$, since in the coordinate system of Ref. [54] the orbital angular momentum is closely aligned with the negative $z$ direction.

[56] M. Burgay, A. Possenti, R. N. Manchester, M. Kramer, M. A. McLaughlin, D. R. Lorimer, I. H. Stairs, B. C. Joshi, A. G. Lyne, F. Camilo, N. D'Amico, P. C. C. Freire, J. M. Sarkissian, A. W. Hotan, and G. B. Hobbs, *Long-Term Variations in the Pulse Emission from PSR J0737-3039B*, Astrophys. J. 624, L113 (2005).

[57] A. Noutsos, G. Desvignes, M. Kramer, N. Wex, P. C. C. Freire, I. H. Stairs, M. A. McLaughlin, R. N. Manchester, A. Possenti, M. Burgay, A. G. Lyne, R. P. Breton, B. B. P. Perera, and R. D. Ferdman, *Understanding and Improving the Timing of PSR J0737-3039B*, Astron. Astrophys. 643, A143 (2020).

[58] B. B. P. Perera, M. A. McLaughlin, M. Kramer, I. H. Stairs, R. D. Ferdman, P. C. C. Freire, A. Possenti, R. P. Breton, R. N. Manchester, M. Burgay, A. G. Lyne, and F. Camilo, *The Evolution of PSR J0737-3039B and a Model for Relativistic Spin Precession*, Astrophys. J. 721, 1193 (2010).

[59] M. Lyutikov and C. Thompson, *Magnetospheric Eclipses in the Double Pulsar System PSR J0737-3039*, Astrophys. J. 634, 1223 (2005).

[60] R. P. Breton, V. M. Kaspi, M. A. McLaughlin, M. Lyutikov, M. Kramer, I. H. Stairs, S. M. Ransom, R. D. Ferdman, F. Camilo, and A. Possenti, *The Double Pulsar Eclipses. I. Phenomenology and Multi-frequency Analysis*, Astrophys. J. 747, 89 (2012).

[61] M. Lyutikov, *On the Nature of Eclipses in Binary Pulsar J0737-3039*, Mon. Not. R. Astron. Soc. 353, 1095 (2004).

[62] At least up to 1PN order for boost-invariant gravity theories [63].

[63] T. Damour, in *Physics of Relativistic Objects in Compact Binaries: From Birth to Coalescence*, Astrophysics and Space Science Library Vol. 359, edited by M. Colpi, P. Casella, V. Gorini, U. Moschella, and A. Possenti (Springer, New York, 2009), pp. 1-41.

[64] C. Granet, H. Z. Zhang, A. R. Forsyth, G. R. Graves, P. Doherty, K. J. Greene, G. L. James, P. Sykes, T. S. Bird, M. W. Sinclair, G. Moorey, and R. N. Manchester, *The Designing, Manufacturing, and Testing of a Dual-Band Feed System for the Parkes Radio Telescope*, IEEE Trans. Antennas Propag. 47, 13 (2005).

[65] L. Staveley-Smith, W. E. Wilson, T. S. Bird, M. J. Disney, R. D. Ekers, K. C. Freeman, R. F. Haynes, M. W. Sinclair, R. A. Vaile, R. L. Webster, and A. E. Wright, *The Parkes 21 cm Multibeam Receiver*, Pub. Astron. Soc. Aust. 13, 243 (1996).

[66] G. Hobbs et al., *An Ultra-wide Bandwidth (704 to 4032 MHz) Receiver for the Parkes Radio Telescope*, Pub. Astron. Soc. Aust. 37, e012 (2020).

[67] A. Hotan, *High-Precision Observations of Relativistic Binary and Millisecond Pulsars*, Ph.D. thesis, Swinburne University of Technology, 2007.

[68] W. van Straten and M. Bailes, *DSPSR: Digital Signal Processing Software for Pulsar Astronomy*, Pub. Astron. Soc. Aust. 28, 1 (2011).

[69] R. N. Manchester et al., *The Parkes Pulsar Timing Array Project*, Pub. Astron. Soc. Aust. 30, e017 (2013).

[70] P. B. Demorest, *Measuring the Gravitational Wave Background Using Precision Pulsar Timing*, Ph.D. thesis, University of California, Berkeley, 2007.

[71] R. DuPlain, S. Ransom, P. Demorest, P. Brandt, J. Ford, and A. L. Shelton, in *Advanced Software and Control for Astronomy II*, SPIE Proceedings Vol. 7019 (SPIE-International Society for Optical Engineering, Bellingham, WA, 2008), p. 70191D.

[72] R. D. Ferdman, *Binary Pulsars: Evolution and Fundamental Physics*, Ph.D. thesis, University of British Columbia, 2008.

[73] http://ascl.net/1509.003.

[74] W. van Straten, P. Demorest, and S. Oslowski, *Pulsar Data Analysis with PSRCHIVE*, Astron. Res. Technol. 9, 237 (2012).

[75] V. R. Grandy, *Wideband Timing of the Double Pulsar (J0737-3039a)*, Bachelor's thesis, University of British Columbia, 2008.

[76] 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).

[77] P. Lazarus, R. Karuppusamy, E. Graikou, R. N. Caballero, D. J. Champion, K. J. Lee, J. P. W. Verbiest, and M. Kramer, *Prospects for High-Precision Pulsar Timing with the New Effelsberg PSRIX Backend*, Mon. Not. R. Astron. Soc. 458, 868 (2016).

[78] C. G. Bassa, G. H. Janssen, R. Karuppusamy, M. Kramer, K. J. Lee, K. Liu, J. McKee, D. Perrodin, M. Purver, S. Sanidas, R. Smits, and B. W. Stappers, *LEAP: The Large European Array for Pulsars*, Mon. Not. R. Astron. Soc. 456, 2196 (2016).

[79] R. Karuppusamy, B. Stappers, and W. van Straten, *Pulsar II: A Wide Band Pulsar Machine for the Westerbork Synthesis Radio Telescope*, Publ. Astron. Soc. Pac. 120, 191 (2008).

[80] G. Desvignes et al., *High-Precision Timing of 42 Millisecond Pulsars with the European Pulsar Timing Array*, Mon. Not. R. Astron. Soc. 458, 3341 (2016).

[81] R. M. Shannon, S. Oslowski, S. Dai, M. Bailes, G. Hobbs, R. N. Manchester, W. van Straten, C. A. Raithel, V. Ravi, L. Toomey, N. D. R. Bhat, S. Burke-Spolaor, W. A. Coles, M. J. Keith, M. Kerr, Y. Levin, J. M. Sarkissian, J.-B. Wang, L. Wen, and X.-J. Zhu, *Limitations in Timing Precision due to Single-Pulse Shape Variability in Millisecond Pulsars*, Mon. Not. R. Astron. Soc. 443, 1463 (2014).

[82] W. van Straten, *Radio Astronomical Polarimetry and Point-Source Calibration*, Astrophys. J. Suppl. Ser. 152, 129 (2004).
[83] A. W. Hotan, M. Bailes, and S. M. Ord, Geodetic Precession in PSR J1141-6545, Astrophys. J. **624**, 906 (2005).

[84] P. B. Demorest et al., Limits on the Stochastic Gravitational Wave Background from the North American Nano-hertz Observatory for Gravitational Waves, Astrophys. J. **762**, 94 (2013).

[85] https://naif.jpl.nasa.gov/pub/naif/JUNO/kernels/spk/de436s.bsp.lbl.

[86] M. Vallisneri et al., Modeling the Uncertainties of Solar System Ephemerides for Robust Gravitational-Wave Searches with Pulsar-Timing Arrays, Astrophys. J. **893**, 112 (2020).

[87] R. N. Caballero et al., Studying the Solar System with the International Pulsar Timing Array, Mon. Not. R. Astron. Soc. **481**, 5501 (2018).

[88] M. J. Keith, W. Coles, R. M. Shannon, G. B. Hobbs, R. N. Manchester, M. Bailes, N. D. R. Bhat, S. Burke-Spolaor, D. J. Champion, A. Chaudhary, A. W. Hotan, J. Khoo, J. Kocz, S. Osowski, V. Ravi, J. E. Reynolds, J. Sarkissian, W. van Straten, and D. R. B. Yardley, Measurement and Correction of Variations in Interstellar Dispersion in High-Precision Pulsar Timing, Mon. Not. R. Astron. Soc. **429**, 2161 (2013).

[89] W. Coles, G. Hobbs, D. J. Champion, R. N. Manchester, and J. P. W. Verbiest, Pulsar Timing Analysis in the Presence of Correlated Noise, Mon. Not. R. Astron. Soc. **418**, 561 (2011).

[90] We note that we have also explored alternative methods, e.g., using the “DMX” method implemented in Tempo, in order to determine the astrometry and its uncertainties. Using simulated timing data with parallaxes of known sizes injected into our datasets, we determined that the method implemented here appears to be the most reliable in recovering the correct astrometry with appropriate uncertainties.

[91] C. E. Rasmussen and C. K. I. Williams, Gaussian Processes for Machine Learning (MIT, Cambridge, MA, 2006).

[92] B. J. Rickett, W. A. Coles, C. F. Nava, M. A. McLaughlin, S. M. Ransom, F. Camilo, R. D. Ferdman, P. C. C. Freire, M. Kramer, A. G. Lyne, and I. H. Stairs, Interstellar Scintillation of the Double Pulsar J0737-3039, Astrophys. J. **787**, 161 (2014).

[93] A. T. Deller, M. Bailes, and S. J. Tingay, Implications of a VLBI Distance to the Double Pulsar J0737-3039A/B, Science **323**, 1327 (2009).

[94] A. T. Deller, S. J. Vigeland, D. L. Kaplan, W. M. Goss, W. F. Brisken, S. Chatterjee, J. M. Cordes, G. H. Janssen, T. J. W. Lazio, L. Petrov, B. W. Stappers, and A. Lyne, Microarcsecond VLBI Pulsar Astrometry with PSRb. I. Two Binary Millisecond Pulsars with White Dwarf Companions, Astrophys. J. **828**, 8 (2016).

[95] R. Blandford and S. A. Teukolsky, Arrival-Time Analysis for a Pulsar in a Binary System, Astrophys. J. **205**, 580 (1976).

[96] Depending on the literature, this is called the argument of periastron. We follow the naming that is more common in pulsar astronomy and, for instance, also used in Ref. [97] (see their Fig. 3.2).

[97] E. Poisson and C. M. Will, Gravity (Cambridge University Press, Cambridge, England, 2014).

[98] R. T. Edwards, G. B. Hobbs, and R. N. Manchester, Tempo2, A New Pulsar Timing Package—II. The Timing

Model and Precision Estimates, Mon. Not. R. Astron. Soc. **372**, 1549 (2006).

[99] G. Schäfer and N. Wex, Second Post-Newtonian Motion of Compact Binaries, Phys. Lett. A **174**, 196 (1993).

[100] G. Schäfer and N. Wex, Second Post-Newtonian Motion of Compact Binaries, Phys. Lett. A **177**, 461 (1993).

[101] N. Wex, The Second Post-Newtonian Motion of Compact Binary-Star Systems with Spin, Classical Quantum Gravity **12**, 983 (1995).

[102] T. Damour and G. Schaefer, The Two-Body Problem in General Relativity, Acad. Sci. Paris C. R. Ser. B **305**, 839 (1987).

[103] The rotational spin frequency parameters are estimated with a standard Tempo analysis. As discussed in Ref. [89], this procedure underestimates the true uncertainty in the presence of unmodeled red timing noise. Comparison of the standard analysis with the results of a full generalized least-squares analysis shows that the true uncertainty in the spin frequency parameters is about an order of magnitude larger than the values quoted in Table IV. We emphasize that the unmodeled red noise has no effect on the orbital parameters as expected because of their much shorter timescale.

[104] Mass determination from PK parameters yields the product $Gm_i$ ($i = A, B$). Consequently, whenever we express the masses in terms of solar masses $M_⊙$, the numerical value is calculated via $Gm_i/(GM)_⊙^N$, where $(GM)_⊙^N = 1.3271244 \times 10^{26} \text{ cm}^3 \text{ s}^{-2}$ is the nominal solar mass parameter as defined by the IAU 2015 Resolution B3 [105].

[105] E. E. Mamajek et al., IAU 2015 Resolution B3 on Recommended Nominal Conversion Constants for Selected Solar and Planetary Properties, arXiv:1510.07674.

[106] T. Damour, P. Jaranowski, and G. Schäfer, Dynamical Invariants for General Relativistic Two-Body Systems at the Third Post-Newtonian Approximation, Phys. Rev. D **62**, 044024 (2000).

[107] R.-M. Memmesheimer, A. Gopakumar, and G. Schäfer, Third Post-Newtonian Accurate Generalized Quasi-Keplerian Parametrization for Compact Binaries in Eccentric Orbits, Phys. Rev. D **70**, 104011 (2004).

[108] C. Königsdörffer and A. Gopakumar, Post-Newtonian Accurate Parametric Solution to the Dynamics of Spinning Compact Binaries in Eccentric Orbits: The Leading Order Spin-Orbit Interaction, Phys. Rev. D **71**, 024039 (2005).

[109] L. Iorio, Prospects for Measuring the Moment of Inertia of Pulsar J0737-3039A, New Astron. **14**, 40 (2009).

[110] 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. R. Astron. Soc. **497**, 3118 (2020).

[111] For the precision needed, we can safely use $g_{S_A}^‖$ instead of $g_{S_A}$ in Ref. [36], since $\delta_A < 3.2^\circ$, while $g_{S_A}^‖ \cos \delta_A$ for the double pulsar (cf. Ref. [33]).

[112] T. Dietrich, M. W. Coughlin, P. T. H. Pang, M. Bulla, J. Heinzel, L. Issa, I. Tews, and S. Antier, Multimessenger Constraints on the Neutron-Star Equation of State and the Hubble Constant, Science **370**, 1450 (2020).

[113] Although the mass of A (1.34 $M_⊙$) is slightly below 1.4 $M_⊙$, within the accuracy needed here, the radius constraints in Ref. [112] for a 1.4 - $M_⊙$ NS can also be applied to pulsar A.
[114] J. M. Lattimer, *Neutron Star Mass and Radius Measurements*, Universe **5**, 159 (2019).

[115] There are other distributions that have been derived for the MOI of pulsar A; see, e.g., Refs. [116–118]. However, for the results of this paper, this does not make any difference.

[116] Y. Lim, J. W. Holt, and R. J. Stahulak, *Predicting the Moment of Inertia of Pulsar J0737–3039A from Bayesian Modeling of the Nuclear Equation of State*, Phys. Rev. C **100**, 035802 (2019).

[117] P. Landry and B. Kumar, *Constraints on the Moment of Inertia of PSR J0737–3039A from GW170817*, Astrophys. J. **868**, L22 (2018).

[118] H. O. Silva, A. M. Holgado, A. Cárdenas-Avendaño, and N. Yunes, *Astrophysical and Theoretical Physics ImPLICATIONS FROM MULTIMESSAGER NEUTRON STAR OBSERVATIONS*, Phys. Rev. Lett. **126**, 181101 (2021).

[119] J. M. Lattimer and B. F. Schutz, *Constraining the Equation of State with Moment of Inertia Measurements*, Astrophys. J. **629**, 979 (2005).

[120] Z. Arzoumanian, K. Joshi, F. Rasio, and S. E. Thorsett, in *Pulsars: Problems and Progress, IAU Colloquium 160*, edited by S. Johnston, M. A. Walker, and M. Bailes (Astronomical Society of the Pacific, San Francisco, 1996), pp. 525–530.

[121] S. M. Kopeikin, *Proper Motion of Binary Pulsars as a Source of Secular Variations of Orbital Parameters*, Astrophys. J. Lett. **467**, L93 (1996).

[122] T. Damour and J. H. Taylor, *On the Orbital Period Change of the Binary Pulsar PSR 1913+16*, Astrophys. J. **366**, 501 (1991).

[123] L. Blanchet, *Gravitational Radiation from Post-Newtonian Sources and Inspiralling Compact Binaries*, Living Rev. Relativity **17**, 2 (2014).

[124] P. C. Peters, *Gravitational Radiation and the Motion of Two Point Masses*, Phys. Rev. **136**, B1224 (1964).

[125] B. R. Iyer and C. M. Will, *Post-Newtonian Gravitational Radiation Reaction for Two-Body Systems: Nonspinning Bodies*, Phys. Rev. D **52**, 6882 (1995).

[126] M. E. Pati and C. M. Will, *Post-Newtonian Gravitational Radiation and Equations of Motion via Direct Integration of the Relaxed Einstein Equations. II. Two-Body Equations of Motion to Second Post-Newtonian Order, and Radiation Reaction to 3.5 Post-Newtonian Order*, Phys. Rev. D **65**, 104008 (2002).

[127] C. Königsdörffer, G. Faye, and G. Schäfer, *Binary Black-Hole Dynamics at the Third-and-a-Half Post-Newtonian Order in the ADM Formalism*, Phys. Rev. D **68**, 044004 (2003).

[128] L. Blanchet and G. Schäfer, *Higher Order Gravitational Radiation Losses in Binary Systems*, Mon. Not. R. Astron. Soc. **239**, 845 (1989).

[129] I. S. Shklovskii, *Possible Causes of the Secular Increase in Pulsar Periods*, Sov. Astron. **13**, 562 (1970).

[130] P. J. McMillan, *The Mass Distribution and Gravitational Potential of the Milky Way*, Mon. Not. R. Astron. Soc. **465**, 76 (2017).

[131] R. Abuter, A. Amorim, M. Bauböck, J. P. Berger, H. Bonnet, W. Brandner, Y. Clénet, V. C. Du Foresto, P. T. de Zeeuw, J. Dexter, G. Duvert et al. (Gravity Collaboration), *A Geometric Distance Measurement to the Galactic Center Black Hole with 0.3% Uncertainty*, Astron. Astrophys. **625**, L10 (2019).

[132] D. M. Eardley, Astrophys. J. **196**, L59 (1975).

[133] C. M. Will and H. W. Zaglauer, *Gravitational Radiation, Close Binary Systems, and the Brans-Dicke Theory of Gravity*, Astrophys. J. **346**, 366 (1989).

[134] N. Wex and R. Rieth, in *Proceedings of the 2nd Gauss Symposium, Conference A: Mathematics and Theoretical Physics, Munich, Germany, 1993*, edited by M. Behara, R. Fritsch, and R. G. Lintz (Walter de Gruyter, Berlin, 1995), pp. 681–694.

[135] See http://www.skatelescope.org.

[136] I. I. Shapiro, Fourth Test of General Relativity, Phys. Rev. Lett. **13**, 789 (1964).

[137] C. M. Will, *Theory and Experiment in Gravitational Physics* (Cambridge University Press, Cambridge, England, 1993).

[138] J. Schneider, *Gravitational Beam Deviation Effects in the Timing Formula of Binary Pulsars*, Astron. Astrophys. **232**, 62 (1990).

[139] D. Lai and R. R. Rafikov, *Effects of Gravitational Lensing in the Double Pulsar System J0737–3039*, Astrophys. J. **621**, L41 (2005).

[140] Simulations are based on the integration of null geodesics in a Schwarzschild spacetime, comoving with pulsar B. Since additive constant terms can be ignored, the use of Eddington-Finkelstein coordinates turns out to be of great advantage in the numerical integration for the propagation time [141].

[141] M. Ali, *Relativistic Propagation Effects in Binary Pulsar Signals for Next Generation Radio Telescopes*, Master's thesis, Friedrich-Wilhelms-Universität Bonn, Bonn, Germany, 2011.

[142] V. M. Kaspi, S. M. Ransom, D. C. Backer, R. Ramachandran, P. Demorest, J. Arons, and A. Spitkovsky, *Green Bank Telescope Observations of the Eclipse of Pulsar “A” in the Double Pulsar Binary PSR J0737–3039*, Astrophys. J. **613**, L137 (2004).

[143] The (intermittent) eclipses of A by the rotating magnetosphere of B lead to a significant reduction in the number of pulses available for the integrated pulse profile. This, in turn, reduces its S/N, leading to significantly fewer viable TOAs, which at the same time have considerably larger errors.

[144] S. A. Klioner and S. Zschocke, *Numerical versus Analytical Accuracy of the Formulae for Light Propagation, Classical Quantum Gravity* **27**, 075015 (2010).

[145] V. Perlick, *Gravitational Lensing from a Spacetime Perspective*, Living Rev. Relativity **7**, 9 (2004).

[146] S. M. Kopeikin and G. Schäfer, *Lorentz Covariant Theory of Light Propagation in Gravitational Fields of Arbitrary-Moving Bodies*, Phys. Rev. D **60**, 124002 (1999).

[147] R. R. Rafikov and D. Lai, *Effects of Gravitational Lensing and Companion Motion on the Binary Pulsar Timing*, Phys. Rev. D **73**, 063003 (2006).

[148] Lensing corrections to the retardation effect as discussed in Ref. [147] are totally negligible. The corresponding postfit residuals are below 0.02 µs, restricted to a narrow range around conjunction where, because of the eclipse, the timing precision is anyway greatly compromised.
[149] M. Kramer, I. H. Stairs, R. N. Manchester, M. A. McLaughlin, A. G. Lyne, R. D. Ferdman, M. Burgay, D. R. Lorimer, A. Possenti, N. D'Amico, J. Sarkissian, B. C. Joshi, P. C. C. Freire, and F. Camilo, *Strong-Field Tests of Gravity with the Double Pulsar*, Ann. Phys. (Berlin) **15**, 34 (2006).

[150] O. V. Doroshenko and S. M. Kopeikin, *Relativistic Effect of Gravitational Deflection of Light in Binary Pulsars*, Mon. Not. R. Astron. Soc. **274**, 1029 (1995).

[151] R. R. Rafikov and D. Lai, *Effects of Pulsar Rotation on Timing Measurements of the Double Pulsar System J0737-3039*, Astrophys. J. **641**, 438 (2006).

[152] Note that the expressions in Refs. [147,151] are based on the lensing equation, which makes the assumption of small angles (cf., e.g., Ref. [153]), and, therefore, they are strictly speaking valid only near conjunction, which also limits their usability in a Tempo implementation.

[153] P. Schneider, J. Ehlers, and E. E. Falco, *Gravitational Lenses* (Springer, New York, 1992).

[154] M. Kramer et al., *The Relativistic Binary Programme on MeerKAT: Science Objectives and First Results*, Mon. Not. R. Astron. Soc. **504**, 2094 (2021).

[155] Periodic variations in the pulsar profile due to a latitudinal shift in the LOS caused by aberration are observed in PSR B1534+12 [156].

[156] I. H. Stairs, S. E. Thorsett, and Z. Arzoumanian, *Measurement of Gravitational Spin-Orbit Coupling in a Binary-Pulsar System*, Phys. Rev. Lett. **93**, 141101 (2004).

[157] The DDS timing model is a modification of the DD model [29], where the Shapiro shape parameter *s* is replaced by the logarithmic Shapiro shape parameter *z*<sub>s</sub> via *s*≡1−exp(−*z*<sub>s</sub>) (see Ref. [149] for details). It is part of the Tempo standard distribution.

[158] W. W. Zhu, I. H. Stairs, P. B. Demorest, D. J. Nice, J. A. Ellis, S. M. Ransom, Z. Arzoumanian, K. Crowter, T. Dolch, R. D. Ferdman, E. Fonseca, M. E. Gonzalez, G. Jones, M. L. Jones, M. T. Lam, L. Levin, M. A. McLaughlin, T. Pennucci, K. Stovall, and J. Swiggum, *Testing Theories of Gravitation Using 21-Year Timing of Pulsar Binary J1713+0747*, Astrophys. J. **809**, 41 (2015).

[159] Z. Arzoumanian, A. Brazier, S. Burke-Spolaor, S. Chamberlin, S. Chatterjee, B. Christy, J. M. Cordes et al. (NANOGrav Collaboration), *The NANOGrav Nine-Year Data Set: Observations, Arrival Time Measurements, and Analysis of 37 Millisecond Pulsars*, Astrophys. J. **813**, 65 (2015).

[160] K. Liu, G. Desvignes, I. Cognard, B. W. Stappers, J. P. W. Verbiest, K. J. Lee, D. J. Champion, M. Kramer, P. C. C. Freire, and R. Karuppusamy, *Measuring Pulse Times of Arrival from Broad-band Pulsar Observations*, Mon. Not. R. Astron. Soc. **443**, 3752 (2014).

[161] T. T. Pennucci, *Frequency-Dependent Template Profiles for High-Precision Pulsar Timing*, Astrophys. J. **871**, 34 (2019).

[162] N. Capitaine, in *Journées Systèmes de Référence Spatio-temporels 2007*, edited by N. Capitaine (2008), p. 61, https://ui.adsabs.harvard.edu/abs/2008jsrs.meet...61C.

[163] Even though the factor $(1 - \kappa)$ is small, such correction is required before using the listed timing results for predictions of the folding parameters with Tempo2, which implements TCB. According to IAU Resolution 2006 B3,

the conversion takes place as $\text{TDB} = \text{TCB} - (1 - \kappa) \times (\text{JD}_{\text{TCB}} - 2443144.5003725) \times 86400 - 6.55 \times 10^{-5}$, where $\text{JD}_{\text{TCB}}$ is the relevant Julian date [164]. Note that this described conversion needs to be applied as given to the quoted epoch of periastron passage ($T_0$) value before it is used with Tempo2.

[164] See https://www.iau.org/static/resolutions/IAU2006_Resol3.pdf.

[165] G. Hobbs, A. G. Lyne, and M. Kramer, *An Analysis of the Timing Irregularities for 366 Pulsars*, Mon. Not. R. Astron. Soc. **402**, 1027 (2010).

[166] R. P. Breton, *Radio Pulsars in Binary Systems*, Ph.D. thesis, McGill University, 2009.

[167] J. M. Lattimer and M. Prakash, *Neutron Star Structure and the Equation of State*, Astrophys. J. **550**, 426 (2001).

[168] B. P. Abbott, R. Abbott, T. D. Abbott, F. Acernese, K. Ackley, C. Adams, T. Adams, P. Addesso, R. X. Adhikari, V. B. Adya et al. (LIGO Scientific Collaboration and Virgo Collaboration), *GW170817: Observation of Gravitational Waves from a Binary Neutron Star Inspiral*, Phys. Rev. Lett. **119**, 161101 (2017).

[169] J. M. Lattimer, in *Proceedings of the Xiamen-CUSTIPEN Workshop on the Equation of State of Dense Neutron-Rich Matter in the Era of Gravitational Wave Astronomy*, American Institute of Physics Conference Series Vol. 2127 (AIP, New York, 2019) p. 020001.

[170] This Doppler-shift rescaling, $m^{\text{measured}} = m^{\text{intrinsic}}/D$, between the pulsar-timing measured masses and the intrinsic ones [29] is due to the same effect which causes the masses directly deduced from gravitational-wave-template fitting in LIGO/Virgo data analysis to differ from the intrinsic, source-frame masses by a redshift factor: $m^{\text{obs}} = (1+z)m^{\text{source}}$ [171].

[171] A. Krolak and B. F. Schutz, *Coalescing Binaries—Probe of the Universe*, Gen. Relativ. Gravit. **19**, 1163 (1987).

[172] J. M. Weisberg and Y. Huang, *Relativistic Measurements from Timing the Binary Pulsar PSR B1913+16*, Astrophys. J. **829**, 55 (2016).

[173] D. J. Nice and J. H. Taylor, PSRs J2019+2425 and J2322+2057 and the Proper Motions of Millisecond Pulsars, Astrophys. J. **441**, 429 (1995).

[174] K. Lazaridis, N. Wex, A. Jessner, M. Kramer, B. W. Stappers, G. H. Janssen, G. Desvignes, M. B. Purver, I. Cognard, G. Theureau, A. G. Lyne, C. A. Jordan, and J. A. Census, *Generic Tests of the Existence of the Gravitational Dipole Radiation and the Variation of the Gravitational Constant*, Mon. Not. R. Astron. Soc. **400**, 805 (2009).

[175] M. J. Reid and A. Brunthaler, *The Proper Motion of Sagittarius A*\*. III. *The Case for a Supermassive Black Hole*, Astrophys. J. **892**, 39 (2020).

[176] R. Abuter et al. (Gravity Collaboration), *Improved GRAVITY Astrometric Accuracy from Modeling Optical Aberrations*, Astron. Astrophys. **647**, A59 (2021).

[177] M. S. Kehl, N. Wex, M. Kramer, and K. Liu, in *Proceedings of the Fourteenth Marcel Grossmann Meeting—MG14*, edited by M. Bianchi, R. T. Jansen, and R. Ruffini (World Scientific, Singapore, 2018), pp. 1860–1865.

[178] B. P. Abbott, R. Abbott, T. D. Abbott, M. R. Abernathy, F. Acernese, K. Ackley, C. Adams, T. Adams, P. Addesso et al. (LIGO Scientific and Virgo Collaborations), *Tests of*
*General Relativity with GW150914*, Phys. Rev. Lett. **116**, 221101 (2016).

[179] We take the combined results of these references, obtained from the IMRPhenomPv2 waveform family.

[180] B. P. Abbott, R. Abbott, T. D. Abbott, S. Abraham, F. Acernese, K. Ackley, C. Adams, R. X. Adhikari, V. B. Adya *et al.* (LIGO Scientific Collaboration and Virgo Collaboration), *Tests of General Relativity with the Binary Black Hole Signals from the LIGO-Virgo Catalog GWTC-1*, Phys. Rev. D **100**, 104036 (2019).

[181] R. Abbott, T. D. Abbott, S. Abraham, F. Acernese, K. Ackley, A. Adams, C. Adams, R. X. Adhikari, V. B. Adya, C. Affeldt *et al.* (LIGO Scientific Collaboration and Virgo Collaboration), *Tests of General Relativity with Binary Black Holes from the Second LIGO-Virgo Gravitational-Wave Transient Catalog*, Phys. Rev. D **103**, 122002 (2021).

[182] S. Khan, S. Husa, M. Hannam, F. Ohme, M. Pürrer, X. J. Forteza, and A. Bohé, *Frequency-Domain Gravitational Waves from Nonprecessing Black-Hole Binaries. II. A Phenomenological Model for the Advanced Detector Era*, Phys. Rev. D **93**, 044007 (2016).

[183] S. Mirshekari and C. M. Will, *Compact Binary Systems in Scalar-Tensor Gravity: Equations of Motion to 2.5 Post-Newtonian Order*, Phys. Rev. D **87**, 084070 (2013).

[184] K. Yagi, L. C. Stein, and N. Yunes, *Challenging the Presence of Scalar Charge and Dipolar Radiation in Binary Pulsars*, Phys. Rev. D **93**, 024010 (2016).

[185] E. Barausse, C. Palenzuela, M. Ponce, and L. Lehner, *Neutron-Star Mergers in Scalar-Tensor Theories of Gravity*, Phys. Rev. D **87**, 081506(R) (2013).

[186] J. Alsing, E. Berti, C. M. Will, and H. Zaglauer, *Gravitational Radiation from Compact Binary Systems in the Massive Brans-Dicke Theory of Gravity*, Phys. Rev. D **85**, 064041 (2012).

[187] Short-range modifications of GR that significantly affect the MOI could, in principle, show up in the Lense-Thirring test of Eq. (50); see discussion in Ref. [110].

[188] I. Ciufolini, A. Paolozzi, E. C. Pavlis, G. Sindoni, J. Ries, R. Matzner, R. Koenig, C. Paris, V. Gurzadyan, and R. Penrose, *An Improved Test of the General Relativistic Effect of Frame-Dragging Using the LARES and LAGEOS Satellites*, Eur. Phys. J. C **79**, 872 (2019).

[189] C. W. F. Everitt *et al.*, *Gravity Probe B: Final Results of a Space Experiment to Test General Relativity*, Phys. Rev. Lett. **106**, 221101 (2011).

[190] A. Ritz and J. Ward, *Weyl Corrections to Holographic Conductivity*, Phys. Rev. D **79**, 066003 (2009).

[191] This is, for instance, the case in the mEIH formalism, where (to leading order) the deflection angle is proportional to $(1 + \gamma_{B0})G_{B0}$ (see e.g., Appendix C in Ref. [192]), which is already tightly constrained by the test of the Shapiro range parameter $r$. More generally, the NLO contribution to $\Delta_A$ is parametrized by $r$ and $s$ according to Eqs. (30) and (31).

[192] T. Damour and G. Esposito-Farèse, *Testing Gravity to Second Post-Newtonian Order: A Field-Theory Approach*, Phys. Rev. D **53**, 5541 (1996).

[193] It follows from the analysis in Ref. [51] that the rotational axis of A is (practically) perpendicular to the orbital plane. See also Sec. VIII C.

[194] An additional test is provided by the PK parameter $\Omega_B^{\text{spin}}$; see Table V.

[195] B. P. Abbott, R. Abbott, T. D. Abbott, F. Acernese, K. Ackley, C. Adams, T. Adams, P. Addesso, Adhikari *et al.* (LIGO Scientific Collaboration and Virgo Collaboration), *Tests of General Relativity with GW170817*, Phys. Rev. Lett. **123**, 011102 (2019).

[196] K. Nordtvedt, *G/G and a Cosmological Acceleration of Gravitationally Compact Bodies*, Phys. Rev. Lett. **65**, 953 (1990).

[197] T. Damour, G. W. Gibbons, and J. H. Taylor, *Limits on the Variability of G Using Binary-Pulsar Data*, Phys. Rev. Lett. **61**, 1151 (1988).

[198] K. Nordtvedt, *Cosmological Time Evolution in the Newtonian Dynamics of Compact Celestial Bodies*, Astrophys. J. **407**, 5 (1993).

[199] There are special cases, like Barker's constant-G theory [200], where $\dot{P}_b^G$ cannot be written in this form.

[200] B. M. Barker, *General Scalar-Tensor Theory of Gravity with Constant g*, Astrophys. J. **219**, 5 (1978).

[201] W. W. Zhu *et al.*, *Tests of Gravitational Symmetries with Pulsar Binary J1713+0747*, Mon. Not. R. Astron. Soc. **482**, 3249 (2019).

[202] A. Genova, E. Mazarico, S. Goossens, F. G. Lemoine, G. A. Neumann, D. E. Smith, and M. T. Zuber, *Solar System Expansion and Strong Equivalence Principle as Seen by the NASA MESSENGER Mission*, Nat. Commun. **9**, 289 (2018).

[203] T. Damour and G. Esposito-Farèse, *Tensor-Multi-scalar Theories of Gravitation*, Classical Quantum Gravity **9**, 2093 (1992).

[204] A. S. Eddington, *The Mathematical Theory of Relativity* (Cambridge University Press, Cambridge, England, 1922).

[205] To leading order, one finds in DEF gravity $B_D \simeq (5/96) \times (\alpha_A - \alpha_B)^2 [1 + \mathcal{O}(\alpha_A, \alpha_B)][1 + \mathcal{O}(e^2)]$ (cf. Refs. [203,206]).

[206] T. Damour and G. Esposito-Farèse, *Gravitational-Wave versus Binary-Pulsar Tests of Strong-Field Gravity*, Phys. Rev. D **58**, 042001 (1998).

[207] B. Bertotti, L. Iess, and P. Tortora, *A Test of General Relativity Using Radio Links with the Cassini Spacecraft*, Nature (London) **425**, 374 (2003).

[208] 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. R. Astron. Soc. **423**, 3328 (2012).

[209] 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 (London) **559**, 73 (2018).

[210] J.-P. Bruneton and G. Esposito-Farèse, *Field-Theoretical Formulations of MOND-like Gravity*, Phys. Rev. D **76**, 124012 (2007).

[211] J. P. W. Verbiest, D. R. Lorimer, and M. A. McLaughlin, *Lutz-Kelker Bias in Pulsar Parallax Measurements*, Mon. Not. R. Astron. Soc. **405**, 564 (2010).

[212] J. P. W. Verbiest, J. M. Weisberg, A. A. Chael, K. J. Lee, and D. R. Lorimer, *On Pulsar Distance Measurements and Their Uncertainties*, Astrophys. J. **755**, 39 (2012).
[213] J. Yao (private communication).

[214] C. R. Purcell, B. M. Gaensler, X. H. Sun, E. Carretti, G. Bernardi, M. Haverkorn, M. J. Kesteven, S. Poppi, D. H. F. M. Schnitzeler, and L. Staveley-Smith, *A Radio-Polarisation and Rotation Measure Study of the Gum Nebula and Its Environment*, Astrophys. J. **804**, 22 (2015).

[215] A. T. Deller, W. M. Goss, W. F. Brisken, S. Chatterjee, J. M. Cordes, G. H. Janssen, Y. Y. Kovalev, T. J. W. Lazio, L. Petrov, B. W. Stappers, and A. Lyne, *Microarcsecond VLBI Pulsar Astrometry with PSR$\pi$ II. Parallax Distances for 57 Pulsars*, Astrophys. J. **875**, 100 (2019).

[216] M. T. Lam, J. M. Cordes, S. Chatterjee, and T. Dolch, *Pulsar Timing Errors from Asynchronous Multi-frequency Sampling of Dispersion Measure Variations*, Astrophys. J. **801**, 130 (2015).

[217] X. P. You, G. Hobbs, W. A. Coles, R. N. Manchester, R. Edwards, M. Bailes, J. Sarkissian, J. P. W. Verbiest, W. van Straten, A. Hotan, S. Ord, F. Jenet, N. D. R. Bhat, and A. Teoh, *Dispersion Measure Variations and Their Effect on Precision Pulsar Timing*, Mon. Not. R. Astron. Soc. **378**, 493 (2007).

[218] Since the DM is the column integral of the electron density along the line-of-sight (z) axis, its 2D spatial spectrum is the same as that of the 3D density with the wave number in z set to zero. Therefore, the 2D spectral exponent of DM (β) is the same as that of the 3D electron-density fluctuations, i.e., -11/3 for a Kolmogorov process. The corresponding 2D structure function has a spectral exponent α = -(β + 2) = +5/3 for Kolmogorov.

[219] One astronomical unit (A.U.) is approximately 1.496 × 10⁸ km.

[220] The fluctuations in DM(t) are those of a 1D cut across the 2D spatial distribution. Therefore, the 1D PSD of DM(t) is the 2D spatial spectrum of DM (scaled by the velocity) with the perpendicular wave number integrated out. This reduction to 1D reduces the Kolmogorov spectral exponent from -11/3 to -8/3.

[221] O. P. Lay, *The Temporal Power Spectrum of Atmospheric Fluctuations due to Water Vapor*, Astron. Astrophys. Suppl. Ser. **122**, 535 (1997).

[222] This value agrees well with a pulsar-based estimate [223]. There is some spread in measured values for $z_{SSB}$ (see, e.g., Table 1 in Ref. [223]), but this has little impact on our conclusions.

[223] J. M. Yao, R. N. Manchester, and N. Wang, *Determination of the Sun's Offset from the Galactic Plane Using Pulsars*, Mon. Not. R. Astron. Soc. **468**, 3289 (2017).

[224] R. N. Manchester, M. Kramer, I. H. Stairs, M. Burgay, F. Camilo, G. B. Hobbs, D. R. Lorimer, A. G. Lyne, M. A. McLaughlin, C. A. McPhee, A. Possenti, J. E. Reynolds, and W. van Straten, *Observations and Modeling of Relativistic Spin Precession in PSR J1141-6545*, Astrophys. J. **710**, 1694 (2010).

[225] V. V. Krishnan, M. Bailes, W. van Straten, E. F. Keane, M. Kramer, N. D. R. Bhat, C. Flynn, and S. Osłowski, *Relativistic Spin Precession in the Binary PSR J1141-6545*, Astrophys. J. **873**, L15 (2019).

[226] M. Burgay, A. Possenti, R. N. Manchester, M. Kramer, M. A. McLaughlin, D. R. Lorimer, I. H. Stairs, B. C. Joshi,

A. G. Lyne, F. Camilo, N. D'Amico, P. C. C. Freire, J. M. Sarkissian, A. W. Hotan, and G. B. Hobbs, *Long-Term Variations in the Pulse Emission from PSR J0737-3039B*, Astrophys. J. **624**, L113 (2005).

[227] K. Akiyama, A. Alberdi, W. Alef, K. Asada, R. Azulay, A.-K. Baczko, D. Ball, M. Baloković, J. Barrett, D. Bintley et al. (Event Horizon Telescope Collaboration), *First M87 Event Horizon Telescope Results. I. The Shadow of the Supermassive Black Hole*, Astrophys. J. **875**, L1 (2019).

[228] R. Abuter, A. Amorim, M. Bauböck, J. P. Berger, H. Bonnet, W. Brandner, V. Cardoso, Y. Clénet, P. T. de Zeeuw, J. Dexter et al. (Gravity Collaboration), *Detection of the Schwarzschild Precession in the Orbit of the Star S2 near the Galactic Centre Massive Black Hole*, Astron. Astrophys. **636**, L5 (2020).

[229] R. Abuter, A. Amorim, M. Bauböck, J. P. Berger, H. Bonnet, W. Brandner, Y. Clénet, V. C. Du Foresto, P. T. de Zeeuw, C. Deen, J. Dexter, G. Duvert et al. (Gravity Collaboration), *Detection of Orbital Motions near the Last Stable Circular Orbit of the Massive Black Hole SgrA**, Astron. Astrophys. **618**, L10 (2018).

[230] N. Yunes, K. Yagi, and F. Pretorius, *Theoretical Physics Implications of the Binary Black-Hole Mergers GW150914 and GW151226*, Phys. Rev. D **94**, 084002 (2016).

[231] T. Baker, D. Psaltis, and C. Skordis, *Linking Tests of Gravity on All Scales: From the Strong-Field Regime to Cosmology*, Astrophys. J. **802**, 63 (2015).

[232] M. Kramer, *Probing Einstein's Universe and Its Physics*, Astron. Geophys. **58**, 3.31(2017).

[233] R. Abuter, A. Amorim, N. Anugu, M. Bauböck, M. Benisty, J. P. Berger, N. Blind, H. Bonnet, W. Brandner, A. Buron, C. Collin, F. Chapron, Y. Clénet, V. C. Du Foresto, P. T. de Zeeuw, C. Deen, F. Delplancke-Ströbele, R. Dembet, J. Dexter, G. Duvert et al. (Gravity Collaboration), *Detection of the Gravitational Redshift in the Orbit of the Star S2 near the Galactic Centre Massive Black Hole*, Astron. Astrophys. **615**, L15 (2018).

[234] T. Do et al., *Relativistic Redshift of the Star S0-2 Orbiting the Galactic Center Supermassive Black Hole*, Science **365**, 664 (2019).

[235] D. Psaltis, L. Medeiros, P. Christian, F. Ozel, K. Akiyama, A. Alberdi, W. Alef, K. Asada, R. Azulay, D. Ball et al. (The EHT Collaboration), *Gravitational Test beyond the First Post-Newtonian Order with the Shadow of the M87 Black Hole*, Phys. Rev. Lett. **125**, 141104 (2020).

[236] J. J. Hermes, M. Kilic, W. R. Brown, D. E. Winget, C. A. Prieto, A. Gianninas, A. S. Mukadam, A. Cabrera-Lavers, and S. J. Kenyon, *Rapid Orbital Decay in the 12.75-Minute Binary White Dwarf J0651+2844*, Astrophys. J. **757**, L21 (2012).

[237] B. J. Rickett, *Radio Propagation through the Turbulent Interstellar Plasma*, Annu. Rev. Astron. Astrophys. **28**, 561 (1990).

[238] W. A. Coles, B. J. Rickett, J. J. Gao, G. Hobbs, and J. P. W. Verbiest, *Scattering of Pulsar Radio Emission by the Interstellar Plasma*, Astrophys. J. **717**, 1206 (2010).

[239] T. Damour, in *Three Hundred Years of Gravitation*, edited by S. Hawking and W. Israel (Cambridge University Press, Cambridge, England, 1987), pp. 128-198.