# Realizing Topological Superconductivity in Tunable Bose-Fermi Mixtures with
Transition Metal Dichalcogenide Heterostructures

Caterina Zerba  Technical University of Munich, TUM School of Natural
Sciences, Physics Department, 85748 Garching, Germany  Munich Center for
Quantum Science and Technology (MCQST), Schellingstr. 4, 80799 München,
Germany  Clemens Kuhlenkamp  Technical University of Munich, TUM School of
Natural Sciences, Physics Department, 85748 Garching, Germany  Munich Center
for Quantum Science and Technology (MCQST), Schellingstr. 4, 80799 München,
Germany  Institute for Quantum Electronics, ETH Zürich, CH-8093 Zürich,
Switzerland  Ataç Imamoğlu  Institute for Quantum Electronics, ETH Zürich,
CH-8093 Zürich, Switzerland  Michael Knap  Technical University of Munich, TUM
School of Natural Sciences, Physics Department, 85748 Garching, Germany
Munich Center for Quantum Science and Technology (MCQST), Schellingstr. 4,
80799 München, Germany

(August 28, 2024)

###### Abstract

Heterostructures of two-dimensional transition metal dichalcogenides (TMDs)
are emerging as a promising platform for investigating exotic correlated
states of matter. Here, we propose to engineer Bose-Fermi mixtures in these
systems by coupling inter-layer excitons to doped charges in a trilayer
structure. Their interactions are determined by the inter-layer trion, whose
spin-selective nature allows excitons to mediate an attractive interaction
between charge carriers of only one spin species. Remarkably, we find that
this causes the system to become unstable to topological p+ip
superconductivity at low temperatures. We then demonstrate a general mechanism
to develop and control this unconventional state by tuning the trion binding
energy using a solid-state Feshbach resonance.

Figure 1: Tunable Bose-Fermi mixtures. (a) We consider a TMD trilayer
heterostructure. The relevant constituents are the charge carriers (holes),
that coherently tunnel between the upper and the middle layer with strength
$t$, an interlayer exciton condensate, and the trion bound state. (b)
Effective attractive interactions between charge carriers are mediated by the
virtual decay of a trion into a hole and a condensate excitation, which then
scatters with another hole to create again a trion. (c) The holes in the trion
form a singlet. When the exciton condensate is created with a certain spin
polarization, for example by circularly polarized light, attractive
interactions are only mediated between holes of opposite spin. (d) The
superconducting state is controlled by a solid-state Feshbach resonance.
Effective interactions are enhanced by tuning the threshold of the open
channel close to the trion bound state in the closed channel by applying an
external electric field $E_{z}\propto\Delta$.

Introduction.—Heterostructures of two-dimensional semiconductors are
particularly promising platforms to realize exotic phases of matter, such as
correlated insulators, Mott-Wigner states, and superconductors [1, 2, 3, 4, 5,
6, 7, 8, 9, 10, 11, 12, 13, 14]. Due to their two-dimensional nature and high
effective electron and hole masses, transition metal dichalcogenides (TMDs) in
particular can host strongly correlated Bose-Fermi mixtures with fermions
realized by doped charges and composite bosons by tightly bound excitons.
Doped holes and inter-layer excitons interact by forming bound molecular
states called trions. The bound-state formation in Bose-Fermi mixtures gives
rise to a particularly rich phase diagram, potentially featuring among others
coexisting liquids, unconventional superfluids, supersolids, and quantum phase
transtions between them [15, 16, 17, 18, 19, 20, 21, 22, 23, 24].

In multi-layer semiconductor systems scattering between excitons and doped
charges can be controlled and strongly enhanced when the energy of the trion
is tuned by an external electric field, which is understood as a solid state
Feshbach resonance [25, 26]. Therefore, two-dimensional materials open new
opportunities to realize and study Bose-Fermi mixtures in challenging low-
temperature regimes.

A natural route to search for unconventional superconducting states in TMDs
[27, 28, 29, 30, 31, 32, 33, 34], such as p-wave superconductivity which hosts
exotic Majorana edge modes that are pertinent for topological quantum
computing [35, 36, 37], are thus provided by Bose-Fermi mixtures, as they
allow for boson-mediated attractive interactions among charge carriers.

However, a controlled realization of p-wave superconductivity remains
challenging, as generically s-wave pairing is favored over higher angular
momentum states.

Here, we propose a general mechanism to develop and enhance unconventional
p-wave superconductivity by tuning Bose-Fermi mixtures of interlayer excitons
and holes in a TMD heterostructure. Superconductivity is achieved by
hybridizing trions and holes, while the pairing is mediated by excitations of
an exciton condensate; as illustrated in Fig. 1. Due to the spin-selective
nature of the trion, only carriers of a single spin species interact
attractively, which effectively restricts our model to spinless charge
carriers. Then, the Pauli principle forces the Cooper pairs to carry odd
angular momentum. Our analysis shows that the superconducting state is
topological with px+$i$py order. Moreover, the superconducting transition
temperature is highly tunable and controllable by an external electric field
via a Feshbach resonance.

Experimental setting.—We consider a trilayer TMD heterostructure as shown in
Fig. 1 (a). Interlayer excitons can be optically injected in the bottom two
layers at high densities ($\simeq 10^{12}$ cm-2) due to their long lifetime
[38, 39]. They are strongly bound and essentially act as composite bosons. We
assume that excitons are condensed [40, 41]. Fermions are introduced in the
form of doped holes, which tunnel coherently between the top and middle layer.
Hexagonal Boron Nitride (hBN) isolates the layers and, to a limited extent,
controls the tunneling of charge carriers. Inter-layer excitons and charge
carriers can bind to a molecular trion state, which allows for controlling the
boson-fermion scattering. Tunable effective interactions between holes are
then mediated by the exciton condensate; see Fig. 1 (b).

We assume that the exciton condensate is optically injected in one valley by
circularly polarized light, and thus through spin-momentum locking in one spin
state, say the $K$-valley with spin-$\uparrow$; Fig. 1 (c). Due to Pauli
blocking, interlayer trion formation is only possible when the holes of the
trion, that are both in the same layer, are in a spin singlet, as the triplet
is generally unbound [42, 43, 44]. For small temperatures compared to the
trion bound state energy, two charge carriers that interact with the same
exciton to form a trion are thus necessarily in the same spin state.
Therefore, even if the exciton polarization is not perfect, no attractive
interactions are mediated between opposite-spin holes, which allows us to
consider spin-polarized fermions [42, 43, 44, 40, 45].

The scattering process between holes and excitons can be controlled by a solid
state Feshbach resonance [25, 26], where open and closed channels are tuned
relatively to each other by an external electric field $E_{z}$; Fig. 1 (d).
Strong interactions are obtained when the open channel and the trion bound
state are approximately in resonance.

The system of holes $\hat{c}$, trions $\hat{m}$ (also referred to as
molecules), and excitons $\hat{x}$ is described by the following effective
Hamiltonian (see supplemental materials [46]):

$\hat{H}=\hat{H}_{\text{holes}}+\hat{H}_{\text{trions}}+\hat{H}_{\text{ex}}+g\sum_{\mathbf{k},\mathbf{k}^{\prime}}[\hat{m}^{\dagger}_{\mathbf{k}}\hat{c}_{\mathbf{k}^{\prime}}\hat{x}_{\mathbf{k}-\mathbf{k}^{\prime}}+\text{h.c.}],$
(1)

where
$\hat{H}_{\text{holes}}=\sum_{\mathbf{k}}\frac{k^{2}}{2m}\hat{c}^{\dagger}_{\mathbf{k}}\hat{c}_{\mathbf{k}}$
is the kinetic energy of holes of mass $m$ and
$\hat{H}_{\text{trions}}=\sum_{\mathbf{k}}\big{(}\frac{k^{2}}{2M}+\Delta\big{)}\hat{m}^{\dagger}_{\mathbf{k}}\hat{m}_{\mathbf{k}}$
is the kinetic energy of trions of mass $M$ with $\Delta=edE_{z}$ being the
energy difference of a charge in the upper and the middle layer separated by
distance $d$ that is tunable by the electric field $E_{z}$. The interlayer
excitons are governed by
$\hat{H}_{\text{ex}}=\sum_{\mathbf{k}}\frac{k^{2}}{2m_{x}}\hat{x}^{\dagger}_{\mathbf{k}}\hat{x}_{\mathbf{k}}+\frac{U}{2}\sum_{\mathbf{k},\mathbf{p},\mathbf{q}}\hat{x}_{\mathbf{k}+\mathbf{q}}^{\dagger}\hat{x}_{\mathbf{p}-\mathbf{q}}^{\dagger}\hat{x}_{\mathbf{k}}\hat{x}_{\mathbf{p}}$,
where $m_{x}$ is the exciton mass and $U$ the strength of the exciton
interactions that depends on the dielectric material between the two lowest
layers and their distance [38]. Since excitons are injected with low momenta
and long-range electron-hole exchange is vanishingly small, due to spatial
separation of the electron and the hole, a quadratic dispersion is justified.
The last term in the Hamiltonian (1) describes the decay of the molecule into
an exciton and a hole of strength $g$ that mediates the effective attractive
interactions.

We consider the regime in which the binding energy of the exciton is larger
than the one of the trion. Then, the exciton’s internal structure does not
influence the properties of the trion. We are thus allowed to describe the
excitons as a Bose gas [47] and assume weak contact interactions. When the
excitons condense [48, 49], we can expand
$\hat{x}_{\mathbf{r}}=\sqrt{n_{0}}+\frac{1}{\sqrt{V}}\sum_{\mathbf{k}\not=0}e^{i\mathbf{k}\mathbf{r}}\hat{x}_{\mathbf{k}}$,
where $n_{0}$ is the macroscopic condensate density. Performing a Bogoliubov
analysis by expanding the purely excitonic part of the Hamiltonian, we obtain

$\displaystyle\hat{H}$
$\displaystyle=\hat{H}_{\text{holes}}+\hat{H}_{\text{trions}}+\sum_{\mathbf{k}\not=0}e_{x}(k)\hat{b}^{\dagger}_{\mathbf{k}}\hat{b}_{\mathbf{k}}$
$\displaystyle+g\sqrt{n_{0}}\sum_{\mathbf{k}}\big{[}\hat{m}^{\dagger}_{\mathbf{k}}\hat{c}_{\mathbf{k}}+\text{h.c.}\big{]}+g\sum_{\mathbf{k}\not=\mathbf{k}^{\prime}}\big{[}\hat{m}^{\dagger}_{\mathbf{k}}\hat{c}_{\mathbf{k}-\mathbf{k}^{\prime}}\hat{b}_{\mathbf{k^{\prime}}}+\text{h.c.}\big{]},$
(2)

where $\hat{b}$ are the condensate excitations with energy
$e_{x}(k)=\sqrt{\frac{k^{2}}{2m_{x}}\big{(}\frac{k^{2}}{2m_{x}}+2Un_{0}\big{)}}$.
The condensate hybridizes the trion and charge bands with strength
$g\sqrt{n_{0}}$. The three-body interaction contains sound excitations of the
condensate $\hat{b}$, which are reminiscent of phonon mediated interactions.

The effective trion-exciton-hole coupling $g$ is fixed by the tunnel coupling
$t$ of the electrons between the two layers as
$g\simeq\sqrt{\frac{2\pi}{\mu_{\text{red}}|E_{T,0}|}}t$, where we estimate the
bare trion binding energy $E_{T,0}\simeq 6$ meV to be lower than the one
reported in Ref. [38] due to the additional hBN layer between the TMDs, and
$\mu_{\text{red}}$ is the reduced mass of the hole and the exciton; see
supplemental materials [46]. In our calculations, we assume the tunneling
between the two layers to be $t\simeq 2$ meV [26]. We typically consider
exciton condensate densities of $n_{0}\sim 10^{12}$ cm-2 and estimate the
interactions $U$ between the excitons using a plate capacitor formula as
$Un_{0}\simeq 10$ meV for $n_{0}\sim 10^{12}$ cm-2 [38]. For the long-range
interactions between excitons, decaying as $\sim r^{-3}$, bosons undergo a BKT
transition [50]. To simplify the situation, we consider a contact
interactions, that, however, capture the linear spectrum at small momenta
excitations. The interactions determine the phase fluctuations and a
condensate description of the problem is justified for scales of the sample
smaller than the typical coherence length and low temperatures [51, 52].

