Searching for gravitational waves via Doppler tracking by future missions to Uranus and Neptune

The past year has seen numerous publications underlining the importance of a space mission to the ice giants in the upcoming decade. Proposed mission plans involve a ∼ 10 yr cruise time to the ice giants. This cruise time can be utilized to search for low-frequency gravitational waves (GWs) by observing the Doppler shift caused by them in the Earth–spacecraft radio link. We calculate the sensitivity of prospective ice giant missions to GWs. Then, adopting a steady-state black hole binary population, we derive a conservative estimate for the detection rate of extreme mass ratio inspirals (EMRIs), supermassive black hole (SMBH), and stellar mass binary black hole (sBBH) mergers. We link the SMBH population to the fraction of quasars f bin resulting from Galaxy mergers that pair SMBHs to a binary. For a total of 10 40-d observations during the cruise of a single spacecraft, O ( f bin ) ∼ 0 . 5 detections of SMBH mergers are likely, if Allan deviation of Cassini-era noise is improved by ∼ 10 2 in the 10 − 5 − 10 − 3 Hz range. For EMRIs the number of detections lies between O (0 . 1) and O (100). Furthermore, ice giant missions combined with the Laser Interferometer Space Antenna (LISA) would improve the localization by an order of magnitude compared to LISA by itself.


INTRODUCTION
Uranus and Neptune are the outermost planets in our solar system, orbiting at roughly 20 and 30 AU from the Sun, respectively.Not surprisingly, they are the least explored planets in the system, having been visited only once by the Voyager II spacecraft in the late 1980s.The past year has seen numerous white papers calling for missions to the ice giants (Simon et al. 2020a;Rymer et al. 2019;Beddingfield et al. 2020;Cartwright et al. 2020;Dahl et al. 2020).The scientific potential of possible missions, along with various mission designs have also been extensively discussed in Hofstadter et al. (2019); Fletcher et al. (2020a,b); Helled & Fortney (2020); Jarmak et al. (2020); Simon et al. (2020b); Kollmann et al. (2020).A consensus is being reached in the community regarding the launch date of a possible mission that would maximize the payload and have a high science yield.The time frame is reported to be around 2029-2030 for Neptune and early 2030s for Uranus, especially if a Jupiter Gravity Assist (JGA) is used to reach the ice giants (Hofstadter et al. 2019).
Considering that they will spend most of their time in interplanetary space (rather than orbiting the planets they are destined for), the science potential of such mission configurations is limited to a fraction of their lifetime.However, if the transponders on the spacecraft were used to detect gravitational waves (GWs) via Doppler tracking during the cruise phase, these missions would also offer a unique opportunity to look for low-frequency GWs.
GWs passing through the space between the transponder and the transmitter/receiver at Earth cause variations in the light travel time ★ E-mail: deniz.soyuer@uzh.chbetween the two, which correspond to a Doppler shift in the frequency of the transmitted/received signal Δ/ 0 , where  0 is the signal carrier frequency.Analysis of the resulting Doppler shift allows to reconstruct the strain due to GWs between Earth and the spacecraft, making the Earth-satellite system an arm of a GW observatory.We expand on the details later and note that this process is elegantly described in Armstrong (2006) and the references therein.
This use of Doppler tracking with interplanetary spacecraft was previously suggested for many space missions; most notably in the Pioneer 11 data analysis Armstrong et al. (1987), the Galileo-Ulysses-Mars Observer coincidence experiment (Anderson et al. 1992;Bertotti et al. 1992), and the Cassini (Comoretto et al. 1992;Bertotti et al. 1999) mission.However, there were no significant GW candidates due to insufficient signal levels (Armstrong 2006).Nevertheless, with increasing capabilities to combat detection noise, we suggest Doppler tracking could be a cheap and efficient way to do science during the cruise phase of prospective ice giant missions.

