# Implementation of space-division multiplexed entanglement-based quantum
cryptography over multicore fiber

Evelyn A. Ortega<EMAIL_ADDRESS>Institute for Quantum Optics and
Quantum Information - IQOQI Vienna, Austrian Academy of Sciences,
Boltzmanngasse 3, 1090 Vienna, Austria  Vienna Center for Quantum Science and
Technology (VCQ), Vienna, Austria  Current address: ICFO-Institut de Ciencies
Fotoniques, The Barcelona Institute of Science and Technology, 08860
Castelldefels (Barcelona), Spain  Jorge Fuenzalida  jorge.fuenzalida@tu-
darmstadt.de Technical University of Darmstadt, Department of Physics,
Institute for Applied Physics, Schloßgartenstraße 7, 64289 Darmstadt, Germany
Krishna Dovzhik  Institute for Quantum Optics and Quantum Information - IQOQI
Vienna, Austrian Academy of Sciences, Boltzmanngasse 3, 1090 Vienna, Austria
Vienna Center for Quantum Science and Technology (VCQ), Vienna, Austria
Current address: Department of Neuroscience and Biomedical Engineering, Aalto
School of Science, Aalto University, Espoo, Finland  Rodrigo F. Shiozaki
Departamento de Física, Universidade Federal de São Carlos, Rodovia Washington
Luís, km 235—SP-310, 13565-905 São Carlos, SP, Brazil  Juan Carlos Alvarado-
Zacarias  CREOL, The University of Central Florida, Orlando, Florida 32816,
USA  Rodrigo Amezcua-Correa  CREOL, The University of Central Florida,
Orlando, Florida 32816, USA  Martin Bohmann  Institute for Quantum Optics and
Quantum Information - IQOQI Vienna, Austrian Academy of Sciences,
Boltzmanngasse 3, 1090 Vienna, Austria  Vienna Center for Quantum Science and
Technology (VCQ), Vienna, Austria  Current address: Quantum Technology
Laboratories GmbH, Clemens-Holzmeister-Straße 6/6, 1100 Vienna, Austria  Sören
Wengerowsky  Institute for Quantum Optics and Quantum Information - IQOQI
Vienna, Austrian Academy of Sciences, Boltzmanngasse 3, 1090 Vienna, Austria
Vienna Center for Quantum Science and Technology (VCQ), Vienna, Austria
Current address: ICFO-Institut de Ciencies Fotoniques, The Barcelona Institute
of Science and Technology, 08860 Castelldefels (Barcelona), Spain  Rupert
Ursin  Institute for Quantum Optics and Quantum Information - IQOQI Vienna,
Austrian Academy of Sciences, Boltzmanngasse 3, 1090 Vienna, Austria  Vienna
Center for Quantum Science and Technology (VCQ), Vienna, Austria  Current
address: Quantum Technology Laboratories GmbH, Clemens-Holzmeister-Straße 6/6,
1100 Vienna, Austria

###### Abstract

Quantum communication implementations require efficient and reliable quantum
channels. Optical fibers have proven to be an ideal candidate for distributing
quantum states. Thus, today’s efforts address overcoming issues towards high
data transmission and long-distance implementations. Here, we experimentally
demonstrate the secret key rate enhancement via space-division multiplexing
using a multicore fiber. Our multiplexing technique exploits the momentum
correlation of photon pairs generated by spontaneous parametric down-
conversion. We distributed polarization-entangled photon pairs into opposite
cores within a 19-core multicore fiber. We estimated the secret key rates in a
configuration with 6 and 12 cores from the entanglement visibility after
transmission through $411$ m long multicore fiber.

##  I Introduction

Quantum key distribution (QKD) is an established quantum communication method
that allows generation of identical secret keys between two remotely separated
parties with unconditional security [1, 2]. The communication protocols to
share the key can be based on entangled states [3, 4] or weak coherent states
[5, 6]. This key needs to be shared by trusted channels such as free-space
links on ground [7, 8], ground-to-space [9, 10] and optical fiber links [11,
12, 13]. Over the years, improvements in optical fiber fabrication have
enabled access to new fibers, which offer flexibility, reliability, and high
fidelity in data transmission [14]. This has allowed implementing quantum key
distribution schemes based on single-photon sources over single-mode fibers
(SMF) [15], few-mode fibers [16, 17], multimode fibers [18, 19], vortex fibers
[20], and multicore fibers (MCF) [21, 22]. Recently, schemes based on
entangled-photon states have been shown to generate a positive key rate over a
long-distance fiber link [23, 24]. These schemes represent a challenge due to
high propagation losses in fiber, polarization mode dispersion, and
stabilization of the state in the propagation.

In order to implement quantum communication schemes over long fiber links, it
is crucial to develop techniques to increase channel transmission capacity.
One of the most widely used techniques for that is multiplexing, by
transmitting more than one signal through the same communication channel. This
has been implemented based on the wavelength (WDM) degree of freedom to
increase the secure key rate [25, 26] and carry out novel entanglement-based
quantum networks [27, 28]. But it is also crucial to explore new domains such
as space-division multiplexing (SDM) to overcome WDM’s capacity [29, 30, 31].
In particular, SDM technology over multicore fiber (MCF) has captured
attention due to the low crosstalk between cores. MCF has been employed as a
quantum communication link to increase the secret key rate [32],
implementation of QKD schemes included integrated into silicon photonic chips
[22], simultaneous transmission of classical and quantum signals over a single
MCF [21], and generation of high-dimensional states [33]. Furthermore, the
multiple fibers available over the same fiber link show an excellent
performance to distribute entangled states based on their spatial correlation
[34, 35].

In this work, we experimentally demonstrate generation of parallel and
independent quantum keys based on SDM over MCF using a single entangled photon
source. The multiplexing is performed through the momentum correlation of the
photon pairs, whereas the secure key rate is extracted from the polarization
measurement outcomes by each correlated channel separately. In particular, we
report the key rate associated with the 19-core configurations into the MCF,
corresponding to a SDM scheme over six and twelve cores. We extrapolate the
achievable key rate for the case of a fast-switching QKD setup as a function
of MCF length to explore its applicability. In addition, we demonstrate a
stable and low quantum bit error rate (QBER) below the threshold limit, along
with a positive key rate for one pair of opposite cores. For at least 24 hours
of continuous operation, polarization changes were observed without any active
stabilization technique. Thus, our results show the first step towards the
performance of SDM on entanglement-based quantum communication protocol and
the development of efficient quantum communication distribution.

##  II Results

###  II.1 Experimental setup

Our experimental setup is depicted in Fig. 1. We generated photon pairs at
telecommunication wavelength from a type-0, non-collinear, quasi-phase-matched
spontaneous parametric down-conversion (SPDC) process. By bidirectionally
pumping a magnesium oxide doped periodically poled lithium niobate (MgO:ppLN)
crystal in a Sagnac configuration, the SPDC source produces polarization-
entangled photon pairs in the Bell state,

$\left|\Phi^{+}\right\rangle=\frac{1}{\sqrt{2}}(\left|H_{s},H_{i}\right\rangle+\left|V_{\mathrm{s}},V_{i}\right\rangle),$
(1)

where H (V) indicates horizontal (vertical) polarization for signal (s) and
idler (i) photons. Exploiting momentum conservation in the SPDC process,
signal and idler photons were collected onto the MCF’s opposite cores, which
were located at the far-field plane of the crystal [34],

$\displaystyle|\boldsymbol{\Phi}\rangle=\sum_{m=1}^{N}g_{m}|c_{m},c_{m}^{\prime}\rangle_{\mathrm{cores}},$
(2)

where $c_{m}$ and $c_{m}^{\prime}$ corresponds to the opposite cores, $g_{m}$
are normalized probability amplitudes, $\sum_{m=1}^{N}|g_{m}|^{2}=1$, and $N$
is the number of opposite cores. The far-field plane contains the transverse
momentum of the signal ($\mathbf{p}_{s}=\hbar\mathbf{q}_{s}$) and idler
($\mathbf{p}_{i}=\hbar\mathbf{q}_{i}$), where $\mathbf{q}$ represents the
transverse wave vector [36, 37]. This plane is accessed by performing a
Fourier transform of the photons source’s plane with a configuration of lenses
as shown in Fig. 1.

The MCF selected for the experiment is $411$ m in length, with 19 cores
distributed in a central core and two hexagonal rings of cores, which we refer
to as the inner ring (6 cores) and the outer ring (12 cores). This MCF
presents very low cross-talk between adjacent cores thanks to a trench-
assisted structure around each core [38]. Additionally, each of these cores
ends into an independent SMF via a fan-in/fan-out (FIFO) device [39]. In order
to couple the polarization-entangled photon pairs into opposite cores of each
MCF ring independently, the emission angle of the SPDC photons was adjusted
via temperature control of the nonlinear crystal [40]. Thus, we can reliably
and quickly switch between a different number of channels (cores), extending
the generation of the simultaneous secure keys without changing optical
components in the SPDC source.

To estimate the quantum key rate from each pair of opposite cores (one at a
time), entangled photons were distributed into the polarization analysis
modules to Alice and Bob through the SMF from the FIFO device. In both
modules, the correlated photons were detected by a superconducting nanowire
single-photon detector (SNSPD) operated at efficiencies of $\sim 80\%$ with
dark count rates of $\sim 100$ counts per second. A time-tagging module
recorded the detection events of each pair of opposite cores. Correlated
photon pairs were identified through coincidence counts with the help of the
temporal cross-correlation functions within $300$ ps time-window.