Effective attractive interactions.—To determine the interactions between doped
charges we diagonalize the free fermionic Hamiltonian [16, 53]:
$\hat{H}_{0}=\hat{H}_{\text{holes}}+\hat{H}_{\text{trions}}+g\sqrt{n_{0}}\sum_{\mathbf{k}}\big{[}\hat{m}^{\dagger}_{\mathbf{k}}\hat{c}_{\mathbf{k}}+\text{h.c.}\big{]}$
to obtain the dressed fermions
$\hat{f}_{\mathbf{k}}=\alpha_{\mathbf{k}}\hat{m}_{\mathbf{k}}-\beta_{\mathbf{k}}\hat{c}_{\mathbf{k}}$
and
$\hat{e}_{\mathbf{k}}=\alpha_{\mathbf{k}}\hat{m}_{\mathbf{k}}+\beta_{\mathbf{k}}\hat{c}_{\mathbf{k}}$,
with eigenenergies
$\lambda_{e,f}(\mathbf{k})=\frac{k^{2}}{4}\frac{M+m}{Mm}+\frac{\Delta}{2}\pm\sqrt{\left(\frac{k^{2}}{4}\frac{M-m}{Mm}-\frac{\Delta}{2}\right)^{2}+g^{2}n_{0}}$,
which define the upper and lower band; see Fig. 2. The gap between the two
bands, $\sqrt{\Delta^{2}+4g^{2}n_{0}}$, is finite for all values of the
electric field and has the minimal value of $2g\sqrt{n_{0}}$ for $\Delta=0$.
Assuming the Fermi energy to be lower than the energy difference between the
two bands and the temperature to be smaller than the gap, we restrict our
analysis the lowest filled band. Then, depending on the sign of the electric
field $\Delta$, charges will mainly have trion ($\Delta<0$) or hole
($\Delta>0$) character. The resulting Hamiltonian of fermions $\hat{f}$ and
condensate excitations $\hat{b}$ is then
$\hat{H}=\sum_{\mathbf{k}}[\lambda_{f}(\mathbf{k})-\mu]\hat{f}_{\mathbf{k}}^{\dagger}\hat{f}_{\mathbf{k}}-g\sum_{\mathbf{k}\not=\mathbf{k}^{\prime}}\big{[}\alpha_{\mathbf{k}}\beta_{\mathbf{k}^{\prime}}\hat{f}^{\dagger}_{\mathbf{k}}\hat{f}_{\mathbf{k^{\prime}}}\hat{b}_{\mathbf{k}-\mathbf{k}^{\prime}}+\text{h.c.}\big{]}$.

Figure 2: Hybridized trion and hole bands. Hybridized bands $\lambda_{e/f}(k)$
(solid lines) and bare bands (dotted lines) as a function of momentum $k$ in
units of the inverse exciton condensate healing length $\xi^{-1}$. The Fermi
energy (horizontal lines) is chosen to be below the upper band. (a) For
$\Delta<0$, the heavier band is filled, implying that doped charges mainly
have trion character. (b) For $\Delta>0$, the lighter band is filled and the
doped charges will mainly have hole character.

We integrate out the condensate interactions using Wegner flow equations [54],
which gives us an effective interaction between charge carries mediated by the
excitons. In contrast to a Schrieffer-Wolff transformation this leads to a
divergence free potential [55, 56, 54, 22]; see supplemental material [46].
The resulting effective Hamiltonian is then:

$\hat{H}^{\text{eff}}=\sum_{\mathbf{k}}[\lambda_{f}(\mathbf{k})-\mu]\hat{f}^{\dagger}_{\mathbf{k}}\hat{f}_{\mathbf{k}}+\sum_{\mathbf{k},\mathbf{q}}V_{k,q}\hat{f}^{\dagger}_{\mathbf{k}-\mathbf{q}}\hat{f}^{\dagger}_{-\mathbf{k}+\mathbf{q}}\hat{f}_{-\mathbf{k}}\hat{f}_{\mathbf{k}},$
(3)

with the interactions
$V_{k,q}=-\frac{M_{\mathbf{k},\mathbf{q}}M^{*}_{-\mathbf{k}+\mathbf{q},\mathbf{q}}e_{x}(\mathbf{q})}{[\lambda_{f}(\mathbf{k})-\lambda_{f}(\mathbf{k}-\mathbf{q})]^{2}+e_{x}(\mathbf{q})^{2}}$
and $M_{k,q}=-g\alpha_{\mathbf{k}}\beta_{\mathbf{k}-\mathbf{q}}$. The
effective interaction is thus attractive and directly proportional to the
hybridization of molecule and hole bands.

Ginzburg-Landau analysis.—To analyze the superconducting order, we identify
the pairing instability by looking for unstable solutions in the temporal
evolution of the superconducting order parameter [57, 58, 59]. This method is
based on the Ginzburg-Landau free energy which to leading-order reads

$\mathcal{F}^{\text{GL}}=\sum_{\mathbf{k},\mathbf{p}}\psi_{\mathbf{k}}^{*}F^{\text{GL}}_{\mathbf{k},\mathbf{p}}\psi_{\mathbf{p}}+\mathcal{O}(\psi^{4}).$
(4)

In the Ginzburg-Landau formalism, the eigenvalues of the matrix
$F^{\text{GL}}$ are proportional to $\propto T-T_{c}$. Hence, the dominant
eigenvalues change their sign when driving the system into the ordered phase.
The associated eigenvectors then determine the structure of the instable mode.
To study the instability, we write the time evolution of $\psi_{\mathbf{k}}$
by means of time-dependent Ginzburg Landau equations [60, 61, 62]

$\frac{d}{dt}{\psi}_{\mathbf{k}}(t)=-\frac{\delta\mathcal{F}^{GL}}{\delta\psi^{*}_{\mathbf{k}}}=-\sum_{{\mathbf{p}}>0}F^{GL}_{\mathbf{k},\mathbf{p}}\psi_{\mathbf{p}}.$
(5)

Performing a mean-field analysis for
$\psi_{\mathbf{k}}=\langle\hat{f}_{\mathbf{k}}\hat{f}_{-\mathbf{k}}\rangle$,
we obtain the corresponding evolution equation in the microscopic formalism

$\frac{d}{dt}\langle\hat{f}_{\mathbf{k}}\hat{f}_{-\mathbf{k}}\rangle=-\frac{i}{\hbar}\sum_{\mathbf{p}>0}\big{[}A_{\mathbf{k},\mathbf{p}}+B_{\mathbf{k},\mathbf{p}}\big{]}\langle\hat{f}_{\mathbf{p}}\hat{f}_{-\mathbf{p}}\rangle,$
(6)

with contributions from the non-interacting part of the Hamiltonian
$A_{\mathbf{k},\mathbf{p}}=2(\lambda_{f}(\mathbf{k})-\mu)$ and the two-body
interaction
$B_{\mathbf{k},\mathbf{p}}=-(2n_{\mathbf{k}}-1)\big{[}V_{\mathbf{k},\mathbf{k}-\mathbf{p}}-V_{\mathbf{k},\mathbf{k}+\mathbf{p}}-V_{\mathbf{-k},\mathbf{-k}-\mathbf{p}}+V_{-\mathbf{k},-\mathbf{k}+\mathbf{p}}\big{]}$.
Here, both $A$ and $B$ are real matrices. Comparing Eq. 5 and Eq. 6, we obtain
the relation $-i\hbar F^{GL}=C=A+B$; the superconducting gap order parameter
is related to $\psi_{\mathbf{k}}$ by
$\Gamma_{\mathbf{k}}=\sum_{\mathbf{p}}V_{\mathbf{k},\mathbf{p}-\mathbf{k}}\psi_{\mathbf{p}}$.
The set of equations can be written in a matrix form as
$i\hbar\frac{d}{dt}|{\psi}\rangle=C|\psi\rangle$. We solve the differential
equation by finding the eigenvectors $|\psi\rangle$ and eigenvalues $E_{\psi}$
of $-i\hbar F^{\text{GL}}$. The critical point is then obtained by searching
for the temperature at which the imaginary part of some of the $E_{\psi}$
turns positive, indicating the onset of the unstable mode. From the associated
eigenvector, we infer the momentum-space structure of the instability. By
expanding the order parameter in its angular momentum contribution we find
numerically that p-wave is the only relevant component; see supplemental
material [46].

The precise p-wave channel that drives the superconducting instability can be
determined separately by minimization of the free energy. To this end, we
compare the free energy for the chiral (px+ipy) and nematic (px or py,
respectively) phases. We find that for the instabilities determined as above,
the free energy is indeed minimized for the topological px+ipy superconductor,
see supplementary materials [46], which hosts chiral Majorana edge modes. The
Ginzburg-Landau analysis identifies the leading instability, and hence
captures the phases near the transition temperature $T_{c}$. At even lower
temperatures, non-linearities may stabilize other competing phases.

Figure 3: Critical temperature. We evaluate the critical temperature $T_{c}$
as a function of the electric field $\Delta$ for (a) different Fermi energies
at condensate density $n_{0}=2$ $10^{12}$ cm-2 and for (b) different
condensate densities and Fermi energy $E_{F}=2$ meV. By tuning the electric
field $\Delta$ we find that superconductivity can be resonantly enhanced.

Tuning topological superconductivity.—We have seen that the attractive
interactions sensitively depend on the electric field via the hybridization of
molecule and hole bands. This suggests that the superconducting order is
controllable by a solid state Feshbach resonance. When adjusting $\Delta$ only
over a couple of meV, we find that the superconducting transition temperature
is tuned by more than an order of magnitude; see Fig. 3. The shape of the
critical temperature as a function of $\Delta$ qualitatively follows the
numerically evaluated effective potential in the vicinity of the Fermi
surface. By tuning $\Delta$, the hybridization between the two bands can be
increased at the relevant energy scales and consequently enhance the critical
temperature $T_{c}$. The largest $T_{c}$ is obtained when the hybridization of
the bands is largest, which arises approximately when $\Delta\sim 3$meV. When
increasing the Fermi energy, we find from our numerical analysis that the
maximum of the critical temperature increases; Fig. 3 (a).This can be
understood by the hybridization between the light hole and heavy trion band.
When increasing the Fermi energy the band becomes more trion-like for
$\Delta\sim|E_{T}|>0$. Hence, the hybridized charge carriers $\hat{f}$ attain
a larger effective mass when increasing the Fermi energy. This in turn leads
to a larger critical temperature.

Varying the condensate density $n_{0}$ at a fixed Fermi energy changes both
the shape of the hybridized bands $\lambda_{e/f}(\mathbf{k})$ and the
condensate healing length $\xi=1/\sqrt{2m_{x}Un_{0}}$. The healing length
determines the length scale at which the condensate effectively responds to a
perturbation. For larger healing length, we thus expect that the effective
attractive interactions are mediated more efficiently. This is confirmed
numerically, where we find that for decreasing condensate density $n_{0}$ and
hence increasing healing length $\xi$, superconductivity is more stable; Fig.
3 (b). This is further supported by evaluating the effective attractive
interaction for scattering near the Fermi surface which for small fields is
given by $V_{k,q}\sim\sqrt{\frac{m_{x}}{Un_{0}q^{2}}}$. The attractive
potential near the Fermi surface thus increases when the condensate density
decreases, as expected from the scaling of the healing length $\xi$.

For small values of the Fermi energy or for small values of the condensate
density, the critical temperature does not follow qualitatively the effective
potential in the vicinity of the Fermi surface anymore but develops additional
structure in the form of a local minimum around zero electric field; Fig. 3.
These features can be attributed to the finite width of the effective
scattering potential, which decreases for $\Delta\sim 0$. Such a local minimum
similarly appears for box-shaped potentials with small widths.