MISSION PLAN
There are many proposed designs for a possible ice giant mission.A notable one consists of a spacecraft separating into two just before engaging in a JGA, after which each payload travels toward their destined planets (see Fig. 1).The mission timeline is projected as: We approximate a trajectory for both spacecraft with above times-tamps, assuming they leave Jupiter's sphere of influence soon after the JGA.Meaning, their trajectories evolve under solar gravity with the initial asymptotic velocities acquired after the JGA.Fig. 1 shows orbital positions of the planets and spacecraft trajectories after the JGA, the Earth-spacecraft distance, and the angle subtended by both spacecraft from the Earth.Although, knowing the exact trajectory is crucial for the actual measurement, it is not important for our sensitivity curve estimation, which will remain qualitatively unaffected by the details of the mission.Both spacecraft enter near-conjunction  > 150 • for ∼2 months at a time, for a total of 9 times for the Uranus spacecraft and 11 times for the Neptune one.As noted in Armstrong (2006), the elongation  is directly related to plasma scintillation noise due to solar winds and irregularities in Earth's ionosphere, where angles  < 150 • are not ideal for observations.For the purposes of our SNR calculations, we take a single Earthspacecraft system that can collect data for ten 40-day observations, evenly spread within 10 years cruising time between 8 and 30 AU.

GW response of a Doppler tracking system
Following Armstrong (2006), we define the fractional frequency fluctuation of a two-way Doppler system, with monochromatic carrier frequency  0 , as  2 () = Δ/ 0 , where  2 is the two-way light time between Earth and the spacecraft, and Δ = ( −  2 ) − ().Then, the frequency fluctuation due to a passing GW can be expressed as where  = k • n is the projection of the unit wavevector k of the wave onto the unit vector connecting the Earth and the spacecraft n, and Ψ is the projection of the GW amplitude onto the Doppler link (Wahlquist 1987).The latter is given by Ψ where h() = h + () e + + h × () e × is the GW amplitude (i.e. the strain) and e +,× are the usual "plus" and the "cross" polarization states of a transverse, traceless plane GW.Compact binary systems are the most promising sources of GWs, as they produce sine waves with the form where  orb is the angular frequency of the binary system and  the monochromatic instantaneous amplitude of the passing GW (Moore et al. 2015).The spectral power response  GW  2 of this fluctuation is a measure of the noiseless sensitivity of a Doppler system to a GW of a particular frequency, and can be calculated analytically for a sinusoidal source by taking the Fourier transform of the fluctuation ỹ2 = F [ GW  2 ] and reading the coefficient of the Dirac delta function Fig. 2 shows  GW  2 of a source located in two extreme directions and the sky-averaged response, for two different light travel times.

Sensitivity of an ice giant mission
The total signal produced by a monochromatic GW source does not only depend on its instantaneous gravitational luminosity, but also on the number of GW cycles passing through the detector.The frequency and amplitude of the GW evolve as the binary shrinks over many orbits.The characteristic strain takes into account the cycles that a binary completes in the proximity of a given GW frequency  For the purposes of this letter it is sufficient to use phenomenological waveform amplitudes like the ones in Ajith et al. (2007), which correctly model the frequency scaling of the inspiral, merger and ringdown phases.The characteristic strain amplitudes read where M = ( 1  2 ) 3/5 /( tot ) 1/5 is the chirp mass and  L the luminosity distance of the binary. and L model the decay of the so called "quasi-normal modes"; oscillations of the shape of the event horizon just after merger, given by L () = ( 0 / 1 ) 2/3 (4( −  1 ) 2 / 2 2 − 1) −1 .The frequencies  0 ,  1 ,  2 and  3 correspond to the merger, ringdown, ringdown decay-width and cut-off frequencies, respectively.Their values are computed as   =  3 (   2 +    +   )/(  tot ), where  = (M/ tot ) 5/3 is the symmetric mass ratio and   ,   and   are conveniently tabulated in Robson et al. (2019).
It is important to note that eq. ( 4) is only valid if the source can complete all of the cycles at a given frequency within one observation time  obs .The limiting frequency for this condition to be true reads For GWs at frequencies  ≤  lim , the number of cycles in a frequency bin is bounded by the observation time.In this case we have The last step required to produce a sensitivity curve for the Doppler system is to equate the spectral power response of the signal  GW  2 and that of the noise   , and subsequently solve for the required instantaneous strain amplitude that leads to an SNR of 1: where   is for 40 days observation time, same as that calculated by Bertotti et al. (1999).If the characteristic strain of a source ℎ  is equal to ℎ  , then that source would have an accumulated SNR = 1.A compilation of noise sources of the Cassini-era observations can be seen in table 2 of Armstrong (2006).We follow Bertotti et al. (1999) to model   () as a triple power-law, corresponding to three frequency regimes that are dominated by different noise sources where  1 <  2 <  3 , and  0 is related to the Allan deviation   via For Cassini, the powers read  1 = −2,  2 = −1/2 and  3 = 2, where the last two model the frequency propagation noise and the onset of thermal noise, respectively.The steep scaling of the low-end does not have a physical origin, but rather it is a conservative estimate due to the lack of sophisticated data analysis at that regime.
Fig. 3 shows the sky-averaged sensitivities of various experiments, as well as those of the prospective ice giant missions.We focus on three cases, where the total Allan deviation is improved by a factor of 3, 30 and 100 w.r.t.Cassini-era values.While the noise scaling at low and high-end is uncertain, the cutoff frequencies at  min = 2/ obs ∼10 −6.2 Hz and  max = 2/ res ∼10Hz are set by observation and resolution time, respectively.Hence, sensitivity of an ice giant Doppler tracking GW detector spans from the high-end of currently operating Pulsar Timing Arrays (PTAs Burke-Spolaor et al. 2019), through the LISA band (Amaro-Seoane et al. 2017), reaching the low-end of ground based detectors (e.g., aLIGO Abbott et al. 2016).