Figure 1: Experimental setup. A MgO:ppLN crystal is pumped in a Sagnac
configuration by a continuous wave laser at $775$ nm (depicted as a blue beam)
to create polarization-entangled photon pairs with the center wavelength of
$1550$ nm (depicted as a red beam). The signal and idler photons were
separated from the pump with a dichroic mirror (DM) and long-pass filter (LPF)
to be distributed through opposite cores of the MCF employing the spatial
correlations of the SPDC photons. These correlations are imaged onto the fiber
end-face by a combination of lenses $\text{L}_{1}$, $\text{L}_{2}$, and
$\text{L}_{3}$. The SPDC photons are emitted in a cone, whose size is adjusted
to fit hexagonal rings of cores (either inner or outer) by controlling the
crystal temperature [40]. The focal lengths of the lenses are
$\text{f}_{1}=200$ mm, $\text{f}_{2}=150$ mm and $\text{f}_{3}=4.51$ mm. The
19 cores are split into individual SMF via a fan-in/fan-out (FIFO) device.
Each pair of opposite cores are connected to a polarization analysis module
(Alice and Bob) consisting of a half-wave plate (HWP) in front of a polarizing
beam splitter (PBS). D1 through D4 refer to detectors (superconducting
nanowire single-photon detectors, SNSPD) used to detect the photons.
Abbreviations: YVO4: yttrium orthovanadate plate, BPF: 1550$\pm$3 nm band-pass
interference filter, PC: polarization controllers.

###  II.2 Estimation of the secure key rate over MCF

We calculated the secure key rate based on the count rates for each pair of
opposite cores [41, 42] and theoretical prediction of its propagation,
assuming an MCF average loss of $0.2$ dB/km per core. A series of polarization
measurements in two mutually unbiased bases allowed us to estimate the secure
key rates on opposite cores in the Horizontal-Vertical (HV)- basis and in the
Diagonal-Antidiagonal (DA) basis. Since our system lacked the ability to have
fast and random basis choices, we estimated the key rate from the polarization
visibility measurement that we would have observed if we had a passive basis
choice on a beam splitter and four detectors, two for each basis. Instead, we
used two motor-driven rotation stages on Alice and Bob to set the half-wave
plates between $0^{\circ}$ (H/V basis) and $22.5^{\circ}$ (D/A basis). They
were followed by a polarizing beam-splitter (PBS) with a fiber-coupler in
order to connect the transmitted and reflected ports to separate detectors,
respectively. For each pair of opposite cores, three pairs for the inner ring
and six pairs for the outer ring, the data was accumulated for $60\,$ s.

Figure 2: Estimated secret key rate as a function of the MCF length. The green
circle represents the overall secret key rate of the inner ring obtained after
a propagation distance of $411$ m in a MCF. A total key rate of $7.3$ kbit/s
was estimated based on the count rates and visibilities measured for a SDM
scheme over 6 cores, i.e., three pairs of opposite cores. The orange square
corresponds to the overall secret key rate of the outer ring obtained after
the same propagation length. In this case, we increased the total key rate up
to $34.5$ kbit/s through an SDM scheme over 12 cores, i.e., six pairs of
cores.

Figure 2 shows the estimated key rate produced for the inner and outer rings
corresponding to green and orange dots, respectively and a theoretical
prediction of how much key rate would be observed under the same conditions as
a function of channel length represented by the solid lines. For the three
pairs of cores of the inner ring, the average coincidence rate was $2287$
counts per second on the H/V basis and $2275$ counts per second on the D/A
basis. For the six pairs of cores of the outer ring, the average coincidence
rate was $7832$ counts per second on the H/V basis and $7770$ counts per
second on the D/A basis. Based on these coincidence count rates, we estimated
an overall secure key rate of $7.3$ kbit/s for the inner ring and $34.5$
kbit/s for the outer ring. Therefore, using one entanglement photon source and
multiplexing in different spatial modes, we achieve up to $4.7$ times higher
key rates in the outer ring than in the inner ring. One would expect a key
rate two times higher in the outer ring than the inner ring due to the number
of cores, the same number of photon pairs generated, and the identical
detection modules. However, we can observe a higher average coincidence rate
in the outer ring, which may be due to the different crystal temperatures of
both configurations, which directly affect SPDC photon emission and spatial
correlations [40]. The spatial correlation increases for a lower crystal
temperature, helping to increase coincidence counts of the polarization-
entangled photon pairs and reducing accidental coincidence counts. Indeed, our
experimental implementation turns out to have different loss levels for the
inner ring of $40.06$ dB and the outer ring of $35.48$ dB. From the
theoretical prediction shown in Fig. 2, the expected key rates that can be
attested to a MCF link have a limiting distance of $180$ km and $190$ km for
the inner and outer rings, respectively. Up to which the key rate would be
positive in the infinite key limit.

Additionally, we performed a study on the system stability by a long-term
measurement of coincidence counts in the H/V and D/A bases for one pair of
opposite cores from the inner ring throughout $24$ hours. The coincidence
counts were collected every $30$ min, alternating the measurement basis for
periods of $60$ s. As evidenced in Fig. 3 (a), QBER (blue dots) remains stable
around $3\%$ for the whole measurement without any active stabilization. It is
worth mentioning that a QBER below $11\%$ ensures a positive secure key rate
[42]. In the other cores of the ring, we estimated the QBER over $60\,$ s,
obtaining in all the combinations QBERs lower than $8.2\%$. The lowest QBER
that can be reached in the inner ring was $3.02\%$ on a H/V basis and $2.42\%$
on a D/A basis. Furthermore, these minima in the outer ring were $2.65\%$ on a
H/V basis and $2.72\%$ on a D/A basis. While Fig. 3 (b) shows the estimated
key rate around $2.3$ kbit/s (green dots) based on the coincidence rate,
considering the error correction efficiency of 1.2 [43].

Figure 3: Estimated a) QBER and b) key rate (KR) over $24$ hours for one pair
of opposite cores from the inner ring. Each data point was calculated from the
polarization visibility measurement in the HV and DA bases over $60$ s,
changing bases every $30$ min.

Hence, employing intrinsic correlations in an SPDC source, it is feasible to
implement a deterministic spatial multiplexing scheme for quantum
communication by extracting independent secure keys of entangled photon pairs.
Note that in our SDM implementation, increasing the number of cores is always
possible by adjusting the crystal temperature, which optimizes the momentum
correlation and allows for avoiding crosstalk. Up to now, MCFs with 37 cores
have been used to implement experimental QKD by employing a classical signal
[32]. This approach adds a third hexagonal ring of 18 cores around our cores
configuration (19 cores) [44]. With such configuration, the third ring would
allow us to connect nine additional core pairs, which would allow obtaining a
secret key rate at least 1.5 times larger than the key rate of the outer ring,
considering the same loss level.

##  III Discussion

We have experimentally shown a space-division multiplexing technique to
enhance the total key rate and, therefore, with promising applications in
quantum communication. Our SDM scheme uses a 19-core MCF wherein we
distributed polarization entangled-photon pairs between 3 pairs of users (6
cores; inner ring) and 6 pairs of users (12 cores; outer ring). The switch
between these two configurations was done by tuning the crystal temperature.
The secure key rate was computed by polarization visibility measurements
between independent pairs of cores. We obtained a coincidence rate for the 6
and 12 cores in the order of $2000$ and $7000$ counts per second,
respectively. With these counts, we computed a key rate of $7.3$ kbit/s
($34.5$ kbit/s) in the inner (outer) ring. In a second analysis, we verified
the estimated QBER and key rate for a pair of cores over 10 hours. Although no
active polarization stabilization techniques were employed, the values
remained stable around $2.2$ kbit/s and a QBER of $3.8\%$.

Our results show a feasible implementation where a 19-core MCF can generate
parallel key rates in independent cores that can be combined, distributed
individually, or grouped for a multi-user quantum network. This scheme can be
easily adapted to a higher number of cores or combined with other degrees of
freedom for enhancing the key rate through a MCF [35]. However, our
experimental implementation can be improved in several aspects, such as the
coupling efficiency of the lens system to implement the far-field plane, which
can be overcome by employing, e.g., a microlens array [45]. Also, the
parameters of the SPDC source, such as the length of the crystal or pump beam
waist, could be optimized to increase the spatial entanglement. An enhancement
in the coupling efficiency could allow the integration of such an SDM scheme
into existing telecommunication infrastructures operating at telecommunication
wavelength. Therefore, our work opens the possibility of employing MCF for
various quantum communication applications.

## Acknowledgments

We acknowledge funding by the European Union’s Horizon 2020 programme grant
agreement No.857156 (OpenQKD) and the Austrian Academy of Sciences. J.C.A.-Z.
and R.A.-C. acknowledge support US. Army Research Office W911NF1710553 and
National Science Foundation ECCS-1711230. E.A.O. acknowledges ANID for the
financial support (Becas de doctorado en el extranjero “Becas Chile”/2016 –
No. 72170402).

##  IV Methods

###  IV.1 Entangled photon pair source