Discussion and Outlook.—In this work, we have demonstrated the potential of
two-dimensional TMD heterostructures as a tunable platform for realizing
correlated Bose-Fermi mixtures of interlayer excitons and doped charges. We
have shown that at low temperatures the system stabilizes p-wave topological
superconductivity. We propose to utilize solid state Feshbach resonances to
modify the effective interaction strength between the constituents.
Experimentally, the superconducting order parameter can be extracted by
measuring the ac conductivity. In the superconducting phase the dissipative
response of the conductivity then only exhibits spectral weight at frequencies
above twice the gap. Within BCS theory, the critical temperature sets the
energy scale for the gap and this places our estimates of these frequencies
below $\sim 200$ GHz, which are in the experimentally accessible microwave
regime. First experimental signatures for the superconducting phase could also
be obtained from the single-particle gap, measured optically via the change of
the line shape of the top layer exciton [6]. Moreover, tunneling spectroscopy
may be able to measure the gapless chiral Majorana edge modes, which is a key
signature of the p+ip topological superconductor.

For future studies it will be interesting to explore competing instabilities
at low temperatures, such as nematic or density-wave phases, for which non-
linearities and the competition of order parameters need to be taken into
account. In this work, we consider the regime of high condensate densities,
such that a possible back-action of the charge carriers on the exciton
condensate is largely negligible. New exciting phenomena appear in degenerate
Bose-Fermi mixtures for low condensate densities, which are interesting to be
explored both theoretically and experimentally [63, 24, 64].

Note added.—Milczewski et al. explored a complementary scheme for inducing
s-wave superconductivity in Bose-Fermi mixtures realized with TMD
heterostructures [65].

Data and Code availability.—Data analysis and simulation codes are available
on Zenodo upon reasonable request [66].

Acknowledgements.—We thank L. Classen, E. Demler, J. Knolle, and R. Schmidt,
for insightful discussions. We acknowledge support from the Deutsche
Forschungsgemeinschaft (DFG, German Research Foundation) under Germany’s
Excellence Strategy–EXC–2111–390814868, DFG grants No. KN1254/1-2, KN1254/2-1,
and TRR 360 - 492547816 and from the European Research Council (ERC) under the
European Unions Horizon 2020 research and innovation programme (Grant
Agreement No. 851161), as well as the Munich Quantum Valley, which is
supported by the Bavarian state government with funds from the Hightech Agenda
Bayern Plus. The work of A.I. was supported by the Swiss National Science
Foundation (SNSF) under Grant Number 200021-204076.

## References

  * Cao  _et al._ [2018a] Y. Cao, V. Fatemi, A. Demir, S. Fang, S. L. Tomarken, J. Y. Luo, J. D. Sanchez-Yamagishi, K. Watanabe, T. Taniguchi, E. Kaxiras, R. C. Ashoori, and P. Jarillo-Herrero, Correlated insulator behaviour at half-filling in magic-angle graphene superlattices, Nature 556, 80 (2018a).
  * Cao  _et al._ [2018b] Y. Cao, V. Fatemi, S. Fang, K. Watanabe, T. Taniguchi, E. Kaxiras, and P. Jarillo-Herrero, Unconventional superconductivity in magic-angle graphene superlattices, Nature 556, 43 (2018b).
  * Lu  _et al._ [2019] X. Lu, P. Stepanov, W. Yang, M. Xie, M. A. Aamir, I. Das, C. Urgell, K. Watanabe, T. Taniguchi, G. Zhang, A. Bachtold, A. H. MacDonald, and D. K. Efetov, Superconductors, orbital magnets and correlated states in magic-angle bilayer graphene, Nature 574, 653 (2019).
  * Sharpe  _et al._ [2019] A. L. Sharpe, E. J. Fox, A. W. Barnard, J. Finney, K. Watanabe, T. Taniguchi, M. A. Kastner, and D. Goldhaber-Gordon, Emergent ferromagnetism near three-quarters filling in twisted bilayer graphene, Science 365, 605 (2019).
  * Yankowitz  _et al._ [2019] M. Yankowitz, S. Chen, H. Polshyn, Y. Zhang, K. Watanabe, T. Taniguchi, D. Graf, A. F. Young, and C. R. Dean, Tuning superconductivity in twisted bilayer graphene, Science 363, 1059 (2019).
  * Shimazaki  _et al._ [2020] Y. Shimazaki, I. Schwartz, K. Watanabe, T. Taniguchi, M. Kroner, and A. Imamoğlu, Strongly correlated electrons and hybrid excitons in a moiré heterostructure, Nature 580, 472 (2020).
  * Xu  _et al._ [2020] Y. Xu, S. Liu, D. A. Rhodes, K. Watanabe, T. Taniguchi, J. Hone, V. Elser, K. F. Mak, and J. Shan, Correlated insulating states at fractional fillings of moiré superlattices, Nature 587, 214 (2020).
  * Regan  _et al._ [2020] E. C. Regan, D. Wang, C. Jin, M. I. B. Utama, B. Gao, X. Wei, S. Zhao, W. Zhao, Z. Zhang, K. Yumigeta, M. Blei, J. D. Carlström, K. Watanabe, T. Taniguchi, S. Tongay, M. Crommie, A. Zettl, and F. Wang, Mott and generalized wigner crystal states in WSe2/WS2 moiré superlattices, Nature 579, 359 (2020).
  * Wang  _et al._ [2020] L. Wang, E.-M. Shih, A. Ghiotto, L. Xian, D. A. Rhodes, C. Tan, M. Claassen, D. M. Kennes, Y. Bai, B. Kim, K. Watanabe, T. Taniguchi, X. Zhu, J. Hone, A. Rubio, A. N. Pasupathy, and C. R. Dean, Correlated electronic phases in twisted bilayer transition metal dichalcogenides, Nature Materials 19, 861 (2020).
  * Tang  _et al._ [2020] Y. Tang, L. Li, T. Li, Y. Xu, S. Liu, K. Barmak, K. Watanabe, T. Taniguchi, A. H. MacDonald, J. Shan, and K. F. Mak, Simulation of hubbard model physics in WSe2/WS2 moiré superlattices, Nature 579, 353 (2020).
  * Shimazaki  _et al._ [2021] Y. Shimazaki, C. Kuhlenkamp, I. Schwartz, T. Smoleński, K. Watanabe, T. Taniguchi, M. Kroner, R. Schmidt, M. Knap, and A. Imamoğlu, Optical signatures of periodic charge distribution in a mott-like correlated insulator state, Phys. Rev. X 11, 021027 (2021).
  * Li  _et al._ [2021a] H. Li, S. Li, E. C. Regan, D. Wang, W. Zhao, S. Kahn, K. Yumigeta, M. Blei, T. Taniguchi, K. Watanabe, S. Tongay, A. Zettl, M. F. Crommie, and F. Wang, Imaging two-dimensional generalized wigner crystals, Nature 597, 650 (2021a).
  * Ma  _et al._ [2021] L. Ma, P. X. Nguyen, Z. Wang, Y. Zeng, K. Watanabe, T. Taniguchi, A. H. MacDonald, K. F. Mak, and J. Shan, Strongly correlated excitonic insulator in atomic double layers, Nature 598, 585 (2021).
  * Huang  _et al._ [2021] X. Huang, T. Wang, S. Miao, C. Wang, Z. Li, Z. Lian, T. Taniguchi, K. Watanabe, S. Okamoto, D. Xiao, S.-F. Shi, and Y.-T. Cui, Correlated insulating states at fractional fillings of the WS2/WSe2 moiré lattice, Nature Physics 17, 715 (2021).
  * Viverit  _et al._ [2000] L. Viverit, C. J. Pethick, and H. Smith, Zero-temperature phase diagram of binary boson-fermion mixtures, Phys. Rev. A 61, 053605 (2000).
  * Powell  _et al._ [2005] S. Powell, S. Sachdev, and H. P. Büchler, Depletion of the bose-einstein condensate in bose-fermi mixtures, Phys. Rev. B 72, 024534 (2005).
  * Büchler and Blatter [2003] H. P. Büchler and G. Blatter, Supersolid versus phase separation in atomic bose-fermi mixtures, Phys. Rev. Lett. 91, 130404 (2003).
  * Enss and Zwerger [2009] T. Enss and W. Zwerger, Superfluidity near phase separation in bose-fermi mixtures, The European Physical Journal B 68, 383 (2009).
  * Fratini and Pieri [2010] E. Fratini and P. Pieri, Pairing and condensation in a resonant bose-fermi mixture, Phys. Rev. A 81, 051605 (2010).
  * Ludwig  _et al._ [2011] D. Ludwig, S. Floerchinger, S. Moroz, and C. Wetterich, Quantum phase transition in bose-fermi mixtures, Phys. Rev. A 84, 033629 (2011).
  * Matuszewski  _et al._ [2012] M. Matuszewski, T. Taylor, and A. V. Kavokin, Exciton supersolidity in hybrid bose-fermi systems, Phys. Rev. Lett. 108, 060401 (2012).
  * Cotleţ  _et al._ [2016] O. Cotleţ, S. Zeytinoǧlu, M. Sigrist, E. Demler, and A. Imamoǧlu, Superconductivity and other collective phenomena in a hybrid bose-fermi mixture formed by a polariton condensate and an electron system in two dimensions, Phys. Rev. B 93, 054510 (2016).
  * Kinnunen  _et al._ [2018] J. J. Kinnunen, Z. Wu, and G. M. Bruun, Induced $p$-wave pairing in bose-fermi mixtures, Phys. Rev. Lett. 121, 253402 (2018).
  * Duda  _et al._ [2023] M. Duda, X.-Y. Chen, A. Schindewolf, R. Bause, J. von Milczewski, R. Schmidt, I. Bloch, and X.-Y. Luo, Transition from a polaronic condensate to a degenerate fermi gas of heteronuclear molecules, Nature Physics 19, 720 (2023).
  * Kuhlenkamp  _et al._ [2022] C. Kuhlenkamp, M. Knap, M. Wagner, R. Schmidt, and A. Imamoğlu, Tunable feshbach resonances and their spectral signatures in bilayer semiconductors, Phys. Rev. Lett. 129, 037401 (2022).
  * Schwartz  _et al._ [2021] I. Schwartz, Y. Shimazaki, C. Kuhlenkamp, K. Watanabe, T. Taniguchi, M. Kroner, and A. Imamoğlu, Electrically tunable feshbach resonances in twisted bilayer semiconductors, Science 374, 336 (2021).
  * Hsu  _et al._ [2017] Y.-T. Hsu, A. Vaezi, M. H. Fischer, and E.-A. Kim, Topological superconductivity in monolayer transition metal dichalcogenides, Nat. Commun. 8, 1 (2017).
  * Margalit  _et al._ [2021] G. Margalit, E. Berg, and Y. Oreg, Theory of multi-orbital topological superconductivity in transition metal dichalcogenides, Annals of Physics 435, 168561 (2021).
  * Lane and Zhu [2022] C. Lane and J.-X. Zhu, Identifying topological superconductivity in two-dimensional transition-metal dichalcogenides, Phys. Rev. Mater. 6, 094001 (2022).
  * Scherer  _et al._ [2022] M. M. Scherer, D. M. Kennes, and L. Classen, Chiral superconductivity with enhanced quantized Hall responses in moiré transition metal dichalcogenides, npj Quantum Mater. 7, 1 (2022).
  * Julku  _et al._ [2022] A. Julku, J. J. Kinnunen, A. Camacho-Guardian, and G. M. Bruun, Light-induced topological superconductivity in transition metal dichalcogenide monolayers, Phys. Rev. B 106, 134510 (2022).
  * Crépel  _et al._ [2023] V. Crépel, D. Guerci, J. Cano, J. H. Pixley, and A. Millis, Topological superconductivity in doped magnetic moiré semiconductors, Phys. Rev. Lett. 131, 056001 (2023).
  * Li  _et al._ [2021b] Y. W. Li, H. J. Zheng, Y. Q. Fang, D. Q. Zhang, Y. J. Chen, C. Chen, A. J. Liang, W. J. Shi, D. Pei, L. X. Xu, S. Liu, J. Pan, D. H. Lu, M. Hashimoto, A. Barinov, S. W. Jung, C. Cacho, M. X. Wang, Y. He, L. Fu, H. J. Zhang, F. Q. Huang, L. X. Yang, Z. K. Liu, and Y. L. Chen, Observation of topological superconductivity in a stoichiometric transition metal dichalcogenide 2M-WS2, Nat. Commun. 12, 1 (2021b).
  * Nayak  _et al._ [2021] A. K. Nayak, A. Steinbok, Y. Roet, J. Koo, G. Margalit, I. Feldman, A. Almoalem, A. Kanigel, G. A. Fiete, B. Yan, Y. Oreg, N. Avraham, and H. Beidenkopf, Evidence of topological boundary modes with topological nodal-point superconductivity, Nature Physics 17, 1413 (2021).
  * Qi and Zhang [2011] X.-L. Qi and S.-C. Zhang, Topological insulators and superconductors, Rev. Mod. Phys. 83, 1057 (2011).
  * Alicea [2012] J. Alicea, New directions in the pursuit of majorana fermions in solid state systems, Reports on Progress in Physics 75, 076501 (2012).
  * Kallin and Berlinsky [2016] C. Kallin and J. Berlinsky, Chiral superconductors, Reports on Progress in Physics 79, 054502 (2016).
  * Jauregui  _et al._ [2019] L. A. Jauregui, A. Y. Joe, K. Pistunova, D. S. Wild, A. A. High, Y. Zhou, G. Scuri, K. De Greve, A. Sushko, C.-H. Yu, T. Taniguchi, K. Watanabe, D. J. Needleman, M. D. Lukin, H. Park, and P. Kim, Electrical control of interlayer exciton dynamics in atomically thin heterostructures, Science 366, 870 (2019).
  * Wilson  _et al._ [2021] N. P. Wilson, W. Yao, J. Shan, and X. Xu, Excitons and emergent quantum phenomena in stacked 2d semiconductors, Nature 599, 383 (2021).
  * Wang  _et al._ [2019] Z. Wang, D. A. Rhodes, K. Watanabe, T. Taniguchi, J. C. Hone, J. Shan, and K. F. Mak, Evidence of high-temperature exciton condensation in two-dimensional atomic double layers, Nature 574, 76 (2019).
  * Troue  _et al._ [2023] M. Troue, J. Figueiredo, L. Sigl, C. Paspalides, M. Katzer, T. Taniguchi, K. Watanabe, M. Selig, A. Knorr, U. Wurstbauer, and A. W. Holleitner, Extended spatial coherence of interlayer excitons in ${\mathrm{mose}}_{2}/{\mathrm{wse}}_{2}$ heterobilayers, Phys. Rev. Lett. 131, 036902 (2023).
  * Witham  _et al._ [2018] O. Witham, R. J. Hunt, and N. D. Drummond, Stability of trions in coupled quantum wells modeled by two-dimensional bilayers, Phys. Rev. B 97, 075424 (2018).
  * Tempelaar and Berkelbach [2019] R. Tempelaar and T. C. Berkelbach, Many-body simulation of two-dimensional electronic spectroscopy of excitons and trions in monolayer transition metal dichalcogenides, Nat. Commun. 10, 1 (2019).
  * Dai and Fu [2023] D. D. Dai and L. Fu, Strong-coupling phases of trions and excitons in electron-hole bilayers at commensurate densities (2023), arXiv:2308.00825 [cond-mat.str-el] .
  * Rivera  _et al._ [2015] P. Rivera, J. R. Schaibley, A. M. Jones, J. S. Ross, S. Wu, G. Aivazian, P. Klement, K. Seyler, G. Clark, N. J. Ghimire, J. Yan, D. G. Mandrus, W. Yao, and X. Xu, Observation of long-lived interlayer excitons in monolayer MoSe2–WSe2 heterostructures, Nat. Commun. 6, 1 (2015).
  * [46] See supplemental materials for details.
  * Amelio  _et al._ [2023] I. Amelio, N. D. Drummond, E. Demler, R. Schmidt, and A. Imamoglu, Polaron spectroscopy of a bilayer excitonic insulator, Phys. Rev. B 107, 155303 (2023).
  * Fogler  _et al._ [2014] M. M. Fogler, L. V. Butov, and K. S. Novoselov, High-temperature superfluidity with indirect excitons in van der Waals heterostructures, Nat. Commun. 5, 1 (2014).
  * Gupta  _et al._ [2020] S. Gupta, A. Kutana, and B. I. Yakobson, Heterobilayers of 2D materials as a platform for excitonic superfluidity, Nat. Commun. 11, 1 (2020).
  * Filinov  _et al._ [2010] A. Filinov, N. V. Prokof’ev, and M. Bonitz, Berezinskii-kosterlitz-thouless transition in two-dimensional dipole systems, Phys. Rev. Lett. 105, 070401 (2010).
  * Petrov  _et al._ [2004] D. S. Petrov, D. M. Gangardt, and G. V. Shlyapnikov, Low-dimensional trapped gases, J. Phys. IV France. 116, 5 (2004).
  * Boudjemâa and Shlyapnikov [2013] A. Boudjemâa and G. V. Shlyapnikov, Two-dimensional dipolar bose gas with the roton-maxon excitation spectrum, Phys. Rev. A 87, 025601 (2013).
  * Marchetti  _et al._ [2008] F. M. Marchetti, C. J. M. Mathy, D. A. Huse, and M. M. Parish, Phase separation and collapse in bose-fermi mixtures with a feshbach resonance, Phys. Rev. B 78, 134517 (2008).
  * Wegner [2000] F. J. Wegner, Flow equations for Hamiltonians, Nucl. Phys. B Proc. Suppl. 90, 141 (2000).
  * Glazek and Wilson [1994] S. D. Glazek and K. G. Wilson, Perturbative renormalization group for hamiltonians, Phys. Rev. D 49, 4214 (1994).
  * Mielke [1997] A. Mielke, Calculating critical temperatures of superconductivity from a renormalized hamiltonian, Europhysics Letters 40, 195 (1997).
  * Abrikosov  _et al._ [1975] A. Abrikosov, L. Gor’kov, and I. Dzyaloshinski,  _Methods of Quantum Field Theory in Statistical Physics_, Dover Books on Physics Series (Dover Publications, 1975).
  * Pekker  _et al._ [2011] D. Pekker, M. Babadi, R. Sensarma, N. Zinner, L. Pollet, M. W. Zwierlein, and E. Demler, Competition between pairing and ferromagnetic instabilities in ultracold fermi gases near feshbach resonances, Phys. Rev. Lett. 106, 050402 (2011).
  * Knap  _et al._ [2016] M. Knap, M. Babadi, G. Refael, I. Martin, and E. Demler, Dynamical cooper pairing in nonequilibrium electron-phonon systems, Phys. Rev. B 94, 214504 (2016).
  * Stephen and Suhl [1964] M. J. Stephen and H. Suhl, Weak time dependence in pure superconductors, Phys. Rev. Lett. 13, 797 (1964).
  * Anderson  _et al._ [1965] P. W. Anderson, N. R. Werthamer, and J. M. Luttinger, An additional equation in the phenomenology of superconductivity: Resistive effects, Phys. Rev. 138, A1157 (1965).
  * Abrahams and Tsuneto [1966] E. Abrahams and T. Tsuneto, Time variation of the ginzburg-landau order parameter, Phys. Rev. 152, 416 (1966).
  * Marco  _et al._ [2019] L. D. Marco, G. Valtolina, K. Matsuda, W. G. Tobias, J. P. Covey, and J. Ye, A degenerate fermi gas of polar molecules, Science 363, 853 (2019).
  * Yan  _et al._ [2023] Z. Z. Yan, Y. Ni, A. Chuang, P. E. Dolgirev, K. Seetharam, E. Demler, C. Robens, and M. Zwierlein, Dissipationless flow in a bose-fermi mixture (2023), arXiv:2304.07663 [cond-mat.quant-gas] .
  * von Milczewski  _et al._ [2023] J. von Milczewski, X. Chen, A. Imamoǧlu, and R. Schmidt, Superconductivity induced by strong electron-exciton coupling in doped atomically thin semiconductor heterostructures (2023).
  * Zerba  _et al._ [2023] C. Zerba, C. Kuhlenkamp, A. Imamoglu, and M. Knap, Realizing Topological Superconductivity in Tunable Bose-Fermi Mixtures with Transition Metal Dichalcogenide Heterostructures (2023).

