# Quantum state preparation of gravitational waves

Fergus Hayes ID<EMAIL_ADDRESS>SUPA, School of Physics and
Astronomy, University of Glasgow, Glasgow G12 8QQ, UK  Sarah Croke ID SUPA,
School of Physics and Astronomy, University of Glasgow, Glasgow G12 8QQ, UK
Chris Messenger ID SUPA, School of Physics and Astronomy, University of
Glasgow, Glasgow G12 8QQ, UK  Fiona C. Speirits ID SUPA, School of Physics and
Astronomy, University of Glasgow, Glasgow G12 8QQ, UK

###### Abstract

We detail a quantum circuit capable of efficiently encoding analytical
approximations to gravitational wave signal waveforms of compact binary
coalescences into the amplitudes of quantum bits using both quantum arithmetic
operations and hybrid classical-quantum generative modelling. The gate cost of
the proposed method is considered and compared to a state preparation routine
for arbitrary amplitudes, where we demonstrate up to a four orders of
magnitude reduction in gate cost when considering the encoding of
gravitational waveforms representative of binary neutron star inspirals
detectable to the Einstein telescope. We demonstrate through a quantum
simulation, that is limited to $28$ qubits, the encoding of a second post-
Newtonian inspiral waveform with a fidelity compared to the desired state of
$0.995$ when using the Grover-Rudolph algorithm, or $0.979$ when using a
trained quantum generative adversarial network with a significant reduction of
required gates.

##  I Introduction

The Laser Interferometer Gravitational-Wave Observatory (LIGO) made the first
detection of a gravitational wave in September of 2015 from a binary black-
hole merger dubbed GW150914 [1], which began the era of gravitational wave
astronomy. Since then, the LIGO-Virgo-KAGRA collaboration has catalogued over
90 additional detections [2, 3, 4] including neutron star-black hole mergers
[5] as well as the first binary neutron star merger event GW170817 [6] that
was detected in conjunction with electromagnetic observations [7, 8] — the
first time that these two observation channels had coincided for a single
event. However, there is still much that gravitational wave astronomy can
offer us: from the discoveries of gravitational waves from sources such as
continuous gravitational waves [9], supernovae [10], the stochastic background
[11] and cosmic strings [12]; to testing our understanding of cosmology [13]
and general relativity [14]. To realize these aims will require future
generations of gravitational wave [15, 16, 17] and electromagnetic detectors
[18] to be constructed, the observations from which will need to be processed
by ever more sophisticated data analysis methods.

Quantum computing as a field has grown independently but alongside the
development of gravitational wave detectors [19]. It is well known that
quantum computers offer the potential of performing some computational tasks
faster than their classical counterparts, such as prime factorisation using
Shor’s algorithm [20] and unstructured searches using Grover’s algorithm [21].
Since the conception of the idea in the 1980’s, quantum computing has now
progressed to the widespread availability of noisy quantum processors.
Repeated, real time error-correction for a single logical qubit has recently
been demonstrated experimentally, providing a first step towards fault
tolerant quantum computing devices [22, 23].

Given the development of this quantum technology and the potential of quantum
algorithms to speed up some computational tasks, it is fitting that
astronomers can now look to see whether quantum computing can help overcome
current and future data analysis challenges. In previous work, it has been
shown that the standard signal detection analysis for gravitational waves,
matched filtering [24], can acquire a quadratic speed-up using Grover’s
algorithm [25], a procedure later built upon in [26]. Bayesian inference of
gravitational wave parameters has also been shown to acquire a polynomial
speed-up [27] with a quantum version of the famous Metropolis-Hastings
algorithm [28].

There are a number of reasons why gravitational wave data analysis is a
particularly interesting case for exploring applications for quantum
computing. Firstly, it is common for data analysis methods to be proposed that
rely upon future technology that is decades in the making due to the
challenges involved in the instrumentation of gravitational wave detections.
Second, the computational expense that gravitational wave data analysis often
faces is not due to the quantity of data but in the size of the solution
space, as opposed to big data problems where techniques such as dimensionality
reduction and data compression are required. Although quantum algorithms with
super-polynomial speed-ups and potential applications to big data problems
have been proposed [29, 30], we must be careful to account for the cost of
loading classical data into quantum states [31]. In the worst case, this can
be sufficient to negate any quantum advantage. Gravitational wave data
analysis problems are not of this flavour. Lastly, gravitational wave signals
are well modelled through general relativity, allowing their waveforms to be
efficiently prepared in quantum states. It is this final aspect that is the
focus of this paper, where we demonstrate how analytical gravitational
waveforms can be encoded into quantum states.

There are different methods for encoding information into quantum states.
Digital encoding involves storing in the computational basis of a string of
qubits, just as information is stored on bit strings in the classical case. A
denser format is angular encoding where each qubit can store up to two real
values as rotation angles about the Bloch sphere [32]. The densest form of
information representation is amplitude encoding, where the components of a
classical vector are stored as the amplitudes of a superposition of quantum
states, the number of which grows exponentially with the number of qubits.
Amplitude encoded data is an essential assumption in a multitude of quantum
machine learning algorithms [33, 34, 35, 36, 37, 38, 39, 40]. The efficient
amplitude encoding of vectors is also a required subroutine of quantum
algorithms for solving linear systems of equations [29, 41, 42]. While
amplitude encoding is the most space efficient, the computational cost of
preparing an arbitrary state scales exponentially with the number of qubits
[43], making it extremely inefficient to perform [31]. It has however been
demonstrated that amplitude encoding can be performed efficiently to prepare
specific states [44, 45].

In this work, we demonstrate by explicit construction that quantum states that
amplitude encode analytic gravitational waveforms for inspiralling binary
systems can be prepared efficiently. This is a crucial first step to study
space-efficient quantum algorithms for gravitational wave data analysis, which
may be implementable on small, near-term quantum devices. We anticipate that
our state preparation algorithm may find application as a subroutine to load
templates to e.g. train a variational quantum circuit to distinguish signal
from noise [46], or in the longer term to construct an oracle for Grover’s
search using amplitude encoded data [25], or for more sophisticated machine
learning techniques [29].

We outline the proposed state preparation routine in Sec. II, where we also
introduce the necessary background from the literature on the subroutines used
in the procedure. The state preparation routine is demonstrated with an
example on a quantum simulator in Sec. III. We relate our state preparation
routine to the challenge of analysing binary neutron star inspiral data faced
during the third generation detector in Sec. IV. The conclusion is provided in
Sec. V.

##  II State preparation algorithms

Our central goal in this paper is to introduce methods for loading
gravitational wave templates into a quantum register, using amplitude
encoding. To achieve this we make use of the fact that analytical expressions
that closely approximate the waveforms of interest may be derived from general
relativity. In this section we begin by outlining our approach to state
preparation, which first prepares the amplitude and then the phase of the
desired coefficients. For each of these steps we draw on existing techniques,
and the remainder of the section is then dedicated to introducing the
background needed for our implementation.

###  II.1 Overview of state preparation procedure

An arbitrary complex vector $\vec{h}$ of length $N$ has components
$\\{\tilde{A}(0)e^{i\Psi(0)},\dots,\tilde{A}(N-1)e^{i\Psi(N-1)}\\}$ where
$\tilde{A}(i)$ and $\Psi(i)$ are real numbers. These components may be
represented as amplitudes of a quantum state $|h\rangle$ given by:

$|h\rangle=\frac{1}{\|\tilde{A}\|}\sum^{2^{n}-1}_{j=0}\tilde{A}(j)e^{i\Psi(j)}|j\rangle,$
(1)

where $|j\rangle$ are computational basis states of
${n=\lceil\log_{2}N\rceil}$ qubits and $\|\cdot\|$ represents the norm of the
vector. The state $|h\rangle$ is then said to be an amplitude encoding of the
normalised vector $\vec{h}$ onto the $n$ qubits. In general preparing an
arbitrary state of this form requires exponential resources [47], however if
these coefficients are given by functions with analytical expressions, there
exist more efficient procedures, which we exploit in this work.

To perform the amplitude encoding of Eq. 1, our strategy will be to divide the
task into two steps. In the first, we construct an operator $\hat{U}_{A}$, the
role of which is to encode the real amplitudes $\tilde{A}(j)$ such that:

$\displaystyle|F\rangle$ $\displaystyle=\hat{U}_{A}|0\rangle^{\otimes n}$
$\displaystyle=\frac{1}{\|\tilde{A}\|}\sum^{2^{n}-1}_{j=0}\tilde{A}(j)|j\rangle.$
(2)

We propose two methods to implement this operator as a quantum circuit: either
by using the Grover-Rudolph algorithm [44], or by using a trained
parameterised quantum circuit [48], each described below.

In the second step, we construct an operator $\hat{U}_{\Psi}$, the role of
which is to prepare the phases $\Psi(j)$. Thus, the action on computational
basis states $|j\rangle$ is such that:

$\hat{U}_{\Psi}|j\rangle=e^{i\Psi(j)}|j\rangle.$ (3)

We show below by explicit construction that if the phase is given analytically
as a function of $j$ this operator may be implemented efficiently as a quantum
circuit. Our procedure is straight-forward and closely related to methods used
in the quantum Fourier transform [49].

Figure 1: Illustration of the Grover-Rudolph algorithm applied to two qubits
to encode $p(f)\propto f^{-7/3}$ (shown as the solid curved line) into the
quantum state amplitudes. Initially no qubits are encoded and therefore $m=0$
and there is only one bin $j=0$. The domain is then divided into two and the
ratio of the left-most shaded region under the distribution to the whole
domain is calculated as $\cos^{2}\zeta_{0,0}$. The $m=0$ qubit is then put
into a superposition conditioned on the corresponding $\zeta_{0,0}$ to produce
$m=1$. The process is repeated for $m=1$ where $j=0,1$ and each region is
split proportionally by $\cos^{2}\zeta_{1,0}$ and $\cos^{2}\zeta_{1,1}$
respectively to acquire $m=2$. Note that a continuous distribution is plotted
for illustrative purposes; the discretised probability distribution is
considered for formulating $\zeta_{m,j}$ in Eq. 5.

###  II.2 Grover-Rudolph algorithm