The SPDC-source is based on $40$ mm long Magnesium Oxide doped periodically
poled Lithium Niobate (MgO:ppLN) crystal with a poling period of $19.2$
$\upmu$m placed in a Sagnac interferometer [46]. The nonlinear crystal was
pumped by a continuous-wave laser at the center wavelength $775$ nm. In the
down-conversion process, one pump photon is converted with low probability to
a pair of signal and idler photons at telecommunication wavelength around
$1550$ nm with the same polarization as the pump photon.

The far-field plane of the crystal was accessed by the lens $\text{L}_{1}$,
where the transverse detection positions of the photon pairs are spatially
distributed in opposite directions. So the photon pairs in the far-field plane
were coupled into the MCF by the 4f-imaging system consisting of the lenses
$\text{L}_{2}$ and $\text{L}_{3}$. This demagnified image in the far-field
plane ensured that the emission cone diameter covered all cores of the
corresponding ring. Additionally, for adjusting the SPDC emission cone
diameter to the inner and outer ring, the crystal temperature is set at
$82.5^{\circ}$C and $82^{\circ}$C, respectively.

###  IV.2 Estimation of key rate

For each pair of opposite cores, the secret key rate has been estimated based
on the assumption that the setup used a fast and random basis choice for each
photon [41, 42]. Considering $H_{2}$ the binary Shannon entropy function
defined as

$H_{2}\left(x\right)=-x\log_{2}\left(x\right)-\left(1-x\right)\log_{2}\left(1-x\right).$
(3)

The secret key rate is given by:

$\begin{split}R&=\frac{1}{2}C_{HV}\left(1-\left(1+f\right)H_{2}\left(Q_{HV}\right)\right)\\\
&\quad+\frac{1}{2}C_{DA}\left(1-(1+f)H_{2}\left(Q_{DA}\right)\right),\end{split}$
(4)

where $C_{HV}$ and $C_{DA}$ are the total number of coincidence counts in the
H/V and D/A bases, respectively. $Q_{HV}$ and $Q_{DA}$ are the QBER in the H/V
and D/A bases, respectively, which can be calculated from the coincidence
counts at the given polarization combination over $60$ s via

$\centering\begin{split}Q_{HV}&=\frac{\textup{error~{}counts}}{\textup{total~{}counts}}\\\
&=\frac{1-V_{HV}}{2}.\end{split}\@add@centering$ (5)

Where $V_{HV}$ corresponds to the visibility in H/V basis,

$V_{HV}=\frac{C_{HH}+C_{VV}-C_{HV}-C_{VH}}{C_{HH}+C_{VV}+C_{HV}+C_{VH}}$ (6)

with $C_{s_{1},s_{2}}$ is the number of coincidence counts for each
polarization setting. The QBER on DA basis is calculated analogously.