##  I Supplemental Material

###  I.1 Effective two-channel model

We obtain the effective Hamiltonian of Eq. (1) starting from a microscopic
model of excitons and doped charges:

$\displaystyle\hat{H}=$
$\displaystyle\sum_{\mathbf{k}}x^{\dagger}_{\mathbf{k}}\frac{k^{2}}{2m_{x}}x_{\mathbf{k}}+\begin{pmatrix}c^{\dagger}_{\mathbf{k},T}\\\
c^{\dagger}_{{\mathbf{k}},M}\end{pmatrix}\begin{pmatrix}\frac{{k}^{2}}{2m}+\Delta&t\\\
t&\frac{k^{2}}{2m}\end{pmatrix}\begin{pmatrix}c_{{\mathbf{k}},T}\\\
c_{{\mathbf{k}},M}\end{pmatrix}$ (7)
$\displaystyle+\frac{\tilde{U}}{V}\sum_{\mathbf{k},\mathbf{k}^{\prime},\mathbf{q}}c_{{\mathbf{k}},M}^{\dagger}c_{{\mathbf{k}}+{\mathbf{q}},M}x^{\dagger}_{{\mathbf{k}^{\prime}}}x_{{\mathbf{k}^{\prime}}-{\mathbf{q}}},$

where the index $\\{T,M\\}$ labels the top and middle layer; the interactions
between inter-layer excitons and holes is modelled by an effective contact
interaction of strength $\tilde{U}$. In the dilute limit, this model leads to
the Hamiltonian Eq. (1) reported in the main text provided the scattering
amplitudes between excitons and electrons are correctly matched.

Close to a Feshbach resonance, when $\Delta\simeq E_{T}^{0}$, the scattering
is dominated by the inter-layer trion bound state. We can thus consider an
effective model for the T-Matrix where scattering between exciton and top
layer electrons is mediated by a virtual molecule. Assuming $g$ is the bare
molecule-exciton-hole coupling, the off-shell T-matrix to second order in $g$
reads $T^{E}(\omega)\sim\frac{g^{2}}{\omega-E^{0}_{T}+\Delta+i\epsilon}$,
which has a pole at the molecular energy.

When considering the microscopic model for excitons and charge carriers in the
same layer, the scattering amplitude in the vicinity of the resonance is
$\frac{2\pi E^{0}_{T}}{\mu_{\text{red}}(\omega-E^{0}_{T}+\Delta+i\epsilon)}$,
where $\mu_{\text{red}}=m\,m_{x}/(m+m_{x})$ is the reduced mass of the hole
and the exciton. This relation is obtained by expanding the logarithm of two-
dimensional T-matrix [25]. The tunneling process can be added perturbatively
to this by introducing a vertex with tunneling strength $t$ and a propagator
that describes tunneling of the upper layer charge carrier to the middle
layer. Then, close to resonance, the tunneling process to the middle layer and
back contributes a factor $\frac{t^{2}}{\Delta^{2}}$. The global T-matrix is
thus $\sim\frac{2\pi
E^{0}_{T}}{\mu_{\text{red}}\Delta^{2}(\omega-E^{0}_{T}+\Delta+i\epsilon)}t^{2}$.
Comparing the scattering amplitudes to second order in the tunneling strength
$t$ yields a simple matching condition:

$g=\sqrt{\frac{2\pi}{\mu_{\text{red}}E_{T}^{0}}}t,$ (8)

relating the scattering matrix in the microscopic and effective model.

###  I.2 Wegner Flow

We present the Wegner flow method that we used to trace out the condensate
excitations and to obtain an effective interaction between the fermions
$\hat{f}$ [54]. Given the full Hamiltonian

$\displaystyle\hat{H}$ $\displaystyle=\hat{H}_{f,0}+\hat{H}_{f,b}$ (9)
$\displaystyle=\sum_{\mathbf{k}}\big{[}\lambda_{f}(\mathbf{k})-\mu\big{]}\hat{f}_{\mathbf{k}}^{\dagger}\hat{f}_{\mathbf{k}}+\sum_{\mathbf{k},\mathbf{q}}\big{[}M_{\mathbf{k},\mathbf{q}}\hat{f}^{\dagger}_{\mathbf{k}}\hat{f}_{\mathbf{k}-\mathbf{q}}\hat{b}_{\mathbf{q}}+\text{h.c.}\big{]},$
(10)

our aim is to obtain a transformation $U(l)$ depending on the continuous
parameter $l$ such that the Hamiltonian $\hat{H}(l)=U^{\dagger}(l)\hat{H}U(l)$
depends more on $\hat{f}$ only as $l$ increases. The generator of the
transformation is the operator $\eta(l)$,