The Grover-Rudolph algorithm [44] is one way to implement the operator
$\hat{U}_{A}$. This algorithm can efficiently prepare a target quantum state
$|\psi\rangle$ of the form:

$|\psi\rangle=\sum^{2^{n}-1}_{j=0}\sqrt{p(j)}|j\rangle,$ (4)

given probability mass function $p$.

The Grover-Rudolph algorithm is an iterative process applied to each
${m=0,\dots,n-1}$ qubits, where each step involves dividing the target
distribution into ${j=0,\dots,2^{m}-1}$ bins. For the $m$th step, each $j$th
bin is divided into two further bins of equal width, and the corresponding
probability amplitude divided between the two new bins is calculated and
stored in an ancillary register. This is illustrated in Fig. 1. The fraction
of the probability residing in the leftmost half of the $j$th bin is equal to:

$\cos^{2}\zeta_{m,j}=\frac{\sum_{i=j2^{n-m}}^{(j+1/2)2^{n-m}}p(i)}{\sum_{l=j2^{n-m}}^{(j+1)2^{n-m}}p(i)}.$
(5)

Consider the operation $\hat{Q}^{(m)}_{\zeta}$ that calculates the angle
$\zeta_{m,j}$, storing this in the computational basis of an ancillary
register $|0\rangle_{\text{a}}$:

$\hat{Q}_{\zeta}^{(m)}\left(\sum_{j=0}^{2^{m}-1}\sqrt{p^{(m)}_{j}}|j\rangle\right)|0\rangle_{\text{a}}=\sum_{j=0}^{2^{m}-1}\sqrt{p^{(m)}_{j}}|j\rangle|\zeta_{m,j}\rangle_{\text{a}},$
(6)

where $p^{m}_{j}$ denotes the probability mass function $p$ integrated across
the $j$th bin given $m$ qubits. Controlled $Y$ rotations
$\hat{R}_{Cy}^{(m+1)}$ are then applied onto qubit $m+1$ such that:

$\displaystyle\hat{R}_{Cy}^{(m+1)}\sqrt{p^{(m)}_{j}}|j\rangle|\zeta_{m,j}\rangle_{\text{a}}|0\rangle_{m+1}$
$\displaystyle=\sqrt{p^{(m)}_{j}}|j\rangle|\zeta_{m,j}\rangle_{\text{a}}(\cos\zeta_{m,j}|0\rangle_{m+1}+\sin\zeta_{m,j}|1\rangle_{m+1})$
$\displaystyle=\sqrt{p^{(m+1)}_{j}}|j\rangle|\zeta_{m,j}\rangle_{\text{a}}.$
(7)

The ancilla register is then cleared by uncomputing $\zeta_{m,j}$ using the
inverse operation $\hat{Q}^{(m)}_{\zeta}$ of Eq. 6. The operations of Eq. 6
and Eq. II.2 can be repeated for $m$ incremented by one, until ${m=n-1}$, and
the state of Eq. 4 is produced.

While this routine allows for efficient state preparation of the $n$ qubit
amplitude states, for practical use we will be interested not just in the
scaling of the cost with $n$, but also the absolute number of gates used. The
computational cost for low $m$ iterations may induce overhead costs that are
much greater than simply using a general state preparation routine for loading
arbitrary states [47]. Therefore we can define a critical value $m_{a}$ where
for $m<m_{a}$ a routine for arbitrary amplitude state preparation is applied
instead that we can denote $R_{y}^{(m<m_{a})}$. This routine can be further
parameterised to reduce computational cost by defining a maximum iteration
$m_{b}$. As $\lim_{m\to\infty}\zeta_{m,j}=\pi/4$ for all $j$ in Eq. 5
(assuming continuity of the probability mass function $p$), then replacing the
operation described in Eq. II.2 by applying Hadamard gates to the qubits for
operations ${m>m_{b}}$ reduces the gate cost at the expense of the accuracy of
the final amplitude states. The overall encoding of the frequency amplitude is
performed by the time-ordered product

$\hat{U}_{A}=\hat{H}^{(m>m_{b})}\left(\prod^{m_{b}-2}_{m=m_{a}}\hat{Q}_{\zeta}^{{\dagger}(m)}\hat{R}_{Cy}^{(m+1)}\hat{Q}_{\zeta}^{(m)}\right)\hat{R}_{y}^{(m<m_{a})}$
(8)

where $\hat{H}^{(m>m_{b})}$ denotes Hadamard gates applied onto qubits
corresponding to $m>m_{b}$.

###  II.3 Quantum generative modelling

Another approach to preparing a state of the form of Eq. 4 is to train a
parameterized quantum circuit $\hat{U}(\phi)$ using a hybrid quantum-classical
generative machine learning model given parameters $\phi$. This is done by
measuring an ensemble of outputs from the parameterized quantum circuit and
utilizing a classical loss function to compare the samples from the measured
distribution $q_{\phi}$ to those from the target distribution $p$ which can be
classically generated. This requires the parameterized quantum circuit to be
chosen so that it is constrained to only result in real amplitudes. We choose
a parameterized circuit of $L$ repeating layers, defined by:

$\displaystyle\hat{U}(\phi)=$
$\displaystyle\prod_{j=0}^{n-1}\hat{R}_{y}^{(j)}(\phi_{L,j})$
$\displaystyle\prod_{i=0}^{L-1}\left(\prod_{k=0}^{n-2}\hat{X}^{(k,k+1)}_{C}\prod_{j=0}^{n-1}\hat{R}_{y}^{(j)}(\phi_{i,j})\right)\hat{H}^{\otimes
n},$ (9)

where parameter $\phi_{i,j}$ is the parameter for layer $i$ and qubit $j$,
$\hat{R}^{(j)}_{y}$ is a $Y$ rotation on qubit $j$ by the given angle, and
$\hat{X}^{(k,k+1)}_{C}$ is an $X$ gate on qubit $k+1$ controlled on qubit $k$.

We explore the use of a hybrid quantum-classical generative adversarial
network to train this parameterised quantum circuit [48]. Generative
adversarial networks consider two competing networks: one, called the
discriminator and dependent on parameters $\omega$, is trained to discern an
ensemble of samples from $q_{\phi}$ from those taken from $p$. The other,
called the generator and dependent on parameters $\phi$, is trained to produce
samples $q_{\phi}$ to trick the discriminator into falsely labeling them as
being products from $p$. As both networks have competing objectives, the
training of both simultaneously leads to ${q_{\phi}\approx p}$ as Nash
equilibrium is reached [50]. Quantum generative adversarial networks use a
quantum parameterized circuit in place of the generator as discussed in Ref.
[51]. The discriminator given parameters $\omega$ outputs $D_{\omega}$, where
${0<D_{\omega}<1/2}$ indicates the samples were drawn from $q_{\theta}$, while
${1/2<D_{\omega}<1}$ indicates the samples were drawn from $p$. The training
involves minimizing the generator loss function given $S$ samples
${\\{x_{1},\dots,x_{S}\\}\sim q_{\phi}}$:

$L_{G}=-\frac{1}{S}\sum_{i=1}^{S}\log D_{\omega}(x_{i}),$ (10)

and maximizing the discriminator loss function given $S$ samples
${\\{x_{1}^{\prime},\dots,x_{S}^{\prime}\\}\sim p}$:

$L_{D}=\frac{1}{S}\sum_{i=1}^{S}[\log
D_{\omega}(x_{i}^{\prime})+\log(1-D_{\omega}(x_{i}))].$ (11)

###  II.4 Phase preparation

To encode the phase information we wish to construct an operator
$\hat{U}_{\Psi}$, as defined in Eq. 3. Given an analytic expression for
$\Psi(j)$, which may be computed efficiently classically, we begin by
evaluating $\Psi^{\prime}(j)=\Psi(j)/2\pi$ and storing it as a binary string
in the computational basis of an ancilla register. We define the operator
implementing this as $\hat{Q}_{\Psi}$:

$\hat{Q}_{\Psi}|j\rangle|0\rangle_{\text{a}}=|j\rangle|\Psi^{\prime}(j)\rangle_{\text{a}}.$
(12)

As this operation can be computed efficiently classically, it follows that
$\hat{Q}_{\Psi}$ may be implemented efficiently as a quantum circuit.

The frequency dependent phase can then be readily produced by $Z$ rotations
$\hat{R}_{z}$ applied to the qubits of the ancillary register. Specifically,
any number $x$ stored in $n$-bits, $p$ of which are precision qubits, has
binary representation:

$x=\sum^{n-1}_{i=0}x_{i}2^{i-p}.$

This may be stored in the computational basis of an $n$-qubit register as the
state

$|x\rangle=|x_{n-1}\rangle\otimes|x_{n-2}\rangle\ldots\otimes|x_{0}\rangle$

Applying single qubit $Z$ rotations to the precision qubits, where
$\hat{R}^{(j)}(2^{j-p+1}\pi)$ represents a rotation by $2^{j-p+1}\pi$ applied
to the $j$-th qubit gives:

$\displaystyle\prod e^{2^{j-p+1}\pi
x_{j}}|x_{n-1}\rangle\otimes|x_{n-2}\rangle\ldots|x_{0}\rangle$
$\displaystyle=e^{2\pi\sum_{j=0}^{p-1}x_{j}2^{j-p}}|x\rangle=e^{2\pi
x}|x\rangle$ (13)

. Thus we need only apply $p_{a}$ single qubit gates (where $p_{a}$ is the
number of precision qubits used in the ancilliary register) to prepare:

$\prod^{p_{a}}_{j=1}\hat{R}_{z}^{(p_{a}-j)}(2^{1-j}\pi)|\Psi^{\prime}(j)\rangle_{\text{a}}=e^{i\Psi(j)}|\Psi^{\prime}(j)\rangle_{\text{a}}.$
(14)

Finally, the ancillary register is cleared by uncomputing the calculation of
$\Psi^{\prime}(j)$ with the inverse operation $\hat{Q}^{{\dagger}}_{\Psi}$
giving

$\hat{U}_{\Psi}=\hat{Q}_{\Psi}^{{\dagger}}\prod^{p_{a}}_{j=1}\hat{R}_{z}^{(p_{a}-j)}(2^{1-j}\pi)\hat{Q}_{\Psi}.$
(15)

###  II.5 Linear piece-wise function