## References

  * Pirandola  _et al._ [2020] S. Pirandola, U. L. Andersen, L. Banchi, M. Berta, D. Bunandar, R. Colbeck, D. Englund, T. Gehring, C. Lupo, C. Ottaviani, _et al._ , Advances in quantum cryptography, Adv. Opt. Photon. 12, 1012 (2020).
  * Xu  _et al._ [2020] F. Xu, X. Ma, Q. Zhang, H.-K. Lo, and J.-W. Pan, Secure quantum key distribution with realistic devices, Rev. Mod. Phys. 92, 025002 (2020).
  * Ursin  _et al._ [2007] R. Ursin, F. Tiefenbacher, T. Schmitt-Manderbach, H. Weier, T. Scheidl, M. Lindenthal, B. Blauensteiner, T. Jennewein, J. Perdigues, P. Trojek, _et al._ , Entanglement-based quantum communication over 144 km, Nature Phys. 3, 481 (2007).
  * Cabrejo-Ponce  _et al._ [2022] M. Cabrejo-Ponce, C. Spiess, A. L. M. Muniz, P. Ancsin, and F. Steinlechner, Ghz-pulsed source of entangled photons for reconfigurable quantum networks, Quantum Sci. Technol. 7, 045022 (2022).
  * Lo  _et al._ [2005] H.-K. Lo, X. Ma, and K. Chen, Decoy state quantum key distribution, Phys. Rev. Lett. 94, 230504 (2005).
  * Wang  _et al._ [2022] S. Wang, Z.-Q. Yin, D.-Y. He, W. Chen, R.-Q. Wang, P. Ye, Y. Zhou, G.-J. Fan-Yuan, F.-X. Wang, Y.-G. Zhu, _et al._ , Twin-field quantum key distribution over 830-km fibre, Nat. Photon. 16, 154 (2022).
  * Scheidl  _et al._ [2009] T. Scheidl, R. Ursin, A. Fedrizzi, S. Ramelow, X.-S. Ma, T. Herbst, R. Prevedel, L. Ratschbacher, J. Kofler, T. Jennewein, and A. Zeilinger, Feasibility of 300 km quantum key distribution with entangled states, New J. Phys. 11, 085002 (2009).
  * Kržič  _et al._ [2023] A. Kržič, S. Sharma, C. Spiess, U. Chandrashekara, S. Töpfer, G. Sauer, L. J. G.-M. del Campo, T. Kopf, S. Petscharnig, T. Grafenauer, R. Lieger, B. Ömer, C. Pacher, R. Berlich, T. Peschel, C. Damm, S. Risse, M. Goy, D. Rieländer, A. Tünnermann, and F. Steinlechner, Towards metropolitan free-space quantum networks, npj Quantum Inf. 9, 95 (2023).
  * Liao  _et al._ [2017] S.-K. Liao, W.-Q. Cai, W.-Y. Liu, L. Zhang, Y. Li, J.-G. Ren, J. Yin, Q. Shen, Y. Cao, Z.-P. Li, _et al._ , Satellite-to-ground quantum key distribution, Nature 549, 43 (2017).
  * Yin  _et al._ [2020] J. Yin, Y.-H. Li, S.-K. Liao, M. Yang, Y. Cao, L. Zhang, J.-G. Ren, W.-Q. Cai, W.-Y. Liu, S.-L. Li, _et al._ , Entanglement-based secure quantum cryptography over 1,120 kilometres, Nature 582, 501 (2020).
  * Boaron  _et al._ [2018] A. Boaron, G. Boso, D. Rusca, C. Vulliez, C. Autebert, M. Caloz, M. Perrenoud, G. Gras, F. Bussières, M.-J. Li, D. Nolan, A. Martin, and H. Zbinden, Secure quantum key distribution over 421 km of optical fiber, Phys. Rev. Lett. 121, 190502 (2018).
  * Shi  _et al._ [2020] Y. Shi, S. M. Thar, H. S. Poh, J. A. Grieve, C. Kurtsiefer, and A. Ling, Stable polarization entanglement based quantum key distribution over a deployed metropolitan fiber, Appl. Phys. Lett. 117, 124002 (2020).
  * Chen  _et al._ [2022] J.-P. Chen, C. Zhang, Y. Liu, C. Jiang, D.-F. Zhao, W.-J. Zhang, F.-X. Chen, H. Li, L.-X. You, Z. Wang, Y. Chen, X.-B. Wang, Q. Zhang, and J.-W. Pan, Quantum key distribution over 658 km fiber with distributed vibration sensing, Phys. Rev. Lett. 128, 180502 (2022).
  * Winzer  _et al._ [2018] P. J. Winzer, D. T. Neilson, and A. R. Chraplyvy, Fiber-optic transmission and networking: the previous 20 and the next 20 years  $invited$ , Opt. Express 26, 24190 (2018).
  * Ribezzo  _et al._ [2023] D. Ribezzo, M. Zahidy, I. Vagniluca, N. Biagi, S. Francesconi, T. Occhipinti, L. K. Oxenløwe, M. Lončarić, I. Cvitić, M. Stipčević, Ž. Pušavec, R. Kaltenbaek, A. Ramšak, F. Cesa, G. Giorgetti, F. Scazza, A. Bassi, P. De Natale, F. S. Cataliotti, M. Inguscio, D. Bacco, and A. Zavatta, Deploying an inter-european quantum network, Adv. Quantum Technol. 6, 2200061 (2023).
  * Wang  _et al._ [2020] B.-X. Wang, Y. Mao, L. Shen, L. Zhang, X.-B. Lan, D. Ge, Y. Gao, J. Li, Y.-L. Tang, S.-B. Tang, J. Zhang, T.-Y. Chen, and J.-W. Pan, Long-distance transmission of quantum key distribution coexisting with classical optical communication over a weakly-coupled few-mode fiber, Opt. Express 28, 12558 (2020).
  * Alarcón  _et al._ [2021] A. Alarcón, J. Argillander, G. Lima, and G. Xavier, Few-mode-fiber technology fine-tunes losses in quantum communication systems, Phys. Rev. Applied 16, 034018 (2021).
  * Zhou  _et al._ [2021] Y. Zhou, B. Braverman, A. Fyffe, R. Zhang, J. Zhao, A. E. Willner, Z. Shi, and R. W. Boyd, High-fidelity spatial mode transmission through a 1-km-long multimode fiber via vectorial time reversal, Nat. Commun. 12, 1866 (2021).
  * Amitonova  _et al._ [2020] L. V. Amitonova, T. B. Tentrup, I. M. Vellekoop, and P. W. Pinkse, Quantum key establishment via a multimode fiber, Opt. Express 28, 5965 (2020).
  * Sit  _et al._ [2018] A. Sit, R. Fickler, F. Alsaiari, F. Bouchard, H. Larocque, P. Gregg, L. Yan, R. W. Boyd, S. Ramachandran, and E. Karimi, Quantum cryptography with structured photons through a vortex fiber, Opt. Lett. 43, 4108 (2018).
  * Dynes  _et al._ [2016] J. F. Dynes, S. J. Kindness, S. W.-B. Tam, A. Plews, A. W. Sharpe, M. Lucamarini, B. Fröhlich, Z. L. Yuan, R. V. Penty, and A. J. Shields, Quantum key distribution over multicore fiber, Opt. Express 24, 8081 (2016).
  * Bacco  _et al._ [2017] D. Bacco, Y. Ding, K. Dalgaard, K. Rottwitt, and L. K. Oxenløwe, Space division multiplexing chip-to-chip quantum key distribution, Sci. Rep. 7, 1 (2017).
  * Wengerowsky  _et al._ [2020] S. Wengerowsky, S. K. Joshi, F. Steinlechner, J. R. Zichi, B. Liu, T. Scheidl, S. M. Dobrovolskiy, R. van der Molen, J. W. N. Los, V. Zwiller, M. A. M. Versteegh, A. Mura, D. Calonico, M. Inguscio, A. Zeilinger, A. Xuereb, and R. Ursin, Passively stable distribution of polarisation entanglement over 192 km of deployed optical fibre, npj Quantum Inf. 6, 5 (2020).
  * Neumann  _et al._ [2022a] S. P. Neumann, A. Buchner, L. Bulla, M. Bohmann, and R. Ursin, Continuous entanglement distribution over a transnational 248 km fiber link, Nat. Commun. 13, 6134 (2022a).
  * Neumann  _et al._ [2022b] S. P. Neumann, M. Selimovic, M. Bohmann, and R. Ursin, Experimental entanglement generation for quantum key distribution beyond 1 gbit/s, Quantum 6, 822 (2022b).
  * Brambila  _et al._ [2023] E. Brambila, R. Gómez, R. Fazili, M. Gräfe, and F. Steinlechner, Ultrabright polarization-entangled photon pair source for frequency-multiplexed quantum communication in free-space, Opt. Express 31, 16107 (2023).
  * Wengerowsky  _et al._ [2018] S. Wengerowsky, S. K. Joshi, F. Steinlechner, H. Hübel, and R. Ursin, An entanglement-based wavelength-multiplexed quantum communication network, Nature 564, 225 (2018).
  * Joshi  _et al._ [2020] S. K. Joshi, D. Aktas, S. Wengerowsky, M. Lončarić, S. P. Neumann, B. Liu, T. Scheidl, G. C. Lorenzo, Ž. Samec, L. Kling, _et al._ , A trusted node–free eight-user metropolitan quantum communication network, Sci. Adv. 6, eaba0959 (2020).
  * Richardson  _et al._ [2013] D. J. Richardson, J. M. Fini, and L. E. Nelson, Space-division multiplexing in optical fibres, Nature Photon. 7, 354 (2013).
  * Puttnam  _et al._ [2021] B. J. Puttnam, G. Rademacher, and R. S. Luís, Space-division multiplexing for optical fiber communications, Optica 8, 1186 (2021).
  * Xavier and Lima [2020] G. B. Xavier and G. Lima, Quantum information processing with space-division multiplexing optical fibres, Commun. Phys. 3, 9 (2020).
  * Bacco  _et al._ [2019] D. Bacco, B. D. Lio, D. Cozzolino, F. D. Ros, X. Guo, Y. Ding, Y. Sasaki, K. Aikawa, S. Miki, H. Terai, T. Yamashita, J. S. Neergaard-Nielsen, M. Galili, K. Rottwitt, U. L. Andersen, T. Morioka, and L. K. Oxenløwe, Boosting the secret key rate in a shared quantum and classical fibre communication system, Commun. Phys. 2, 140 (2019).
  * Cañas  _et al._ [2017] G. Cañas, N. Vera, J. Cariñe, P. González, J. Cardenas, P. W. R. Connolly, A. Przysiezna, E. S. Gómez, M. Figueroa, G. Vallone, P. Villoresi, T. F. da Silva, G. B. Xavier, and G. Lima, High-dimensional decoy-state quantum key distribution over multicore telecommunication fibers, Phys. Rev. A 96, 022317 (2017).
  * Ortega  _et al._ [2021] E. A. Ortega, K. Dovzhik, J. Fuenzalida, S. Wengerowsky, J. C. Alvarado-Zacarias, R. F. Shiozaki, R. Amezcua-Correa, M. Bohmann, and R. Ursin, Experimental space-division multiplexed polarization-entanglement distribution through 12 paths of a multicore fiber, PRX Quantum 2, 040356 (2021).
  * Achatz  _et al._ [2023] L. Achatz, L. Bulla, S. Ecker, E. A. Ortega, M. Bartokos, J. C. Alvarado-Zacarias, R. Amezcua-Correa, M. Bohmann, R. Ursin, and M. Huber, Simultaneous transmission of hyper-entanglement in three degrees of freedom through a multicore fiber, npj Quantum Inf. 9, 45 (2023).
  * Walborn  _et al._ [2010] S. Walborn, C. Monken, S. Pádua, and P. Souto Ribeiro, Spatial correlations in parametric down-conversion, Phys. Rep. 495, 87 (2010).
  * Achatz  _et al._ [2022] L. Achatz, E. A. Ortega, K. Dovzhik, R. F. Shiozaki, J. Fuenzalida, S. Wengerowsky, M. Bohmann, and R. Ursin, Certifying position-momentum entanglement at telecommunication wavelengths, Physica Scripta 97, 015101 (2022).
  * Takenaga  _et al._ [2011] K. Takenaga, Y. Arakawa, S. Tanigawa, N. Guan, S. Matsuo, K. Saitoh, and M. Koshiba, Reduction of crosstalk by trench-assisted multi-core fiber, in  _2011 Optical Fiber Communication Conference and Exposition and the National Fiber Optic Engineers Conference_ (2011) p. OWJ4.
  * Alvarado-Zacarias  _et al._ [2019] J. C. Alvarado-Zacarias, J. E. Antonio-Lopez, M. S. Habib, S. Gausmann, N. Wang, D. Cruz-Delgado, G. Li, A. Schülzgen, A. Amezcua-Correa, L.-A. Demontmorillon, P. Sillard, and R. Amezcua-Correa, Low-loss 19 core fan-in/fan-out device using reduced-cladding graded index fibers, in  _Optical Fiber Communication Conference (OFC) 2019_ (Optical Society of America, 2019) p. Th3D.2.
  * Ortega  _et al._ [2023] E. A. Ortega, J. Fuenzalida, M. Selimovic, K. Dovzhik, L. Achatz, S. Wengerowsky, R. F. Shiozaki, S. P. Neumann, M. Bohmann, and R. Ursin, Spatial and spectral characterization of photon pairs at telecommunication-wavelength from type-0 spontaneous parametric downconversion, J. Opt. Soc. Am. B 40, 165 (2023).
  * Ma  _et al._ [2007] X. Ma, C.-H. F. Fung, and H.-K. Lo, Quantum key distribution with entangled photon sources, Phys. Rev. A 76, 012307 (2007).
  * Neumann  _et al._ [2021] S. P. Neumann, T. Scheidl, M. Selimovic, M. Pivoluska, B. Liu, M. Bohmann, and R. Ursin, Model for optimizing quantum key distribution with continuous-wave pumped entangled-photon sources, Phys. Rev. A 104, 022406 (2021).
  * Elkouss  _et al._ [2009] D. Elkouss, A. Leverrier, R. Alleaume, and J. J. Boutros, Efficient reconciliation protocol for discrete-variable quantum key distribution, in  _2009 IEEE International Symposium on Information Theory_ (2009) pp. 1879–1883.
  * Sasaki  _et al._ [2017] Y. Sasaki, K. Takenaga, K. Aikawa, Y. Miyamoto, and T. Morioka, Single-mode 37-core fiber with a cladding diameter of 248 $\mu$m, in  _2017 Optical Fiber Communications Conference and Exhibition (OFC)_ (IEEE, 2017) p. Th1H.2.
  * Dietrich  _et al._ [2017] P.-I. Dietrich, R. J. Harris, M. Blaicher, M. K. Corrigan, T. J. Morris, W. Freude, A. Quirrenbach, and C. Koos, Printed freeform lens arrays on multi-core fibers for highly efficient coupling in astrophotonic systems, Opt. Express 25, 18288 (2017).
  * Kim  _et al._ [2006] T. Kim, M. Fiorentino, and F. N. C. Wong, Phase-stable source of polarization-entangled photons using a polarization sagnac interferometer, Phys. Rev. A 73, 012316 (2006).


# First Doppler Limits on Binary Planets and Exomoons in the HR 8799 System

Andrew Vanderburg Department of Physics and Kavli Institute for Astrophysics
and Space Research, Massachusetts Institute of Technology, 77 Massachusetts
Avenue, Cambridge, MA 02139, USA  Joseph E. Rodriguez Department of Physics
and Astronomy, Michigan State University, East Lansing, MI 48824, USA Andrew
Vanderburg<EMAIL_ADDRESS>
###### Abstract