$\frac{d}{dl}\hat{H}(l)=\big{[}\eta(l),\hat{H}(l)]$ (11)

and $\big{[}\frac{d}{dl}U(l)\big{]}U^{\dagger}(l)=\eta(l)$. Following [54], we
consider $\eta(l)=[\hat{H}_{f,0}(l),\hat{H}_{f,x}(l)]$ as an ansatz for the
continuous transformation

$\eta(l)=\sum_{\mathbf{k},\mathbf{q}}\left[M_{\mathbf{k},\mathbf{q}}(l)\Omega(\mathbf{k},\mathbf{q})\hat{f}^{\dagger}_{\mathbf{k}}\hat{f}_{\mathbf{k}-\mathbf{q}}\hat{b}_{\mathbf{q}}-\text{h.c.}\right].$
(12)

where
$\Omega(\mathbf{k},\mathbf{q})=\lambda_{f}(\mathbf{k})-\lambda_{f}(\mathbf{k}-\mathbf{q})-e_{x}(\mathbf{q})$.
The flow in Eq. (11) has a contribution from the commutator
$\big{[}\eta(l),\hat{H}_{f,0}]$, which yields a three body term and defines
the flow of the coefficient $M_{\mathbf{k},\mathbf{q}}(l)$ as

$\frac{d}{dl}M_{\mathbf{k},\mathbf{q}}(l)=-\Omega(\mathbf{k},\mathbf{q})^{2}M_{\mathbf{k},\mathbf{q}}(l)$
(13)

Solving this equation gives
$M_{\mathbf{k},\mathbf{q}}(l)=e^{-\Omega(\mathbf{k},\mathbf{q})^{2}l}M_{\mathbf{k},\mathbf{q}}(l=0)$,
and hence $M_{\mathbf{k},\mathbf{q}}(l\to\infty)=0$. For $l\to\infty$ the
three-body term disappears as we required from the flow. However, the
commutator $[\eta(l),\hat{H}_{f,x}]$ introduces complicated many-body terms
for the $\hat{f}$ particles. In particular, they generate a fermion-fermion
interaction $\hat{H}_{f,f}=\sum_{\mathbf{k},\mathbf{k^{\prime}},\mathbf{q}}$
$V_{\mathbf{k},\mathbf{k^{\prime}},\mathbf{q}}\hat{f}_{\mathbf{k}-\mathbf{q}}^{\dagger}\hat{f}_{\mathbf{k^{\prime}}+\mathbf{q}}^{\dagger}\hat{f}_{\mathbf{k^{\prime}}}\hat{f}_{\mathbf{k}}$.
To analyze superconducting instabilities, we focus on the terms in the
potential with $\mathbf{k}=-\mathbf{k}^{\prime}$, which are the most relevant
for the Cooper pair formation. From Eq. (11) we then obtain

$\frac{d}{dl}V_{k,-k,q}=M_{\mathbf{k},\mathbf{q}}(l)M_{-\mathbf{k}+\mathbf{q},\mathbf{q}}(l)\big{[}\Omega(\mathbf{k},\mathbf{q})+\Omega(-\mathbf{k}+\mathbf{q},\mathbf{q})\big{]}.$
(14)

By explicitly integrating the above equation we obtain the effective
interaction

$\int_{0}^{\infty}\frac{d}{dl}V_{k,-k,q}=\frac{-M_{\mathbf{k},\mathbf{q}}(0)M^{*}_{-\mathbf{k}+\mathbf{q},\mathbf{q}}(0)e_{x}(\mathbf{q})}{(\lambda_{f}(\mathbf{k})-\lambda_{f}(\mathbf{k}-\mathbf{q}))^{2}+e_{x}(\mathbf{q})^{2}}.$
(15)

The commutator of Eq. (11) gives rise to additional terms as well, including
density-density terms, that are not relevant for our superconductivity
analysis. Other terms are of higher order in the superfluid density, that is
assumed to be small, or are of higher order in the interaction $g$ and, thus
can be neglected as well [54].

###  I.3  Angular momentum expansion

To identify the contribution of each angular momentum to the superconducting
order parameter, we consider an expansion of the eigenvectors of the matrix
$C$ as
$\psi_{\mathbf{k}}=\sum_{n}\psi^{c}_{k,n}\cos{((2n+1)\theta_{k})}+\psi^{s}_{k,n}\sin{((2n+1)\theta_{k})}$.
Using Eq. (6) we can then write a set of equations for the components
$\psi^{c}_{k,n}$ and $\psi^{s}_{k,n}$

$\frac{d\psi^{c}_{n,k}}{dt}=-\frac{i}{\hbar}\sum_{m}\int
d\theta_{k}\cos{(2n+1)\theta_{k}}A_{\mathbf{k}}\big{[}\cos{((2m+1)\theta_{k})}\\\
\psi^{c}_{m,k}+\sin{((2m+1)\theta_{k})}\psi^{s}_{m,k}\big{]}+\sum_{p,m}\int
d\theta_{k}d\theta_{p}\cos{((2n+1)\theta_{k})}\\\
B_{\mathbf{k},\mathbf{p}}\big{[}\cos{((2m+1)\theta_{p})}\psi^{c}_{m,p}+\sin{((2m+1)\theta_{p})}\psi^{s}_{m,p}\big{]}$
(16) $\frac{d\psi^{s}_{n,k}}{dt}=-\frac{i}{\hbar}\sum_{m}\int
d\theta_{k}\sin{((2n+1)\theta_{k})}A_{\mathbf{k}}\big{[}\cos{((2m+1)\theta_{k})}\\\
\psi^{c}_{m,k}+\sin{((2m+1)\theta_{k})}\psi^{s}_{m,k}\big{]}+\sum_{p,m}\int
d\theta_{k}d\theta_{p}\sin{((2n+1)\theta_{k})}\\\
B_{\mathbf{k},\mathbf{p}}\big{[}\cos{((2m+1)\theta_{p})}\psi^{c}_{m,p}+\sin{((2m+1)\theta_{p})}\psi^{s}_{m,p}\big{]}$
(17)

Consequently, we define a projected matrix $\tilde{C}=\tilde{A}+\tilde{B}$.
The instability can then be identified by exponentially growing solutions of
the order parameter as a function of time as discussed in the main text. In
this new definition, $\tilde{A}$ and $\tilde{B}$ are tensors defined by
indices: $\tilde{A}_{k,n,p,m}^{c,c}=\psi_{k,p}\int
d\theta_{k}\cos{((2n+1)\theta_{k})}A_{\mathbf{k}}\cos{((2m+1)\theta_{k})}$,
$\tilde{A}_{k,n,p,m}^{c,s}=\psi_{k,p}\int
d\theta_{k}\cos{((2n+1)\theta_{k})}A_{\mathbf{k}}\sin{((2m+1)\theta_{k})}$ and
$\tilde{A}_{k,n,p,m}^{s,s}$ and $\tilde{A}_{k,n,p,m}^{s,c}$ are given
correspondingly by considering sine instead of cosine in
$\tilde{A}_{k,n,p,m}^{c,c}$ and $\tilde{A}_{k,n,p,m}^{c,s}$. Similarly,
$\tilde{B}^{c,c}_{k,n,p,m}=\int
d\theta_{k}d\theta_{p}\cos{((2n+1)\theta_{k})}B_{\mathbf{k},\mathbf{p}}\cos{((2m+1)\theta_{p})}$,
$\tilde{B}^{c,s}_{k,n,p,m}=\int
d\theta_{k}d\theta_{p}\cos{((2n+1)\theta_{k})}B_{\mathbf{k},\mathbf{p}}\sin{((2m+1)\theta_{p})}$,
and $\tilde{B}_{k,n,p,m}^{s,s}$ and $\tilde{B}_{k,n,p,m}^{s,c}$ are given
correspondingly by considering sine instead of cosine and vice versa in
$\tilde{B}_{k,n,p,m}^{c,c}$ and $\tilde{B}_{k,n,p,m}^{c,s}$.
We always find a two-fold degeneracy of the eigenvalues of the matrix
$\tilde{C}$. We interpret this as due to an equal-weight contribution of the
sine and cosine to the eigenvector $|\psi\rangle$. In the degenerate
2-dimensional eigenspace, the basis that diagonalizes $\tilde{C}$ can be
chosen in such a way that we have an eigenvector with only sine contributions
and an eigenvector with only cosine contributions. This is expected, since the
Hamiltonian is symmetric under rotations. This result implies that critical
temperature for the sine and cosine eigenvalue is the same.

###  I.4 Topological superconductivity

Figure 4: Analysis of the relative phase $\varphi$ of the angular momentum
contributions to the order parameter. (a) Free energy as a function of the
phase $\varphi$, computed for $T=1.5K$, $\Delta=2$ meV, $E_{F}=1$ meV.
Different magnitude of the instabilities are compared by changing the
parameter r, $max_{k}(\psi_{k}$)=r$\sqrt{\lambda}$, where $\lambda\propto$
Tc-T and $\varrho\propto\frac{1}{r^{2}}$ is the fourth order term in the
Ginzburg-Landau formalism, see text. (b) Difference of Free energy at
$\varphi=0$ and $\varphi=\frac{\pi}{2}$ as a function of the temperature for
different magnitudes of the order parameter shows that $\varphi=\frac{\pi}{2}$
is the solution with lowest free energy.

Here, we discuss the symmetry of the order parameter, which leads to the
conclusion that the state is a topological $p$+$ip$ superconductor.

We observed that the only relevant contribution to the eigenvector of the
maximally unstable solution $|\psi\rangle$ is given by the p-wave
contribution. The instability is generated by two eigenvectors of the matrix
$\tilde{C}$: one with contribution only in the cosine degrees of freedom and
the other one only in the sine degrees of freedom. Furthermore, by Ginzburg-
Landau analysis, we observe that the critical temperature is the same for the
two directions. Because of this, the Ginzburg- Landau free energy in Eq. 4 can
be written as
$\mathcal{F}^{GL}=\sum\psi_{x}^{*}\lambda\psi_{x}+\psi_{y}^{*}\lambda\psi_{y}$,
where we defined $\psi_{x}$ as the eigenstate with only cosine contributes,
$\psi_{y}$ the eigenstate with only sine contribution, or equivalently by
considering any unitary transformation of $(\psi_{x},\psi_{y})$. Among these
possible linear combinations, we find also the topological p-wave
superconductor $(\psi_{x}+i\psi_{y},\psi_{x}-i\psi_{y})$ and the nematic one
$(\psi_{x}+\psi_{y},\psi_{x}-\psi_{y})$; up to a rotation in momentum space,
the generic form of the instability can be written as
$(\psi_{x}+e^{i\varphi}\psi_{y},\psi_{x}-e^{-i\varphi}\psi_{y})$. To determine
the nature of the p-wave superconductivity, we would have to consider further
order terms in the Ginzburg-Landau free energy and find $\varphi$ such that
the Ginzburg-Landau free energy is minimized. The higher-order terms can be
computed perturbatively. However, in order to estimate the type of instability
a simpler strategy is sufficient. We determine the order parameter up to a
constant. To this end, we introduce a mean-field Hamiltonian that describes a
system, that is unstable under the previous determined superconducting order
parameter, and for which the only free parameter is the phase $\varphi$. By
considering several different magnitudes of the order parameter, we show that
the instability directions are the only relevant information needed to
determine the topological nature of the superconductor, see Fig. 4.

First, we show that if the only relevant contribution to the instability is
p-wave, the same holds for the order parameter. Indeed we have
$\Gamma_{\mathbf{k}}=\sum_{\mathbf{p}}V_{\mathbf{k},\mathbf{k}-\mathbf{p}}\psi_{\mathbf{p}}=\sum_{\mathbf{p}>0}V_{\mathbf{k},\mathbf{k}-\mathbf{p}}\psi_{\mathbf{p}}-V_{\mathbf{k},\mathbf{k}+\mathbf{p}}\psi_{\mathbf{p}}$.
Using the fact that the potential is invariant under the transformation
$V_{\mathbf{k},\mathbf{q}}\to V_{-\mathbf{k},-\mathbf{q}}$, we can write
$\Gamma_{\mathbf{k}}=\frac{1}{2}\sum_{\mathbf{p}>0}V_{\mathbf{k},\mathbf{k}-\mathbf{p}}\psi_{\mathbf{p}}-V_{\mathbf{k},\mathbf{k}+\mathbf{p}}\psi_{\mathbf{p}}+V_{\mathbf{-k},\mathbf{-k}+\mathbf{p}}\psi_{\mathbf{p}}-V_{\mathbf{-k},\mathbf{-k}-\mathbf{p}}\psi_{\mathbf{p}}$.
When $\psi_{\mathbf{k}}$ is an eigenvector of the matrix $C$, then we have
from Eq (17)