$\Qcircuit@C=.3em@R=.02em@!{\lstick{|0\rangle_{l}}&\multigate{1}{\text{Label}}\ctrl{1}\qw\ctrl{1}\ctrl{1}\qw\qw\\\
\lstick{|0\rangle_{\text{c}}}\ghost{\text{Label}}\gate{X_{1}}\multigate{2}{\text{Mult}}\gate{X_{1}^{{\dagger}}}\gate{X_{0}}\multigate{1}{\text{Add}}\qw\\\
\lstick{|0\rangle_{o}}\qw\qw\ghost{\text{Mult}}\qw\qw\ghost{\text{Add}}\rstick{|f(x)\rangle_{o}}\qw\\\
\lstick{|x\rangle_{x}}\ctrl{-2}\qw\ghost{\text{Mult}}\qw\qw\qw\qw\\\ }$
Figure 2: Quantum circuit of the piece-wise linear function operation
$\hat{Q}_{f}$ approximating $f(x)\approx A_{1}^{l}x+A_{0}^{l}$ for values of
$x$ within each of the $2^{n_{l}}$ sub-domains of $f(x)$. The label gate
correlates the sub-domains to the label register $|l\rangle_{l}$ before
loading in the coefficients $A_{1}^{l}$ for each sub-domain into the
coefficient register using controlled operations $X_{1}$. The input register
$|x\rangle_{x}$ is then multiplied by the coefficient register with the
arithmetic multiplication operation and the result is saved in the output
register. The coefficient register is cleared before the zeroth order
coefficients are loaded with operation $X_{0}$ before they are added to the
output register with an arithmetic addition operation.

For both Eq. 6 of the Grover-Rudolph algorithm and Eq. 12 of the phase
preparation, we require a quantum circuit that is able to efficiently input
and evaluate functions $f(x)$ for a given $x$. Any function that is computable
efficiently classically may be directly implemented as a quantum circuit,
performing the function evaluation in the computational basis and writing the
result in an ancilla register:

$|x\rangle|0\rangle\rightarrow|x\rangle|f(x)\rangle$

One way to implement function evaluation in practice is by a piece-wise
function approximation, described in Ref. [52], which has previously been
proposed as a way to perform the Grover-Rudolph algorithm efficiently in Ref.
[53]. This involves dividing the domain into $2^{n_{l}}$ sub-domains and
approximating the function in each sub-domain by a polynomial of a given
order. The Remez Algorithm [54] is applied to determine the polynomial
coefficients to the given order that minimizes the $L_{\infty}$ error across
the given domain. This classical pre-processing step is performed with a cost
of $O(2^{n_{l}})$. The domains of the function are then correlated with a
label register $|0\rangle^{\otimes n_{l}}_{l}$ using the label gate described
in Ref. [52]. Here we consider the simplest case where the function $f(x)$ is
approximated to be linear within the given sub-domain over $x$, requiring the
zeroth and first order polynomial coefficients $A_{0}^{l}$ and $A_{1}^{l}$
such that $f(x)\approx A_{1}^{l}x+A_{0}^{l}$ for $x$ in the given sub-domain.
The $x$ argument is given to the input register $|x\rangle_{x}^{\otimes
n_{x}}$ of size $n_{x}$ and the coefficients for each sub-domain are loaded
into the coefficient register $|0\rangle^{\otimes n_{c}}_{c}$, where $n_{c}$
qubits are used. This requires the output register in which the outcome $f(x)$
is stored in $n_{\text{o}}$ qubits. Generally, when approximating a linear
piece-wise function, we introduce the operation $\hat{Q}_{f}$, performing the
operations:

$\displaystyle|x\rangle_{x}|0\rangle_{o}|A_{1}^{l}\rangle_{c}|l\rangle_{l}$
$\displaystyle\xrightarrow{\text{Mult}}|x\rangle_{x}|A_{1}^{l}x\rangle_{o}|A_{1}^{l}\rangle_{c}|l\rangle_{l}$
$\displaystyle\xrightarrow{\text{Load}}|x\rangle_{x}|A_{1}^{l}x\rangle_{o}|A_{0}^{l}\rangle_{c}|l\rangle_{l}$
$\displaystyle\xrightarrow{\text{Add}}|x\rangle_{x}|A_{1}^{l}x+A_{0}^{l}\rangle_{o}|A_{0}^{l}\rangle_{c}|l\rangle_{l}.$

The quantum circuit used to perform this action is shown in Fig. 2, where the
zeroth and first order coefficients are loaded and unloaded using gates
$X_{0,1}$ and $X_{0,1}^{{\dagger}}$ respectively. The gate and space cost of
applying the piece-wise linear function depends on the choice of adder and
multiplier routine. The coefficients of a given domain are loaded and unloaded
into the coefficient register using controlled $X$ gates conditioned on the
label register.

##  III Encoding gravitational wave inspirals

In this section we detail how the techniques outlined in the previous section
are applied to encode the waveforms of gravitational wave binary system
inspirals. An analytical expression of the frequency dependent amplitude of an
inspiral waveform is obtained using the stationary-phase approximation (see
Appendix A), resulting in the expression [55]:

$\tilde{A}_{N}(f)=\frac{Q\mathcal{M}^{5/6}}{D}f^{-7/6}.$ (16)

where $Q$ is dependent on the geometry of the detector and source system,
$\mathcal{M}$ is the chirp mass, and $D$ is the luminosity distance to the
source. Similarly, the frequency dependent phase to a second post-Newtonian
order waveform is of the form:

$\Psi_{\text{2PN}}(f)=2\pi
ft_{c}-\phi_{c}-\frac{\pi}{4}+\frac{3}{128}(\pi\mathcal{M}f)^{-5/3}\left[1+\frac{20}{9}\left(\frac{743}{336}+\frac{11}{4}\eta\right)(\pi
Mf)^{2/3}-4(4\pi-\beta)(\pi Mf)\right.\\\
\left.+10\left(\frac{3058673}{1016064}+\frac{5429}{1008}\eta+\frac{617}{144}\eta^{2}-\sigma\right)(\pi
Mf)^{4/3}\right].$ (17)

where $t_{c}$ is the time of coalescence, $\phi_{c}$ is the phase at
coalescence, $M$ is the total mass, $\eta$ is the reduced mass, $\sigma$ is
the spin-spin parameter, and $\beta$ is the spin-orbit parameter. For further
details on these parameters, we refer the reader to Appendix B.

Consider encoding the frequency dependent waveform over a frequency range with
lower frequency cut-off $f_{\text{min}}$ and upper frequency $f_{\text{max}}$,
discretised into ${N=2^{n}}$ frequency bins of widths ${\Delta
f=(f_{\text{max}}-f_{\text{min}})/2^{n}}$. This gives the number of integer
qubits to be ${n_{\text{int}}=\lceil\log_{2}(2^{n}\Delta f)\rceil}$ and
therefore the number of precision bits as ${p=\lceil\log_{2}T\rceil}$, where
$T$ is the waveform’s temporal duration ${T=\Delta f^{-1}}$. The $j$th
computational basis state has a probability amplitude equal to the integrated
real classical frequency amplitude within the frequency range
${f_{\text{min}}+[j\Delta f,(j+1)\Delta f)}$.

While the Newtonian waveform amplitude depends on the detector response,
chirp-mass and source distance, only the frequency dependence is required in
the encoding due to the normalisation, and therefore these terms are ignored
such that ${\tilde{A}(f)\propto f^{-7/6}}$.

To demonstrate the state-preparation routine of Sec. II, we prepare a waveform
with a Newtonian amplitude of Eq. 16 and a second post-Newtonian phase of Eq.
17, given a spinless (${\beta=\sigma=0}$), near equal-mass system of
${m_{1}=35\,M_{\odot}}$ and ${m_{2}=30\,M_{\odot}}$ using IBM’s quantum
simulator and qiskit software [56]. This simulation is carried out on Python
code that is publicly available on Github [57]. The simulation is run on a
commercial computer with 7.9 gigabytes of available memory, allowing for
simulations of up to 28 qubits when allowing complex numbers to be stored on
16 bytes. This waveform is sampled within the frequency interval
${f_{\text{min}}=40}\,$Hz and ${f_{\text{max}}=168}\,$Hz with ${\Delta
f=2}\,$Hz, requiring ${n=6}$ qubits. A value of $\Delta T=0.02\,$s is chosen.

After the simulation, the output waveform state $|\psi_{\text{out}}\rangle$ is
compared to the classically computed target waveform state
$|\psi_{\text{targ}}\rangle$ and the fidelity $\mathcal{F}$ between the two
states is calculated, defined as
$|\langle\psi_{\text{out}}|\psi_{\text{targ}}\rangle|^{2}$, such that
${1-\sqrt{\mathcal{F}}}$ is synonymous with the mismatch associated between
templates as defined in [58] where the noise is flat across all frequencies.

(a)

(b)

Figure 3: The top panel shows the loss over training iterations for the
quantum generative adversarial network configuration for both the generator
shown as the black solid line as defined in Eq. 10, and the discriminator in
Eq. 11 as the grey solid line. The dialectic between both networks is
represented by the loss, as they exhibit opposing gradients about an
equilibrium of ${\sim}0.7$. The bottom panel shows the resulting generated
quantum state mismatch with the target state amplitude of $\tilde{A}(f)\propto
f^{-7/6}$ over iterations. $(a)$ shows the training of a parameterized quantum
circuit of $L=12$ with $78$ trainable parameters, which reaches a minimum
mismatch of $8.57{\times}10^{-3}$. $(b)$ shows the training of a parameterized
quantum circuit of $L=20$ with $126$ trainable parameters, which reaches a
minimum mismatch of $8.36{\times}10^{-3}$.

(a)

(b)

Figure 4: A scatter of simulated frequency state amplitude compared to the
target state amplitudes of $\tilde{A}(f)\propto f^{-7/6}$ shown as the solid
black line in the top panel. The bottom panel shows the relative difference
between the output and the target frequency states. $(a)$ shows the resulting
amplitudes from simulating the parameterized quantum circuit of Eq. II.3 after
the training illustrated in Fig. 3(a). The resulting state fidelity is $0.983$
with a mismatch of $8.6{\times}10^{-3}$. Discrepancies between the simulated
and target state are mainly due to the limited size of the parameterized
quantum circuit in this case. The number of controlled NOT gates used in the
simulation is $100$. $(b)$ shows the resulting amplitudes from simulating the
circuit described in Fig. 5 resulting in a state fidelity of $0.999$ and
mismatch of $4.1{\times}10^{-4}$. Discrepancies between the two are caused by
the limited number of qubits available for the arithmetic operations involved
in the Grover-Rudolph algorithm when invoking the piece-wise linear function
operation as well as the omission of second-order terms. The number of
controlled NOT gates used is $23,796$.