We place the first constraints on binary planets and exomoons from Doppler
monitoring of directly imaged exoplanets. We model radial velocity
observations of HR 8799 b, c, and d from Ruffio et al. (2021) and determine
upper limits on the $m\sin{i}$ of short-period binary planets and satellites.
At 95% confidence, we rule out companions orbiting the three planets more
massive than $m\sin{i}=2$ MJ with orbital periods shorter than 5 days. We
achieve our tightest constraints on moons orbiting HR 8799 c, where with 95%
confidence we rule out out edge-on Jupiter-mass companions in periods shorter
than 5 days and edge-on half-Jupiter-mass moons in periods shorter than 1 day.
These radial velocity observations come from spectra with resolution 20 times
lower than typical radial velocity instruments and were taken using a
spectrograph that was designed before the first directly imaged exoplanet was
discovered. Similar datasets from new and upcoming instruments will probe
significantly lower exomoon masses.

planetary systems, planets and satellites: detection, stars: individual (HR
8799)

††journal: ApJL††facilities: Keck I (OSIRIS)††software: matplotlib (Hunter,
2007), numpy (Oliphant, 2006), scipy (Virtanen et al., 2020), EXOFASTv2
(Eastman et al., 2019), edmcmc (Vanderburg, 2021)

##  1 Introduction

Three decades of exoplanet discovery have resulted in the detection of over
4500 planets outside the solar
system111https://exoplanetarchive.ipac.caltech.edu/cgi-bin/TblView/nph-
tblView?app=ExoTbls&config=PS&constraint=default_flag=1, accessed 18
September, 2021. (Akeson et al., 2013) and the knowledge that planets
outnumber stars in our galaxy (e.g. Swift et al., 2013; Dressing &
Charbonneau, 2015). However, despite the prevalence of planets in our galaxy
and the ubiquity of moons orbiting the planets in our solar system,
astronomers have not yet securely detected any exomoons and have only a
handful of unconfirmed candidates (e.g. Teachey & Kipping 2018; Lazzoni et al.
2020; Limbach et al. 2021, but see also Kreidberg et al. 2019). Given the
importance of Earth’s moon on our planet’s spin dynamics (Li & Batygin, 2014)
and the potential habitability of icy moons in the outer solar system (e.g.
Reynolds et al., 1983), continuing to search for moons beyond our solar system
is a worthwhile endeavor.

There are two likely reasons for the lack of secure exomoon detections. First,
large exomoons and binary planets (defined here as objects orbiting planets
with enough mass that the system’s center of mass is outside the primary
planet, Stern & Levison 2002) are apparently rare. The Kepler mission would
have be highly sensitive to any nearly equal-sized binary-planet companions,
but detected none (Lewis et al., 2015; Teachey et al., 2018). In a number of
cases, upper limits on the presence of exomoons and binary-planet companions
have been placed using transit photometry that rule out moons larger than
Earth, and in the best cases rule out moons almost as small as Ganymede
(Kipping et al., 2014, 2015). The rarity of large exomoons can be understood
by the fact that moons forming circumplanetary disks are expected to have
masses smaller than $10^{-4}$ times that the planetary host (Canup & Ward,
2006), while the gravitational capture of giant moons or binary planets is
likely to be a rare dynamical event (Ochiai et al., 2014). Second, our
observations are not yet sensitive to small moons. There are over 200 moons
orbiting the 8 major Solar System planets and many more orbiting smaller
bodies, but none of these moons would be detectable from afar with our current
technology. Many different exomoon detection methods have been proposed (e.g.
Sartoretti & Schneider, 1999; Hook, 2005; Cabrera & Schneider, 2007; Lewis et
al., 2008; Kipping, 2009; Agol et al., 2015; Heller, 2016; Sengupta & Marley,
2016; Forgan, 2017; Hwang et al., 2018; Vanderburg et al., 2018; Lazzoni et
al., 2020) but so far only transit photometry (Kipping et al., 2015; Teachey
et al., 2018) and direct imaging (Lazzoni et al., 2020) have yielded
constraints on the presence of exomoons. Future space observations from the
Roman Space Telescope microlensing survey (Spergel et al., 2013; Penny et al.,
2019) are expected to either detect or constrain planets/moons orbiting low-
mass brown dwarfs or free-floating planets, although unambiguous detections
will be challenging (Hwang et al., 2018).

In this work, we place the first limits on exomoons and binary planets from
radial velocity measurements of directly imaged exoplanets. This method has
been mentioned occasionally in literature since 2005 (Hook, 2005; Cabrera &
Schneider, 2007; Heller, 2016; Lillo-Box et al., 2018) and was recently
discussed in detail by Vanderburg et al. (2018) who showed that existing and
forthcoming instruments could yield strong constraints on moons orbiting
bright directly imaged planets. We apply this method to three such planets in
the HR 8799 multi-planet system (Marois et al., 2008, 2010). HR 8799 hosts
four directly imaged super-Jovian exoplanets named HR 8799 e, HR 8799 d, HR
8799 c, and HR 8799 b (in order of increasing orbital semimajor axis). HR 8799
c, d, and e all have masses of $7.2_{-0.7}^{+0.6}$MJ, while the outermost
planet HR 8799 b has a slightly lower mass of $5.8\pm 0.5$MJ (Wang et al.,
2018b). Recently, Ruffio et al. (2021) observed HR 8799 b, c, and d and
reported repeated planetary radial velocity measurements with a typical
precision of 1 $\rm km\,s^{-1}$. Here, we show that these radial velocity
observations are sensitive to massive exomoons in short-period orbits and use
the time series to place upper limits on the presence of exomoon and binary
companions to the HR 8799 planets. This Letter is organized as follows:
Section 2 describes the observations and our radial velocity analysis, Section
3 describes the constraints we are able to place, and Section 4 discusses the
implications of our upper limits and the prospects for improved constraints or
exomoon detections with future datasets.

##  2 Observations and Analysis

###  2.1 Keck OSIRIS Radial Velocities

We use the radial velocity observations of HR 8799 b, c, and d reported by
Ruffio et al. (2021). These observations were made over the course of a decade
(between UT 2010 July 11 and UT 2020 August 03) with the OH-Suppressing
Infrared Integral Field Spectrograph (OSIRIS, Larkin et al. 2006; Mieda et al.
2014) on the Keck telescopes. OSIRIS is a moderate-resolution
($\lambda/\Delta\lambda\approx 4000$) integral field spectrograph capable of
measuring the spectra of multiple objects across a small field of view. While
almost all observations were taken with the smallest plate scale (with
0$\farcs$02 pixels and a 0$\farcs$32$\times$1$\farcs$28 field of view), one
was taken with a slightly larger plate scale (0$\farcs$035 pixels and a
0$\farcs$56$\times$2$\farcs$24 field of view). Most of the observations
reported were made in the K-band, some were derived from both H- and K-bands,
and others were made entirely from the H-band observations (see Table 1).
Ruffio et al. (2021) extracted spectra and measured radial velocities using
least-squares template fitting, while simultaneously fitting for stellar
contamination and modeling telluric absorption as a sum of vectors from a
principal component analysis of many spectroscopic observations. The radial
velocity calibration was performed by measuring the location of OH- sky
emission lines. In total, Ruffio et al. (2021) reported 11 observations of HR
8799 b, 8 observations of HR 8799 c, and 7 observations of HR 8799 d. We list
the these radial velocity measurements in Table 1 and show them in Figure 1.

Two other groups have reported radial velocity measurements of the HR 8799
planets using different instruments on the Keck telescopes. Wang et al.
(2018a) used the NIRSPEC near infrared spectrograph behind the Keck adaptive
optics system to measure the radial velocity of HR 8799 c and Wang et al.
(2021a) used NIRSPEC with the the Keck Planet Imager and characterizer (KPIC)
front-end to measure the radial velocity of all four planets in the HR 8799
system. Both groups achieved radial velocity uncertainties as small as 1-2
km/s, at the level of the OSIRIS observations used here, but only reported one
measurement per planet. Because the different instruments used to collect the
observations and the different data analysis methods likely introduce offsets
between these measurements and the OSIRIS datasets, we do not include them in
our analysis.

Figure 1: Radial velocity observations of the HR 8799 planets: HR 8799 b
(top), HR 8799 c (middle), and HR 8799 d (bottom). The points show the radial
velocity measurements from Ruffio et al. (2021) with their uncertainties
inflated by adding the jitter term from our analysis in Section 2.2 in
quadrature. The transparent rectangular regions show the best-fit mean radial
velocities from our analysis for each planet and their 1$\sigma$
uncertainties. Table 1: Radial Velocities of HR 8799 b, c, and d