$2(\lambda_{f}(\mathbf{k})-\mu)\psi_{\mathbf{k}}-(2n_{\mathbf{k}}-1)\sum_{\mathbf{p}>0}\big{[}V_{\mathbf{k},\mathbf{k}-\mathbf{p}}+\\\
-V_{\mathbf{k},\mathbf{k}+\mathbf{p}}-V_{\mathbf{-k},\mathbf{-k}-\mathbf{p}}+V_{\mathbf{-k},\mathbf{-k}+\mathbf{p}}\big{]}\psi_{\mathbf{p}}=E_{\psi}\psi_{\mathbf{k}}$
(18)

Hence,

$2(\lambda_{f}(\mathbf{k})-\mu)\psi_{\mathbf{k}}-2(2n_{\mathbf{k}}-1)\Gamma_{\mathbf{k}}=E_{\psi}\psi_{\mathbf{k}}.$
(19)

from which we obtain an expression for $\Gamma_{\mathbf{k}}$

$\Gamma_{\mathbf{k}}=\frac{2(\lambda_{f}(\mathbf{k})-\mu)-E_{\psi}}{2(2n_{\mathbf{k}}-1)}\psi_{\mathbf{k}}$
(20)

where the only dependence from $\theta_{k}$ is coming from
$\psi_{\mathbf{k}}$. This implies that to the eigenvector of the maximal
unstable solution $\ket{\psi^{c}_{k}}$, which contributes only to the cosine
degrees of freedom, corresponds to an order parameter, which also contributes
only on the cosine degrees of freedom. The same is true for the corresponding
sine degrees of freedom. Notice that $\Gamma_{\mathbf{k}}$ has been determined
up to a constant which is the same for all momenta $\mathbf{k}$.

To determine whether the p-wave superconducting phase is topological or
trivial, we discard the contributes of the potential that are not relevant for
the formation of superconductivity and consider the mean-field model

$H=\sum_{\mathbf{k}}\xi(\mathbf{k})\hat{f}^{\dagger}_{\mathbf{k}}\hat{f}_{\mathbf{k}}+\sum_{\mathbf{k}}\bigg{(}\big{[}\Gamma_{k}\cos(\theta_{k})+\\\
+e^{i\varphi}\Gamma_{k}\sin{\theta_{k}}\big{]}\hat{f}^{\dagger}_{\mathbf{k}}\hat{f}^{\dagger}_{-\mathbf{k}}+\text{h.c.}\bigg{)},$
(21)

where we defined $\xi(\mathbf{k})=\lambda_{f}(\mathbf{k})-\mu$. When tuning
the phase $\varphi$ different p-wave superconductors can be realized,
including the topological chiral and the non-topological nematic ones.
Specifically, phases of $\varphi=0,\pi$ correspond to nematic order in which
the only contribution to the superconductor is $p_{x}$ or $p_{y}$ and the gap
is not rotational invariant leading to gapless points on the Fermi surface. By
contrast, phases of $\varphi=\frac{\pi}{2},\frac{3\pi}{2}$ correspond to
chiral topological phases in which the gap is rotational invariant and the
order parameter is fully gapped.

For this model and finite temperatures we can evaluate the free energy,
$F=E-TS$; we compute the entropy as
$S=-k_{B}[n_{\mathbf{k}}\ln{n_{\mathbf{k}}}+(1-n_{\mathbf{k}})\ln({1-n_{\mathbf{k}}})]$
using the occupation $n_{\mathbf{k}}$ of the Bogoliubov excitations,
$n_{\mathbf{k}}=\frac{1}{\exp[\beta E_{\mathbf{k}}]+1}$, $\beta=k_{B}T$,
$E_{\mathbf{k}}=\sqrt{{\xi(\mathbf{k})^{2}+|\Gamma_{k}|^{2}(1+\cos(\varphi)\sin(2\theta_{k}))}}$,
and $E=\sum_{\mathbf{k}}E_{\mathbf{k}}n_{\mathbf{k}}+E_{GS}$. We observe that
the dependence on $k$ of the order parameter guarantees that
$\varphi=\pm\frac{\pi}{2}$ are minima, while $\varphi=0,\pi$ are maxima., see
Fig. 4 (a). We test this result against different magnitudes of the
instability, effectively exploring two order of magnitudes of values for the
prefactor of the fourth order term in the Ginzburg-Landau theory
$\varrho\propto\frac{1}{r^{2}}$, defined by
$\mathcal{F}=\lambda|\psi|^{2}+\frac{1}{2}\varrho|\psi|^{4}$. In Fig. 4 (b) we
report the difference between the free energy at $\varphi=0$ and
$\varphi=\frac{\pi}{2}$ for different values of the temperature $T\sim T_{c}$,
showing that the same holds in a finite range of temperature in the vicinity
of $T_{c}$.


# ChatGPT scores a bad birdie in counting gravitational-wave chirps

Floor S. Broekgaarden,1
1Center for Astrophysics | Harvard $\&$ Smithsonian, 60 Garden St., Cambridge, MA 02138, USA
E-mail<EMAIL_ADDRESS>
(April 1 2023)

###### Abstract

How many gravitational-wave observations from compact object mergers have we
seen to date? This seemingly simple question has a surprisingly complex answer
that even ChatGPT struggles to answer. To shed light on this, we present a
database with the literature’s answers to this question. We find values
spanning 67–100 for the number of detections from double compact object
mergers to date, emphasizing that the exact number of detections is uncertain
and depends on the chosen data-analysis pipeline and underlying assumptions.
We also review the number of gravitational-wave detections expected in the
coming decades with future observing runs, finding values up to millions of
detections per year in the era of Cosmic Explorer and Einstein Telescope. We
present a publicly available code to visualize the detection numbers,
highlighting the exponential growth in gravitational-wave observations in the
coming decades and the exciting prospects of gravitational-wave astrophysics.
See http://www.broekgaarden.nl/floor/wordpress/elementor-967/. We plan to keep
this database up-to-date and welcome comments and suggestions for additional
references.

###### keywords:

Aprils Fools – Visualizations – Do not count on ChatGPT for GWs just yet –
thesis procrastination

††pubyear: 2023††pagerange: ChatGPT scores a bad birdie in counting
gravitational-wave chirps–ChatGPT scores a bad birdie in counting
gravitational-wave chirps

##  1 Introduction

Since the first detection of gravitational waves from a binary black hole
merger in 2015 (Abbott et al., 2016a), GW interferometers from advanced LIGO
(Aasi et al., 2015; Abbott et al., 2016b), advanced Virgo (Acernese et al.,
2015, 2019), and advanced KAGRA (Somiya, 2012; Aso et al., 2013; Akutsu et
al., 2021)111Henceforth referred to as the LVK collaboration have rapidly
increased the number of observed double-compact object mergers. To share this
excitement, numerous papers and talks on gravitational waves often will
mention the current number of detections from double-compact object mergers:
“To date, a total of $X$ gravitational-wave sources have been
observed”222Often followed by cheesy sentences including things like ‘a
revolution in astronomy’ ‘opening a new window into the Universe’, ‘golden
era’, ‘new frontier’, indeed even ChatGPT loves these, see 2.. However, as we
will discuss in more detail below, this number $X$ of confident detections
from double-compact object mergers depends on the used data-analysis pipelines
and the underlying assumptions such as the detection or signal-to-noise
threshold. Moreover, for future observation runs and next-generation
detectors, different expected detection numbers can especially arise from
different assumptions about the underlying population properties of double
compact object mergers in our Universe.

In this research note, we shed light on the number of detections in existing
and future GW observation runs as reported in the literature. We present an
overview of the reported number of detections in the literature from past
observing runs O1, O2, and O3 in Section 2 and for future observing runs in
Section 3. A publicly available code to create visualizations to showcase the
(growing) number of detections is provided at need-to-make-this, which
includes videos made for outreach and presentations. We hope this note
provides a fun tool for astronomers interested in this topic. We plan to
update this regularly in the future and welcome additional sources and
references from the community.

##  2 Current Detections from O1, O2, and O3

At the time of writing LVK have shared the data and detections from observing
runs O1, O2, O3a and O3b in the Gravitational-wave Transient Catalogs (GWTC)
GWTC-1 (Abbott et al., 2019), GWTC-2 (Abbott et al., 2021b), GWTC-2.1 (Abbott
et al., 2021a), and GWTC-3 (The LIGO Scientific Collaboration et al., 2021a).
O1 ran from September 12th, 2015 to January 19th, 2016, O2 ran from November
30th, 2016 to August 25th, 2017, O3a ran from 1 April 2019 until 1 October
2019, and O3b ran between 1 November 2019 and 27 March 2020. Based on the data
in these finished runs, several catalogs have been presented with different
numbers of events, based on different underlying assumptions for the detection
pipelines and thresholds. We list the main catalog options considered in the
literature below, and summarize them in Table 1 and Figure 1

Table 1: Summary of Detection numbers. Acronyms are: FAR (false alarm rate),
GWTC (Gravitational Wave Transient Catalog), GWOSC (Gravitational Wave Open
Science Center).
References:
(a) GWOSC Dataset (Collaboration et al., 2022) events with
$p_{\rm{astro}}>0.5$.
(b): Table 1 from LVK population paper (The LIGO Scientific Collaboration et
al., 2021b).
(c) All detections with $p_{\rm{astro}}>0.5$ or $\rm{FAR}<1/100\,\rm{yr}^{-1}$
from Table 3 in (Nitz et al., 2021a).
(d) Union of GWOSC and (Nitz et al., 2021a) catalog with all additional events
found with $p_{\rm{astro}}\geq 0.5$ in the independent catalogs by Venumadhav
et al. (2019, 2020); Zackay et al. (2019b); Nitz et al. (2020, 2021b, 2021a).
(e) Using https://chat.openai.com/chat, see Figure 2.  Label | Total | BHBH | BHNS | NSNS | Ref
---|---|---|---|---|---
Detections from O1, O2, O3 |  |  |  |  |
GWOSC (GWTC–1 $+$ GWTC–2.1 $+$ GWTC–3): | 90 | 85 | 3 | 2 | (a)
GWTC–3 population events with FAR $<0.25\,\rm{yr}^{-1}$ | 67 | 63 | 2 | 2 | (b)
GWTC–3 population events with FAR $<1\,\rm{yr}^{-1}$ | 76 | 69 | 4 | 2 | (b)
4-OGC Catalog (Nitz et al., 2021a) with $p_{\rm{astro}}\geq 0.5$ | 94 | 90 | 2 | 2 | (c)
Optimistic (all possible events with $p_{\rm{astro}}\geq 0.5$) | 101 | 95 | 4 | 2 | (d)
ChatGPT | 50 | 44 | 1 | 5 | (e)
Figure 1: Cumulative histogram of the number of detected mergers from binary
black holes (BH-BH; black), black hole-neutron stars (BHNS; teal), and binary
neutron stars (NS-NS; orange) for different data analysis pipelines in the
literature from Table 1. The exact number of detections depends on the
underlying assumptions such as the chosen threshold of assigning
gravitational-wave events to be of (likely) astrophysical origin. The number
of events by the end of O3 for different catalogs is quoted in the middle
right. An animated version of this figure is available at
http://www.broekgaarden.nl/floor/wordpress/elementor-967/.

###  2.1 GWOSC catalog (GWTC–1, GWTC–2.1, and GWTC–3)