###  III.1 Amplitude preparation

$\Qcircuit@C=.0em@R=.0em@!{\lstick{|0\rangle_{l}}&\qw\multigate{7}{Q_{\zeta}^{(4)}}\qw\multigate{7}{Q_{\zeta}^{{\dagger}(4)}}\qw\\\
\lstick{|0\rangle_{c}}\qw\ghost{Q_{\zeta}^{(4)}}\qw\ghost{Q_{\zeta}^{{\dagger}(4)}}\qw\\\
\lstick{|0\rangle_{\text{a}}}\qw\ghost{Q_{\zeta}^{(4)}}\ctrl{6}\ghost{Q_{\zeta}^{{\dagger}(4)}}\qw\\\
\ustick{0}\multigate{4}{R_{Y}^{(m<5)}}\ghost{Q_{\zeta}^{(4}}\qw\ghost{Q_{\zeta}^{{\dagger}(4)}}\qw\\\
\ustick{1}\ghost{R_{Y}^{(m<5)}}\ghost{Q_{\zeta}^{(4)}}\qw\ghost{Q_{\zeta}^{{\dagger}(4)}}\qw\\\
\ustick{2}\ghost{R_{Y}^{(m<5)}}\ghost{Q_{\zeta}^{(4)}}\qw\ghost{Q_{\zeta}^{{\dagger}(4)}}\qw\\\
\ustick{3}\ghost{R_{Y}^{(m<5)}}\ghost{Q_{\zeta}^{(4)}}\qw\ghost{Q_{\zeta}^{{\dagger}(4)}}\qw\\\
\ustick{4}\ghost{R_{Y}^{(m<5)}}\ghost{Q_{\zeta}^{(4)}}\qw\ghost{Q_{\zeta}^{{\dagger}(4)}}\qw\\\
\ustick{5}\qw\qw\gate{R_{Y}^{(5)}}\qw\qw\\\ }$
Figure 5: The circuit of the amplitude preparation step using the Grover-
Rudolph algorithm where the steps for $m=0,1,2,3,4$ are performed by
controlled rotations of angles $\zeta_{m,j}$ for all $j=0,\dots,2^{m}-1$
denoted by operation $R_{y}^{(m<5)}$, while the operations for $m=5$ are
performed by the operations
$\hat{Q}^{{\dagger}(m-1)}_{\zeta}\hat{R}_{y}^{(m)}\hat{Q}_{\zeta}^{(m-1)}$.
$\Qcircuit@C=.3em@R=.2em@!{\lstick{|0\rangle_{l}}&\multigate{3}{Q_{\Psi}}\qw\multigate{3}{Q_{\Psi}^{{\dagger}}}\qw\\\
\lstick{|0\rangle_{c}}\ghost{Q_{\Psi}}\qw\ghost{Q_{\Psi}^{{\dagger}}}\qw\\\
\lstick{|0\rangle_{\text{a}}}\ghost{Q_{\Psi}}\gate{R_{z}}\ghost{Q_{\Psi}^{{\dagger}}}\qw\\\
\lstick{|F\rangle}\ghost{Q_{\Psi}}\qw\ghost{Q_{\Psi}^{{\dagger}}}\rstick{|h\rangle}\qw}$
Figure 6: Quantum circuit for the phase preparation step where the function
$\Psi^{\prime}$ is encoded into the computational basis of the ancillary
register $|0\rangle_{\text{a}}$ through the operation $\hat{Q}_{\Psi}$,
implemented through the linear piece-wise function using coefficient and label
registers $|0\rangle_{c}$ and $|0\rangle_{l}$. The phase $e^{i\Psi(j)}$ on
computational basis state $|j\rangle$ is produced through $Z$ rotations on the
precision bits of the binary string stored in the ancillary register. Finally,
the ancillary, label and coefficient registers are cleared by the inverse
operation $\hat{Q}_{\Psi}^{{\dagger}}$.

We demonstrate the preparation of the real amplitudes using both the Grover-
Rudolph algorithm method and a parameterized quantum circuit trained using a
generative adversarial network for this case.

####  III.1.1 Quantum generative adversarial network

First the quantum generative adversarial network described in Sec. II.3 is
applied for the case of $n=6$. This is done by assuming the parameterized
quantum circuit of Eq. II.3. The parameterized circuit is trained using
qiskit’s quantum generative adversarial network implementation with PyTorch
[59]. The classical learning is performed using the Adam optimizer with a
learning rate of $0.01$, and first and second momentum parameters of $0.7$ and
$0.999$. The training was run for $1,500$ iterations, taking $10,000$ samples
from the generated quantum state and the target distribution to pass to the
classical discriminator. The results of the training are displayed in Fig. 3
for a network size of ${L=12}$ and ${L=20}$ in Fig. 3(a) and Fig. 3(b)
respectively. For both cases, the loss function of both the generator and
discriminator are plotted in the top panel (black and grey respectively) over
training iterations, and can be seen to oscillate about an equilibrium of
${\sim}0.7$ as the decrease of loss of one network necessarily requires the
increase of loss of the other. The bottom panels display the corresponding
mismatch between states generated by the parameterised quantum circuit and the
target state at different points of the training. The decrease in the mismatch
between generated and target state corresponds to the loss functions of both
generator and discriminator converging to the equilibrium values. Both
networks converge to a similar mismatch of $8.57\times 10^{-3}$ and
$8.36\times 10^{-3}$ for the case of ${L=12}$ and ${L=20}$ respectively. While
not plotted, further training proved to result in a mismatch increase for both
cases. The mismatch of Fig. 3(a) has a more gradual decline than Fig. 3(b),
where the convergence to a desired solution occurs more suddenly and generally
after more training iterations. These two training cases exemplify a general
trend found that while increasing the network size can improve the mismatch,
the training becomes more sporadic and less stable. The work required to
stabilize this training for larger network sizes is beyond the scope of this
paper.

The resulting amplitudes after applying the trained parameterized quantum
circuit with $L=20$ is shown in Fig. 4(a) scattered in black and compared to
the target state shown as the solid black line, and the relative difference
between the two amplitudes scattered in the bottom panel. The fidelity between
the output state and the target state of $0.983$ is achieved, corresponding to
a mismatch of $8.36\times 10^{-3}$. The simulation is performed using $100$
controlled NOT gates.

####  III.1.2 Grover-Rudolph algorithm

The amplitude preparation step using the Grover-Rudolph algorithm described in
Sec. II.2 is applied to obtain the state described in Eq. 2. The circuit
performing the operation $\hat{U}_{A}$ is shown in Fig. 5. The operation of
Eq. 6 is performed using simple controlled $Y$ rotation and $X$ gates for
$m<5$, denoted $\hat{R}_{Y}^{(m<5)}$, and a piece-wise linear function
approximation applied to higher values. The size of the ancilla and
coefficient registers is chosen so that $n_{c}=9$ and $n_{\text{o}}=9$, and
the domain of the function is divided into $16$ sub-domains using $n_{l}=4$
qubits where the boundaries of each sub-domain are spread uniformly between
the domain of the function.

Figure 7: The top panel shows a scatter of the simulated values of
$\Psi^{\prime}$ stored in the computational basis in the ancillary register
for each frequency register state compared to the target function, plotted as
the solid line. The boundaries of the piece-wise linear function domains are
shown in the dashed vertical lines. The deviations between the values stored
in the ancillary register and the target function are scattered in the bottom
panel. The simulated states show a maximum deviation from the target function
of $|\Delta\Psi|<0.04$. The number of controlled NOT gates used in the
simulation is $9,464$.

(a)

(b)

Figure 8: The top panel shows the simulated real and imaginary amplitudes of
the frequency register states, scattered in black and grey respectively. The
simulated state amplitudes are compared to the target normalised gravitational
waveform, plotted in black and grey solid lines. The deviations between the
simulated and the target waveform are scattered in the bottom panel. $(a)$
shows the result after applying the trained parameterized quantum circuit to
approximate the amplitude preparation routine as described in Sec. III.1.1.
The simulated waveform has a fidelity of $0.979$ and corresponding mismatch of
$1.0{\times}10^{-2}$ with the target state. The number of controlled NOT gates
used in the simulation is $19,028$. $(b)$ gives the same result but after the
phase preparation routine described in Sec. III.2 after performing the Grover-
Rudolph amplitude preparation routine described in Sec. III.1.2 where the
simulated waveform has a fidelity of $0.995$ and corresponding mismatch of
$2.4\times 10^{-3}$ with the target state. The number of controlled NOT gates
used is $42,724$.

After the application of the circuit shown in Fig. 5, the amplitude of the
output state $|F\rangle$ is plotted in the top panel of Fig. 4(b) as black
dots and compared to the target state shown as the solid black line. The
fidelity between the output state and the target state is calculated to be
$0.999$ corresponding to a mismatch of $4.1{\times}10^{-4}$, while the
relative difference in amplitude between each output and target frequency
state is plotted in the bottom panel. Deviations from the target state are due
to the limited number of ancillary qubits to store the rotation angles
$\zeta_{m,j}$, load the polynomial coefficients and define the sub-domains of
the function, as well as the omission of higher order terms from the linear
function approximation. The simulation is performed using $23,796$ controlled
NOT gates.

###  III.2 Phase preparation

For the phase preparation step, the operation $\hat{Q}_{\Psi}$ of Eq. 12 is
again applied by the piece-wise linear function approximation described in
Sec. II.5. The sizes of the coefficient and ancillary registers are set to
$n_{c}=8$ and $n_{\text{a}}=10$ given the piece-wise function coefficients,
while the label register size remains unchanged with $n_{l}=4$. The resulting
joint state across both the frequency and ancillary register corresponds to a
superposition in which each term has a single binary string representing
$\Psi^{\prime}$ stored in the ancilla registers, correlated to the
corresponding frequency state in the frequency register. The $\Psi^{\prime}$
values for each of the binary strings are scattered in the top panel of Fig. 7
across their corresponding frequency states and compared to the target
function $\Psi_{\text{PN}}/2\pi$, plotted as the solid black line. The
resulting joint state amplitudes of Fig. 7 clearly follow a smooth continuous
function comparable to the target function. The bottom panel shows the
relative difference between the target function $\Psi_{\text{PN}}/2\pi$ and
those stored in the ancillary register, which deviate only such that
$|\Delta\Psi|<0.04$. The boundaries of the piece-wise polynomial function are
shown as vertical dotted lines and are uniformly spread across the space. This
step is simulated using $9,464$ controlled NOT gates.

The $\hat{R}_{z}$ rotations of Eq. 15 are applied to the ancilla register that
$\Psi^{\prime}$ in the computational basis to produce the frequency dependent
phase of Eq. 14. The ancilla register is then uncomputed of $\Psi^{\prime}$
with the inverse piece-wise linear function, costing an additional $9,464$
controlled NOT gates. With enough qubits to adequately account for the
required precision of the multiplication operation, the uncomputing of
$\Psi^{\prime}$ leaves the ancillary and label register in a singular state of
$|0\rangle_{l}|0\rangle_{\text{a}}$, which can now be discarded from the
circuit. Fig. 8(a) and Fig. 8(b) depicts the amplitudes of the final simulated
states compared to the target state for when the trained parameterised quantum
circuit of Fig. 3(b) and the Grover-Rudolph algorithm are used for amplitude
preparation steps respectively. The real and imaginary parts of the target
state are plotted as solid lines in their respective colours, with the
difference between the output and target state displayed in the bottom panel.
When the Grover-Rudolph algorithm is employed, the resulting state fidelity is
$\mathcal{F}=0.995$ with an associated mismatch of $2.4{\times}10^{-3}$ using
$42,724$ controlled NOT gates, while when the trained parameterised quantum
circuit is used, the state fidelity is $\mathcal{F}=0.979$ with a mismatch of
$1.0{\times}10^{-2}$ using $19,028$ controlled NOT gates.

##  IV Case study: third generation gravitational wave detectors

Figure 9: Ratio of gate cost over duration when considering an arbitrary state
preparation routine against the upper bound gate cost (see Appendix E) of
preparing a insprial waveform using the method described in Sec. II with
${n_{c}=16}$, ${n_{l}=4}$ and ${n_{a}=n+n_{c}}$ using the Grover-Rudolph
amplitude preparation routine plotted in black, and using a parameterised
quantum circuit with ${L=100}$ in grey. Long duration waveforms that are
required for third generation detectors, which probe the signals at low
frequencies, can be produced by a factor of up to four orders of magnitude
less gates.

We consider the case of waveforms required for the application of data
analysis methods given future detector data, with higher sensitivities than
the detectors that are currently in operation. Detectors such as the Einstein
Telescope will probe gravitational wave emissions at frequencies of as low as
$1{-}7$ Hz, where sources such as binary neutron star pairs can be apparent in
the sensitive band for time periods of orders of days, a four orders of
magnitude increase compared to current detectors [60]. The longer observation
period allows for the orbital evolution of these systems to be studied in
detail [15]. However, the increased sensitivity at lower frequencies for third
generation detectors provide unique data analysis challenges and increases
computational demand [61]. The longer duration signals require longer
waveforms to perform coherent analyses, which is required to maximise
sensitivity when considering a multi-detector network. If a sampling rate of
$4096$ Hz is considered, the memory requirements to store a single binary
neutron star inspiral waveform can extend into the gigabytes. The majority of
the matched filter signal-to-noise ratio accumulated by the Einstein telescope
for binary neutron star inspirals will occur in the frequency region below
$10$ Hz [62]. The matched filtering algorithm used to perform signal detection
requires applying the fast Fourier transform to each of the $M$ templates
which has a computational cost that scales as $O(MN\log N)$ [24]. Therefore,
given the size of the template bank $M$ must also grow to account for these
low frequency signals, where generally the number of templates in the bank
scales proportional to $f_{\text{min}}^{-8/3}$ to maintain a constant mismatch
between potential signals and templates [63], longer duration signals will
pose a computational challenge to the signal detection process [61]. This is
further confounded by the need for a higher dimensionality template bank to
account for spin and tidal deformability effects that increase the signal
strength of higher-order terms of the waveform [64].

By considering the method outlined in this work, the longest binary inspiral
gravitational waveforms detectable to third generation detector, that require
gigabytes of classical memory to store, can be stored on less than 32 qubits.
Fig. 9 compares the gate cost of producing waveforms using this method when
assuming ${n_{c}=16}$, ${n_{l}=4}$ and ${n_{a}=n+n_{c}}$ to the cost of using
the arbitrary state preparation routine. The black line shows the relative
speed-up in terms of an arbitrary amplitude state preparation routine that
requires $2^{n}$ controlled NOT gates [65] when compared to the upper bound of
controlled NOT gates of the analysis when using the Grover-Rudolph amplitude
preparation routine specified in Eq. 35 of Appendix E. While the overhead gate
cost ensures that such a routine is inefficient for waveforms of duration less
than $10^{3}$ s, the method demonstrates a clear reduction in gate cost which
can reach two orders of magnitude when waveforms of durations of up to
${\sim}10$ days are considered. Waveforms of such durations are required for
analysing low mass binary inspiral systems that are detectable by third
generation detectors.

We note that the computational cost of preparing a waveform can be reduced
from that described in Fig. 9 at the expense of fidelity of the resulting
amplitude encoded waveform by reducing $n_{c}$, $n_{l}$, $n_{a}$ or $m_{b}$.
Conversely the fidelity can be increased by increasing $n_{c}$ or $n_{l}$ with
higher gate cost. The solid grey line of Fig. 9 depicts the relative speed-up
when a trained parameterized quantum circuit like that of Eq. II.3 with
${L=100}$ layers of number of controlled NOT gates given in Eq. 36 when
compared to the $2^{n}$ controlled NOT gates of the arbitrary state
preparation routine. The parameterised quantum circuit gains an order of
magnitude reduction in controlled NOT gates when compared to the Grover-
Rudolph algorithm amplitude preparation routine.

While the application of quantum algorithms to directly address the data
analysis problems of third generation detectors is a task left for future
work, we illustrate the power of the state preparation routine by showing the
longest duration binary neutron star inspiral observable to third generation
detectors can be amplitude encoded onto less than 32 qubits with the aid of 68
ancillary qubits. This can be performed with an exponential speed-up over an
arbitrary state preparation routine, a necessary condition to be met to
qualify for any potential advantage of many quantum algorithms [31].

##  V Conclusion

We provide an efficient routine for preparing amplitude encoded
representations of the inspiral phase waveform of a gravitational wave signal
emitted by the coalescence of a compact binary system. The routine is based on
the use of either the Grover-Rudolph algorithm or a generative hybrid
classical-quantum machine learning method, and the ability to efficiently
evaluate arithmetic functions through a piece-wise linear approximation. We
have demonstrated the application of this routine by simulating the
preparation of a spinless $35~{}M_{\odot}$ and $30~{}M_{\odot}$ binary black-
hole merger waveform in the frequency domain onto 6 qubits. The resulting
state obtains a fidelity of $0.995$ to the desired state when using the
Grover-Rudolph algorithm method, which could further increase by increasing
the number of available ancillary qubits, or assuming higher order polynomial
approximations to evaluate the functions.

A fidelity of $0.979$ was demonstrated when a generative adversarial network
was used instead of the Grover-Rudolph algorithm, leading to a significant
reduction in the gate count. The resulting fidelity may be increased by
training a larger parameterized circuit, however we demonstrate that training
larger circuits requires more careful training regiments. The consideration of
other generative modelling schemes with more stable training may be
beneficial, such as adopting a normalizing flow model [66]. A further gate
count reduction may be achieved by approximating the function evaluation step
using a trained parameterised quantum circuit, rather than the current piece-
wise function approximation that suffers from the high gate count required for
arithmetic operations.

Amplitude encoding allows for exponentially larger waveforms to be encoded
onto the given number of qubits in comparison to other encoding formats, which
may allow for waveforms of sizes comparable to those used for classical
analyses on qubit numbers obtainable on quantum processors in the near future.
We investigate this case when considering binary neutron star inspiral
waveforms of lengths necessary for the analysis of signals from third
generation detectors such as the Einstein Telescope. Such detectors will
encounter computational challenges when probing lower frequencies with higher
sensitivity detectors. We demonstrate that the state preparation method
provides a greater speed-up over arbitrary state preparation methods the
longer the waveform duration is, suggesting quantum advantages could be sought
in this low frequency cut-off regime where the classical computational cost is
greatest.

While we mainly pose the use of this routine for preparing states that
represent gravitational waveforms, the same routine can be applied to prepare
other functions of the form in Eq. 18, including the merger and ringdown
waveforms of gravitational wave binary merger events [67]. The adaption of
this routine to perform state preparation of a combined inspiral, merger and
ringdown waveform is left as future work.

Our work demonstrates that efficient encoding of gravitational waveforms into
quantum states with a number of qubits that scales only logarithmically with
the length of the waveform is possible. This extremely space-efficient
encoding offers the tantalising possibility of near-term applications of
quantum computational devices to gravitational wave astronomy. While the
waveforms may be efficiently prepared, we note that this is not the case for
encoding the detector data itself. This will require an arbritrary state
preparation routine to encode and induce a computational cost that scales as
$O(2^{n})$, which is unavoidable. Nonetheless, efficient encoding of template
waveforms allows us to explore data analysis algorithms, with a potential
speed-up in the training phase of machine learning approaches, or as a basis
to construct measurements discriminating between signals in a measurement-
based approach. We expect that our state preparation routine will find
immediate applications in variational-based learning approaches to data
analysis, as well as in the longer term to more sophisticated quantum
protocols.

##  VI Acknowledgements

We are grateful for the support from the EPSRC under Grant No. EP/T001062/1.

##  Appendix A Stationary phase approximation

The waveforms that are represented analytically in the positive frequency-
domain can be written:

$\tilde{h}(f)\propto\tilde{A}(f)e^{i\Psi(f)},$ (18)

where $\tilde{h}(f)$ is the Fourier transform of real-valued $h(t)$:

$h(t)\propto\Re\left\\{A(t)e^{-i\Phi(t)}\right\\},$ (19)

such that $\tilde{h}(f)=\int h(t)e^{2\pi ift}dt$, and $A$ and $\tilde{A}$ are
real functions. Commonly, $\tilde{h}(f)$ can be approximated from an
analytical expression of Eq. 19 using the stationary phase approximation [68]
such that

$\tilde{A}(f)\mathrel{\vbox{
\offinterlineskip\halign{\hfil$#$\cr\propto\cr\kern
2.0pt\cr\sim\cr\kern-2.0pt\cr}}}A(t_{0})\left(\frac{d^{2}\Phi(t_{0})}{dt^{2}}\right)^{1/2}$
(20)

and

$\Psi(f)\approx 2\pi ft_{0}-\phi(f)-\frac{\pi}{4},$ (21)

where $t_{0}$ is the time at which $d\Phi(t)/dt|_{t_{0}}=2\pi f$, at which
point $A(t)$ varies slowly. The function $\phi(f)$ can be determined from
$\Phi(t)$ given the frequency dependence on time through
${\Phi(t)=\int^{t}fdt^{\prime}}$.

##  Appendix B Gravitational waves from compact binary inspirals

Much of the following background is based on the work presented in Refs. [69,
55]. The gravitational waveform is entirely determined through general
relativity. However, solutions to the non-linear Einstein field equations can
only be found numerically, and analytical forms of the waveforms can only be
approximated. For compact binary systems, analytical waveforms are commonly
determined using a post-Newtonian approximation in the near-field of the
orbital system, while a post-Minkowskian approximation is made to the
surrounding field [70]. Solutions are expanded about $(v/c)^{2}$, where ${v\ll
c}$ is the orbital velocity. Taking only the leading-order terms leads to a
‘Newtonian’ waveform which describes gravitational radiation produced solely
from the change in the mass-quadrupole moment of the binary system, with
frequency twice that of the orbital frequency [71]. However, such a solution
becomes inaccurate in systems with high mass ratios, component spins, or as
the orbital velocity increases for systems with shorter separation where
higher-order modes are excited. This requires higher-order terms to be
considered, leading to ‘post-Newtonian’ waveforms. Analytical waveforms are
also often formed through fitting phenomenological models to sets of waveforms
determined through numerical relativity [67]. For the remainder of this
section, geometric units can be assumed such that ${G=c=1}$.

Gravitational radiation amplitudes depend on two independent polarization
states denoted ‘$+$’ and ‘$\times$’, separated by an angle of $\pi/4$ from one
another. For a binary system of compact objects, these polarization amplitudes
depend on the inclination $\iota$ of the system’s angular momentum vector with
respect to the line-of-sight, such that

$\displaystyle h_{+}=$
$\displaystyle\frac{h_{0}}{2}(1+\cos^{2}\iota)\cos\Phi(t),$ (22)
$\displaystyle h_{\times}=$ $\displaystyle h_{0}\cos\iota\sin\Phi(t).$ (23)

The common amplitude $h_{0}$ is equal to

$h_{0}\approx\frac{\mathcal{M}^{5/3}(2\pi f)^{2/3}}{D}.$ (24)

Here $\mathcal{M}$ is the chirp-mass, related to the total mass
$M=m_{1}+m_{2}$ and reduced mass $\eta=m_{1}m_{2}/(m_{1}+m_{2})$ of the system
of component masses $m_{1,2}$ by $\mathcal{M}=\eta^{3/5}M^{2/5}$, and $D$ is
the distance to the source. To induce a time-dependent strain on a detector,
the polarization amplitudes couple to the detector’s response to each
polarization $F_{+,\times}$ and combine linearly:

$h(t)=F_{+}h_{+}+F_{\times}h_{\times}=Qh_{0}\cos(\Phi(t)-\Phi_{0}).$ (25)

This can be rewritten as

$h(t)=Qh_{0}\cos(\Phi(t)-\Phi_{0}),$ (26)

where ${Q=(A_{+}^{2}+A_{\times}^{2})^{1/2}}$, and
${\Phi_{0}=\arctan(A_{\times}/A_{+})}$, given
${A_{+}=F_{+}(1+\cos^{2}\iota)/2}$ and ${A_{\times}=F_{\times}\cos\iota}$. The
expression in Eq. 26 is in the form of Eq. 19.

The spin of the binary system is parameterised by the spin-orbit parameter:

$\beta=\frac{1}{12}\sum^{2}_{i=1}\boldsymbol{L}\cdot\boldsymbol{\chi}_{i}[113(m_{i}/M)^{2}+75\eta],$
(27)

and spin-spin parameter

$\sigma=\frac{\eta}{48}\left(721(\boldsymbol{L}\cdot\boldsymbol{\chi}_{1}\boldsymbol{L}\cdot\boldsymbol{\chi}_{2})-247(\boldsymbol{\chi}_{1}\cdot\boldsymbol{\chi}_{2})\right).$
(28)

Here, $\boldsymbol{\chi}_{1,2}=\boldsymbol{S}_{1,2}/m_{1,2}^{2}$ where
$\boldsymbol{S}_{1,2}$ is the spin angular momentum of bodies $1$ and $2$,
which is compared against the orbital angular momentum of the system
$\boldsymbol{L}$.

The time-frequency dependence to second post-Newtonian order is

$\frac{df}{dt}=\frac{96}{5\pi\mathcal{M}^{2}}(\pi\mathcal{M}f)^{11/3}\left[1-\left(\frac{743}{336}+\frac{11}{4}\eta\right)(\pi
Mf)^{2/3}\right.\\\ \left.+(4\pi-\beta)(\pi
Mf)+\left(\frac{34103}{18144}+\frac{13661}{2016}\eta+\frac{59}{18}\eta^{2}+\sigma\right)(\pi
Mf)^{4/3}\right].$ (29)

From Eq. B, the frequency dependence on time can be calculated, giving the
phase as $\Phi(t)=\int^{t}fdt^{\prime}$ from which Eq. 17 is calculated from
Eq. 21. The time-frequency relation can be inverted allowing for
$\phi(f)=\Phi(t(f))$ to be determined, resulting in the expression for the
frequency dependent amplitude of Eq. 16 given Eq. 20.

##  Appendix C Fixed-point binary encoding

Our routine to prepare the state of Eq. 1 requires binary strings to be
encoded in the computational basis using signed-magnitude (and unsigned-
magnitude) and two’s-complement representations, which we detail here.

We represent signed binary numbers in a bit-strings of $n$ bits with the
ordering:

$x=\underbrace{x_{n-1}}_{\text{Sign bit}}\underbrace{x_{n-2}\dots
x_{n-n_{\text{int}}-1}}_{n_{\text{int}}\text{ integer
bits}}\underbrace{x_{p-1}\dots x_{0}}_{p\text{ precision bits}},$

so that the sign of the number is represented in the leading order bit, the
integer part in the following $n_{\text{int}}$ bits, and the fraction part in
the final $p$ precision bits. Note that while $n_{\text{int}}$ and $p$ must
follow the condition that $n=n_{\text{int}}+p+1$, the values are not
restricted to positive numbers, allowing for the precision $2^{-p}\geq 1$ or
the upper bound to $x$ that the bit-string can represent to be
$2^{n_{\text{int}}}\leq 1$.

###  C.1 Unsigned-magnitude representation

When only the magnitude of $x$ is required, the sign bit is dropped and the
representation of the number is simply:

$x=\sum^{n-1}_{i=0}x_{i}2^{i-p}.$

###  C.2 Signed-magnitude representation

Positive and negative numbers of equal magnitude are related by flipping the
leading order bit (the sign bit) so that:

$x=(-1)^{x_{n-1}}\sum^{n-2}_{i=0}x_{i}2^{i-p}.$

Note that this representation includes both a positive and negative value for
zero.

###  C.3 Two’s-complement representation

Negative numbers are related to positive values of equal magnitude by a bit-
flip of all bits and incrementing the result by the least significant bit
(equivalent to the addition of $2^{-p}$):

$x=-x_{n-1}2^{n_{\text{int}}}+\sum^{n-2}_{i=0}x_{i}2^{i-p}.$

##  Appendix D Fixed-point multiplication through the quantum Fourier
transform

The multiplication of fixed-point numbers can be carried out without ancillary
qubits through the use of the quantum Fourier transform. The multiplicand is
stored in state $|a\rangle$ and the multiplier in $|b\rangle$. A third,
initially empty, output register $|0\rangle_{o}$ is required to store the
product such that

$\displaystyle|a\rangle|b\rangle|0\rangle_{o}$
$\displaystyle\xrightarrow{\text{Mult}}|a\rangle|b\rangle|ab\rangle_{o}.$

Like classical multiplication algorithms, the multiplication is performed by a
series of additions of $|a\rangle$ onto $|0\rangle_{o}$ of a number depending
on $|b\rangle$. However, here we perform the additions to $|0\rangle_{o}$ by
first performing the quantum Fourier transform across the register. The
addition of $|a\rangle$ onto the output register while in the conjugate basis
space is then a series of controlled rotation operations. The details of how
these additions are performed is given in Ref. [72], which is expanded to
multiplication by further controlling the rotational addition operations on
the multiplier $|b\rangle$. Similar to [52], we simplify the multiplication by
asserting that the multiplicand is positive, and the result is processed
conditioned on the most significant qubit; if the multiplier is negative (the
most significant qubit is $|1\rangle$), the output is bit-flipped and
incremented by 1 to correspond with negative values in the two’s complement
notation. As the multiplication product is encoded in the conjugate basis in
two’s-complement representation, the following addition operation of Fig. 2 is
applied independently of the product’s sign.

##  Appendix E Gate cost

An upper bound to the number of controlled $X$ gates can be assigned to each
of the operations used in the the state preparation routine. This upper bound
could be reduced by optimising the implementation of each of the operations or
by introducing additional ancillary qubits. However, we provide the upper
bound when assuming the form of each operation as used to perform the quantum
simulation described in Sec. III. For the multiplication of a register of
length $n_{1}$ on one of $n_{2}$ where the result is stored on a register of
length $n_{1}+n_{2}$, the upper bound to the number of controlled $X$ gates
is:

$\displaystyle C_{\text{Mult}}(n_{1},n_{2})=$ $\displaystyle
8(n_{1}+n_{2})n_{2}(n_{1}-1)$ $\displaystyle+20(n_{1}+n_{2})^{2}$
$\displaystyle-13(n_{1}+n_{2}).$ (30)

Similarly the addition of registers of length $n_{1}$ and $n_{2}$, where the
result is stored on the second register requires:

$C_{\text{Add}}(n_{1},n_{2})=2n_{1}n_{2}+2n_{2}(n_{2}-1).$ (31)

The label operation when utilizing the quantum Fourier transform is bound by:

$C_{\text{Label}}(n,n_{l})=2n_{l}(n_{l}-1)+2^{n_{l}+1}(6n+n_{l}+1).$ (32)

The inputting of the piece-wise function coefficients requires:

$C_{\text{X}}(n_{c},n_{l})=n_{c}2^{n_{l}}C_{C^{\otimes n}X}(n_{l}),$ (33)

where $C_{C^{\otimes n}X}$ is the number of controlled $X$ gates required to
perform a $n$-controlled $X$ operation. Without the use of ancillary qubits
$C_{C^{\otimes n}X}(n_{l})=2n_{l}^{2}-6n_{l}+5$ [73], otherwise $C_{C^{\otimes
n}X}(n_{l})=20(n_{l}-2)$ with the use of $n-2$ ancillary qubits [74] or
$12(n_{l}-1)+1$ with $n-1$ ancillary qubits [19].

The linear piece-wise function (as shown in Fig. 2) is then constrained by:

$\displaystyle C_{\text{LPF}}(n,n_{c},n_{l})=$ $\displaystyle
C_{\text{Mult}}(n_{c},n)+C_{\text{Add}}(n_{c},n_{c}+n)$
$\displaystyle+C_{\text{Label}}(n,n_{l})+3C_{\text{X}}(n_{c},n_{l})$ (34)

Then the upper bound on the combined amplitude and phase preparation steps is
placed when assuming $m_{a}=1$ and $m_{b}=n$:

$\displaystyle C_{\text{GW}}(n,n_{c},n_{l})=$ $\displaystyle
2\sum_{m=1}^{n}C_{\text{LPF}}(m,n_{c},n_{l})+2(n+n_{c}(n-1)).$ (35)

If a trained parameterized quantum circuit is used, the amplitude preparation
step can be reduced to $L(n-1)$ controlled $X$ gates, reducing the total count
to:

$\displaystyle C_{\text{GW}}(n,n_{c},n_{l})=$ $\displaystyle
2C_{\text{LPF}}(n,n_{c},n_{l})+L(n-1).$ (36)

## References

  * [1] Benjamin P Abbott, Richard Abbott, TD Abbott, MR Abernathy, Fausto Acernese, Kendall Ackley, Carl Adams, Thomas Adams, Paolo Addesso, RX Adhikari, et al.  Observation of gravitational waves from a binary black hole merger.  Physical review letters, 116(6):061102, 2016.
  * [2] BP Abbott, Richard Abbott, TDea Abbott, S Abraham, F Acernese, K Ackley, C Adams, RX Adhikari, VB Adya, Christoph Affeldt, et al.  Gwtc-1: a gravitational-wave transient catalog of compact binary mergers observed by ligo and virgo during the first and second observing runs.  Physical Review X, 9(3):031040, 2019.
  * [3] R Abbott, TD Abbott, S Abraham, F Acernese, K Ackley, A Adams, C Adams, RX Adhikari, VB Adya, Christoph Affeldt, et al.  Gwtc-2: compact binary coalescences observed by ligo and virgo during the first half of the third observing run.  Physical Review X, 11(2):021053, 2021.
  * [4] R Abbott, TD Abbott, F Acernese, K Ackley, C Adams, N Adhikari, RX Adhikari, VB Adya, C Affeldt, D Agarwal, et al.  Gwtc-3: compact binary coalescences observed by ligo and virgo during the second part of the third observing run.  arXiv preprint arXiv:2111.03606, 2021.
  * [5] R Abbott, TD Abbott, S Abraham, F Acernese, K Ackley, A Adams, C Adams, RX Adhikari, VB Adya, Christoph Affeldt, et al.  Observation of gravitational waves from two neutron star–black hole coalescences.  The Astrophysical journal letters, 915(1):L5, 2021.
  * [6] Benjamin P Abbott, Rich Abbott, TD Abbott, Fausto Acernese, Kendall Ackley, Carl Adams, Thomas Adams, Paolo Addesso, RX Adhikari, Vaishali B Adya, et al.  Gw170817: observation of gravitational waves from a binary neutron star inspiral.  Physical review letters, 119(16):161101, 2017.
  * [7] Benjamin P Abbott, Robert Abbott, TD Abbott, F Acernese, K Ackley, C Adams, T Adams, P Addesso, RX Adhikari, VB Adya, et al.  Gravitational waves and gamma-rays from a binary neutron star merger: Gw170817 and grb 170817a.  The Astrophysical Journal Letters, 848(2):L13, 2017.
  * [8] B. P. Abbott et. al.  Multi-messenger observations of a binary neutron star merger.  The Astrophysical Journal Letters, 848(2):L12, oct 2017.
  * [9] R Abbott, TD Abbott, S Abraham, F Acernese, K Ackley, A Adams, C Adams, RX Adhikari, VB Adya, C Affeldt, et al.  All-sky search in early o3 ligo data for continuous gravitational-wave signals from unknown neutron stars in binary systems.  Physical Review D, 103(6):064017, 2021.
  * [10] BP Abbott, R Abbott, TD Abbott, S Abraham, Fausto Acernese, K Ackley, C Adams, VB Adya, C Affeldt, M Agathos, et al.  Optically targeted search for gravitational waves emitted by core-collapse supernovae during the first and second observing runs of advanced ligo and advanced virgo.  Physical Review D, 101(8):084002, 2020.
  * [11] R Abbott, TD Abbott, S Abraham, F Acernese, K Ackley, A Adams, C Adams, RX Adhikari, VB Adya, C Affeldt, et al.  Search for anisotropic gravitational-wave backgrounds using data from advanced ligo and advanced virgo’s first three observing runs.  Physical Review D, 104(2):022005, 2021.
  * [12] R Abbott, TD Abbott, S Abraham, F Acernese, K Ackley, A Adams, C Adams, RX Adhikari, VB Adya, C Affeldt, et al.  Constraints on cosmic strings using data from the third advanced ligo–virgo observing run.  Physical review letters, 126(24):241102, 2021.
  * [13] BP Abbott, R Abbott, TD Abbott, S Abraham, F Acernese, K Ackley, C Adams, RX Adhikari, VB Adya, Christoph Affeldt, et al.  A gravitational-wave measurement of the hubble constant following the second observing run of advanced ligo and virgo.  The Astrophysical Journal, 909(2):218, 2021.
  * [14] R Abbott, TD Abbott, S Abraham, F Acernese, K Ackley, A Adams, C Adams, RX Adhikari, VB Adya, C Affeldt, et al.  Tests of general relativity with binary black holes from the second ligo-virgo gravitational-wave transient catalog.  Physical review D, 103(12):122002, 2021.
  * [15] Michele Maggiore, Chris Van Den Broeck, Nicola Bartolo, Enis Belgacem, Daniele Bertacca, Marie Anne Bizouard, Marica Branchesi, Sebastien Clesse, Stefano Foffa, Juan García-Bellido, et al.  Science case for the einstein telescope.  Journal of Cosmology and Astroparticle Physics, 2020(03):050, 2020\.
  * [16] David Reitze, Rana X. Adhikari, Stefan Ballmer, Barry Barish, Lisa Barsotti, GariLynn Billingsley, Duncan A. Brown, Yanbei Chen, Dennis Coyne, Robert Eisenstein, Matthew Evans, Peter Fritschel, Evan D. Hall, Albert Lazzarini, Geoffrey Lovelace, Jocelyn Read, B. S. Sathyaprakash, David Shoemaker, Joshua Smith, Calum Torrie, Salvatore Vitale, Rainer Weiss, Christopher Wipf, and Michael Zucker.  Cosmic Explorer: The U.S. Contribution to Gravitational-Wave Astronomy beyond LIGO.  In Bulletin of the American Astronomical Society, volume 51, page 35, September 2019.
  * [17] Emanuele Berti, Vitor Cardoso, and Clifford M Will.  Gravitational-wave spectroscopy of massive black holes with the space interferometer lisa.  Physical Review D, 73(6):064030, 2006.
  * [18] LORENZO Amati, Paul O’Brien, D Götz, E Bozzo, C Tenzer, F Frontera, Giancarlo Ghirlanda, CLAUDIO Labanti, JP Osborne, G Stratta, et al.  The theseus space mission concept: science case, design and expected performances.  Advances in Space Research, 62(1):191–244, 2018.
  * [19] Michael A Nielsen and Isaac Chuang.  Quantum computation and quantum information, 2002.
  * [20] Peter W Shor.  Polynomial-time algorithms for prime factorization and discrete logarithms on a quantum computer.  SIAM review, 41(2):303–332, 1999.
  * [21] Lov K Grover.  Quantum mechanics helps in searching for a needle in a haystack.  Physical review letters, 79(2):325, 1997.
  * [22] Ciaran Ryan-Anderson, JG Bohnet, Kenny Lee, Daniel Gresh, Aaron Hankin, JP Gaebler, David Francois, Alexander Chernoguzov, Dominic Lucchetti, NC Brown, et al.  Realization of real-time fault-tolerant quantum error correction.  Physical Review X, 11(4):041058, 2021.
  * [23] Sebastian Krinner, Nathan Lacroix, Ants Remm, Agustin Di Paolo, Elie Genois, Catherine Leroux, Christoph Hellings, Stefania Lazar, Francois Swiadek, Johannes Herrmann, et al.  Realizing repeated quantum error correction in a distance-three surface code.  Nature, 605(7911):669–674, 2022.
  * [24] Bruce Allen, Warren G Anderson, Patrick R Brady, Duncan A Brown, and Jolien DE Creighton.  Findchirp: An algorithm for detection of gravitational waves from inspiraling compact binaries.  Physical Review D, 85(12):122006, 2012.
  * [25] Sijia Gao, Fergus Hayes, Sarah Croke, Chris Messenger, and John Veitch.  Quantum algorithm for gravitational-wave matched filtering.  Physical Review Research, 4(2):023006, 2022.
  * [26] Koichi Miyamoto, Gonzalo Morrás, Takahiro S Yamamoto, Sachiko Kuroyanagi, and Savvas Nesseris.  Gravitational wave matched filtering by quantum monte carlo integration and quantum amplitude amplification.  arXiv preprint arXiv:2205.05966, 2022.
  * [27] Gabriel Escrig, Roberto Campos, Pablo AM Casares, and MA Martin-Delgado.  Parameter estimation of gravitational waves with a quantum metropolis algorithm.  arXiv preprint arXiv:2208.05506, 2022.
  * [28] W Keith Hastings.  Monte carlo sampling methods using markov chains and their applications.  Biometrika, 1970.
  * [29] Aram W. Harrow, Avinatan Hassidim, and Seth Lloyd.  Quantum algorithm for linear systems of equations.  Phys. Rev. Lett., 103:150502, Oct 2009.
  * [30] Seth Lloyd, Masoud Mohseni, and Patrick Rebentrost.  Quantum principal component analysis.  Nature Physics, 10(9):631–633, 2014.
  * [31] Scott Aaronson.  Read the fine print.  Nature Physics, 11(4):291–293, 2015.
  * [32] Ryan LaRose and Brian Coyle.  Robust data encodings for quantum classifiers.  Physical Review A, 102(3):032420, 2020.
  * [33] Iordanis Kerenidis and Anupam Prakash.  Quantum recommendation systems.  arXiv preprint arXiv:1603.08675, 2016.
  * [34] Kwok Ho Wan, Oscar Dahlsten, Hlér Kristjánsson, Robert Gardner, and MS Kim.  Quantum generalisation of feedforward neural networks.  npj Quantum information, 3(1):1–8, 2017.
  * [35] Jonathan Romero, Jonathan P Olson, and Alan Aspuru-Guzik.  Quantum autoencoders for efficient compression of quantum data.  Quantum Science and Technology, 2(4):045001, 2017.
  * [36] Edward Farhi and Hartmut Neven.  Classification with quantum neural networks on near term processors.  arXiv preprint arXiv:1802.06002, 2018.
  * [37] Maria Schuld, Alex Bocharov, Krysta M Svore, and Nathan Wiebe.  Circuit-centric quantum classifiers.  Physical Review A, 101(3):032308, 2020.
  * [38] Dmytro Bondarenko and Polina Feldmann.  Quantum autoencoders to denoise quantum data.  Physical review letters, 124(13):130502, 2020.
  * [39] XiMing Wang, YueChi Ma, Min-Hsiu Hsieh, and Man-Hong Yung.  Quantum speedup in adaptive boosting of binary classification.  Science China Physics, Mechanics & Astronomy, 64(2):1–10, 2021\.
  * [40] Xiao-Ming Zhang, Weicheng Kong, Muhammad Usman Farooq, Man-Hong Yung, Guoping Guo, and Xin Wang.  Generic detection-based error mitigation using quantum autoencoders.  Physical Review A, 103(4):L040403, 2021.
  * [41] B David Clader, Bryan C Jacobs, and Chad R Sprouse.  Preconditioned quantum linear system algorithm.  Physical review letters, 110(25):250504, 2013.
  * [42] Andrew M Childs, Robin Kothari, and Rolando D Somma.  Quantum algorithm for systems of linear equations with exponentially improved dependence on precision.  SIAM Journal on Computing, 46(6):1920–1950, 2017.
  * [43] Mikko Mottonen, Juha J Vartiainen, Ville Bergholm, and Martti M Salomaa.  Transformation of quantum states using uniformly controlled rotations.  arXiv preprint quant-ph/0407010, 2004.
  * [44] Lov Grover and Terry Rudolph.  Creating superpositions that correspond to efficiently integrable probability distributions.  arXiv preprint quant-ph/0208112, 2002.
  * [45] Jian Zhao, Yu-Chun Wu, Guang-Can Guo, and Guo-Ping Guo.  State preparation based on quantum phase estimation.  arXiv preprint arXiv:1912.05335, 2019.
  * [46] Marco Cerezo, Andrew Arrasmith, Ryan Babbush, Simon C Benjamin, Suguru Endo, Keisuke Fujii, Jarrod R McClean, Kosuke Mitarai, Xiao Yuan, Lukasz Cincio, et al.  Variational quantum algorithms.  Nature Reviews Physics, 3(9):625–644, 2021.
  * [47] Raban Iten, Roger Colbeck, Ivan Kukuljan, Jonathan Home, and Matthias Christandl.  Quantum circuits for isometries.  Physical Review A, 93(3):032318, 2016.
  * [48] Marcello Benedetti, Erika Lloyd, Stefan Sack, and Mattia Fiorentini.  Parameterized quantum circuits as machine learning models.  Quantum Science and Technology, 4(4):043001, 2019.
  * [49] Stephen Barnett.  Quantum information, volume 16.  Oxford University Press, 2009.
  * [50] Ian Goodfellow, Jean Pouget-Abadie, Mehdi Mirza, Bing Xu, David Warde-Farley, Sherjil Ozair, Aaron Courville, and Yoshua Bengio.  Generative adversarial networks.  Communications of the ACM, 63(11):139–144, 2020.
  * [51] Christa Zoufal, Aurélien Lucchi, and Stefan Woerner.  Quantum generative adversarial networks for learning and loading random distributions.  npj Quantum Information, 5(1):103, 2019.
  * [52] Thomas Häner, Martin Roetteler, and Krysta M Svore.  Optimizing quantum circuits for arithmetic.  arXiv preprint arXiv:1805.12445, 2018.
  * [53] Ethan T Wang.  Efficient quantum state preparation for the cauchy distribution based on piecewise arithmetic.  IEEE Transactions on Quantum Engineering, 3:1–9, 2021.
  * [54] Eugene Y Remez.  Sur la détermination des polynômes d’approximation de degré donnée.  Comm. Soc. Math. Kharkov, 10(196):41–63, 1934.
  * [55] Eric Poisson and Clifford M Will.  Gravitational waves from inspiraling compact binaries: Parameter estimation using second-post-newtonian waveforms.  Physical Review D, 52(2):848, 1995.
  * [56] Md Sajid Anis et al.  Qiskit: An open-source framework for quantum computing, 2021.
  * [57] Fergus Hayes.  Quantum state preparation of gravitational waves.  https://github.com/Fergus-Hayes/Quantum_state_preparation_of_gravitational_waves, 2023\.
  * [58] Benjamin J Owen.  Search templates for gravitational waves from inspiraling binaries: Choice of template spacing.  Physical Review D, 53(12):6749, 1996.
  * [59] Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, et al.  Pytorch: An imperative style, high-performance deep learning library.  Advances in neural information processing systems, 32, 2019.
  * [60] Man Leong Chan, Chris Messenger, Ik Siong Heng, and Martin Hendry.  Binary neutron star mergers and third generation detectors: Localization and early warning.  Physical Review D, 97(12):123014, 2018.
  * [61] Leone Bosi and Edward K Porter.  Data analysis challenges for the einstein telescope.  General Relativity and Gravitation, 43:519–535, 2011.
  * [62] Amber K Lenon, Duncan A Brown, and Alexander H Nitz.  Eccentric binary neutron star search prospects for cosmic explorer.  Physical Review D, 104(6):063011, 2021.
  * [63] Benjamin J Owen and Bangalore Suryanarayana Sathyaprakash.  Matched filtering of gravitational waves from inspiraling compact binaries: Computational cost and template placement.  Physical Review D, 60(2):022002, 1999.
  * [64] Xisco Jiménez Forteza, Tiziano Abdelsalhin, Paolo Pani, and Leonardo Gualtieri.  Impact of high-order tidal terms on binary neutron-star waveforms.  Physical Review D, 98(12):124014, 2018.
  * [65] Vivek V Shende, Stephen S Bullock, and Igor L Markov.  Synthesis of quantum logic circuits.  In Proceedings of the 2005 Asia and South Pacific Design Automation Conference, pages 272–275, 2005.
  * [66] Danilo Rezende and Shakir Mohamed.  Variational inference with normalizing flows.  In International conference on machine learning, pages 1530–1538. PMLR, 2015.
  * [67] Parameswaran Ajith, Stanislav Babak, Yanbei Chen, Martin Hewitson, Badri Krishnan, AM Sintes, John T Whelan, B Brügmann, P Diener, Nils Dorband, et al.  Template bank for gravitational waveforms from coalescing binary black holes: Nonspinning binaries.  Physical Review D, 77(10):104017, 2008.
  * [68] Kip S Thorne, S Hawking, and W Israel.  300 years of gravitation, 1987.
  * [69] Curt Cutler and Eanna E Flanagan.  Gravitational waves from merging compact binaries: How accurately can one extract the binary’s parameters from the inspiral waveform?  Physical Review D, 49(6):2658, 1994.
  * [70] Luc Blanchet, Thibault Damour, Bala R Iyer, Clifford M Will, and Alan G Wiseman.  Gravitational-radiation damping of compact binary systems to second post-newtonian order.  Physical Review Letters, 74(18):3515, 1995.
  * [71] Lee Samuel Finn and David F Chernoff.  Observing binary inspiral in gravitational radiation: One interferometer.  Physical Review D, 47(6):2198, 1993.
  * [72] Lidia Ruiz-Perez and Juan Carlos Garcia-Escartin.  Quantum arithmetic with the quantum fourier transform.  Quantum Information Processing, 16(6):1–14, 2017.
  * [73] Mehdi Saeedi and Massoud Pedram.  Linear-depth quantum circuits for n-qubit toffoli gates with no ancilla.  Physical Review A, 87(6):062318, 2013.
  * [74] Adriano Barenco, Charles H Bennett, Richard Cleve, David P DiVincenzo, Norman Margolus, Peter Shor, Tycho Sleator, John A Smolin, and Harald Weinfurter.  Elementary gates for quantum computation.  Physical review A, 52(5):3457, 1995.