Planet | Time | RV | $\sigma_{\rm RV}$ | Pixel | Filter/
---|---|---|---|---|---
Name | BJD | $\rm{km}\,\rm{s}^{-1}$ | $\rm{km}\,\rm{s}^{-1}$ | Scale(″) | Band
HR 8799 b | 2455389.05 | -10.6 | 1.1 | 0.02 | K
HR 8799 b | 2455390.04 | -12.4 | 1.1 | 0.02 | K
HR 8799 b | 2455391.12 | -3.7 | 2.2 | 0.02 | H
HR 8799 b | 2456499.05 | -9.5 | 1.0 | 0.02 | K
HR 8799 b | 2456500.0 | -8.4 | 1.3 | 0.02 | K
HR 8799 b | 2456501.01 | -9.1 | 1.9 | 0.02 | K
HR 8799 b | 2457698.82 | -8.2 | 2.2 | 0.02 | K
HR 8799 b | 2457699.82 | -9.7 | 1.5 | 0.02 | K
HR 8799 b | 2457700.84 | -3.2 | 2.1 | 0.02 | K
HR 8799 b | 2458322.02 | -6.2 | 1.3 | 0.035 | K
HR 8799 b | 2459064.96 | -9.9 | 2.3 | 0.02 | K
HR 8799 c | 2455393.05 | -11.8 | 0.7 | 0.02 | K
HR 8799 c | 2455504.81 | -11.0 | 0.7 | 0.02 | H/K
HR 8799 c | 2455766.04 | -10.0 | 1.3 | 0.02 | K
HR 8799 c | 2455767.06 | -10.8 | 2.1 | 0.02 | H/K
HR 8799 c | 2455768.07 | -11.0 | 2.2 | 0.02 | H/K
HR 8799 c | 2456500.08 | -15.8 | 3.5 | 0.02 | K
HR 8799 c | 2458060.75 | -10.5 | 0.9 | 0.02 | H/K
HR 8799 c | 2459060.0 | -10.6 | 1.3 | 0.02 | K
HR 8799 d | 2457224.01 | -13.8 | 2.8 | 0.02 | K
HR 8799 d | 2457227.02 | -18.3 | 2.7 | 0.02 | K
HR 8799 d | 2457263.04 | -13.6 | 3.5 | 0.02 | K
HR 8799 d | 2459060.09 | -11.4 | 1.6 | 0.02 | K
HR 8799 d | 2459061.06 | -12.9 | 1.9 | 0.02 | K
HR 8799 d | 2459062.04 | -10.8 | 1.1 | 0.02 | K
HR 8799 d | 2459065.05 | -10.1 | 1.8 | 0.02 | K

These observations are as reported by Ruffio et al. (2021).

###  2.2 MCMC Radial Velocity Analysis

We searched the radial velocity time series of all three planets for evidence
of orbiting companions. We calculated Lomb-Scargle periodograms (Lomb, 1976;
Scargle, 1982) of the three time series, and saw no clear peaks at periods
between 0.1 and 10 days. A boot-strap analysis yielded false alarm
probabilities for the highest peaks in the power spectra for HR 8799 b, c, and
d of 68%, 49%, and 15%, respectively. We concluded that there is no evidence
for orbiting companions in the OSIRIS radial velocity observations.

We therefore determined upper limits on the presence of putative companions to
the HR 8799 planets using Markov Chain Monte Carlo (MCMC) explorations of
plausible short-period orbits. Radial velocity observations are unable to
constrain the true mass of planetary companions because of a well-known
degeneracy with the companion’s orbital inclination, $i$. We therefore place
constraints on $m\sin{i}$ rather than the true masses of orbiting binary
planets or moon. In our MCMC fits, we modeled each planets’ radial velocity
time series as the sum of a single Keplerian signal describing a moon/binary-
planet companion (restricted to circular orbits), a constant radial velocity
offset, and a linear velocity trend in time (to model each planets’ orbital
motion). We fit for a white noise “jitter” term, $j$, which we added in
quadrature to the error bars, $\sigma$, reported by Ruffio et al. (2021) to
account for any additional sources of radial velocity uncertainty. That is,
the effective uncertainty, $\sigma_{\rm new},$ of each observation was
inflated to $\sigma_{\rm new}\equiv\sqrt{\sigma^{2}+j^{2}}$. We imposed a
Gaussian prior on the masses of the three planets with central values and
widths for HR 8799 b, c, and d of $M_{b}=5.8\pm 0.5$, $M_{c}=7.2\pm 0.7$, and
$M_{d}=7.2\pm 0.7$ (Wang et al., 2018b). The radial velocity datasets from
Ruffio et al. (2021) are sparsely sampled, so we were forced to restrict the
putative companion orbital periods to be between 0 and 10 days. This
restricted parameter space is still interesting though; more than a quarter of
the moons in the solar system, comprising nearly half the total moon mass,
orbit in periods shorter than 10 days222https://ssd.jpl.nasa.gov/sats/elem/.
We forced the $m\sin{i}$ of the putative companion to be nonnegative with a
uniform prior over the interval $[0,\infty)$. Finally, we restricted the
orbital phase of conjunction to be between 0 and 1.

All in all, our model included seven parameters: the mass of each host planet,
the $m\sin{i}$, orbital period, and phase of conjunction for the putative
moon/binary-planet companion, a radial velocity offset, a linear radial
velocity trend, and an excess radial velocity “jitter” term. We explored this
seven-dimensional parameter space using edmcmc (Vanderburg, 2021), a Python
implementation333https://github.com/avanderburg/edmcmc of the differential
evolution MCMC sampler of Ter Braak (2006). For each planet, we ran 50
parallel chains for 400,000 links, discarding the first 20,000 links as burn-
in. We assessed convergence using the Gelman-Rubin diagnostic (Gelman & Rubin,
1992), which was below 1.02 for all parameters. As an independent check on our
analysis, we ran a similar MCMC fit for HR 8799 b using the EXOFASTv2 software
(Eastman et al., 2019) and found consistent results.

Figure 2: Limits on the presence of orbiting companions to the HR 8799
planets: HR 8799 b (top), HR 8799 c (middle), and HR 8799 d (bottom). The
background shading shows the posterior probability distribution of allowed
moon/binary planet $m\sin{i}$ as a function of orbital period. The dark,
medium, and light purple lines show the $1\sigma$, $2\sigma$, and $3\sigma$
(68%, 95%, and 99.7% confidence) contours as a function of orbital period. The
grey lines show the boundary between a binary-planet system and a planet/moon
system based on the definition of Stern & Levison (2002) – objects below the
grey lines are moons, while objects above the grey lines are binary-planet
companions. The spiky vertical structures in the probability distribution and
contours at longer orbital periods are due to the dataset’s sparse sampling.

##  3 Results

Once our MCMC analysis was complete, we marginalized over most of our free
parameters to yield constraints on the presence of exomoons and binary-planet
companions as a function of moon $m\sin{i}$ and orbital period. Our
constraints are summarized in Figure 2, and our $1\sigma$, $2\sigma$, and
$3\sigma$ limits are provided as Data Behind the Figure.

We find that for all three planets, we are able to rule out most companions
with $m\sin{i}$ greater than 2MJ in periods shorter than 5 days. In
particular, our 95% confidence limits for companions in periods shorter than 5
days are 1.1 MJ, 0.82 MJ, and 1.7 MJ for HR 8799 b, c, and d respectively.
Between 5 and 10 days, our results are less constraining: we rule out
companions more massive than $m\sin{i}$ = 2.4 MJ, 1.8 MJ, and 2.3 MJ for the
three planets at 95% confidence. Our constraints are stronger for short-period
companions for two reasons: 1) the radial velocity signal scales with period
as $P^{-1/3}$ so longer-period companions induce smaller signals and 2) the
sparse sampling of the RV observations causes poor phase coverage at longer
orbital periods, making it more difficult to rule out all companions. The
effects of the sparse radial velocity sampling is visible in Figure 2 as
vertical spikes in both the posterior probability distribution and our
$1\sigma$, $2\sigma$, and $3\sigma$ limits at specific orbital periods. The
spikes become more prominent at periods longer than 5 days, reflecting the
typical sampling of the dataset, which have observations on a few consecutive
nights interspersed between years-long gaps without data.

To determine whether our limits predominantly rule out exomoons versus binary
planets, we calculated the boundary between these two classes of objects in
$m\sin{i}$/period space (assuming edge-on inclinations). Stern & Levison
(2002) define a moon as an orbiting companion with a small enough mass and
semimajor axis that the moon/planet system center of mass remains within the
radius of the primary planet. A binary planet is any companion massive enough
that the system center of mass lies outside the primary planet. Using this
rule, we calculate the mass boundary, $M_{\rm boundary}$ above which a
companion is a binary planet rather than a moon:

$M_{\rm boundary}=\frac{M_{p}\,R_{p}}{a-R_{p}}$ (1)

where $M_{p}$, $R_{p}$, and $a$ are the planet mass, planet radius, and
companion semimajor axis, respectively. We plot this curve in
$m\sin{i}$/period space in Figure 2 and find that our radial velocity
constraints are sufficient to rule out almost all edge-on binary planets and
some massive edge-on exomoons with periods shorter than 1-2 days. At orbital
periods longer than about 2 days, our observations are only sensitive to
binary planets.

A natural by-product of our MCMC analysis is that we measure the mean velocity
of each planet over the span of the observations. Since we also fit for a
linear velocity trend, we calculate the velocity of each planet at the mean
time of all of its observation timestamps. We measure mean velocities of
$-8.4\pm 1.7$ $\rm km\,s^{-1}$, $-11.1\pm 1.1$ $\rm km\,s^{-1}$, and $-13.0\pm
1.8$ $\rm km\,s^{-1}$ for HR 8799 b, c, and d, respectively. These velocities
are largely consistent with the values ($-9.1\pm 0.4$ $\rm km\,s^{-1}$,
$-11.1\pm 0.4$ $\rm km\,s^{-1}$, and $-11.6\pm 0.8$ $\rm km\,s^{-1}$ for HR
8799 b, c, and d) calculated by Ruffio et al. (2021), but have larger error
bars (due to our use of a jitter term to inflate the individual error bars)
and are slightly more extreme (likely because including the jitter term gives
slightly more relative credence to apparently outlying datapoints). The more
extreme velocities we measure for the planets match the predicted velocities
from orbital fits ($-8.6^{+0.5}_{-0.6}$ $\rm km\,s^{-1}$,
$-10.5^{+0.5}_{-0.6}$ $\rm km\,s^{-1}$, and $-13.2^{+0.6}_{-0.7}$ $\rm
km\,s^{-1}$ for HR 8799 b, c, and d, Wang et al. 2021a) slightly better than
the Ruffio et al. (2021) averages, indicating that this type of jitter term
may be a useful model of the dataset’s noise properties.