EVENT RATE ESTIMATES
Expected sources of GWs in this range are the late inspiral and mergers of SBHBs (Klein et al. 2016a) at the low-end, to the EMRIs (Amaro-Seoane et al. 2007) at mid-frequencies, to pre-merger inspirals of sBBHs that will merge weeks to years later in the LIGO band (Sesana 2016;Gerosa et al. 2019).Fig. 3 shows the characteristic strains of various merger events with generic masses and distances.
To estimate event rates for each type of source, we consider binaries with mass ratio  =  2 / 1 ≤ 1, and zero orbital eccentricity such that they emit GWs at twice the orbital frequency  = 2 orb .We assume that the number of binaries per orbital frequency is governed by a steady-state continuity equation with orbital frequency evolution dominated by GW emission,  G .The first assumption is valid when the binaries are formed at  orb lower than the detection frequencies (Christian & Loeb 2017), and the second for binaries at late inspiral, both well motivated here.Thus, the number density of binaries is where dR/d is the volumetric merger rate per total binary mass for a specified population, and where we have further assumed that the binary mass and mass ratios are fixed over the observation time.
The total number of these binaries that a detector with sensitivity   may detect above a threshold SNR   , is found by integrating over volume  and including only sources above the detection threshold: with where d/d is the angle integrated cosmological volume element in a flat universe (Hogg 1999) and H is the Heaviside function.We have assumed the population can be modeled with a single representative mass ratio , which is a variable of the model.The detection probability D is built from the ℎ  of eq.s (4) and ( 7) and the the characteristic noise, ℎ  , explained surrounding eq. ( 8).We note that D undercounts the SNR of a source as it compares ℎ  and ℎ  at a single frequency rather than integrating over all source frequencies.

Stellar mass binary black hole mergers
For sBBHs that will eventually merge in the LIGO band, we can reliably tie the merger rate to the measured LIGO merger rate.As only the most massive sBBHs are expected to be observable at low frequencies and sensitivities of interest here (Gerosa et al. 2019) we estimate the differential merger rate as a fraction of the LIGO inferred merger rate  R LIGO ( − * ).For binary parameters representative of the massive sBBHs detected by LIGO,  * = 60 and  = 1 (Abbott et al. 2019), eq. (11) shows that even for a 10 2 improvement in   , we would expect 0.02 detectable binaries with   = 1.This is consistent with assuming D (, , , ) results in a cut of sources beyond a maximum distance  max and minimum frequency  min : For comparison, with LISA SNR cuts of   = 2, 5, 8, one expects to detect 200, 14, and 3 sBBHs respectively, in line with recent estimates (Sesana 2016;Gerosa et al. 2019).SMBH Merger @ z = 5,  2015)).Orange, blue and pink curves represent sensitivities of the ice giants mission, with an Allan deviation improvement of 3, 30 and 100 times compared to the Cassini-era, respectively.We sample 10, 40-day observations between  2 = 8000s and 30000s.The red curves correspond to various equal mass GW sources.The green and black curves are the  = 2 (solid) and  = 4 (dashed) harmonics of a eccentric EMRIs at  = 0.3, 0.1.Bottom panels: Number of detectable SBHBs per log mass, redshift, and observed GW frequency.Dashed lines represent all SBHBs predicted by our model, while the solid lines represent the number of SBHBs detectable with the specified detectors  SBHB via eq.( 11), with a SNR threshold of   = 1 (partially justified due to undercounting by D).A luminosity distance of 600 Mpc is denoted in the redshift plot for clarity.