LVK published several event catalogs based on observing runs O1, O2, and O3.
In GWTC-1 the combined detections from O1 and O2 are presented with 3 binary
black hole (BHBH) mergers in O1, one detection of a binary neutron star (NSNS)
in O2, and 7 more BHBH mergers in O2. In GWTC-2 LVK reported a total of 39
compact object coalescence detections from the O3a run. In GWTC-2.1, the LVK
collaboration reported an updated, deeper, list of candidate events for the
same observation period O3a that removed three GW events from GWTC-2 and added
eight additional events, due to improved calibration and subtraction of excess
noise. This gives a total of 44 events with an astrophysical probability of at
least 0.5 in O3a based on GWTC-2.1, including one NSNS. In GWTC-3 the LVK
collaboration presented the detections from O3b which included 35 detections
with an astrophysical probability $>0.5$, including three black hole–neutron
stars. All in all, this gives a total of 90 events, with 85 BHBHs, 3 BHNSs,
and 2 NSNSs that have an astrophysical probability $>0.5$. The combined GWTC
catalog is presented and kept up-to-date as the GWTC by the Gravitational Wave
Open Science Center (GWOSC) available via Zenodo (Collaboration et al.,
2022)333As listed on the GWOSC website https://www.gw-
openscience.org/eventapi/html/GWTC/.. Noteworthy, BHNS merger GW200105 is not
included in this GWOSC catalog, as it has an astrophysical probability of
$<0.5$. On the other hand, GW200105 is included in the paper representing the
first two detections of a BHNS by the LVK collaboration (Abbott et al., 2021c)
and the GWTC-3 population inference catalogs (see below) which uses a false-
alarm-rate (FAR) instead of the astrophysical probability as a threshold for
assigning events to a catalog. Larger numbers of BHNS detections in the future
will help characterize the properties of the BHNS population and might address
whether GW200105 is likely of astrophysical origin. We also note that in the
GWOSC catalog (and all other catalogs mentioned below), GW190814 with a
secondary mass of $m_{\rm{2}}=2.6^{+0.1}_{-0.1}\,\rm{M}_{\odot}$ is classified
as a BHBH, but the possibility exists that this is instead a BHNS source
(Abbott et al., 2020; The LIGO Scientific Collaboration et al., 2021b).

###  2.2 GWTC-3 population catalog with FAR $<0.25\,\rm{yr}^{-1}$

In The LIGO Scientific Collaboration et al. (2021b), LVK present the
population properties of compact object mergers based on the data from O1, O2,
and O3. The authors use an updated set of events, based on a chosen false-
alarm threshold of $<0.25\,\rm{yr}^{-1}$, to assure high purity, particularly
among mergers with a neutron star. This omits several candidates of moderate
significance including GW190403$\\_$051519 and GW200220$\\_$061928 from the
GWTC catalogs (GWTC-1, GWTC-2, GWTC-2.1, and GWTC-3). At this threshold LVK
expects approximately one event to not be of astrophysical origin. They
present the population catalog in their Table 1, consisting of 63 BHBH, 2
BHNS, and 2 NSNS mergers.

###  2.3 GWTC-3 population catalog events with FAR $<1\,\rm{yr}^{-1}$

Same as above, but now including all events with FAR $<1\,\rm{yr}^{-1}$, which
is a slightly less strict threshold. At this threshold LVK expects
approximately 4.6 events to not be of astrophysical origin. Using this LVK
presents the population catalog in Table 1 with 76 events consisting of 69
BHBH, 4 BHNS, and 2 NSNS mergers (The LIGO Scientific Collaboration et al.,
2021b). This population catalog is used by the authors to analyze the redshift
evolution of BHBH mergers.

###  2.4 4-OGC: independent catalog events with $p_{\rm{astro}}>0.5$

In Nitz et al. (2021a), the authors present an independent data analysis of
the O1, O2, and O3 data with slightly different data analysis pipelines.
Compared to earlier catalogs (1-OGC, 2-OGC, and 3-OGC), the authors find 7 new
unreported BHBH events compared to GWTC-3, and also find that three earlier
BHBH events now pass the threshold for a detection, based on updated
assumptions. We take all the detections with $p_{\rm{astro}}>0.5$ or an
inverse false alarm rate $>100\,\rm{years}$ from Table 3 in (Nitz et al.,
2021a), which includes a total of 94 observed events with 90 BHBH, 2 BHNS, and
2 NSNS mergers.

###  2.5 Optimistic catalog events (union of all claimed events with
$p_{\rm{astro}}>0.5$)

Several gravitational-wave searches independent of the LVK collaboration
analysis have been performed. These searches typically use slightly different
data analysis pipelines and methods to calculate the astrophysical probability
($p_{\rm{astro}}$) or threshold or a significant source. Many of the
gravitational-wave events reported by LVK are also detected in these
pipelines, but a few additional events not reported in the GWTC catalogs of
the LVK collaboration have also been found. For the purpose of visualization,
we decide to make an ‘optimistic’ catalog that consists of all events from the
GWOSC and 4-OGC catalogs but also includes additional events from independent
catalogs that have an astrophysical origin of $>0.5$. This includes in O1 one
BHBH event444Where we note that GW151205 with $p_{\rm{astro}}=0.53$ was
reported by (Nitz et al., 2020) but was found to have $p_{\rm{astro}}<0.5$ in
the updated 4-OGC catalog Nitz et al. (2021a) and is omitted by us from the
optimistic catalog.: GW151216 with astrophysical probability
$p_{\rm{astro}}\sim 0.71$ Venumadhav et al. (2019); Zackay et al. (2019b). In
O2 three additional BHBH events with likely astrophysical origin
($p_{\rm{astro}}\geq 0.98$) have been reported by Venumadhav et al. (2020):
GW170121, GW170304, and GW170727. These events were found independently with
$p_{\rm{astro}}\geq 0.7$ in Nitz et al. (2020). Three other BHBH events,
GW170425, GW170202, and GW170403, with astrophysical probabilities between
$0.5$–$0.8$ were also reported in Venumadhav et al. (2020) and confirmed, but
with astrophysical probabilities $<0.25$ by 2-OGC in (Nitz et al., 2020), but
then updated to astrophysical probabilities of 0.37, 0.86, and 0.72 in Nitz et
al. (2021a). In Nitz et al. (2021a), the authors presented with 4-OGC seven
additional events with $p_{\rm{astro}}\geq 0.5$ which we also add to the
catalog. We also include GW170817A, which has been reporterd to have an
astrophysical probability of 0.86 (Zackay et al., 2019a). The same study found
the potential gravitational-wave source GWC170402, but we do not include it,
as the event cannot fully be described by a binary black hole waveform. We
also do not include events like GW170502 (Udall et al., 2020), which are
considered marginal triggers in all literature as far as we are aware.

We combine all additional events with the GWOSC and 4-OGC catalogs to create a
new ‘Optimistic’ catalog that presents a very optimistic catalog of
detections. This catalog has a total of 101 events with 95 BHBH, 4 BHNS, and 2
NSNS.

###  2.6 ChatGPT