Another by-product of our MCMC analysis is an estimate of the radial velocity
“jitter,” or excess radial velocity scatter in each dataset. We measure an
excess radial velocity uncertainty of $2.1^{+1.3}_{-1.0}$ $\rm km\,s^{-1}$,
$0.7^{+1.0}_{-0.5}$ $\rm km\,s^{-1}$, and $1.3^{+1.9}_{-0.9}$ $\rm km\,s^{-1}$
for the HR 8799 b, c, and d observations respectively. For both HR 8799 c and
d the most probable value of the jitter value is 0, but the HR 8799 b dataset
shows a preference for non-zero jitter. This likely because HR 8799 b is the
faintest of the HR 8799 planets, and therefore is the most susceptible to low-
level residual systematic errors in the spectral extraction or the telluric
and starlight correction.

##  4 Discussion

It is not particularly surprising that the HR 8799 planets do not have
detectable short-period Jovian moons or binary-planet companions. Among the
thousands of known exoplanets, there are still no securely detected exomoons
or binary planets, even though companions the mass and size of Jupiter should
have been detectable around most Kepler target stars (Lewis et al., 2015).
Nevertheless, our constraints ruling out massive edge-on exomoons and binary
planets around HR 8799 are novel. Almost all published limits on the presence
of exomoons (e.g. Kipping et al., 2014, 2015) come from transiting exoplanet
systems, which due to selection effects orbit close to their host stars
(within 1-2 AU). At these short orbital periods, the planets’ Hill spheres are
smaller, limiting the parameter space in which moons can survive. Indeed, only
three of the hundreds of known moons in the solar system are found within 2 AU
of the Sun. So far, only Lazzoni et al. (2020) have published constraints on
moons orbiting planets at similar orbital radii to the HR 8799 planets, and
their direct imaging constraints only probe moons orbiting far from their host
planets. Our constraints are first probes of moons and binary planets in a new
parameter space (planet semimajor axis $\gtrsim$ 10 AU and moon orbital period
$\lesssim$ 10 days).

One limitation of the radial velocity method of exomoon (and exoplanet)
detection is that it is most sensitive to companions with close to edge-on
orbital inclinations. The radial velocity signal induced by a companion is
proportional to the sine of its orbital inclination, and generally it is not
possible to measure a companion’s true mass without an independent inclination
measurement. It is therefore possible that massive exomoons or binary planets
are hiding around the HR 8799 planets in nearly face-on orbits. The HR 8799
planets orbit their host star with orbital inclinations of about 20 degrees;
if moons or binary planets in the system orbit in the same plane, their radial
velocity signals would be attenuated by about a factor of 3. However, such
massive moons and binary planets likely must have formed by dynamical capture,
and would not necessarily orbit in the same plane as the planetary orbits.

Another promising way to detect exomoons around directly imaged planets is by
measuring the astrometric wobble of the planet due to an orbiting satellite
(Agol et al., 2015). Recently, the GRAVITY instrument on the Very Large
Telescope Interferometer (VLTI) has achieved extremely precise astrometric
measurements (uncertainties ranging from 20-150 microarcseconds) of directly
imaged planets (Gravity Collaboration et al., 2020; Lagrange et al., 2020;
Lacour et al., 2021), including HR 8799 e (Gravity Collaboration et al.,
2019). The astrometric signal induced by a Jupiter-mass companion in a 10 day
orbit around one of the HR 8799 planets would be about 60 microarcseconds,
roughly 3 times larger than the smallest single-point GRAVITY astrometric
uncertainties. Since the radial velocity semiamplitude of a Jupiter-mass
companion in a 10 day orbit around one of the HR 8799 planets ($\approx
2.4~{}\rm km\,s^{-1}$) is also about 3 times larger than the best smallest
single-point OSIRIS radial velocity uncertainty measurement, astrometry and
radial velocities are likely to be similarly sensitive to companions at 10 day
orbital periods. These techniques are highly complementary because a
companion’s radial velocity amplitude decreases with orbital period (radial
velocity semiamplitude $k\propto P^{-1/3}$), while the astromeric amplitude
increases with orbital period (astrometric amplitude $\theta\propto P^{2/3}$).
So while radial velocities will likely provide the tightest constraints on
close-in exomoons, astrometry could be a more sensitive probe of long-period
companions. Owing to its extremely high angular resolution, GRAVITY can
spatially resolve the Hill spheres of many nearby directly imaged exoplanets
(Wang et al., 2021b), including those around HR 8799. It therefore may also be
possible to use GRAVITY to directly detect luminous orbiting exomoons and
binary-planet companions on wide orbits around their primary planet. Detecting
exomoons orbiting directly imaged planets using the transit method (e.g.
Cabrera & Schneider, 2007) will be very challenging given the need for high-
duty cycle observations on 10 meter class telescopes (required to resolve the
individual planets). However, targeted observations following up on known
signals from either radial velocities or astrometry could be feasible.

Perhaps the most interesting result of this work is not the limits on exomoons
and binary planets themselves, but the fact that it was possible to constrain
the presence of companions in the first place. Detecting exoplanets around
stars with the radial velocity method requires many high signal-to-noise
observations with stabilized high-resolution spectrographs, but the
observational requirements for detecting exomoons are dramatically lessened
due to the low mass of the host planet. This is why it was possible for us to
place meaningful constraints on the presence of exomoons using a moderate
resolution spectrograph ($\lambda/\Delta\lambda\approx 4000$) that was
designed (Larkin et al., 2003) before the discovery of the first directly
imaged exoplanet (Chauvin et al., 2004). A new generation of specially
designed instruments for taking high resolution spectra of directly imaged
planets will dramatically improve the detection limits. Doppler uncertainties
scale with the full width at half maximum of spectral lines to the 1.5 power
(Lovis & Fischer, 2010), so increasing the spectral resolution from
$\lambda/\Delta\lambda\approx 4000$ (with OSIRIS) to
$\lambda/\Delta\lambda\approx 40000$ (with an instrument like KPIC) with
signal-to-noise held constant could decrease radial velocity uncertainties by
a factor of 30. Upgrades to existing instruments like KPIC/Keck (Jovanovic et
al., 2020), CRIRES+/SPHERE/HiRISE/VLT (Dorn et al., 2014; Vigan et al., 2018),
and IRD/SCExAO/REACH/Subaru (Kotani et al., 2020) and new instruments on 30
meter class telescopes like GMTNIRS/GMT (Jaffe et al., 2016) and HIRES/ELT
(Marconi et al., 2021) will yield high signal-to-noise spectra, potentially
enabling Doppler detections of rocky exomoons.

We thank Jason Eastman for helpful discussions. We thank the anonymous referee
for a quick and constructive report and for the idea that GRAVITY’s high
angular resolution could allow it to actually resolve luminous moons and
binary-planet companions. This research has made use of NASA’s Astrophysics
Data System, the NASA Exoplanet Archive, which is operated by the California
Institute of Technology, under contract with the National Aeronautics and
Space Administration under the Exoplanet Exploration Program, and the SIMBAD
database, operated at CDS, Strasbourg, France.