Extreme mass ratio inspirals
We follow a similar procedure for EMRIs, where instead of a measured merger rate, we rely on a range of rates predicted in the literature for various EMRI production channels (see Chen & Han 2018).This varies from ∼ 10 −9 − 10 −6 /gal/yr.Assuming a characteristic EMRI consisting of a 30 BH inspiraling to a 10 6  SMBH, we find (see Fig. 3), such EMRIs fall above the intermediate (blue) sensitivity curve for  0.1.The total number of galaxies within  = 0.1 is approximately 10 −2 gal/Mpc 3 × 4(420Mpc) 3 /3 = 3 × 10 6 , so we find that the expected number of detectable EMRIs ranges from 3 × 10 −3 to 3 over the lifetime of the proposed mission.Considering the highest sensitivity (pink) mission, the expected number of EMRI events approaches ∼ 0.08 to 80 over the mission lifetime.

Supermassive black hole binary mergers
To compute a differential merger rate for SBHBs, we assume a plausible merger channel, namely that a fraction  bin of quasars result from galaxy mergers that pair SMBHs into a binary at the nucleus of the newly formed galaxy (Haiman et al. 2009;D'Orazio & Loeb 2019).In other words, if  bin = 1, then all quasars counted in the observationally determined quasar luminosty function (QLF) have a binary at some stage in its evolution, (albeit not necessarily at the stage we want for the frequency to be high enough to detect).The number of SBHBs is traced by the quasar population which is itself traced by the observationally determined QLF: d 2   /(dd where we relate the quasar bolometric luminosity  to , assuming the binary accretes at some fraction of the Eddington rate  Edd , with  Edd = 4   /  the Eddington luminosity.We assume an average value of  Edd = 0.1 for bright quasars (Shankar et al. 2013) that are traced by the QLF from (Hopkins et al. 2007).We then need to know the distribution of these SBHBs in emitted GW frequency and characteristic strain.The number of binaries per frequency in a steady-state population, assuming circular orbits, can be estimated by assuming binaries are driven together by GW emission.Under this assumption, the fraction of binaries at frequency  is approximately the residence time of the binary /  GW ∝  −8/3 (Sesana et al. 2005;Christian & Loeb 2017;D'Orazio & Di Stefano 2020) divided by the binary lifetime.As fiducial parameters we choose a quasar lifetime of   = 10 7 yr (Martini 2004), and  = 0.3, typical of major galactic mergers (Volonteri et al. 2003).Using eqs.( 15) and ( 16) in eq. ( 11), and integrating yields estimates for the detectable SBHB population scalable in terms of the  bin and the   .Bottom panels of Fig. 3 show our estimates for the corresponding ice giant Doppler ranging missions and for a nominal LISA mission, for  bin = 1.We see that our most optimistic detector could observe of order a few SBHB inspiral/mergers at redshifts of ∼ 0.1 − 1.0 and frequencies of 10 −4.5 − 10 −3.5 Hz.Interestingly, our model predicts that such a detection is most likely for an inspiralling, lighter SBHB at  ≤ 10 6  .This is because the quasar (and hence SMBH) population is largest at these smaller masses.While our model relies on the uncertain means by which SBHBs are brought together and merge, we note that by including only a merger channel that traces the quasars, we tie our estimate to a rele- vant observable quantity while also conservatively underestimating high redshift mergers.Our estimates are in approximate agreement with the low redshift predictions from more general SBHB population estimates (Wyithe & Loeb 2003;Klein et al. 2016b).