We also asked ChatGPT (https://chat.openai.com/chat) what the number of
gravitational-wave detections by the time of the publication of GWTC-3 was.
The responses are shown in Figure 2. Despite being aware of GWTC-3 and the
third observing run, ChatGPT vastly underestimates the number of GW detections
to be with certainty 50 total events, consisting of 44 BHBH, 5 NSNS and 1
BHNS. ChatGPT clearly also uses different thresholds for detections, as
rephrasing a few times the question resulted in answers by ChatGPT for the
total number of events to be 50, 59, or 81. Interestingly, despite that
ChatGPT’s knowledge cutoff date of September 2021 is before GWTC-3 was
published, it claims to know the GWTC-3 detections with much certainty,
despite the fact that the claimed events are a combination of real events,
triggers from the public trigger data, and made-up events. See Figure 3, with
a few funny responses. All in all, it seems clear that ChatGPT performs under
par when counting gravitational-wave events and that the work from GW
scientists will be required to understand the current number of GW detections.
We hope this research note helps in this goal.

Figure 2: Answer by ChatGPT to the question “How many and what type of
confident gravitational wave sources have been detected by the time GWTC-3 was
published?”. ChatGPT answers this with 50 total events of which 44 BHBH, 5
NSNS, and 1 BHNS mergers Figure 3: Answers by ChatGPT to questions about the
number of gravitational-wave detections.

##  3 Future Detections

###  3.1 Expected Detections for O4 and O5

The next upcoming LVK observing runs are O4 in 2023 and O5 around 2027. O4 is
expected to start on May 24, 2023, and it is planned to have an 18-calendar
month observing run, and LIGO detector sensitivity for binary neutron stars of
about $160$–$190\,\rm{Mpc}$555See
https://observing.docs.ligo.org/plan/index.html for the most up-to-date
information.. O5 will start around 2027 and is planned to run for about two
years with a target LIGO sensitivity of $240$–$325\,\rm{Mpc}$. Two recent
examples of expected detection numbers are the following. First, LVK reported
the prospects for the number of detections in O4 and O5 (Abbott et al.,
2018)666Where we use the updated information for O4 and O5 from
https://emfollow.docs.ligo.org/.. Assuming among other things that the mass
distributions follows the “Power law $+$ Dip $+$ Break” model from (The LIGO
Scientific Collaboration et al., 2021b), they find that they expect to detect
in O4 about 11 NSNS, 2 BHNS, and 77 BHBH mergers777Where we have scaled down
the median numbers in the summary Table from https://emfollow.docs.ligo.org/
by a factor 3.375 and rounded the value to go from public alert count to an
approximate number of confident detections, see their footnote.. For O5 the
expected numbers are about 53 NSNS, 9 BHNS, and 258 BHBH mergers. For both O4
and O5, the expected numbers have uncertainties that typically are of order
unity or twice unity, reflecting both the merger rate uncertainty and the
poison counting uncertainty (see also the discussion in Petrov et al. 2022).

It is important to realize that the exact number of BHBH, BHNS, and NSNS
mergers that will be detected from O4, and O5 is highly uncertain. This number
will depend on many factors including the actual sensitivity and observing
duty cycle of the instruments, and the Poisson statistics of the events (i.e.
the probability of a given number of events occurring in a fixed interval of
time). Another important uncertainty is the underlying population properties.
One example is the mass distributions of the double compact object merger
populations. Since the LVK detector sensitivity scales strongly with the
source masses, the number of detected events is currently heavily biased
towards more massive events, making the frequency and properties of lower-mass
black holes more uncertain (The LIGO Scientific Collaboration et al., 2021b).
This makes it challenging to estimate the number of detections for future
observing runs that probe larger distances and therefore see more of the
underlying astrophysical population. Another important uncertainty is the
redshift evolution of the rate and masses of the compact object merger events.
From theory and observations, it is thought that the rate and properties of
the double-compact objects change at higher redshifts (e.g. van Son et al.,
2021; The LIGO Scientific Collaboration et al., 2021b; Callister & Farr,
2023). Moreover, different formation channels with different properties might
be at play (Ng et al., 2021; Zevin et al., 2021). This can all impact the
number of detections in future detectors.

Given the many uncertain factors that can impact the number of detected events
in future detectors, we decide to take a very simplistic approach for the
purpose of the visualizations presented here. We caveat to the reader that
these visualizations are to sketch a potential landscape of the number of
expected gravitational-wave detections and that the actual expected number is
still highly uncertain. We instead simply use the GWOSC catalog as a proxy for
the ratio between BHBH, BHNS, and NSNS detections as $85:3:2$, and assume the
total number of 50 BHBH mergers to be somewhat representative for a year-long
observing run at O3 sensitivity. To obtain the expected number of detections
for O4 and O5, we then scale these numbers to the larger detection volume by
using the O4 and O5 binary neutron star ranges and multiplying them by the
detector observing times of 1.5 years and 2 years, respectively, and using a
duty cycle of 0.7 in both O4 and O5. This gives the very very rough estimates
for O4 and O5 of

$N_{\rm{inO4}}=50\times(175/120)^{3}*1.5=217\rm{\ total\ BHBH\ detections}$
(1) $N_{\rm{inO5}}=50\times(282.5/120)^{3}*2=913\rm{\ total\ BHBH\
detections},$ (2)

with $175\,\rm{Mpc}$ and $282.5\,\rm{Mpc}$ being the middle of the estimated
NSNS distance for O4 and O5, respectively, and $120\,\rm{Mpc}$ as proxy for
the O3 NSNS distance888See https://observing.docs.ligo.org/plan/. For the
uncertainty, we simply use a factor of three in both directions. For the BHNS
and NSNS rates we simply scale all numbers with the ratio of BHBH, BHNS, and
NSNS detections as $85:3:2$.

###  3.2 Expected Detections for Voyager, Cosmic Explorer and Einstein
Telescope

The exact sensitivity, and hence the expected number of detections that will
be made with Voyager is uncertain. In Adhikari et al. (2020), the authors
mention that Voyager might detect 100 times more sources compared to advanced
LIGO/Virgo/KAGRA. As a proxy, we decide to use a factor of 50 times more
sources compared to our O4 estimate (which comes down to about a factor 100
compared to O3)

The number of detections expected with Cosmic Explorer (CE) and Einstein
Telescope (ET) is also still highly uncertain. In recent work, Iacovelli et
al. (2022) present the most up-to-date expectations for the number of
detections. In their Tables 3, 4, and 5 they quote the expected number of
BHBH, BHNS, and NSNS detections with an ET or ET and two CE-like
interferometer networks. We summarize these rates in our Table 3 and use the
ET$+$2CE rates as a proxy for the CE/ET rates. A source of uncertainty
includes the uncertain local merger rate, which the authors use as
normalization in the underlying merger rate population for their calculation.
In their Section 5.4 they discuss the impact on the detection number by
changing these bounds. We take the most optimistic and pessimistic numbers
from these ranges for our uncertainty bounds in the visualization plots, which
we also quote in Table 3. Earlier work has also been performed by authors
including (Borhanian & Sathyaprakash, 2022; Ronchini et al., 2022), whose
thorough analysis of the prospects for CE and ET are aligned with the bounds
by Iacovelli et al. (2022) (but these authors use slightly different models or
local merger rates). In the near future, we plan to further update our catalog
with these and other results (e.g., Petrov et al. 2022).

We show a summary of the possible expected observations in the coming decades
in Figure 4. We emphasize again that these are very rough estimates for the
detection number and uncertainty, and provide publicly available code that
allows the user to change these assumptions. For outreach purposes it can be
fun to show a video with a linear y-axis to emphasize the exponential growth.
We provide this video online and show a snapshot of the video in Figure 5.

Figure 4: Number of detected mergers expected from binary black holes (BH-BH; black), black hole-neutron stars (BHNS; teal), and binary neutron stars (NS-NS; orange) systems for different data analysis pipelines and assumptions in the literature. The exact number of detections depends on the data analysis pipeline and is sensitive to the chosen threshold of assigning gravitational-wave events to be of (likely) astrophysical origin. See http://www.broekgaarden.nl/floor/wordpress/elementor-967/ for an animated version. Figure 5: Snapshots of a movie that visualizes the growth of the gravitational-wave catalog in the coming decades on a linear scale. See http://www.broekgaarden.nl/floor/wordpress/elementor-967/ for the full animation. Table 2: Summary of Detection numbers for O4 and O5. Acronyms are: FAR (false alarm rate), GWTC (Gravitational Wave Transient Catalog), GWOSC (Gravitational Wave Open Science Center). References: (a) taken from Abbott et al. (2018) using their updated table from https://emfollow.docs.ligo.org/ and scaling down by a factor 3.375 to convert to confident detections (representing SNR$\gtrsim 12$). (b) Numbers are taken from Tables 3, 4, and 5 from Iacovelli et al. (2022).  Label | Total | BHBH | BHNS | NSNS | Ref
---|---|---|---|---|---
Detections for O4 |  |  |  |  |
Abbott et al. (2018) PBL model: | 90 | 77 | 2 | 11 | (a)
Iacovelli et al. (2022): | 90 | 86 | 3 | 4 | (b)
Detections for O5 |  |  |  |  |
Abbott et al. (2018) PBL model: | 320 | 258 | 9 | 53 | (a)
Table 3: Summary of Detection numbers per year for CE/ET era. Acronyms are: FAR (false alarm rate), GWTC (Gravitational Wave Transient Catalog), GWOSC (Gravitational Wave Open Science Center). References: (a) Taken from Tables 3, 4, and 5 from Iacovelli et al. (2022), where we use the ‘Analysed’ events that have SNR$\gtrsim 12$. (b) taken from the lowest (pessimistic) and highest (optimistic) bounds in Section from Iacovelli et al. (2022).  Label | Total | BHBH | BHNS | NSNS | Ref
---|---|---|---|---|---
Detections for ET/CE |  |  |  |  |
Iacovelli et al. (2022) (ET): | 66 993 | 48456 | 9384 | 9153 | (a)
Iacovelli et al. (2022) (ET$+$2CE): | 138 656 | 69610 | 29462 | 39584 | (a)
Iacovelli et al. (2022) Pessimistic: | 33810 | $\bf 3.1\times 10^{4}$ | $\bf 1.9\times 10^{3}$ | $\bf 9.1\times 10^{2}$ | (b)
Iacovelli et al. (2022) Optimistic: | 794000 | $\bf 1.1\times 10^{5}$ | $\bf 3.4\times 10^{4}$ | $\bf 6.5\times 10^{5}$ | (b)

##  4 Discussion: beyond numbers

Importantly, future detectors such as O4, O5, Voyager, Einstein Telescope and
Cosmic Voyager will transform astronomy much beyond just increasing the number
of detections, which has been the focus of this paper. Indeed, these next-
generation detectors provide data with ever-increasing precision and out to
larger distances, pushing the reach of GW astronomy to the edge of the
observable Universe; revolutionizing our view of the cosmos by allowing
detections of gravitational-wave sources possibly beyond redshifts of ten.
This will allow us to address the redshift evolution of these sources, and
probe compact objects and their progenitor (stars) out to extremely high
redshifts. See Reitze et al. 2019 and Maggiore et al. 2020 and the references
therein for more details. Moreover, the future detectors will not just see
more events, but also many more events with higher precision (SNR), which will
allow scientists to infer the source properties with more precision (see e.g.,
the discussion and results presented in Abbott et al. (2021a); Borhanian &
Sathyaprakash (2022); Iacovelli et al. (2022); Ronchini et al. (2022). All in
all, there is no doubt that this will be some exciting decades to come with
many new discoveries and surprises along the way!

## Acknowledgements

FSB is grateful for the thesis writing procrastination (aka I have to hand in
a Ph.D thesis in seven days from today) that allowed her to write this work.
FSB is grateful for the many suggestions and feedback received for this work
from a list of people including Carl-Johan Haster, Christopher Berry, Karan
Jani, Alex Nitz, Sohen Mondhai, Eric Thrane, Ilya Mandel, Edo Berger, Yvette
Cendes, Maryam Hussaini, Lieke van Son. This work (and much more serious work)
was performed under a NASA FINESST scholarship 80NSSC22K1601, awarded to the
President and Fellows of Harvard College. FSB is also grateful for the support
on visualizing research received through the Harvard Horizons Scholarship
program.

## Data Availability

All data is publicly available, see
http://www.broekgaarden.nl/floor/wordpress/elementor-967/ for more details and
links.

## References

  * Aasi et al. (2015) Aasi J., et al., 2015, Class. Quant. Grav., 32, 074001
  * Abbott et al. (2016a) Abbott B. P., et al., 2016a, Phys. Rev. Lett., 116, 061102
  * Abbott et al. (2016b) Abbott B. P., et al., 2016b, Phys. Rev. Lett., 116, 131103
  * Abbott et al. (2018) Abbott B. P., et al., 2018, Living Rev. Rel., 21, 3
  * Abbott et al. (2019) Abbott B. P., et al., 2019, Physical Review X, 9, 031040
  * Abbott et al. (2020) Abbott R., et al., 2020, ApJ, 896, L44
  * Abbott et al. (2021a) Abbott R., Abbott T. D., Abraham S., Acernese F., et al., 2021a, arXiv e-prints, p. arXiv:2108.01045
  * Abbott et al. (2021b) Abbott R., Abbott T. D., Abraham S., Acernese F., et al., 2021b, Physical Review X, 11, 021053
  * Abbott et al. (2021c) Abbott R., et al., 2021c, ApJ, 915, L5
  * Acernese et al. (2015) Acernese F., et al., 2015, Classical and Quantum Gravity, 32, 024001
  * Acernese et al. (2019) Acernese F., et al., 2019, Phys. Rev. Lett., 123, 231108
  * Adhikari et al. (2020) Adhikari R. X., et al., 2020, Classical and Quantum Gravity, 37, 165003
  * Akutsu et al. (2021) Akutsu T., et al., 2021, Progress of Theoretical and Experimental Physics, 2021, 05A102
  * Aso et al. (2013) Aso Y., Michimura Y., Somiya K., Ando M., Miyakawa O., Sekiguchi T., Tatsumi D., Yamamoto H., 2013, Phys. Rev. D, 88, 043007
  * Borhanian & Sathyaprakash (2022) Borhanian S., Sathyaprakash B. S., 2022, arXiv e-prints, p. arXiv:2202.11048
  * Callister & Farr (2023) Callister T. A., Farr W. M., 2023, arXiv e-prints, p. arXiv:2302.07289
  * Collaboration et al. (2022) Collaboration L. S., Collaboration V., Collaboration K., 2022, GWOSC Event Portal Snapshots, doi:10.5281/zenodo.7249086, https://doi.org/10.5281/zenodo.7249086
  * Iacovelli et al. (2022) Iacovelli F., Mancarella M., Foffa S., Maggiore M., 2022, ApJ, 941, 208
  * Maggiore et al. (2020) Maggiore M., et al., 2020, J. Cosmology Astropart. Phys., 2020, 050
  * Ng et al. (2021) Ng K. K. Y., Vitale S., Farr W. M., Rodriguez C. L., 2021, ApJ, 913, L5
  * Nitz et al. (2020) Nitz A. H., et al., 2020, ApJ, 891, 123
  * Nitz et al. (2021a) Nitz A. H., Kumar S., Wang Y.-F., Kastha S., Wu S., Schäfer M., Dhurkunde R., Capano C. D., 2021a, arXiv e-prints, p. arXiv:2112.06878
  * Nitz et al. (2021b) Nitz A. H., Capano C. D., Kumar S., Wang Y.-F., Kastha S., Schäfer M., Dhurkunde R., Cabero M., 2021b, ApJ, 922, 76
  * Petrov et al. (2022) Petrov P., et al., 2022, ApJ, 924, 54
  * Reitze et al. (2019) Reitze D., et al., 2019, in Bulletin of the American Astronomical Society. p. 35 (arXiv:1907.04833), doi:10.48550/arXiv.1907.04833
  * Ronchini et al. (2022) Ronchini S., et al., 2022, A&A, 665, A97
  * Somiya (2012) Somiya K., 2012, Classical and Quantum Gravity, 29, 124007
  * The LIGO Scientific Collaboration et al. (2021a) The LIGO Scientific Collaboration et al., 2021a, arXiv e-prints, p. arXiv:2111.03606
  * The LIGO Scientific Collaboration et al. (2021b) The LIGO Scientific Collaboration et al., 2021b, arXiv e-prints, p. arXiv:2111.03634
  * Udall et al. (2020) Udall R., Jani K., Lange J., O’Shaughnessy R., Clark J., Cadonati L., Shoemaker D., Holley-Bockelmann K., 2020, ApJ, 900, 80
  * van Son et al. (2021) van Son L. A. C., et al., 2021, arXiv e-prints, p. arXiv:2110.01634
  * Venumadhav et al. (2019) Venumadhav T., Zackay B., Roulet J., Dai L., Zaldarriaga M., 2019, Phys. Rev. D, 100, 023011
  * Venumadhav et al. (2020) Venumadhav T., Zackay B., Roulet J., Dai L., Zaldarriaga M., 2020, Phys. Rev. D, 101, 083030
  * Zackay et al. (2019a) Zackay B., Dai L., Venumadhav T., Roulet J., Zaldarriaga M., 2019a, arXiv e-prints, p. arXiv:1910.09528
  * Zackay et al. (2019b) Zackay B., Venumadhav T., Dai L., Roulet J., Zaldarriaga M., 2019b, Phys. Rev. D, 100, 023007
  * Zevin et al. (2021) Zevin M., et al., 2021, ApJ, 910, 152

  *[GW]: gravitational wave
  *[BHBH]: binary black hole
  *[NSNS]: binary neutron star
  *[BHBHs]: binary black hole
  *[BHNSs]: black hole–neutron star
  *[NSNSs]: binary neutron star
  *[BHNS]: black hole–neutron star