## References

  * Agol et al. (2015) Agol, E., Jansen, T., Lacy, B., Robinson, T. D., & Meadows, V. 2015, ApJ, 812, 5, doi: 10.1088/0004-637X/812/1/5
  * Akeson et al. (2013) Akeson, R. L., Chen, X., Ciardi, D., et al. 2013, PASP, 125, 989, doi: 10.1086/672273
  * Cabrera & Schneider (2007) Cabrera, J., & Schneider, J. 2007, A&A, 464, 1133, doi: 10.1051/0004-6361:20066111
  * Canup & Ward (2006) Canup, R. M., & Ward, W. R. 2006, Nature, 441, 834, doi: 10.1038/nature04860
  * Chauvin et al. (2004) Chauvin, G., Lagrange, A. M., Dumas, C., et al. 2004, A&A, 425, L29, doi: 10.1051/0004-6361:200400056
  * Dorn et al. (2014) Dorn, R. J., Anglada-Escude, G., Baade, D., et al. 2014, The Messenger, 156, 7
  * Dressing & Charbonneau (2015) Dressing, C. D., & Charbonneau, D. 2015, ApJ, 807, 45, doi: 10.1088/0004-637X/807/1/45
  * Eastman et al. (2019) Eastman, J. D., Rodriguez, J. E., Agol, E., et al. 2019, arXiv e-prints, arXiv:1907.09480.  https://arxiv.org/abs/1907.09480
  * Forgan (2017) Forgan, D. H. 2017, MNRAS, 470, 416, doi: 10.1093/mnras/stx1217
  * Gelman & Rubin (1992) Gelman, A., & Rubin, D. B. 1992, Statistical science, 457
  * Gravity Collaboration et al. (2019) Gravity Collaboration, Lacour, S., Nowak, M., et al. 2019, A&A, 623, L11, doi: 10.1051/0004-6361/201935253
  * Gravity Collaboration et al. (2020) Gravity Collaboration, Nowak, M., Lacour, S., et al. 2020, A&A, 633, A110, doi: 10.1051/0004-6361/201936898
  * Heller (2016) Heller, R. 2016, A&A, 588, A34, doi: 10.1051/0004-6361/201527496
  * Hook (2005) Hook, I. 2005, The science case for the European Extremely Large Telescope : the next step in mankind’s quest for the Universe
  * Hunter (2007) Hunter, J. D. 2007, Computing in Science and Engineering, 9, 90, doi: 10.1109/MCSE.2007.55
  * Hwang et al. (2018) Hwang, K. H., Udalski, A., Bond, I. A., et al. 2018, AJ, 155, 259, doi: 10.3847/1538-3881/aac2cb
  * Jaffe et al. (2016) Jaffe, D. T., Barnes, S., Brooks, C., et al. 2016, in Proc. SPIE, Vol. 9908, Ground-based and Airborne Instrumentation for Astronomy VI, 990821, doi: 10.1117/12.2232994
  * Jovanovic et al. (2020) Jovanovic, N., Calvin, B., Porter, M., et al. 2020, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 11447, Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, 114474U, doi: 10.1117/12.2563107
  * Kipping (2009) Kipping, D. M. 2009, MNRAS, 392, 181, doi: 10.1111/j.1365-2966.2008.13999.x
  * Kipping et al. (2014) Kipping, D. M., Nesvorný, D., Buchhave, L. A., et al. 2014, ApJ, 784, 28, doi: 10.1088/0004-637X/784/1/28
  * Kipping et al. (2015) Kipping, D. M., Schmitt, A. R., Huang, X., et al. 2015, ApJ, 813, 14, doi: 10.1088/0004-637X/813/1/14
  * Kotani et al. (2020) Kotani, T., Kawahara, H., Ishizuka, M., et al. 2020, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 11448, Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, 1144878, doi: 10.1117/12.2561755
  * Kreidberg et al. (2019) Kreidberg, L., Luger, R., & Bedell, M. 2019, ApJ, 877, L15, doi: 10.3847/2041-8213/ab20c8
  * Lacour et al. (2021) Lacour, S., Wang, J. J., Rodet, L., et al. 2021, arXiv e-prints, arXiv:2109.10671.  https://arxiv.org/abs/2109.10671
  * Lagrange et al. (2020) Lagrange, A. M., Rubini, P., Nowak, M., et al. 2020, A&A, 642, A18, doi: 10.1051/0004-6361/202038823
  * Larkin et al. (2006) Larkin, J., Barczys, M., Krabbe, A., et al. 2006, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 6269, Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, ed. I. S. McLean & M. Iye, 62691A, doi: 10.1117/12.672061
  * Larkin et al. (2003) Larkin, J. E., Quirrenbach, A., Krabbe, A., et al. 2003, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 4841, Instrument Design and Performance for Optical/Infrared Ground-based Telescopes, ed. M. Iye & A. F. M. Moorwood, 1600–1610, doi: 10.1117/12.461775
  * Lazzoni et al. (2020) Lazzoni, C., Zurlo, A., Desidera, S., et al. 2020, A&A, 641, A131, doi: 10.1051/0004-6361/201937290
  * Lewis et al. (2015) Lewis, K. M., Ochiai, H., Nagasawa, M., & Ida, S. 2015, ApJ, 805, 27, doi: 10.1088/0004-637X/805/1/27
  * Lewis et al. (2008) Lewis, K. M., Sackett, P. D., & Mardling, R. A. 2008, ApJ, 685, L153, doi: 10.1086/592743
  * Li & Batygin (2014) Li, G., & Batygin, K. 2014, ApJ, 790, 69, doi: 10.1088/0004-637X/790/1/69
  * Lillo-Box et al. (2018) Lillo-Box, J., Kipping, D., Rebollido, I., et al. 2018, ArXiv e-prints.  https://arxiv.org/abs/1803.04010
  * Limbach et al. (2021) Limbach, M. A., Vos, J. M., Winn, J. N., et al. 2021, ApJ, 918, L25, doi: 10.3847/2041-8213/ac1e2d
  * Lomb (1976) Lomb, N. R. 1976, Ap&SS, 39, 447, doi: 10.1007/BF00648343
  * Lovis & Fischer (2010) Lovis, C., & Fischer, D. 2010, Radial Velocity Techniques for Exoplanets, ed. S. Seager, 27–53
  * Marconi et al. (2021) Marconi, A., Abreu, M., Adibekyan, V., et al. 2021, The Messenger, 182, 27, doi: 10.18727/0722-6691/5219
  * Marois et al. (2008) Marois, C., Macintosh, B., Barman, T., et al. 2008, Science, 322, 1348, doi: 10.1126/science.1166585
  * Marois et al. (2010) Marois, C., Zuckerman, B., Konopacky, Q. M., Macintosh, B., & Barman, T. 2010, Nature, 468, 1080, doi: 10.1038/nature09684
  * Mieda et al. (2014) Mieda, E., Wright, S. A., Larkin, J. E., et al. 2014, PASP, 126, 250, doi: 10.1086/675784
  * Ochiai et al. (2014) Ochiai, H., Nagasawa, M., & Ida, S. 2014, ApJ, 790, 92, doi: 10.1088/0004-637X/790/2/92
  * Oliphant (2006) Oliphant, T. E. 2006, A guide to NumPy, USA: Trelgol Publishing
  * Penny et al. (2019) Penny, M. T., Gaudi, B. S., Kerins, E., et al. 2019, ApJS, 241, 3, doi: 10.3847/1538-4365/aafb69
  * Reynolds et al. (1983) Reynolds, R. T., Squyres, S. W., Colburn, D. S., & McKay, C. P. 1983, Icarus, 56, 246, doi: 10.1016/0019-1035(83)90037-4
  * Ruffio et al. (2021) Ruffio, J.-B., Konopacky, Q. M., Barman, T., et al. 2021, arXiv e-prints, arXiv:2109.07614.  https://arxiv.org/abs/2109.07614
  * Sartoretti & Schneider (1999) Sartoretti, P., & Schneider, J. 1999, A&AS, 134, 553, doi: 10.1051/aas:1999148
  * Scargle (1982) Scargle, J. D. 1982, ApJ, 263, 835, doi: 10.1086/160554
  * Sengupta & Marley (2016) Sengupta, S., & Marley, M. S. 2016, ApJ, 824, 76, doi: 10.3847/0004-637X/824/2/76
  * Spergel et al. (2013) Spergel, D., Gehrels, N., Breckinridge, J., et al. 2013, arXiv e-prints, arXiv:1305.5422.  https://arxiv.org/abs/1305.5422
  * Stern & Levison (2002) Stern, S. A., & Levison, H. F. 2002, Highlights of Astronomy, 12, 205
  * Swift et al. (2013) Swift, J. J., Johnson, J. A., Morton, T. D., et al. 2013, ApJ, 764, 105, doi: 10.1088/0004-637X/764/1/105
  * Teachey & Kipping (2018) Teachey, A., & Kipping, D. M. 2018, Science Advances, 4, eaav1784, doi: 10.1126/sciadv.aav1784
  * Teachey et al. (2018) Teachey, A., Kipping, D. M., & Schmitt, A. R. 2018, AJ, 155, 36, doi: 10.3847/1538-3881/aa93f2
  * Ter Braak (2006) Ter Braak, C. J. 2006, Statistics and Computing, 16, 239
  * Vanderburg (2021) Vanderburg, A. 2021, avanderburg/edmcmc: v1.0.0, v1.0.0, Zenodo, doi: 10.5281/zenodo.5599854
  * Vanderburg et al. (2018) Vanderburg, A., Rappaport, S. A., & Mayo, A. W. 2018, AJ, 156, 184, doi: 10.3847/1538-3881/aae0fc
  * Vigan et al. (2018) Vigan, A., Otten, G. P. P. L., Muslimov, E., et al. 2018, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 10702, Ground-based and Airborne Instrumentation for Astronomy VII, ed. C. J. Evans, L. Simard, & H. Takami, 1070236, doi: 10.1117/12.2313681
  * Virtanen et al. (2020) Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nature Methods, 17, 261, doi: 10.1038/s41592-019-0686-2
  * Wang et al. (2018a) Wang, J., Mawet, D., Fortney, J. J., et al. 2018a, AJ, 156, 272, doi: 10.3847/1538-3881/aae47b
  * Wang et al. (2018b) Wang, J. J., Graham, J. R., Dawson, R., et al. 2018b, AJ, 156, 192, doi: 10.3847/1538-3881/aae150
  * Wang et al. (2021a) Wang, J. J., Ruffio, J.-B., Morris, E., et al. 2021a, arXiv e-prints, arXiv:2107.06949.  https://arxiv.org/abs/2107.06949
  * Wang et al. (2021b) Wang, J. J., Vigan, A., Lacour, S., et al. 2021b, AJ, 161, 148, doi: 10.3847/1538-3881/abdb2d