DISCUSSION AND CONCLUSION
In this letter, we have calculated the feasibility of detecting GWs from BH mergers with a prospective ice giants mission by constructing a Doppler tracking sensitivity curve and implementing a population synthesis model to estimate the merger detection rate.As shown in Fig. 4, the Allan deviation,   , is a crucial factor in determining whether the mission will be successful as a GW observatory.Within our conservative merger population model, we find that an improvement in   of at least 2 orders of magnitude over the Cassini parameters is required to have a reasonable chance of detecting a few GW events by SBHBs, and a few to tens of EMRIs, depending on the accuracy of the merger rate estimates.Their is also some uncertainty in  bin , some estimates put it around  bin ∼ 0.25 (Charisi et al. 2016).
It is hard to predict how much the total Allan deviation is likely to improve over the years, since it depends on several different technologies (as discussed thoroughly in Armstrong 2006).However, a few qualitative arguments suggest that a factor 10 2 might not be improbable.A hypothetical ice giant mission would be scheduled for the 2030s, meaning three decades of technological development from the Cassini-era.Moreover, the expertise gained by the recent success of LIGO and the LISA pathfinder mission is likely to have significant crossover to a Doppler tracking system.
In the mission scenario detailed in Section 2, two spacecraft will form a ∼90 • angle for most of the cruise time (see Fig. 1).While our calculations are for a single spacecraft, adding a second independent tracking system would allow to reduce the noise via cross-correlation of the signals.Furthermore, if the signal analysis issues with low frequencies mentioned in Bertotti et al. (1999) were resolved, we would expect significantly more detections around 10 −5 − 10 −4 Hz, where the detection rate per frequency bin peaks.With a significant improvement in the total Allan deviation, a Doppler tracking experiment might become as capable as LISA at such low frequencies, and bridge the gap between mHz detectors and PTAs.
It is a major advantage that an ice giant mission would be concurrent with LISA, since the two experiments would be likely to detect the same signal with completely independent systems, reducing systematic noise and increasing the precision of parameter estimation.This, combined with the cost efficiency of an ice giant Doppler tracking experiment, demonstrates that a future ice giant mission would be a great opportunity for low-frequency GW astronomy.

•
Feb. 2031: Space Launch System (SLS) departure from Earth.• Dec. 2032: Separation of the spacecraft and subsequent JGA.• Apr.2042: Arrival of the first spacecraft at Uranus.• Sep.2044: Arrival of the second spacecraft at Neptune.

Figure 1 .
Figure 1.Top panel: Orbital positions of Earth (red), Jupiter (orange), Uranus (cyan) and Neptune (blue) centered around the Sun, plotted until 01/Jan/2050.Markers represent different timestamps, where full ones indicate the presence of a spacecraft.Location data have been acquired from the JPL HORIZONS System using the Astroquery tool Ginsburg et al. (2019).Mission details given in: https://github.com/ice-giants/papers/raw/master/presentation/IGs2020_missiondesign_elliott.pdf.Center panel: Angle between the Uranus spacecraft (1), Earth (E) and the Neptune spacecraft (2).Relevant timestamps of the mission are shown with vertical lines.Bottom panel: The distances between spacecraft and the Earth.

Figure 3 .
Figure 3. Top panel: Sensitivity of various experiments expressed in log ℎ  vs. log  (Ulysses and Cassini with data fromBertotti et al. (1999), LISA curve fromRobson et al. (2019), aLIGO and the SKA-era PTA fromMoore et al. (2015)).Orange, blue and pink curves represent sensitivities of the ice giants mission, with an Allan deviation improvement of 3, 30 and 100 times compared to the Cassini-era, respectively.We sample 10, 40-day observations between  2 = 8000s and 30000s.The red curves correspond to various equal mass GW sources.The green and black curves are the  = 2 (solid) and  = 4 (dashed) harmonics of a eccentric EMRIs at  = 0.3, 0.1.Bottom panels: Number of detectable SBHBs per log mass, redshift, and observed GW frequency.Dashed lines represent all SBHBs predicted by our model, while the solid lines represent the number of SBHBs detectable with the specified detectors  SBHB via eq.(11), with a SNR threshold of   = 1 (partially justified due to undercounting by D).A luminosity distance of 600 Mpc is denoted in the redshift plot for clarity.

Figure 4 .
Figure 4. Expected detections of SMBH binary mergers vs. Allan deviation for 400 days of total observation.Black curves represent SNRs of 1,3 and 5. Vertical colored lines are the ice giant missions with Allan deviation improvements of 3, 30 and 100 with respect to the Cassini-era measurements.
).Assuming this fraction of the quasars facilitates a SBHB over the quasar lifetime   , the differential merger rate for SBHBs reads