# Exact and approximate continuous-variable
gate decompositions

Timjan Kalajdzievski  Xanadu, Toronto, ON, M5G 2C8, Canada  Nicolás Quesada
Xanadu, Toronto, ON, M5G 2C8, Canada

###### Abstract

We gather and examine in detail gate decomposition techniques for continuous-
variable quantum computers and also introduce some new techniques which expand
on these methods. Both exact and approximate decomposition methods are studied
and gate counts are compared for some common operations. While each having
distinct advantages, we find that exact decompositions have lower gate counts
whereas approximate techniques can cover decompositions for all continuous-
variable operations but require significant circuit depth for a modest
precision.

##  1 Introduction

For a quantum computer to be able to execute an algorithm, often defined by a
product of high-level unitary transformations [1, 2, 3, 4, 5, 6], the
operations of the algorithm must be translated into the small set of logical
gates directly built into the hardware. These built-in operations are the
quantum computing analog to the logical gates which are hardwired into the
transistors of a classical computer. In order to program a quantum computer to
execute an algorithm, a combination of the built-in elementary gates must be
found that can reproduce the operation of the algorithm. Certain sets of gates
exist such that any arbitrary unitary operation can be expressed as a finite
product of gates from the set, to any desired precision [7, 8, 9, 10, 11, 12].
Due to their ability to express any given transformation, these are referred
to as universal gates sets. Thus, in order to program a quantum computer, a
method is required to decompose any given unitary in terms of elements of the
universal gate sets. Ideally, a decomposition method will allow for the
implementation of an algorithm with high precision and as few gates as
possible. Gate decomposition methods and techniques have been broadly studied,
especially in the discrete variable (qubit) domain [7, 10, 13, 14, 15].

The qubit model of quantum computing uses two-level quantum systems that act
as the quantum analog of bits on a classical computer. The logical gates
acting on qubits can be represented by $2\times 2$ unitary matrices, and the
process of translating general operations into the logical gates of this
model, or gate decomposition, has been well studied. For example, the Solovay-
Kitaev theorem [15, 5], stated informally tells us that any single-qubit
operation can be approximated to high precision using a sequence of only a few
logical gates. Other results have also been discovered which strengthen
decomposition techniques in qubit systems by allowing for less logical gates
to be needed, or the resulting decomposition to be of a greater precision.
These techniques are used for single-qubit operations [16, 17, 18, 19] as well
as general multi-qubit operations [6].

In this manuscript we focus on gate decomposition techniques for continuous-
variable (CV) quantum computers, which have not been as extensively studied.
In the CV model of quantum computing, registers are infinite-dimensional
quantum systems called qumodes (the infinite-dimensional generalization of
qubits) and physically correspond to quantum harmonic oscillators. The logic
gates in this model are unitary operators which act on the infinite-
dimensional Hilbert space [11, 20, 21, 22, 1, 2, 3, 4]. CV quantum computing
using quadrature operators to form logic gates was first proposed by Lloyd and
Braunstein [11]. The logic gates in that case were formed using Hamiltonians
which were polynomial in the quadrature operators. These gates allow one to
express analog computation using quantum-mechanical operations, and moreover,
by constructing a universal set one can construct any gate generated by a
Hamiltonian that is a polynomial of finite order in the quadrature operations,
like the ones typically found in the simulation of bosonic systems [23, 2].
The problem of finding optimal ways to decompose an arbitrary gate into a
product of gates from a given universal set is precisely the focus of this
manuscript. The study of this problem also finds application when continuous
degrees of freedom are used to encode error-correctable qubits, since the
logical operations of the logical qubits need to be ultimately executed as
continuous-variable gates [24, 25, 26, 27].

We have structured our manuscript as follows: in Sec. 2 we set up some basic
notations. Then in Sec. 3 we examine a comprehensive method for decomposing
arbitrary CV transformations using commutator identities introduced by Sefi
and van Loock in Ref. [13]. This method is approximate, meaning that the final
decomposition only implements the desired operation up to a certain error. The
error is decreased by repeating the sequence of gates in the decomposition
[28, 13, 29]. Exact decompositions which do not rely on Trotterization (see
next section) are also known but are not universally applicable [1, 13, 3]. We
first study exact decompositions of quadratic Hamiltonians in Sec. 4.1 using
their symplectic representation [30, 31]. Then, in Sec. 4.2 we study exact
decomposition of higher than quadratic polynomials of mutually commuting
operators [14, 32]. These methods are further extended and generalized in Sec.
4.3. In Sec. 5 we examine another approximate method which relies on using
squeezing and displacement operations to assemble a single mode cubic-phase
gate starting from a Kerr gate [33]. Finally we provide a comparison of the
methods as well as some discussion in Sec. 6.

##  2 Background and notation

Registers in the CV model are quantum harmonic oscillators with corresponding
creation and annihilation operators $\hat{a}^{\dagger}_{j}$ and $\hat{a}_{j}$,
where the subscript refers to the mode they act upon. For a system with $M$
modes, these operators satisfy the bosonic commutation relations
$[\hat{a}_{j},\hat{a}^{\dagger}_{k}]=\delta_{jk}$, and
$[\hat{a}_{j},\hat{a}_{k}]=[\hat{a}^{\dagger}_{j},\hat{a}^{\dagger}_{k}]=0$.
An equivalent operator description of a bosonic system uses the quadrature
operators $\hat{x}$ and $\hat{p}$, which can be represented in terms of the
annihilation and creation operators as

$\displaystyle\hat{x}_{j}=\sqrt{\frac{\hbar}{2}}\left(\hat{a}^{\dagger}_{j}+\hat{a}_{j}\right),\quad\hat{p}_{j}=i\sqrt{\frac{\hbar}{2}}\left(\hat{a}^{\dagger}_{j}-\hat{a}_{j}\right),$
(1)

with commutators $[\hat{x}_{j},\hat{x}_{k}]=[\hat{p}_{j},\hat{p}_{k}]=0$ and
$[\hat{x}_{j},\hat{p}_{k}]=i\hbar\delta_{j,k}$. For convenience we group the
operators as follows,

$\displaystyle\bm{\hat{x}}$
$\displaystyle=(\hat{x}_{1},\ldots,\hat{x}_{M})^{T},\quad\bm{\hat{p}}=(\hat{p}_{1},\ldots,\hat{p}_{M})^{T},$
(2) $\displaystyle\bm{\hat{r}}$ $\displaystyle=\begin{pmatrix}\hat{\bm{x}}\\\
\hat{\bm{p}}\end{pmatrix}=(\hat{x}_{1},\ldots,\hat{x}_{M},\hat{p}_{1},\ldots,\hat{p}_{M})^{T}.$
(3)

This notation allows us to write the commutation relations as
$[\hat{r}_{j},\hat{r}_{k}]=i\hbar\Omega_{j,k}$ where

$\displaystyle\bm{\Omega}=\begin{pmatrix}0&\ \mathbb{I}_{M}\\\
-\mathbb{I}_{M}&\ 0\end{pmatrix}$ (4)

is the so-called symplectic form, and $\mathbb{I}_{M}$ is the $M$-dimensional
identity matrix.

As mentioned previously, a universal gate set is a collection of gates such
that any arbitrary unitary operation can be expressed as a finite sequence of
gates from the set, to any desired precision. The CV decomposition techniques
examined in later sections focus on the universal set specified by the gates

$\displaystyle R_{j}(\theta)$
$\displaystyle=\exp\left(i\tfrac{\theta}{2\hbar}\left[\hat{x}_{j}^{2}+\hat{p}_{j}^{2}\right]\right),$
(5a) $\displaystyle Z_{j}(t_{1})$
$\displaystyle=\exp\left(i\frac{t_{1}}{\hbar}\hat{x}_{j}\right),$ (5b)
$\displaystyle P_{j}(t_{2})$
$\displaystyle=\exp\left(i\frac{t_{2}}{2\hbar}\hat{x}^{2}_{j}\right),$ (5c)
$\displaystyle V_{j}(t_{3})$
$\displaystyle=\exp\left(i\frac{t_{3}}{3\hbar}\hat{x}^{3}_{j}\right),$ (5d)
$\displaystyle CZ(\tau)$
$\displaystyle=\exp\left(i\frac{\tau}{\hbar}\hat{x}_{j}\hat{x}_{k}\right),$
(5e)

where $\theta$, $t_{1}$, $t_{2}$, $t_{3}$, and $\tau$ are real parameters.
This particular universal set is often chosen for mathematical convenience.
Gates generated by a Hamiltonian that is at most quadratic are referred to as
Gaussian operations. In the universal set above the rotation operator (5a),
pure-momentum displacement (5b), quadratic-phase (5c), and the multi-mode
controlled-phase gate (5e) are all Gaussian. These operations are implemented
experimentally with linear optics.

The cubic phase gate (5d), which is the only non-Gaussian operator in the
universal set, provides the non-linearity needed to allow for decompositions
of arbitrarily higher order. This gate can be realized either by using optical
nonlinearities [33] or measurements [27, 34, 35, 36, 37, 38, 39]. For example,
in [27] a photon counting measurement is used to introduce the higher order
non-linearity. In this case the full implementation involves a displaced two-
mode squeezed state for which $R^{\dagger}\hat{n}R$ (photon counting in a
rotated basis) is measured on one arm. The desired cubic operation is then
collapsed onto the second unmeasured mode. In [36] repeat-until-success photon
subtractions and Gaussian operations are used, and in [39] quadrature
detection for feed-forward manipulation of parameters produces the necessary
nonlinear interaction. The cubic phase gate can also be realized in the
microwave domain using the Josephson nonlinearity [40].

A rotation operation by $\theta=\frac{\pi}{2}$ is referred to as the Fourier
transform gate [41, 42]. This operation maps between the quadrature operators
$\hat{x}$ and $\hat{p}$ and is a cheap gate to implement experimentally.

$\displaystyle\mathcal{F}^{\dagger}\hat{x}\mathcal{F}=-\hat{p},\quad\mathcal{F}^{\dagger}\hat{p}\mathcal{F}=\hat{x}.$
(6)

An equivalent universal set is one where the non-Gaussian gate $V_{j}(t_{3})$
is replaced by the Kerr gate
$\exp\left(it_{3}(\hat{a}^{\dagger}\hat{a})^{2}\right)$. In general, the non-
Gaussian gates in the set can be replaced with any other and the set will
retain universality. The reason why the cubic phase gate is chosen in the
universal set here, as opposed to a Kerr gate, is due to mathematical
convenience in the exact decompositions. As well, at the moment it seems as
though gates with a mixture of $\hat{x}$ and $\hat{p}$ in the same mode cause
problems for these same methods. This is because the polynomials used to
create the desired powers in the decompositions will no longer be able to be
constructed with unitary conjugation. Finally, note that in Sec. 5 we examine
in detail a method derived in Ref. [33] where a driven-Kerr Hamiltonian is
used to synthesize a cubic phase gate by unitary conjugation with squeezing
and displacement gates.

One can also remove one of the Gaussian gates by showing that it can be
decomposed using other elements in the set [13]. In both of these cases the
choice is often made such that the gates in the set are not costly to
implement experimentally. If the decomposition of one gate in the set includes
multiple uses of another harder to implement gate, then it might be beneficial
to still include it. For example, the quadratic phase gate can be decomposed
in terms of other gates, but its decomposition includes multiple uses of the
cubic phase gate. This means that it can be taken out of the universal set,
but as the cubic phase gate is much more costly to implement, it is common to
still include it. In general, Gaussian gates may be decomposed exactly in
terms of squeezing, displacements and beamsplitters. We demonstrate how to
achieve these decompositions in section 4.1.

Throughout this manuscript we express an arbitrary unitary operator as
$U=e^{itH}$ where the Hermitian operator $H=\sum_{j=1}^{N}H_{j}$ is the
generator of the gate. When decomposing arbitrary gates into a universal set,
it is often necessary to approximate this sum of operators in the exponent as
a product of exponential operators. More specifically, for $H=\hat{A}+\hat{B}$
where $\hat{A}$ and $\hat{B}$ are Hermitian operators, the Zassenhaus formula
[43] states that

$e^{it(\hat{A}+\hat{B})}=e^{it\hat{A}}e^{it\hat{B}}e^{\frac{t^{2}}{2}[\hat{A},\hat{B}]}e^{\frac{-it^{3}}{6}\left(2[\hat{B},[\hat{A},\hat{B}]]+[\hat{A},[\hat{A},\hat{B}]]\right)}\cdots$
(7)

When $[\hat{A},\hat{B}]=0$ the product ends immediately after the first two
operations. In general, however, it is possible that this product never
terminates, resulting in a decomposition that is no longer finite. In this
case, the product can be truncated at a designated stage in the expansion and
the remaining commutators can be neglected. This strategy is referred to as a
Trotter-Suzuki approximation [44, 45], which in general can be written as

$e^{itH}=\left(\prod_{j=1}^{N}e^{i\frac{t}{K}H_{j}}\right)^{K}+O(t^{2}/K),$
(8)

where $H=\sum_{j=1}^{N}H_{j}$. This approximation requires
$K=O(1/\varepsilon)$ gates to achieve precision $\varepsilon$ for fixed $t$.
Note that the subscript $j$ on $H_{j}$ is only an index and does not refer to
a mode, as each $H_{j}$ may contain any number of modes. In general, the
Trotter-Suzuki approximation can be useful in combination with a variety of
other decomposition techniques and approximations as we show in the next
section.

##  3 Approximate decompositions using commutators

In this section we examine an approximate decomposition technique detailed in
Ref. [13]. This method expresses sums and products of the quadrature operators
in terms of commutators and then approximates the exponentials of these
commutators as repeated products of their arguments. More specifically, given
a general unitary $U=e^{it\hat{H}}$ where $\hat{H}$ is a sum of products of
quadrature operators, each term in the sum can be represented in terms of
commutators with the following identities:

$\displaystyle\hat{x}^{m}$
$\displaystyle=-\frac{2}{3(m-1)(2\hbar)^{2}}[\hat{x}^{m-1},[\hat{x}^{3},\hat{p}^{2}]],$
(9) $\displaystyle\hat{x}^{m}\hat{p}^{n}+\hat{p}^{n}\hat{x}^{m}$
$\displaystyle=-\frac{4i}{(n+1)(m+1)(2\hbar)}[\hat{x}^{m+1},\hat{p}^{n+1}]-\frac{1}{(n+1)(2\hbar)^{2}}\sum_{k=1}^{n-1}[\hat{p}^{n-k},[\hat{x}^{m},\hat{p}^{k}]],$
(10)

and for two-mode terms

$\displaystyle\hat{p}_{1}^{m}\hat{p}_{2}^{n}=-\frac{1}{(m+1)(n+1)(2\hbar)^{2}}[\hat{p}_{2}^{n+1},[\hat{p}_{1}^{m+1},\hat{x}_{1}\hat{x}_{2}]].$
(11)

Once the exponent of the general unitary operator is written in terms of a sum
of commutators, it can be separated with the Trotter-Suzuki technique as in
Eq. (8), and then expanded in the following way: given two Hermitian operators
$\hat{A}$ and $\hat{B}$, it holds that [46]

$\displaystyle e^{t^{2}[\hat{A},\hat{B}]}$
$\displaystyle=\left(e^{i\frac{t}{K}\hat{B}}e^{i\frac{t}{K}\hat{A}}e^{-i\frac{t}{K}\hat{B}}e^{-i\frac{t}{K}\hat{A}}\right)^{K^{2}}+O(t^{4}/K).$
(12)

Note that for fixed $t$, $K=O(1/\varepsilon)$ gates are required to achieve an
error of $\varepsilon$ in this approximation, but the resulting circuit will
have a depth of $O(1/\varepsilon^{2})$. It is possible then that very large
circuits may be required to achieve a desired precision. A similar expression
can also be used in the case of a nested commutator in the form
$e^{t^{3}[\hat{B},[\hat{B},\hat{A}]]}$, given by

$\displaystyle e^{it^{3}[\hat{B},[\hat{B},\hat{A}]]}$
$\displaystyle=\left(e^{i\frac{t}{K}\hat{B}}e^{\frac{t^{2}}{K^{2}}[\hat{B},\hat{A}]}e^{-i\frac{t}{K}\hat{B}}e^{-\frac{t^{2}}{K^{2}}[\hat{B},\hat{A}]}\right)^{K^{2}}+O(t^{4}/K).$
(13)

To illustrate the use of this decomposition technique, consider the
decomposition of the unitary
$e^{it\left(\hat{x}^{2}\hat{p}+\hat{p}\hat{x}^{2}\right)}$. First, using the
equality from Eq. (10) we have that

$e^{it\left(\hat{x}^{2}\hat{p}+\hat{p}\hat{x}^{2}\right)}=e^{\frac{t}{3\hbar}[\hat{x}^{3},\hat{p}^{2}]}.$
(14)

Next, using Eq. (12) with $\hat{A}=\hat{x}^{3}$ and $\hat{B}=\hat{p}^{2}$
leads to

$\displaystyle e^{\frac{t}{3\hbar}[\hat{x}^{3},\hat{p}^{2}]}\approx$
$\displaystyle\left(e^{i\frac{\sqrt{\frac{t}{3\hbar}}}{K}\hat{p}^{2}}e^{i\frac{\sqrt{\frac{t}{3\hbar}}}{K}\hat{x}^{3}}e^{-i\frac{\sqrt{\frac{t}{3\hbar}}}{K}\hat{p}^{2}}e^{-i\frac{\sqrt{\frac{t}{3\hbar}}}{K}\hat{x}^{3}}\right)^{K^{2}}+O\left[\left(\frac{t}{3\hbar}\right)^{2}/K\right],$
(15)

where each of the resulting gates on the right-hand side are part of the
universal set in Eq. (5) up to Fourier transforms. In order to obtain a
precision of $O(1/K)$, this product of gates must be repeated $O(K^{2})$
times. For example, if $t=1$ and the desired precision of the decomposition is
$10^{-3}$, then the product of gates on the right-hand side needs to be
repeated approximately $10^{5}$ times.

Another operator we can decompose to demonstrate this technique is the Kerr
operation. We also examine the decomposition of this operation with a
different method in section 5. The Kerr operation in terms of quadrature
operators can be written as
$e^{it\left(\hat{x}^{4}+\hat{x}^{2}\hat{p}^{2}+\hat{p}^{2}\hat{x}^{2}+\hat{p}^{4}-\hbar\hat{x}^{2}-\hbar\hat{p}^{2}\right)}$.
After splitting the terms in the exponent using Trotter-Suzuki as in Eq. (8),
and neglecting Fourier transforms, we are left to decompose operators which
look like $e^{it\hat{x}^{4}}$ and
$e^{it(\hat{x}^{2}\hat{p}^{2}+\hat{p}^{2}\hat{x}^{2})}$. These operators can
be written in terms of commutators with Eq. (9) and (10) respectively to get

$\displaystyle
e^{it\hat{x}^{4}}=e^{-\frac{it}{18\hbar^{2}}[\hat{x}^{3},[\hat{x}^{3},\hat{p}^{2}]]},$
(16)

and

$\displaystyle
e^{it(\hat{x}^{2}\hat{p}^{2}+\hat{p}^{2}\hat{x}^{2})}=e^{\frac{2t}{9\hbar}[\hat{x}^{3},\hat{p}^{3}]}.$
(17)

Expanding with Eq. (12) and (13), we get

$\displaystyle
e^{-\frac{it}{18\hbar^{2}}[\hat{x}^{3},[\hat{x}^{3},\hat{p}^{2}]]}$ (18)
$\displaystyle\ \
\approx\left(e^{-i\frac{\left(\frac{t}{18\hbar^{2}}\right)^{1/3}}{K}\hat{x}^{3}}e^{\frac{\left(\frac{t}{18\hbar^{2}}\right)^{2/3}}{K^{2}}[\hat{x}^{3},\hat{p}^{2}]}e^{i\frac{\left(\frac{t}{18\hbar^{2}}\right)^{1/3}}{K}\hat{x}^{3}}e^{-\frac{\left(\frac{t}{18\hbar^{2}}\right)^{2/3}}{K^{2}}[\hat{x}^{3},\hat{p}^{2}]}\right)^{K^{2}}+O\left[\left(\frac{t}{18\hbar^{2}}\right)^{4/3}/K\right],$

and

$\displaystyle e^{\frac{2t}{9\hbar}[\hat{x}^{3},\hat{p}^{3}]}\approx$
$\displaystyle\left(e^{i\frac{\sqrt{\frac{2t}{9\hbar}}}{K}\hat{p}^{3}}e^{i\frac{\sqrt{\frac{2t}{9\hbar}}}{K}\hat{x}^{3}}e^{-i\frac{\sqrt{\frac{2t}{9\hbar}}}{K}\hat{p}^{3}}e^{-i\frac{\sqrt{\frac{2t}{9\hbar}}}{K}\hat{x}^{3}}\right)^{K^{2}}+O\left[\left(\frac{2t}{9\hbar}\right)^{2}/K\right],$
(19)

where the second and fourth operators on the right-hand side of Eq. (18) can
be expanded further, similar to the expansion in Eq. (15). Once again the
product of gates on the right-hand sides of each equation must be repeated
$O(K^{2})$ times to obtain a precision of $O(1/K)$. In this case, because the
original operator needs to be expanded first with a Trotter-Suzuki
approximation, the final sequence of gates needs to again be repeated to
increase the precision as in Eq. (8). For derivations of some of these
identities, as well as discussion on the necessity of good approximations, see
the supplementary materials section of [13]

To potentially avoid the cost of repeating gates for greater precision, exact
decomposition methods can be used as we show in the next section.

##  4 Exact decompositions

###  4.1 Gaussian operations

In this section we provide an exact decomposition for an arbitrary Gaussian
operator in terms of displacement, interferometer and squeezing gates [47].
Before developing the decompositions we provide definitions for these unitary
operations. We start with the single mode displacement operator,

$\displaystyle D(\alpha)=\exp(\alpha\hat{a}^{\dagger}-\alpha^{*}\hat{a}).$
(20)

which takes as argument the complex number $\alpha$. The annihilation and
quadrature operators transform as

$\displaystyle D^{\dagger}(\alpha)\hat{a}D(\alpha)=\hat{a}+\alpha,\quad
D^{\dagger}(\alpha)\hat{x}D(\alpha)=\hat{x}+\sqrt{2\hbar}\Re(\alpha),\quad
D^{\dagger}(\alpha)\hat{p}D(\alpha)=\hat{p}+\sqrt{2\hbar}\Im(\alpha),$ (21)

where $\alpha=\Re(\alpha)+i\Im(\alpha)$. Note in particular that the pure-
momentum displacement in Eq. (5b) is simply a displacement operator by a
purely imaginary amount, $Z(t_{1})=D(it_{1}/\sqrt{2\hbar})$. We also define a
multimode displacement operator as follows

$\displaystyle\mathcal{W}(\bm{r})=\bigotimes_{j=1}^{M}D_{j}\left(\tfrac{1}{\sqrt{2\hbar}}(x_{j}+ip_{j})\right),\text{
where }\bm{r}=(x_{1},\ldots,x_{M},p_{1}\ldots,p_{M})^{T}.$ (22)

By inserting the prefactor $\tfrac{1}{\sqrt{2\hbar}}$ in the argument of the
single-mode displacement operators we make the following relation
$\hbar$-independent

$\displaystyle\mathcal{W}^{\dagger}(\bm{r})\bm{\hat{r}}\mathcal{W}(\bm{r})=\bm{\hat{r}}+\bm{r}.$
(23)

We continue with the single-mode squeezing operator defined as

$\displaystyle\mathcal{T}(s)=\exp\left[\tfrac{s}{2}\left(\hat{a}^{\dagger
2}-\hat{a}^{2}\right)\right],$ (24)

which transforms the annihilation and quadrature operators as follows

$\displaystyle\mathcal{T}^{\dagger}(s)\hat{a}\mathcal{T}(s)=\cosh(s)\hat{a}+\sinh(s)\hat{a}^{\dagger},\quad\mathcal{T}^{\dagger}(s)\hat{x}\mathcal{T}(s)=e^{s}\hat{x},\quad\mathcal{T}^{\dagger}(s)\hat{p}\mathcal{T}(s)=e^{-s}\hat{p}.$
(25)

Finally, a multi-mode interferometer operator $\mathcal{U}$, parametrized by a
unitary matrix $\bm{U}$ acts as follows

$\displaystyle\mathcal{U}(\bm{U})^{\dagger}a_{i}\mathcal{U}(\bm{U})$
$\displaystyle=\sum_{i=1}^{M}U_{il}a_{l},$ (26)
$\displaystyle\mathcal{U}(\bm{U})^{\dagger}\hat{\bm{r}}\ \mathcal{U}(\bm{U})$
$\displaystyle=\bm{O}_{\bm{U}}\hat{\bm{r}},\quad\bm{O}_{\bm{U}}\equiv\left(\begin{array}[]{c|c}\Re\left(\bm{U}\right)&-\Im\left(\bm{U}\right)\\\
\hline\cr\Im\left(\bm{U}\right)&\Re\left(\bm{U}\right)\end{array}\right).$
(29)

Note that the interferometer operators respect the group structure of the
unitary group. This last property makes it easy to decompose an arbitrary
$M$-mode interferometer $\mathcal{U}(\bm{U})$ into a product of
interferometers acting on at most two modes by simply decomposing the
associated matrix $\bm{U}$ into a product of unitary matrices where each
unitary matrix couples at most two modes [48, 49, 50, 51].

With this notation we are ready to tackle the decomposition of the most
general operation generated by a quadratic Hamiltonian. We essentially follow
the elegant arguments exposed by Serafini [52] and write this operator as

$\displaystyle\mathcal{Q}(\bm{H},\bar{\bm{r}})=\exp\left[\tfrac{i}{\hbar}\left(\tfrac{1}{2}\hat{\bm{r}}^{T}\bm{H}\hat{\bm{r}}+\hat{\bm{r}}^{T}\bar{\bm{r}}\right)\right],$
(30)

where we assume $\bm{H}\in\mathbb{R}^{2M\times 2M}$ is a symmetric matrix and
$\bar{\bm{r}}\in\mathbb{R}^{2M}$. To obtain a decomposition we first eliminate
the linear part in the exponential by inserting resolutions of the identity in
terms of $\mathcal{W}$ as follows

$\displaystyle\mathcal{Q}(\bm{H},\bar{\bm{r}})$
$\displaystyle=\mathcal{W}(\bm{r})\mathcal{W}^{\dagger}(\bm{r})\mathcal{Q}(\bm{H},\bar{\bm{r}})\mathcal{W}(\bm{r})\mathcal{W}^{\dagger}(\bm{r})=\mathcal{W}(\bm{r})\exp\left[\tfrac{i}{\hbar}\left(\tfrac{1}{2}(\hat{\bm{r}}+\bm{r})^{T}\bm{H}(\hat{\bm{r}}+\bm{r})+(\hat{\bm{r}}+\bm{r})^{T}\bar{\bm{r}}\right)\right]\mathcal{W}^{\dagger}(\bm{r})$
$\displaystyle=e^{\frac{i}{2\hbar}\bm{r}^{T}\bm{H}\bm{r}+\frac{i}{\hbar}\bm{r}^{T}\bar{\bm{r}}}\mathcal{W}(\bm{r})\exp\left[\tfrac{i}{\hbar}\left(\tfrac{1}{2}\hat{\bm{r}}^{T}\bm{H}\hat{\bm{r}}+\hat{\bm{r}}^{T}\\{\bar{\bm{r}}+\bm{H}\bm{r}\\}\right)\right]\mathcal{W}^{\dagger}(\bm{r}).$
(31)

Note that the first term in the right-hand side of the last equation is an
inconsequential global phase and that if we pick as displacement
$\bm{r}=-\bm{H}^{-1}\bar{\bm{r}}$ then the term inside the curly braces in the
last equation vanishes and we are only left with an operator with a purely
quadratic generator,

$\displaystyle\mathcal{S}(\bm{H})=\exp\left[\tfrac{i}{2\hbar}\hat{\bm{r}}^{T}\bm{H}\hat{\bm{r}}\right].$
(32)

To decompose this gate we study its action on the canonical operators to find
(cf. Sec. 3.2.2 of Ref. [52])

$\displaystyle\mathcal{S}(\bm{H})^{\dagger}\hat{\bm{r}}\mathcal{S}(\bm{H})$
$\displaystyle=\bm{S}\hat{\bm{r}},$ (33) $\displaystyle\bm{S}$
$\displaystyle=e^{-\bm{\Omega}\bm{H}}.$ (34)

Note that $\mathcal{S}(\bm{H})$ is an operator acting on the Hilbert space of
$M$ modes while $S=e^{-\bm{\Omega}\bm{H}}$ is $2M\times 2M$ matrix. It is not
hard to show that the latter is a member of the real symplectic group [31],
i.e., that it satisfies

$\displaystyle\bm{S}\bm{\Omega}\bm{S}^{T}=\bm{\Omega},$ (35)

where $\bm{\Omega}$ is the symplectic form in Eq. (4). To continue with the
decomposition of the operator $\hat{S}(\bm{H})$ we will now turn our attention
to the singular-value decomposition of the matrix $\bm{S}$,

$\displaystyle\bm{S}=\bm{O}^{(1)}\bm{Z}\bm{O}^{(2)}.$ (36)

Since $\bm{S}$ is symplectic, one can construct the orthogonal matrices
$\bm{O}^{(1)},\bm{O}^{(2)}$ to also be symplectic and satisfy (cf. Appendix B
of Ref. [52])

$\displaystyle\bm{O}^{(i)}=\left(\begin{array}[]{c|c}\Re\left(\bm{U}^{(i)}\right)&-\Im\left(\bm{U}^{(i)}\right)\\\
\hline\cr\Im\left(\bm{U}^{(i)}\right)&\Re\left(\bm{U}^{(i)}\right)\end{array}\right),$
(39)

where $\bm{U}^{(i)}$ is a complex unitary matrix of size $M\times M$. This
form is precisely the same form appearing in Eq. (29). Invoking the membership
of $\bm{S}$ in the symplectic group one can also show that the singular values
of $\bm{S}$ come in pairs of the form $(z,1/z)$ and thus one can write the
diagonal symplectic matrix $\bm{Z}$ as

$\displaystyle\bm{Z}=\text{diag}(z_{1},\ldots,z_{M},1/z_{1},\ldots,1/z_{M}).$
(40)

By comparing Eq. (39) and Eq. (40) with Eq. (29) and Eq. (24) respectively one
identifies [30]

$\displaystyle\mathcal{S}(\bm{H})=\mathcal{U}(\bm{U}^{(1)})\left[\bigotimes_{j=1}^{M}\mathcal{T}_{j}(\log
z_{j})\right]\mathcal{U}(\bm{U}^{(2)}),$ (41)

which concludes the decomposition, since one can verify that the operators on
either of the equations above transform in exactly the same way the
quadratures $\bm{\hat{r}}$.

As a first example consider the decomposition of the single-mode ($M=1$) CV-
phase gate [53]

$\displaystyle
P(s)=\exp\left(i\frac{s}{2\hbar}\hat{x}^{2}\right)=\exp\left(\frac{i}{2\hbar}\hat{\bm{r}}\begin{pmatrix}s&0\\\
0&0\end{pmatrix}\hat{\bm{r}}^{T}\right),$ (42)

where without loss of generality we take $s\geq 0$, since
${P}(-s)={P}^{\dagger}(s)$. From this definition we easily find

$\displaystyle e^{-\bm{\Omega}\bm{H}}$ $\displaystyle=\begin{pmatrix}1&0\\\
s&1\end{pmatrix}=\begin{pmatrix}\cos\theta&-\sin\theta\\\
\sin\theta&\cos\theta\end{pmatrix}\begin{pmatrix}e^{r}&0\\\
0&e^{-r}\end{pmatrix}\begin{pmatrix}\sin\theta&\cos\theta\\\
-\cos\theta&\sin\theta\end{pmatrix},$ (43) $\displaystyle\frac{s}{2}$
$\displaystyle=\sinh
r,\quad\cos\theta=\frac{1}{\sqrt{1+e^{2r}}},\quad\sin\theta=\frac{e^{r}}{\sqrt{1+e^{2r}}},$
(44)

and then we can write

$\displaystyle P(s)=R(\theta)\mathcal{T}(r){R}(\theta-\pi/2),$ (45)

where we used the fact that a single-mode interferometer
$\mathcal{U}(e^{i\theta})$ is identical to a rotation gate ${R}(\theta)$.

As a second example let us consider the two-mode controlled-phase gate,

$\displaystyle
CZ(s)=\exp\left(i\tfrac{s}{\hbar}\hat{x}_{1}\hat{x}_{2}\right)=\exp\left[\frac{i}{2\hbar}\hat{\bm{r}}\left(\begin{smallmatrix}0&s&0&0\\\
s&0&0&0\\\ 0&0&0&0\\\ 0&0&0&0\end{smallmatrix}\right)\hat{\bm{r}}^{T}\right],$
(46) $\displaystyle
e^{-\bm{\Omega}\bm{H}}=\left(\begin{array}[]{cccc}1&0&0&0\\\ 0&1&0&0\\\
0&s&1&0\\\ s&0&0&1\\\ \end{array}\right)$ (51)
$\displaystyle=\left(\begin{array}[]{cccc}\cos\theta&0&0&-\sin\theta\\\
0&\cos\theta&-\sin\theta&0\\\ 0&\sin\theta&\cos\theta&0\\\
\sin\theta&0&0&\cos\theta\\\
\end{array}\right)\left(\begin{array}[]{cccc}e^{r}&0&0&0\\\ 0&e^{r}&0&0\\\
0&0&e^{-r}&0\\\ 0&0&0&e^{-r}\\\
\end{array}\right)\left(\begin{array}[]{cccc}\sin\theta&0&0&\cos\theta\\\
0&\sin\theta&\cos\theta&0\\\ 0&-\cos\theta&\sin\theta&0\\\
-\cos\theta&0&0&\sin\theta\\\ \end{array}\right),$ (64)

where $r$ and $\theta$ have the same functional dependence on the parameter
$s$ as in Eq. (44) and we assumed without loss of generality that $s\geq 0$.
Having the singular-value decomposition of the symplectic matrix we can easily
write

$\displaystyle
CZ(s)=BS(\theta)\left[\mathcal{T}_{1}(r)\otimes\mathcal{T}_{2}(r)\right]BS(\theta-\pi/2),$
(65)

where $BS(\theta)$ is the (symmetric) beamsplitter

$\displaystyle
BS(\theta)=\mathcal{U}\left(\left[\begin{smallmatrix}\cos\theta&i\sin\theta\\\
i\sin\theta&\cos\theta\end{smallmatrix}\right]\right)=\exp\left(i\theta\left[\hat{a}_{1}^{\dagger}\hat{a}_{2}+\hat{a}_{2}^{\dagger}\hat{a}_{1}\right]\right).$
(66)

###  4.2 Non-Gaussian Operators

In this section we study exact decomposition techniques for non-Gaussian
operators. The techniques apply broadly to gates that have generators made of
monomials that mutually commute. These techniques were first introduced by one
of the authors in Ref. [14] and provides a recipe for exact decompositions of
multi-mode gates $e^{itH}$, where the operator $H$ is of the form

$\displaystyle H=\left(\prod_{j=1}^{N-1}\hat{x}_{j}\right)\hat{x}_{N}^{n},$
(67)

or

$\displaystyle H=\hat{x}_{1}^{n_{1}}\hat{x}_{2}^{n_{2}},$ (68)

for $n$, $n_{1}$, and $n_{2}$ positive integers. As well as single-mode gates
$e^{itH}$ with

$H=\hat{x}^{N}.$ (69)

These gates can also be extended to include momentum quadrature operators
$\hat{p}_{j}$ by Fourier transforms acting on individual modes. Note that the
method works for any product where at most one operator has an exponent $n>1$,
so the label of the modes in Eq. (67) is arbitrary. It is also required that
$N$ is divisible by either 2 or 3, and in the multi-mode case, the product
$nN$ is divisible by 2 or 3, as well as $n_{1}$ and $n_{2}$ divisible by 2.
The intuition behind these restrictions is that the method relies on factoring
the overall power $N$ into two lower powers and being able to decompose
operators with exponents to those two lower powers. Because of the recursive
nature of some of the single mode gates which are decomposed with this method,
only decompositions of certain powers are known [32]. More specifically, all
single mode gates with power divisible by 2 or 3 are known, and as a result,
the overall powers of any decomposition must also be divisible by 2 or 3.
Although more restricted than the approximate commutator expansion method, the
set of gates for which these exact decompositions can be obtained encompasses
a large class of operators arising in several CV quantum algorithms and
simulations of bosonic systems [54, 3, 2, 23, 13, 1, 55].

The method relies on using the following three identities:
(i) unitary conjugation

$U^{\dagger}e^{itH}U=e^{itU^{\dagger}HU},$ (70)

(ii) a lemma to the Baker-Campbell-Hausdorff (BCH) formula

$\displaystyle e^{\hat{A}}\hat{B}e^{-\hat{A}}=\hat{B}+[\hat{A},\hat{B}]+$
$\displaystyle\frac{1}{2!}[\hat{A},[\hat{A},\hat{B}]]+\ldots,$ (71)

and (iii), the identity

$\displaystyle e^{i3\alpha^{2}t\hat{p}_{k}\hat{x}_{j}^{2}}=\>$ $\displaystyle
e^{i\frac{\alpha}{\hbar}\hat{x}_{j}\hat{x}_{k}}e^{i\frac{t}{\hbar}\hat{p}_{k}^{3}}e^{-i\alpha\hat{x}_{j}\hat{x}_{k}}e^{-i\frac{t}{\hbar}\hat{p}_{k}^{3}}e^{-i\frac{\alpha}{\hbar}\hat{x}_{j}\hat{x}_{k}}e^{i\frac{t}{\hbar}\hat{p}_{k}^{3}}e^{i\alpha\hat{x}_{j}\hat{x}_{k}}e^{-i\frac{t}{\hbar}\hat{p}_{k}^{3}}e^{i\frac{3\alpha^{3}t}{4\hbar}\hat{x}_{j}^{3}},$
(72)

with $\alpha$ and $t$ real parameters. The identity Eq. (72) comes from
repeated usage of both (i) and (ii), the details of which are shown in
Appendix A of Ref. [32], as well as a few derivations of other similar
identities in the appendix of that reference. A few examples can now be
studied to provide intuition on how these identities can be used.

Assume an exact decomposition is needed for the unitary operator
$e^{i\alpha\hat{x}_{j}\hat{x}_{k}\hat{x}_{l}}$. First, express the operator
$\hat{x}_{j}\hat{x}_{k}\hat{x}_{l}$ as a linear combination of polynomials of
degree three in the quadrature operators $\hat{x}_{j}$, $\hat{x}_{k}$, and
$\hat{x}_{l}$. This is done using the identity

$\displaystyle\hat{x}_{j}\hat{x}_{k}\hat{x}_{l}=$
$\displaystyle\tfrac{1}{6}[(\hat{x}_{j}+\hat{x}_{k}+\hat{x}_{l})^{3}-(\hat{x}_{j}+\hat{x}_{k})^{3}-(\hat{x}_{j}+\hat{x}_{l})^{3}-(\hat{x}_{k}+\hat{x}_{l})^{3}+\hat{x}_{j}^{3}+\hat{x}_{k}^{3}+\hat{x}_{l}^{3}],$
(73)

which implies the identity

$\displaystyle e^{i\alpha\hat{x}_{j}\hat{x}_{k}\hat{x}_{l}}=\>$ $\displaystyle
e^{\frac{i\alpha}{6}\left(\hat{x}_{j}+\hat{x}_{k}+\hat{x}_{l}\right)^{3}}e^{\frac{-i\alpha}{6}\left(\hat{x}_{j}+\hat{x}_{k}\right)^{3}}e^{\frac{-i\alpha}{6}\left(\hat{x}_{j}+\hat{x}_{l}\right)^{3}}e^{\frac{-i\alpha}{6}\left(\hat{x}_{k}+\hat{x}_{l}\right)^{3}}e^{\frac{i\alpha}{6}\hat{x}^{3}_{j}}e^{\frac{i\alpha}{6}\hat{x}^{3}_{k}}e^{\frac{i\alpha}{6}\hat{x}^{3}_{l}},$
(74)

since all the terms in the exponent commute. The right-hand side of this
equation still includes operators which are not contained within the universal
set. To decompose these, the identities in Eq. (70) and Eq. (71) can be used
with $U=e^{2i\hat{p}_{j}\hat{x}_{k}}$ as the unitary of conjugation, to arrive
at the following decompositions

$\displaystyle e^{i\alpha\left(\hat{x}_{j}+\hat{x}_{k}\right)^{3}}$
$\displaystyle=e^{i\hat{p}_{j}\hat{x}_{k}/\hbar}e^{i\alpha\hat{x}_{j}^{3}}e^{-i\hat{p}_{j}\hat{x}_{k}/\hbar},$
(75) $\displaystyle
e^{i\alpha\left(\hat{x}_{j}+\hat{x}_{k}+\hat{x}_{l}\right)^{3}}$
$\displaystyle=e^{i\hat{p}_{j}\hat{x}_{l}/\hbar}e^{i\alpha(\hat{x}_{j}+\hat{x}_{k})^{3}}e^{-i\hat{p}_{j}\hat{x}_{l}/\hbar}.$
(76)

To summarize this example, the exponent of the original operator
$\hat{x}_{j}\hat{x}_{k}\hat{x}_{l}$ is expressed as a linear combination of
polynomials of operators, and as such the full operator
$e^{i\alpha\hat{x}_{j}\hat{x}_{k}\hat{x}_{l}}$ can be written in terms of a
product of gates, each of which is decomposed exactly using two of the three
identities mentioned above.

To illustrate the method for a higher order single-mode gate, assume the
operator we wish to decompose is $e^{i\alpha\hat{x}_{j}^{4}}$. Following the
strategy used for the previous example, express the operator $\hat{x}_{j}^{4}$
as a linear combination of degree-four polynomials. More specifically, the
identity

$\hat{x}_{j}^{4}=(\hat{x}_{j}^{2}+\hat{x}_{k})^{2}-\hat{x}_{k}^{2}-2\hat{x}_{j}^{2}\hat{x}_{k},$
(77)

implies the relation

$e^{i\alpha\hat{x}_{j}^{4}}=e^{i\alpha(\hat{x}_{j}^{2}+\hat{x}_{k})^{2}}e^{-i\alpha\hat{x}_{k}^{2}}e^{-2i\alpha\hat{x}_{j}^{2}\hat{x}_{k}}.$
(78)

Here, the operator $e^{-i\alpha\hat{x}_{k}^{2}}$ is contained in the universal
set, while Eq. (72) gives an exact decomposition for
$e^{-i\alpha\hat{x}_{j}^{2}\hat{x}_{k}}$ up to a Fourier transform. Similar to
the previous example, the remaining term can be decomposed using unitary
conjugation:

$e^{i\alpha(\hat{x}_{j}^{2}+\hat{x}_{k})^{2}}=e^{i\hat{p}_{k}\hat{x}_{j}^{2}/\hbar}e^{i\alpha\hat{x}_{k}^{2}}e^{-i\hat{p}_{k}\hat{x}_{j}^{2}/\hbar},$
(79)

leading to a full decomposition for the target gate
$e^{i\alpha\hat{x}_{j}^{4}}$. Note that an additional ancillary mode $k$ was
required in this decomposition.

To extend the method to higher order single and multi-mode gates a more
general form of Eq. (72) is needed

$\displaystyle e^{2i\alpha^{2}\hat{p}_{k}\hat{x}_{j}^{N}}=\>$ $\displaystyle
e^{i\frac{\alpha}{\hbar}\hat{x}_{j}^{N-2}\hat{x}_{k}}e^{-i\alpha\hat{x}_{j}^{2}\hat{p}_{k}^{2}}e^{-2i\frac{\alpha}{\hbar}\hat{x}_{j}^{N-2}\hat{x}_{k}}e^{i\alpha\hat{x}_{j}^{2}\hat{p}_{k}^{2}}e^{i\frac{\alpha^{3}}{\hbar}\hat{x}_{j}^{2(N-1)}},$
(80)

where decompositions of single mode gates of lower order in $N$ are required.
Although, in this more general case the same strategy is used. Namely, the
target gate is expressed in terms of a linear combination of polynomials and
the resulting gates are decomposed using unitary conjugation or another lower
order decomposition. The derivation of this identity can be found in the
appendix of Ref. [32].

The same procedure involving creating polynomials and increasing overall order
can also be used to characterize decompositions for the group of operators
with $H$ given by Eq. (68). Although in this case the orders of each power
must be increased individually through the build up of polynomials.

Note that this method can provide exact decompositions for a wide variety of
useful operations [54, 3, 2, 23, 13, 1, 55], but the total gate count of the
decomposition scales exponentially as the overall power in each mode is
increased.

###  4.3 Generalizing exact decomposition methods

By allowing for additional constants to be added to the polynomials, the
scaling of the method in the previous section can be significantly improved
and the method extended to allow for the decomposition of a larger category of
operators. This is shown in this section with the use of a general formula for
representing products as sums from Ref. [56].

As was discussed in the previous section, one of the types of Hamiltonians
which can be decomposed exactly is of the form $e^{it\hat{H}}$, where the
operator $H$ is of the form

$\displaystyle
H=\left(\prod_{j=1}^{n-1}\hat{x}_{j}\right)\hat{x}_{n}^{s_{n}}.$ (81)

The method for decomposing these Hamiltonians worked for any product where at
most one operator has an exponent $s_{n}>1$. It also required that $n$ is
divisible by either 2 or 3, and that the product $ns_{n}$ must also be
divisible by 2 or 3.

A lemma from Ref. [56] can be used to find decompositions which no longer have
the first restriction that at most only one operator has an exponent
$s_{n}>1$. In other words, the position quadrature operators in the product
can be taken to any arbitrary powers, as long as the overall sum of the powers
is divisible by 2 or 3.

The following equation is valid for real numbers [56], but because each of the
$\hat{x}_{i}$s commute with each other it is also valid for a product of
operators.

$\displaystyle\hat{x}^{s_{1}}_{1}\hat{x}^{s_{2}}_{2}\cdots\hat{x}^{s_{n}}_{n}=\frac{1}{s!}\sum_{v_{1}=0}^{s_{1}}\cdots\sum_{v_{n}=0}^{s_{n}}(-1)^{\sum_{i=1}^{n}v_{i}}\binom{s_{1}}{v_{1}}\cdots\binom{s_{n}}{v_{n}}\left(\sum_{i=1}^{n}h_{i}\hat{x}_{i}\right)^{s},$
(82)

where $s=s_{1}+s_{2}+\cdots+s_{n}$ and $h_{i}=s_{i}/2-v_{i}$. When at least
one of the $s_{i}$ are odd repeated terms in the sum may be grouped together
reducing the expression to

$\displaystyle\hat{x}^{s_{1}}_{1}\hat{x}^{s_{2}}_{2}\cdots\hat{x}^{s_{n}}_{n}=\frac{2}{s!}\sum_{v_{1}=0}^{(s_{1}-1)/2}\sum_{v_{2}=0}^{s_{2}}\cdots\sum_{v_{n}=0}^{s_{n}}(-1)^{\sum_{i=1}^{n}v_{i}}\binom{s_{1}}{v_{1}}\cdots\binom{s_{n}}{v_{n}}\left(\sum_{i=1}^{n}h_{i}\hat{x}_{i}\right)^{s},$
(83)

where without loss of generality it is assumed that $s_{1}$ is odd. Following
the procedure in the last section, each term on the right-hand side may be
created through unitary conjugation, where each operation in the unitary
conjugation can be decomposed exactly. The challenge then is to be able to
implement gates with the correct coefficients which depend on each $s_{i}$. We
note that, to the best of our knowledge, the idea of using the identity in Eq.
(82) to decompose CV gates has not been presented in the literature before.

As an example, expand the product $\hat{x}_{1}\hat{x}^{2}_{2}\hat{x}^{3}_{3}$
using Eq. (83). Note that this product is not covered by the $H$ in Eq. (67),
but has an overall power $s$ which is even.

$\displaystyle\hat{x}_{1}\hat{x}^{2}_{2}\hat{x}^{3}_{3}$
$\displaystyle=\frac{1}{360}\sum_{v_{1}=0}^{0}\sum_{v_{2}=0}^{2}\sum_{v_{3}=0}^{3}(-1)^{\sum_{i=1}^{n}v_{i}}\binom{1}{v_{1}}\binom{2}{v_{2}}\binom{3}{v_{3}}\left(h_{1}\hat{x}_{1}+h_{2}\hat{x}_{2}+h_{3}\hat{x}_{3}\right)^{6}$
(84)
$\displaystyle=\frac{1}{360}\left(\hat{x}_{2}+\frac{1}{2}\hat{x}_{1}+\frac{3}{2}\hat{x}_{3}\right)^{6}-\frac{1}{120}\left(\hat{x}_{2}+\frac{1}{2}\hat{x}_{1}+\frac{1}{2}\hat{x}_{3}\right)^{6}+\frac{1}{120}\left(\hat{x}_{2}+\frac{1}{2}\hat{x}_{1}-\frac{1}{2}\hat{x}_{3}\right)^{6}$
$\displaystyle\;\;\;-\frac{1}{360}\left(\hat{x}_{2}+\frac{1}{2}\hat{x}_{1}-\frac{3}{2}\hat{x}_{3}\right)^{6}-\frac{1}{11520}\left(\hat{x}_{1}+3\hat{x}_{3}\right)^{6}+\frac{1}{3840}\left(\hat{x}_{1}+\hat{x}_{3}\right)^{6}-\frac{1}{3840}\left(\hat{x}_{1}-\hat{x}_{3}\right)^{6}$
$\displaystyle\;\;\;+\frac{1}{11520}\left(\hat{x}_{1}-3\hat{x}_{3}\right)^{6}+\frac{1}{360}\left(\hat{x}_{2}-\frac{1}{2}\hat{x}_{1}-\frac{3}{2}\hat{x}_{3}\right)^{6}-\frac{1}{120}\left(\hat{x}_{2}-\frac{1}{2}\hat{x}_{1}-\frac{1}{2}\hat{x}_{3}\right)^{6}$
$\displaystyle\;\;\;+\frac{1}{120}\left(\hat{x}_{2}-\frac{1}{2}\hat{x}_{1}+\frac{1}{2}\hat{x}_{3}\right)^{6}-\frac{1}{360}\left(\hat{x}_{2}-\frac{1}{2}\hat{x}_{1}+\frac{3}{2}\hat{x}_{3}\right)^{6}.$

Each of the resulting terms can be made with unitary conjugation, where the
central operator will have the single mode with unit coefficient in the
brackets of each term.

This relation along with the identities used in the last section (Eq. (70),
(71), and (72) ) can also be used to find decompositions for operators that
were already covered in the last section, in order to compare gate counts. For
the operator where $H=\hat{x}_{1}\hat{x}_{2}\hat{x}_{3}$ the total gate count
using the previous method is 17 gates from the universal set, whereas when
using the updated method described here it is 20 gates. For a higher order
operation like $H=\hat{x}_{1}\hat{x}_{2}\hat{x}_{3}\hat{x}_{4}$ the total gate
count from the universal set was 440 but with this method is 280. Other
examples comparing the two exact methods are given in Table 1. In some cases
it seems like the method in this section requires more unitary conjugation but
less usage of single mode gates. So as the overall order increases the single
mode gates become harder to decompose but the number of operators needed for
the unitary conjugation here do not increase as fast.

A potential problem with this method, as well as the exact decomposition
methods of previous sections, are terms with a mix of $\hat{x}$ and $\hat{p}$
operators which act on the same mode. To deal with these terms the commutator
expansion method can be used, which was described in section 3, or one can use
another approximate method introduced in the next section.

Target gate | Commutator approx. | Exact decomposition | Exact decomposition
---|---|---|---
| ($10^{-3}$ precision) |  | using Eq. (82, 83)
$e^{it\hat{x}^{4}}$ |  $1.8\times 10^{4}$ gates | 29 gates | -
$e^{it\hat{x}_{j}^{2}\hat{x}_{k}^{2}}$ |  $2.8\times 10^{4}$ gates | 119 gates | 279 gates
$e^{it\hat{x}_{j}\hat{x}_{k}^{3}}$ |  $2.9\times 10^{8}$ gates | 269 gates | 93 gates
$e^{it\hat{x}_{j}\hat{x}_{k}\hat{x}_{l}}$ |  $4.2\times 10^{8}$ gates | 17 gates | 20 gates
$e^{it\hat{x}_{j}^{2}\hat{x}_{k}\hat{x}_{l}}$ |  $1.4\times 10^{9}$ gates | 249 gates | 198 gates
$e^{it\hat{x}_{j}\hat{x}_{k}\hat{x}_{l}\hat{x}_{m}}$ |  $6.9\times 10^{13}$ gates | 440 gates | 280 gates
$e^{it\hat{x}^{6}}$ |  $1.2\times 10^{13}$ gates | 809 gates | -
$e^{it\hat{x}_{j}^{2}\hat{x}_{k}^{4}}$ |  $2.4\times 10^{13}$ gates | 3320 gates | 12165 gates
Table 1: Gate counts for decompositions of some common operations, using the
standard commutator approximation as well as the exact decomposition methods
demonstrated in this work. The final column uses the generalization for multi-
mode exact decompositions using Eq. (82, 83). The gate counts neglect any
Fourier transforms used as they are inexpensive to implement experimentally.

##  5 Asymptotically exact decompositions with squeezing and displacement

In the last two sections we discussed the exact decompositions of Gaussian
operations and a subset of non-Gaussian operations. To be precise, for the
latter group we focused on non-Gaussian gates that have generators of
monomials that mutually commute. In this section we discuss some recent ideas
on how to decompose more complicated gates in which one mixes in the
generators of position and momentum of the same mode. The method was
introduced by Yanagimoto et al. [33]. Both squeezing and displacement are
operations which can be implemented with linear optics. As a result these
operations are Gaussian, and are not costly to implement. In fact, the
displacement operator is equivalent to a Gaussian operator which is already
part of the universal set in Eq. (5). The squeezing operation was given in Eq.
(24), and the displacement operation was given in Eq. (5b).

In order to use these operations to form a gate decomposition method, we can
examine their effect in unitary conjugation with the quadrature operators.
Squeezing has the effect of multiplying $\hat{x}$ and $\hat{p}$ operators by a
constant

$\displaystyle\mathcal{T}^{\dagger}(\log\lambda)\hat{x}_{i}\mathcal{T}(\log\lambda)=\lambda\hat{x}_{i},$
(85)
$\displaystyle\mathcal{T}^{\dagger}(\log\lambda)\hat{p}_{i}\mathcal{T}(\log\lambda)=\lambda^{-1}\hat{p}_{i},$
(86)

and displacement has the effect of adding a constant as in Eq. (21) (which
will we will specify by $y$) [57]. Applying both of these to the annihilation
and creation operators maps them to a new effective operation

$\hat{a}^{\dagger}\rightarrow\hat{a}^{\dagger}_{\text{eff}}=\mathcal{T}^{\dagger}(\log\lambda)D^{\dagger}(y)\hat{a}^{\dagger}D(y)\mathcal{T}(\log\lambda)=\frac{\lambda}{\sqrt{2\hbar}}\hat{x}-i\frac{\lambda^{-1}}{\sqrt{2\hbar}}\hat{p}+y.$
(87)

With this mapping and the correct choice of constants $\lambda$ and $y$, terms
with a greater overall power in $\hat{x}$ or $\hat{p}$ can be picked out while
the contribution of other terms can be ignored. To illustrate this consider
transforming the Kerr operation into a cubic phase gate, as in [33]. Starting
with the driven Kerr Hamiltonian

$\displaystyle
H_{\text{Kerr}}=-\frac{\chi}{2}\hat{a}^{\dagger^{2}}\hat{a}^{2}+\delta\hat{a}^{\dagger}\hat{a}+\beta(\hat{a}+\hat{a}^{\dagger}),$
(88)

where $\delta$ and $\beta$ are the detuning and drive constants respectively,
and $\chi$ is associated with the self phase modulation of a propagating pulse
in a Kerr medium [33]. We can then apply squeezing and displacement to find
the effective Hamiltonian

$\displaystyle H_{\text{Kerr}}^{\text{eff}}=$
$\displaystyle\mathcal{T}^{\dagger}(\log\lambda)D^{\dagger}(y)H_{\text{Kerr}}D(y)\mathcal{T}(\log\lambda)$
(89) $\displaystyle=$
$\displaystyle-\frac{\chi}{8}(\lambda^{4}\hat{x}^{4}+\hat{p}\hat{x}^{2}\hat{p}+\hat{x}\hat{p}^{2}\hat{x}+\lambda^{-4}\hat{p}^{4})-\frac{\sqrt{\hbar}}{\sqrt{2}}\chi\lambda^{3}y\hat{x}^{3}-\frac{\sqrt{\hbar}}{\sqrt{2}}\chi\lambda^{-1}y\hat{p}\hat{x}\hat{p}$
(90) $\displaystyle+\frac{\hbar}{2}\lambda^{2}(-3\chi
y^{2}+\chi+\delta)\hat{x}^{2}+\frac{\hbar}{2}\lambda^{-2}(-\chi
y^{2}+\chi+\delta)\hat{p}^{2}+\sqrt{2}\hbar^{3/2}\lambda(-\chi y^{3}+\chi
y+\delta y+\beta)\hat{x}.$

Some terms will then cancel by choosing $\delta=3\chi y^{2}-\chi$ and
$\beta=-2\chi y^{3}$ to get

$\displaystyle
H_{\text{Kerr}}=-\frac{\chi\lambda^{3}y}{\sqrt{2}}\left(\hat{x}^{3}+\frac{\lambda
y^{-1}}{4\sqrt{2}}\hat{x}^{4}-\sqrt{2\hbar}\lambda^{-5}y\hat{p}^{2}+\mathcal{O}(\lambda^{-4})\right).$
(91)

With the correct choice of the displacement constant $y$, namely
$y\sim\lambda^{3}$, the effective Hamiltonian takes on the form of a cubic
operation to leading order in $\lambda$. Although this decomposition can be
done with only squeezing and displacement, and a tuning of their respective
parameters, we note that squeezing to arbitrary levels can be difficult to
implement [30].

##  6 Discussion

In this work we outlined approximate and exact gate decomposition methods for
CV quantum computers. For approximate methods we reviewed a commutator
expansion method which expressed quadrature operators in terms of commutators,
as well as another decomposition procedure which used squeezing and
displacement operations in unitary conjugation to approximate gates. We also
examined exact decomposition methods both for Gaussian and non-Gaussian
operations. For Gaussian operators we showed how to express the target gate in
terms of displacement, interferometer and squeezing gates. In order to
decompose non-Gaussian operations exactly we used a technique where the target
gate was expressed in terms of a linear combination of polynomials and then
decomposed using unitary conjugation, a lemma to the Baker-Campbell-Hausdorff
formula (Eq. (70) and (71)) and other known lower order exact decompositions.
This procedure was generalized using Eq. (82, 83) to allow for any product of
quadrature operators to be expressed as a linear combination of polynomials of
the same operators. To illustrate these methods we examined examples of
applying them to some simpler operations. The decompositions were all made
given the chosen universal set in Eq. (5).

The commutator expansion method provided a comprehensive CV decomposition
method that allows for the decomposition of any operator into gates from a
universal set. Although it allowed for the decomposition of any operator, the
recursive structure of the method causes exponential build up in the number of
gates required as the order of the operator to be decomposed increases.
Combined with the need to repeat the final sequence of gates to increase the
precision of the decomposition, the overall circuit depth may become very
large. As an example we examined the decomposition of the gate
$e^{it\left(\hat{x}^{2}\hat{p}+\hat{p}\hat{x}^{2}\right)}$ in Eq. (15), as
well as the Kerr operation. For these examples the method provided
decompositions requiring small sequences of gates from the universal set, but
in order to achieve a modest precision, the product of gates needed to be
repeated many times.

While more limited in the scope of what can be decomposed for non-Gaussian
gates, exact decompositions can in some cases significantly lower the gate
count needed for a decomposition when compared to approximate methods. For
example, the gate used as an example of applying the commutator expansion
method ($e^{it\left(\hat{x}^{2}\hat{p}+\hat{p}\hat{x}^{2}\right)}$) cannot be
decomposed exactly with the methods shown in this work, due to the mixture of
$\hat{x}$ and $\hat{p}$ quadrature operators of the same mode. Although
another multi-mode operation, namely
$e^{it\hat{x}_{j}\hat{x}_{k}\hat{x}_{l}}$, can be decomposed exactly with only
17 gates from the universal set, whereas it would take $4.2\times 10^{8}$
gates using the commutator expansion method for a precision of $10^{-3}$. The
set of gates for which exact decompositions can be obtained still encompass a
large class of useful CV operators [54, 3, 2, 23, 13, 1, 55], and can be
expanded by using Eq. (82, 83). This expanded method for multi-mode operations
can also be used to decompose operations that could already be decomposed with
the regular exact method to compare gate counts. Some example gate counts can
be seen in Table 1. In some cases it seems as though the expanded method
requires more unitary conjugation but less usage of single-mode gates and as a
result can require fewer gates. It is still unclear whether exact
decompositions can be generalized to include all CV operations.

The final method examined in this work uses squeezing and displacement
operators in unitary conjugation to approximate an operator into another
contained within the universal set. While no longer being exact, this method
allowed for the precision of the decomposition to be increased by simply
tuning the parameters of the final sequence of gates instead of having to
repeat them. To demonstrate this method we examined how the Kerr operator can
be approximated in terms of a cubic gate, following the procedure in [33].

By putting together all these results in a single place, it is our hope that
this manuscript becomes an easy to use reference for researchers in the area,
and moreover, that it inspires further work in finding more efficient and
compact decompositions of CV gates.

## Acknowledgments

T.K. thanks C. Weedbrook for valuable discussions. T.K. and N.Q. thank N.
Killoran and L. Madsen for comments on the manuscript.

## References

  * Lau et al. [2017] Hoi-Kwan Lau, Raphael Pooser, George Siopsis, and Christian Weedbrook.  Quantum machine learning over infinite dimensions.  _Phys. Rev. Lett_ , 118:080501, 2017.  doi: 10.1103/PhysRevLett.118.080501.
  * Kalajdzievski et al. [2018] Timjan Kalajdzievski, Christian Weedbrook, and Patrick Rebentrost.  Continuous-variable gate decomposition for the Bose-Hubbard model.  _Phys. Rev. A_ , 97(6):062311, 2018.  doi: 10.1103/PhysRevA.97.062311.
  * Arrazola et al. [201908] Juan Miguel Arrazola, Timjan Kalajdzievski, Christian Weedbrook, and Seth Lloyd.  Quantum algorithm for non-homogeneous linear partial differential equations.  _Phys. Rev. A_ , 100:032306, 201908.  doi: 10.1103/PhysRevA.100.032306.
  * Sefi et al. [2013] Seckin Sefi, Vishal Vaibhav, and Peter van Loock.  A measurement-induced optical kerr interaction.  _Phys. Rev. A_ , 88:012303, 2013.  doi: 10.1103/PhysRevA.88.012303.
  * Dawson and Nielsen [2006] Christopher M. Dawson and Michael A. Nielsen.  The Solovay-Kitaev algorithm.  _Quantum Inf. Comput._ , 6:1, 2006.  doi: 10.5555/2011679.2011685.
  * Amy et al. [2013] Matthew Amy, Dmitri Maslov, Michele Mosca, and Martin Roetteler.  A meet-in-the-middle algorithm for fast synthesis of depth-optimal quantum circuits.  _IEEE Trans. Comput. Aided Des. Integr. Circuits Syst._ , 32(6):818–830, 2013.  doi: 10.1109/TCAD.2013.2244643.
  * Lloyd [1995] Seth Lloyd.  Almost any quantum logic gate is universal.  _Phys. Rev. Lett._ , 75(2):346, 1995.  doi: 10.1103/PhysRevLett.75.346.
  * DiVincenzo [1995] David P DiVincenzo.  Two-bit gates are universal for quantum computation.  _Phys. Rev. A_ , 51(2):1015, 1995.  doi: 10.1103/PhysRevA.51.1015.
  * Barenco et al. [1995] 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.  _Phys. Rev. A_ , 52(5):3457, 1995.  doi: 10.1103/PhysRevA.52.3457.
  * Lloyd [1996] Seth Lloyd.  Universal quantum simulators.  _Science_ , 23:1073, 1996.  doi: 10.1126/science.273.5278.1073.
  * Lloyd and Braunstein [1999] Seth Lloyd and Samuel L. Braunstein.  Quantum computation over continuous variables.  _Phys. Rev. Lett_ , 82:1784, 1999.  doi: 10.1103/PhysRevLett.82.1784.
  * Bravyi and Kitaev [2005] Sergey Bravyi and Alexei Kitaev.  Universal quantum computation with ideal Clifford gates and noisy ancillas.  _Phys. Rev. A_ , 71(2):022316, 2005.  doi: 10.1103/PhysRevA.71.022316.
  * Sefi and van Loock [2011] Seckin Sefi and Peter van Loock.  How to decompose arbitrary continuous-variable quantum operations.  _Phys. Rev. Lett._ , 107:170501, 2011.  doi: 10.1103/PhysRevLett.107.170501.
  * Kalajdzievski and Arrazola [2019] Timjan Kalajdzievski and Juan Miguel Arrazola.  Exact gate decompositions for photonic quantum computing.  _Phys. Rev. A_ , 99:022341, 2019.  doi: 10.1103/PhysRevA.99.022341.
  * Kitaev [1997] A Yu Kitaev.  Quantum computations: algorithms and error correction.  _Russ. Math. Surv._ , 52(6):1191–1249, 1997.  doi: 10.1070/RM1997v052n06ABEH002155.
  * Kliuchnikov et al. [2013] Vadym Kliuchnikov, Dmitri Maslov, and Michele Mosca.  Asymptotically optimal approximation of single qubit unitaries by Clifford and T circuits using a constant number of ancillary qubits.  _Phys. Rev. Lett._ , 110(19):190502, 2013.  doi: 10.1103/PhysRevLett.110.190502.
  * Kliuchnikov et al. [2014] Vadym Kliuchnikov, Alex Bocharov, and Krysta M Svore.  Asymptotically optimal topological quantum compiling.  _Phys. Rev. Lett._ , 112(14):140504, 2014.  doi: 10.1103/PhysRevLett.112.140504.
  * Kliuchnikov and Yard [2015] Vadym Kliuchnikov and Jon Yard.  A framework for exact synthesis.  _arXiv:1504.04350_ , 2015.
  * Bocharov et al. [2015] Alex Bocharov, Martin Roetteler, and Krysta M Svore.  Efficient synthesis of universal repeat-until-success quantum circuits.  _Phys. Rev. Lett._ , 114(8):080502, 2015.  doi: 10.1103/PhysRevLett.114.080502.
  * Menicucci et al. [2006] Nicolas C. Menicucci, Peter van Loock, Mile Gu, Christian Weedbrook, Timothy C. Ralph, and Michael A. Nielsen.  niversal quantum computation with continuous-variable cluster states.  _Phys. Rev. Lett_ , 97:110501, 2006.  doi: 10.1103/PhysRevLett.97.110501.
  * Gu et al. [2009a] Mile Gu, Christian Weedbrook, Nicolas C. Menicucci, Timothy C. Ralph, and Peter van Loock.  Quantum computing with continuous-variable clusters.  _Phys. Rev. A_ , 79:062318, 2009a.  doi: 10.1103/PhysRevA.79.062318.
  * Weedbrook et al. [2012] Christian Weedbrook, Stefano Pirandola, Raúl García-Patrón, Nicolas J. Cerf, Timothy C. Ralph, Jeffrey H. Shapiro, and Seth Lloyd.  Gaussian quantum information.  _Rev. Mod. Phys._ , 84:621, 2012.  doi: 10.1103/RevModPhys.84.621.
  * Sowinski et al. [2012] Tomasz Sowinski, Omjyoti Dutta, Philipp Hauke, Luca Tagliacozzo, and Maciej Lewenstein.  Dipolar molecules in optical lattices.  _Phys. Rev. Lett._ , 108:115301, 2012.  doi: 10.1103/PhysRevLett.108.115301.
  * Myers and Ralph [2011] CR Myers and TC Ralph.  Coherent state topological cluster state production.  _New J. Phys._ , 13(11):115015, 2011.  doi: 10.1088/1367-2630/13/11/115015.
  * Ralph et al. [2003] Timothy C Ralph, Alexei Gilchrist, Gerard J Milburn, William J Munro, and Scott Glancy.  Quantum computation with optical coherent states.  _Phys. Rev. A_ , 68(4):042319, 2003.  doi: 10.1103/PhysRevA.68.042319.
  * Pantaleoni et al. [2020] Giacomo Pantaleoni, Ben Q Baragiola, and Nicolas C Menicucci.  Modular bosonic subsystem codes.  _Phys. Rev. Lett._ , 125(4):040501, 2020.  doi: 10.1103/PhysRevLett.125.040501.
  * Gottesman et al. [2001] Daniel Gottesman, Alexei Kitaev, and John Preskill.  Encoding a qubit in an oscillator.  _Phys. Rev. A_ , 64(1):012310, 2001.  doi: 10.1103/PhysRevA.64.012310.
  * Hatano and Suzuki [2005] Naomichi Hatano and Masuo Suzuki.  Finding exponential product formulas of higher orders.  In A. Das and B.K. Chakrabarti, editors, _Quantum Annealing and Other Optimization Methods_ , pages 37–68. Springer, Berlin, 2005.
  * Wiebe et al. [2010] Nathan Wiebe, Dominic W. Berry, Peter Hoyer, and Barry C. Sanders.  Higher order decompositions of ordered operator exponentials.  _J. Phys. A: Math. Theor._ , 43:065203, 2010.  doi: 10.1088/1751-8113/43/6/065203.
  * Braunstein [2005] Samuel L Braunstein.  Squeezing as an irreducible resource.  _Phys. Rev. A_ , 71(5):055801, 2005.  doi: 10.1103/PhysRevA.71.055801.
  * Dutta et al. [1995] Biswadeb Dutta, N Mukunda, R Simon, et al.  The real symplectic groups in quantum mechanics and optics.  _Pramana_ , 45(6):471–497, 1995.  doi: 10.1007/BF02848172.
  * Kalajdzievski [2020] Timjan Kalajdzievski.  _Exact Gate Decompositions For Photonic Quantum Computers_.  PhD thesis, York University, 2020.  URL https://yorkspace.library.yorku.ca/xmlui/handle/10315/37435.
  * Yanagimoto et al. [2020] Ryotatsu Yanagimoto, Tatsuhiro Onodera, Edwin Ng, Logan G. Wright, Peter L. McMahon, and Hideo Mabuchi.  Engineering a Kerr-based deterministic cubic phase gate via gaussian operation.  _Phys. Rev. Lett._ , 124:240503, 2020.  doi: 10.1103/PhysRevLett.124.240503.
  * Yukawa et al. [2013] Mitsuyoshi Yukawa, Kazunori Miyata, Hidehiro Yonezawa, Petr Marek, Radim Filip, and Akira Furusawa.  Emulating quantum cubic nonlinearity.  _Phys. Rev. A_ , 88(5):053816, 2013.  doi: 10.1103/PhysRevA.88.053816.
  * Gu et al. [2009b] Mile Gu, Christian Weedbrook, Nicolas C Menicucci, Timothy C Ralph, and Peter van Loock.  Quantum computing with continuous-variable clusters.  _Phys. Rev. A_ , 79(6):062318, 2009b.  doi: 10.1103/PhysRevA.79.062318.
  * Marshall et al. [2015] Kevin Marshall, Raphael Pooser, George Siopsis, and Christian Weedbrook.  Repeat-until-success cubic phase gate for universal continuous-variable quantum computation.  _Phys. Rev. A_ , 91(3):032321, 2015.  doi: 10.1103/PhysRevA.91.032321.
  * Sabapathy and Weedbrook [2018] Krishna Kumar Sabapathy and Christian Weedbrook.  ON states as resource units for universal quantum computation with photonic architectures.  _Phys. Rev. A_ , 97(6):062315, 2018.  doi: 10.1103/PhysRevA.97.062315.
  * Sabapathy et al. [2019] Krishna Kumar Sabapathy, Haoyu Qi, Josh Izaac, and Christian Weedbrook.  Production of photonic universal quantum gates enhanced by machine learning.  _Phys. Rev. A_ , 100(1):012326, 2019.  doi: 10.1103/PhysRevA.100.012326.
  * Marek et al. [2018] Petr Marek, Radim Filip, Hisashi Ogawa, Atsushi Sakaguchi, Shuntaro Takeda, Jun ichi Yoshikawa, and Akira Furusawa.  General implementation of arbitrary nonlinear quadrature phase gates.  _Phys. Rev. A_ , 97:022329, 2018.  doi: 10.1103/PhysRevA.97.022329.
  * Hillmann et al. [2020] Timo Hillmann, Fernando Quijandría, Göran Johansson, Alessandro Ferraro, Simone Gasparinetti, and Giulia Ferrini.  Universal gate set for continuous-variable quantum computation with microwave circuits.  _Phys. Rev. Lett._ , 125(16):160501, 2020.  doi: 10.1103/PhysRevLett.125.160501.
  * Weinstein et al. [2001] Yaakov S. Weinstein, Seth Lloyd, and David G. Cory.  Implementation of the quantum Fourier transform.  _Phys. Rev. Lett._ , 86:1889, 2001.  doi: 10.1103/PhysRevLett.86.1889.
  * Patra and Braunstein [2011] Manas K. Patra and Samuel L. Braunstein.  Quantum Fourier transform, heisenberg groups and quasiprobability distributions.  _New J. Phys._ , 13:063013, 2011.  doi: 10.1088/1367-2630/13/6/063013.
  * Magnus [1954] Wilhelm Magnus.  On the exponential solution of differential equations for a linear operator.  _Commun. Pure Appl. Math._ , 7(4):649–673, 1954\.  doi: 10.1002/cpa.3160070404.
  * Trotter [1959] Hale F Trotter.  On the product of semi-groups of operators.  _Proc. Am. Math. Soc._ , 10(4):545–551, 1959\.  doi: 10.2307/2033649.
  * Suzuki [1976] Masuo Suzuki.  Generalized Trotter’s formula and systematic approximants of exponential operators and inner derivations with applications to many-body problems.  _Commun. Math. Phys._ , 51:183, 1976.  doi: 10.1007/BF01609348.
  * Childs et al. [2018] Andrew M Childs, Dmitri Maslov, Yunseong Nam, Neil J Ross, and Yuan Su.  Toward the first quantum simulation with quantum speedup.  _Proc. Natl. Acad. Sci. U.S.A._ , 115:9456–9461, 2018.  doi: 10.1073/pnas.1801723115.
  * Barnett and Radmore [2002] Stephen Barnett and Paul M Radmore.  _Methods in theoretical quantum optics_ , volume 15.  Oxford University Press, 2002.
  * Reck et al. [1994] Michael Reck, Anton Zeilinger, Herbert J Bernstein, and Philip Bertani.  Experimental realization of any discrete unitary operator.  _Phys. Rev. Lett._ , 73(1):58, 1994.  doi: 10.1103/PhysRevLett.73.58.
  * Clements et al. [2016] William R Clements, Peter C Humphreys, Benjamin J Metcalf, W Steven Kolthammer, and Ian A Walmsley.  Optimal design for universal multiport interferometers.  _Optica_ , 3(12):1460–1465, 2016.  doi: 10.1364/OPTICA.3.001460.
  * de Guise et al. [2018] Hubert de Guise, Olivia Di Matteo, and Luis L Sánchez-Soto.  Simple factorization of unitary transformations.  _Phys. Rev. A_ , 97(2):022328, 2018.  doi: 10.1103/PhysRevA.97.022328.
  * Su et al. [2019] Daiqin Su, Ish Dhand, Lukas G Helt, Zachary Vernon, and Kamil Brádler.  Hybrid spatiotemporal architectures for universal linear optics.  _Phys. Rev. A_ , 99(6):062301, 2019.  doi: 10.1103/PhysRevA.99.062301.
  * Serafini [2017] Alessio Serafini.  _Quantum continuous variables: a primer of theoretical methods_.  CRC press, 2017.
  * Fiurášek [2003] Jaromír Fiurášek.  Unitary-gate synthesis for continuous-variable systems.  _Phys. Rev. A_ , 68(2):022304, 2003.  doi: 10.1103/PhysRevA.68.022304.
  * Sparrow et al. [2018] Chris Sparrow, Enrique Martín-López, Nicola Maraviglia, Alex Neville, Christopher Harrold, Jacques Carolan, Yogesh N Joglekar, Toshikazu Hashimoto, Nobuyuki Matsuda, Jeremy L O’Brien, et al.  Simulating the vibrational quantum dynamics of molecules using photonics.  _Nature_ , 557(7707):660, 2018.  doi: 10.1038/s41586-018-0152-9.
  * Rebentrost et al. [2018] Patrick Rebentrost, Brajesh Gupt, and Thomas R. Bromley.  Photonic quantum algorithm for Monte Carlo integration.  _arXiv:1809.02579_ , 2018.
  * Kan [2008] Raymond Kan.  From moments of sum to moments of product.  _J. Multivar. Anal._ , 99:542, 2008.  doi: 10.1016/j.jmva.2007.01.013.
  * Killoran et al. [2019] Nathan Killoran, Josh Izaac, Nicolás Quesada, Ville Bergholm, Matthew Amy, and Christian Weedbrook.  Strawberry fields: A software platform for photonic quantum computing.  _Quantum_ , 3:129, 2019.  doi: 10.22331/q-2019-03-11-129.


s}_{t}\,\mathbb{C}^{1,\nabla}_{t}(\mathbb{r})\cdot\nabla\times\mathbb{D}^{1,\sigma}_{t}(\mathbb{r})\,,$
(24) $\displaystyle\mathcal{E}_{\rm pair}(\mathbb{r})$
$\displaystyle=\sum_{q=p,n}g_{q}(\mathbb{r})\,\Big{[}\tilde{D}_{q}^{1,1}(\mathbb{r})\,\tilde{D}_{q}^{1,1}(\mathbb{r})$
$\displaystyle\phantom{\sum_{q=p,n}g_{q}(\mathbb{r})\,\Big{[}}+\tilde{C}_{q}^{1,1}(\mathbb{r})\,\tilde{C}_{q}^{1,1}(\mathbb{r})\Big{]}\,,$
(25)

where the functions $g_{q}(\mathbb{r})$ are defined in Eq. (9b). This notation
unambiguously separates the time-even contribution
$\propto\tilde{D}_{q}^{1,1}(\mathbb{r})\,\tilde{D}_{q}^{1,1}(\mathbb{r})$ to
the pairing EDF from the time-odd contribution
$\propto\tilde{C}_{q}^{1,1}(\mathbb{r})\,\tilde{C}_{q}^{1,1}(\mathbb{r})$.

##  Appendix C Explanation of the supplementary material

We provide as supplementary material the file
Mass_Table_BSkG2.dat, which contains the calculated ground state properties of
all nuclei with $8\leq Z\leq 110$ lying between the proton and neutron drip
lines. Its contents follows the same conventions as the supplementary material
of Ref. Scamps21 , but we repeat the contents of the columns in Table. 5 for
convenience. For clarity we mention that (i) we defined $M(N,Z)$ in Eq. (2) as
the atomic mass, while columns 3 and 4 of the file list mass _excesses_ and
(ii) we report only the largest value among the three Belyaev MOIs
$\mathcal{I}^{\rm B}_{x/y/z}$ in column 16.

Column | Quantity | Units | Explanation
---|---|---|---
1 | Z | $-$ | Proton number
2 | N | $-$ | Neutron number
3 | $M_{\rm exp}$ | MeV | Experimental atomic mass excess
4 | $M_{\rm th}$ | MeV | BSkG1 atomic mass excess
5 | $\Delta M$ | MeV | $M_{\rm exp}-M_{\rm th}$
6 | $E_{\rm tot}$ | MeV | Total energy, Eq. (1)
7 | $\beta_{20}$ | $-$ | Quadrupole deformation
8 | $\beta_{22}$ | $-$
9 | $\beta$ | $-$
10 | $E_{\rm rot}$ | MeV | Rotational correction
11 | $\langle\Delta\rangle_{n}$ | MeV | Average neutron gap
12 | $\langle\Delta\rangle_{p}$ | MeV | Average proton gap
13 | $r_{\rm BSkG2}$ | fm | Calculated rms charge radius
14 | $r_{\rm exp}$ | fm | Experimental rms charge radius
15 | $\Delta r$ | fm | $r_{\rm exp}-r_{\rm BSkG2}$
16 | $\mathcal{I}^{B}$ |  $\hbar^{2}$ MeV-1 | Calculated Belyaev MOI.
17 | par(p) | $-$ | Parity of proton qp. excitation
18 | par(n) | $-$ | Parity of neutron qp. excitation
Table 5: Contents of the Mass_Table_BSkG2.dat file.

## References

  * (1) M. Arnould and S. Goriely, _Astronuclear Physics: A tale of the atomic nuclei in the skies_ , Prog. Part. Nucl. Phys. 112, 103766 (2020).
  * (2) H. Hergert, _A Guided Tour of ab initio Nuclear Many-Body Theory_ , Frontiers in Physics 8, 379 (2020).
  * (3) M. Bender, P.-H. Heenen, and P.-G. Reinhard, _Self-consistent mean-field models for nuclear structure_ , Rev. Mod. Phys. 75, 121 (2003).
  * (4) S. Goriely, N. Chamel, and J. M. Pearson, _Further explorations of Skyrme-Hartree-Fock-Bogoliubov mass formulas. XVI. Inclusion of self-energy effects in pairing._ , Phys. Rev. C 93, 034337 (2016).
  * (5) M. Dutra _et al._ , _Skyrme interaction and nuclear matter constraints_ , Phy. Rev. C 85, 035201 (2012).
  * (6) M. Dutra _et al._ , _Relativistic mean-field hadronic models under nuclear matter constraints_ , Phys. Rev. C 90, 055203 (2014).
  * (7) S. Goriely, S. Hilaire, M. Girod, and S. Péru, _First Gogny-Hartree-Fock-Bogoliubov Nuclear Mass Model_ , Phys. Rev. Lett. 102, 242501 (2009).
  * (8) D. Peña-Arteaga, S. Goriely, and N. Chamel, _Relativistic mean-field mass models_ , Eur. Phys. J. A 52, 320 (2016).
  * (9) S. Goriely, N. Chamel, and J. M. Pearson, _Further explorations of Skyrme-Hartree-Fock-Bogoliubov mass formulas. XIII. The 2012 atomic mass evaluation and the symmetry coefficient_ , Phys. Rev. C 88, 061302(R) (2013).
  * (10) P. Möller, A. J. Sierk, T. Ichikawa, and H. Sagawa, _Nuclear ground-state masses and deformations: FRDM(2012)_ , At. Data Nucl. Data Tables, 109, 1-204 (2016).
  * (11) G. Scamps, S. Goriely, E. Olsen, M. Bender, and W. Ryssens, _Skyrme-Hartree-Fock-Bogoliubov mass models on a 3D mesh: Effect of triaxial shape_ , Eur. Phys. J. A 57, 333 (2021).
  * (12) W. Ryssens, _Symmetry breaking in nuclear mean-field models_ , Ph.D. thesis, Université Libre de Bruxelles, (2016).
  * (13) W. Ryssens, P.-H. Heenen, and M. Bender, _Numerical accuracy of mean-field calculations in coordinate space_ , Phys. Rev. C 92, 064318 (2015).
  * (14) L. M. Robledo and G. F. Bertsch, _Global systematics of octupole excitations in even-even nuclei_ , Phys. Rev. C 84, 054302 (2011).
  * (15) M. Chen, T. Li, J. Dobaczewski and W. Nazarewicz, _Microscopic origin of reflection-asymmetric nuclear shapes_ , Phys. Rev. C 103, 034303 (2021).
  * (16) T. Duguet, P. Bonche, P.-H. Heenen, and J. Meyer, _Pairing correlations. I. Description of odd nuclei in mean-field theories_ , Phys. Rev. C 65, 014310 (2001).
  * (17) V. Hellemans, P.-H. Heenen, and M. Bender, _Tensor part of the Skyrme energy density functional. III. Time-odd terms at high spin_ , Phys. Rev. C 85, 014326 (2012).
  * (18) J. A. Maruhn, P.-G. Reinhard, P. D. Stevenson, and M. R. Strayer, _Spin-excitation mechanisms in Skyrme-force time-dependent Hartree-Fock calculations_ , Phys. Rev. C 74, 027601 (2006).
  * (19) M. Bender, J. Dobaczewski, J. Engel, and W. Nazarewicz, _Gamow-Teller strength and the spin-isospin coupling constants of the Skyrme energy functional_ , Phys. Rev. C 65, 054322 (2002).
  * (20) M. Kortelainen, N. Hinohara, and W. Nazarewicz, _Multipole modes in deformed nuclei within the finite amplitude method_ , Phys. Rev. C 92, 051302 (2015).
  * (21) P. Bonche, J. Dobaczewski, H. Flocard, P.-H. Heenen, and J. Meyer, _Analysis of the generator coordinate method in a study of shape isomerism in 194Hg_, Nucl. Phys. A 510, 466 (1990).
  * (22) N. Hinohara, T. Nakatsukasa, M. Matsuo, and K. Matsuyanagi, _Effects of time-odd components in mean field on large amplitude collective dynamics_ , Prog. Theoret. Phys. 115, 567, (2006)
  * (23) P. Klüpfel, _Skyrme’s Interaction Beyond the Mean-Field: The DGCM+GOA Hamiltonian of Nuclear Quadrupole Motion_ , Thesis, Friedrich-Alexander-Universität Erlangen-Nürnberg (2008).
  * (24) K. Petrík and M. Kortelainen, _Thouless-Valatin rotational moment of inertia from linear response theory_ , Phys. Rev. C 97, 034321 (2018).
  * (25) K. Washiyama, N. Hinohara, and T. Nakatsukasa, _Finite-amplitude method for collective inertia in spontaneous fission_ , Phys. Rev. C 103, 014306 (2021).
  * (26) P. Bonche, H. Flocard, and P.-H. Heenen, _Self-consistent calculation of nuclear rotations: The complete yrast line of 24Mg_, Nucl. Phys. A 467, 115-135 (1987).
  * (27) J. Dobaczewski and J. Dudek, _Time-odd rotating superdeformed nuclei_ , Phys. Rev. C 52, 1827 (1995).
  * (28) L. M. Robledo, R. N. Bernard, and G. F. Bertsch, _Spin constraints on nuclear energy density functionals_ , Phys. Rev. C 89, 021303 (2014).
  * (29) P. L. Sassarini, J. Dobaczewski, J. Bonnard, and R. F. G. Ruiz, _Global analysis of electromagnetic moments in odd near doubly magic nuclei_ , arxiv.org/abs/2111.04675.
  * (30) D. Peña Arteaga, M. Grasso, E. Khan, and P. Ring, _Nuclear structure in strong magnetic fields: Nuclei in the crust of a magnetar_ , Phys. Rev. C 84, 045806 (2011).
  * (31) J. Margueron, J. Navarro, and N. Van Giai, _Instabilities of infinite matter with effective Skyrme-type interactions_ , Phys. Rev. C 66, 014303 (2002).
  * (32) N. Chamel, S. Goriely and J. M. Pearson, _Further explorations of Skyrme-Hartree-Fock-Bogoliubov mass formulas. XI. Stabilizing neutron stars against a ferromagnetic collapse_ , Phys. Rev. C 80, 065804 (2009).
  * (33) N. Chamel and S. Goriely, _Spin and spin-isospin instabilities in asymmetric nuclear matter at zero and finite temperatures using Skyrme functionals_ , Phys. Rev. C 82, 045804 (2010).
  * (34) L.-G. Cao, G. Colò, and H. Sagawa, _Spin and spin-isospin instabilities and Landau parameters of Skyrme interactions with tensor correlations_ , Phys. Rev. C 81, 044302 (2010).
  * (35) A. Pastore, D. Davesne, K. Bennaceur, J. Meyer, and V. Hellemans, _Fitting Skyrme functionals using linear response theory_ , Phys. Scr. T154, 014014 (2013).
  * (36) A. Pastore, D. Tarpanov, D. Davesne, and J. Navarro, _Spurious finite-size instabilities in nuclear energy density functionals: Spin channel_ , Phys. Rev. C 92, 024305 (2015).
  * (37) K. J. Pototzky, J. Erler, P.-G. Reinhard, and V. O. Nesterenko, _Properties of odd nuclei and the impact of time-odd mean fields: A systematic Skyrme-Hartree-Fock analysis_ , Eur. Phys. J. A 46, 299-313 (2010).
  * (38) A. V. Afanasjev, and H. Abusara, _Time-odd mean fields in covariant density functional theory: Nonrotating systems_ , Phys. Rev. C 81, 014309 (2010).
  * (39) S. Perez-Martin and L. M. Robledo, _Microscopic justification of the equal filling approximation_ , Phys. Rev. C 78, 014304 (2008).
  * (40) N. Schunck _et al._ , _One-quasiparticle states in the nuclear energy density functional theory_ , Phys. Rev. C 81, 024316 (2010).
  * (41) S. Goriely, M. Samyn, and J. M. Pearson, _Further explorations of Skyrme-Hartree-Fock-Bogoliubov mass formulas. VII. Simultaneous fits to masses and fission barriers_ , Phys. Rev. C 75, 064312 (2007).
  * (42) W. Ryssens, G. Scamps, S. Goriely, and M. Bender, _Skyrme-Hartree-Fock-Bogoliubov mass models on a 3D mesh: fission barriers with BSkG2_ , in preparation.
  * (43) D. Lunney, J. M. Pearson, and C. Thibault, _Recent trends in the determination of nuclear masses._ , Rev. Mod. Phys. 75, 1021-1082 (2003).
  * (44) M. Wang et al., _The AME 2020 atomic mass evaluation (II). Tables, graphs and references_ , Chin. Phys. C 45, 030003 (2021).
  * (45) S. J. Krieger, P. Bonche, H. Flocard, P. Quentin, and M. S. Weiss, _An improved pairing interaction for mean field calculations using Skyrme potentials_ , Nucl. Phys. A 517, 275 (1990).
  * (46) S. Goriely and J. M. Pearson, _Further explorations of Skyrme-Hartree-Fock-Bogoliubov mass formulas. VIII. Role of Coulomb exchange_ , Phys. Rev. C 77, 031301 (2008).
  * (47) M. Bender, K. Rutz, P.-G. Reinhard, and J. A. Maruhn, _Consequences of the center-of-mass correction in nuclear mean-field models_ , Eur. Phys. J. A 7, 467-478 (2000).
  * (48) F. Tondeur, S. Goriely, J. M. Pearson, and M. Onsi, _Towards a Hartree-Fock mass formula_ , Phys. Rev. C 62, 024308 (2000).
  * (49) S. Goriely, N. Chamel, and J. M. Pearson, _Further explorations of Skyrme-Hartree-Fock-Bogoliubov mass formulas. XIII. The 2012 atomic mass evaluation and the symmetry coefficient_ , Phys. Rev. C 88, 024308 (2013).
  * (50) S. Goriely, M. Samyn, M. Bender, and J. M. Pearson, _Further explorations of Skyrme-Hartree-Fock-Bogoliubov mass formulas. II. Role of the effective mass_ , Phys. Rev. C 68, 054325 (2003).
  * (51) W. Ryssens, M. Bender and P.-H. Heenen, _Iterative approaches to the self-consistent nuclear energy density functional problem: Heavy ball dynamics and potential preconditioning_ , Eur. Phys. J. A 55, 93 (2019).
  * (52) M. M. Sharma, G. Lalazissis, J. König, and P. Ring, _Isospin Dependence of the Spin-Orbit Force and Effective Nuclear Potentials_ , Phys. Rev. Lett. 74, 3744-3747 (1995).
  * (53) P.-G. Reinhard and H. Flocard, _Nuclear effective forces and isotope shifts_ , Nucl. Phys. A 584, 467-488 (1995).
  * (54) T. Lesinski, M. Bender, K. Bennaceur, T. Duguet, and J. Meyer, _Tensor part of the Skyrme energy density functional: Spherical nuclei_ , Phys. Rev. C 76, 014312 (2007).
  * (55) P. Klüpfel, P.-G. Reinhard, T. J. Bürvenich, and J. A. Maruhn, _Variations on a theme by Skyrme: A systematic study of adjustments of model parameters_ , Phys. Rev. C 79, 034310 (2009).
  * (56) M. Kortelainen _et al._ , _Nuclear energy density optimization_ , Phys. Rev. C 82, 024313 (2010).
  * (57) M. Kortelainen _et al._ , _Nuclear energy density optimization: Shell structure_ , Phys. Rev. C 89, 054314 (2014).
  * (58) S. Goriely, _Further explorations of Skyrme–Hartree–Fock–Bogoliubov mass formulas. XV: The spin-orbit coupling_ , Nucl. Phys. A 933, 68-81 (2015).
  * (59) D. Davesne, A. Pastore, and J. Navarro, _Linear response of homogeneous nuclear matter with energy density functionals_ , Phys. Rep. 563, 1 (2015).
  * (60) V. Hellemans _et al._ , _Spurious finite-size instabilities in nuclear energy density functionals_ , Phys. Rev. C 88, 064323 (2013).
  * (61) D. Davesne, A. Pastore, and J. Navarro, _Linear response theory with finite-range interactions_ , Prog. Part. Nucl. Phys. 120, 103870 (2021).
  * (62) M. Martini, A. De Pace, and K. Bennaceur, _Spurious finite-size instabilities with Gogny-type interactions_ , Eur. Phys. J. A 55, 150 (2019).
  * (63) C. Gonzalez-Boquera, M. Centelles, X. Viñas, and L. M. Robledo, _Finite-size instabilities in finite-range forces_ , Phys. Rev. C 103, 064314 (2021).
  * (64) E. Perlińska, S. G. Rohoziński, J. Dobaczewski, and W. Nazarewicz, _Local density approximation for proton-neutron pairing correlations: Formalism_ , Phys. Rev. C 69, 014316 (2004).
  * (65) S. T. Belyaev, _Concerning the calculation of the nuclear moment of inertia_ , Nucl. Phys. 24, 322-325 (1961).
  * (66) D. Baye and P.-H. Heenen, _Generalised meshes for quantum mechanical problems_ , J. Phys. A. Math. Gen. 19, 2041 (1986).
  * (67) J. Dobaczewski, J. Dudek, S. G. Rohoziński, and T. R. Werner, _Point symmetries in the Hartree-Fock approach. I. Densities, shapes, and currents_ , Phys. Rev. C 62, 014310 (2000).
  * (68) E. P. Wigner, _Normal form of antiunitary operators._ , J. Math. Phys. 1, 409-413 (1960).
  * (69) B. Banerjee, P. Ring, and H. J. Mang, _On the character of the Hartree-Fock-Bogoliubov solutions in a rotating frame_ , Nucl. Phys. A 221, 564-572 (1974).
  * (70) P. Ring and P. Schuck, The Nuclear Many-Body Problem, (Springer Verlag, 1980).
  * (71) G. Bertsch, J. Dobaczewski, W. Nazarewicz, and J. Pei, _Hartree-Fock-Bogoliubov theory of polarized Fermi systems_ , Phys. Rev. A 79, 043602 (2009).
  * (72) M. Bender, N. Schunck, J.-P. Ebran, and T. Duguet, Chapter 3: Single-Reference and Multi-Reference Formulation, In _Energy Density Functional Methods for Atomic Nuclei_ , Nicolas Schunck [ed.], IoP Publishing Ltd (2019).
  * (73) P.-H. Heenen, P. Bonche, and H. Flocard, _Microscopic study of superdeformation in 193Hg_, Nucl. Phys. A 588, 490-500 (1995).
  * (74) J. Dobaczewski _et al._ , _Solution of the Skyrme-Hartree-Fock-Bogolyubov equations in the Cartesian deformed harmonic-oscillator basis. (VI) hfodd (v2.40h): A new version of the program_ , Comp. Phys. Comm. 180, 2361-2391 (2009).
  * (75) R. N. Perez, N. Schunck, R.-D. Lasseri, C. Zhang, and J. Sarich, _Axially deformed solution of the Skyrme-Hartree-Fock-Bogolyubov equations using the transformed harmonic oscillator basis (III) hfbtho (v3.00): A new version of the program_ , Comp. Phys. Comm. 220, 363-375 (2017).
  * (76) N. D. Woods, M. C. Payne, and P. J. Hasnip, _Computing the self-consistent field in Kohn–Sham density functional theory_ , J. Phys. Condens. Matter 31, 453001 (2019).
  * (77) J. L. Egido, J. Lessing, V. Martin, and L. M. Robledo, _On the solution of the Hartree-Fock-Bogoliubov equations by the conjugate gradient method_ , Nucl. Phys. A 594, 70-86 (1995).
  * (78) L. M. Robledo and G. F. Bertsch, _Application of the gradient method to Hartree-Fock-Bogoliubov theory_ , Phys. Rev. C 84, 014312 (2011).
  * (79) B. Gall, P. Bonche, J. Dobaczewski, H. Flocard, and P.-H. Heenen, _Superdeformed rotational bands in the mercury region. A cranked Skyrme-Hartree-Fock-Bogoliubov study_ , Z. Phys. A 348, 183-197 (1994).
  * (80) W. Ryssens, M. Bender, K. Bennaceur, P.-H. Heenen, and J. Meyer, _Impact of the surface energy coefficient on the deformation properties of atomic nuclei as predicted by Skyrme energy density functionals_ , Phys. Rev. C 99, 044315 (2019).
  * (81) C. J. Gallagher and S. A. Moszkowski, _Coupling of angular momenta in odd-odd nuclei_ , Phys. Rev. 111, 1282-1290 (1958).
  * (82) J. P. Boisson, R. Piepenbring, and W. Ogle, _The effective neutron-proton interaction in rare-earth nuclei_ , Phys. Rep. 26, 99-147 (1978).
  * (83) D. J. Thouless and J. G. Valatin, _Time-dependent Hartree-Fock equations and rotational states of nuclei_ , Nucl. Phys. 31, 211-230 (1962).
  * (84) N. Hinohara, _Collective inertia of the Nambu-Goldstone mode from linear response theory_ , Phys. Rev. C 92, 034321 (2015).
  * (85) J. Libert, M. Girod, and J. P. Delaroche, _Microscopic descriptions of superdeformed bands with the Gogny force: Configuration mixing calculations in the A $\sim$ 190 mass region_, Phys. Rev. C 60, 054301 (1999).
  * (86) A. Bohr and B. R. Mottelson, _Single-Particle and Collective Aspects of Nuclear Rotation_ , Phys. Scr. 24, 71-76 (1981).
  * (87) J. Dudek, _Nuclear superdeformation at high spins_ , Prog. Part. Nucl. Phys. 28, 131-185 (1992).
  * (88) B. Singh, R. Zywina and R. B. Firestone, _Table of Superdeformed Nuclear bands and Fission Isomers_ , Nucl. Data Sheets 97, 241-592 (2002).
  * (89) S. Goriely, N. Samyn, and J. M. Pearson, _Further explorations of Skyrme-Hartree-Fock-Bogoliubov mass formulas; VI: Weakened pairing_ , Nucl. Phys. A 773, 279 (2006).
  * (90) M. Bender, K. Rutz, P.-G. Reinhard, and J. A. Maruhn, _Pairing gaps from nuclear mean-field models_ , Eur. Phys. J. A, 8, 59-75 (2000).
  * (91) T. Duguet, P. Bonche, P.-H. Heenen, and J. Meyer, _Pairing correlations. II. Microscopic analysis of odd-even mass staggering in nuclei_ , Phys. Rev. C 65, 014311 (2001).
  * (92) R. Capote _et al_., _RIPL - Reference Input Parameter Library for Calculation of Nuclear Reactions and Nuclear Data Evaluation_ , Nucl. Data Sheets 110, 3107-3214 (2009).
  * (93) E. Chabanat, P. Bonche, P. Haensel, J. Meyer and R. Schaeffer, _A Skyrme parametrization from subnuclear to neutron star densities._ , Nucl. Phys. A 627, 710 (1997).
  * (94) G. Colò, N. V. Giai, J. Meyer, K. Bennaceur, and P. Bonche, _Microscopic determination of the nuclear incompressibility within the nonrelativistic framework_ , Phys. Rev. C 70, 024307 (2004).
  * (95) L. G. Cao, U. Lombardo, and P. Schuck, _Screening effects in superfluid nuclear and neutron matter within Brueckner theory_ , Phys. Rev. C 74, 064301 (2006).
  * (96) W. Zuo, A. Lejeune, U. Lombardo, and J.-F. Mathiot, _Interplay of three-body interactions in the EOS of nuclear matter_ , Nucl. Phys. A 706, 418 (2002).
  * (97) I. Angeli and K. P. Marinova, _Table of experimental nuclear ground state charge radii: An update_ , At. Data Nucl. Data Tables 99, 69 (2013).
  * (98) R. Jodon, M. Bender, K. Bennaceur, and J. Meyer, _Constraining the surface properties of effective Skyrme interactions_ , Phys. Rev. C 94, 024335 (2016).
  * (99) Ph. Da Costa, M. Bender, K. Bennaceur, J. Meyer, and W. Ryssens, _in preparation_.
  * (100) A. Mamdouh, J. M. Pearson, M. Rayet, and F. Tondeur, _Fission barriers of neutron-rich and superheavy nuclei calculated with the ETFSI method_ , Nuc. Phys. A 679, 337-358 (2001).
  * (101) S. A. Giuliani, G. Martínez-Pinedo, and L. M. Robledo, _Fission properties of superheavy nuclei for r-process calculations_ , Phys. Rev. C 97, 034323 (2018).
  * (102) E. Chabanat, P. Bonche, P. Haensel, J. Meyer, and R. Schaeffer, _A Skyrme parametrization from subnuclear to neutron star densities Part II. Nuclei far from stabilities_ , Nucl. Phys. A635, 231 (1998).
  * (103) W. Satuła, _Wigner energy, odd-even mass staggering and the time-odd mean-fields_ , AIP Conference Proceedings 141, 141-150 (1999).
  * (104) D. Vretenar, A. V. Afanasjev, G. A. Lalazissis, and P. Ring, _Relativistic Hartree-Bogoliubov theory: Static and dynamic aspects of exotic nuclear structure_ , Phys. Rep. 409, 101-259 (2005).
  * (105) A. Sulaksono, P.-G. Reinhard, T. J. Bürvenich, P. O. Hess, and J. A. Maruhn, _From self-consistent covariant Effective Field Theories to their Galilean-invariant counterparts_ , Phys. Rev. Lett. 98, 262501 (2007).
  * (106) N. Paar, T. Nikšić, D. Vretenar, and P. Ring, _Quasiparticle random phase approximation based on the relativistic Hartree-Bogoliubov model. II. Nuclear spin and isospin excitations_ , Phys. Rev. C 69, 054303 (2004).
  * (107) K. Rutz, M. Bender, P.-G. Reinhard, and J. A. Maruhn, _Pairing gap and polarisation effects_ , Phys. Lett. B 468, 1-6 (1999).
  * (108) L. M. Robledo, R. Bernard, and G. F. Bertsch, _Pairing gaps in the Hartree-Fock-Bogoliubov theory with the Gogny D1S interaction_ , Phys. Rev. C 86, 064313 (2012).
  * (109) A. S. Jensen, P. G. Hansen, and B. Jonson, _New mass relations and two- and four-nucleon correlations_ , Nucl. Phys. A 431, 393-418 (1984).
  * (110) W. A. Friedman and G. F. Bertsch, _Neutron-proton pairing reexamined_ , Phys. Rev. C 76, 057301 (2007).
  * (111) Z. Wu, S. A. Changizi, and C. Qi, _Empirical residual neutron-proton interaction in odd-odd nuclei_ , Phys. Rev. C 93, 034334 (2016).
  * (112) J. Y. Zeng, T. H. Jin, and Z. J. Zhao, _Blocking effect and odd-even differences in the moments of inertia of rare-earth nuclei_ , Phys. Rev. C 50, 1388 (1994).
  * (113) A. V. Afanasjev, J. König, P. Ring, L. M. Robledo, and J. L. Egido, _Moments of inertia of nuclei in the rare earth region: A relativistic versus nonrelativistic investigation_ , Phys. Rev. C 62, 054306 (2000).
  * (114) J. M. Pearson, Y. Aboussir, A. K. Dutta, R. C. Nayak, M. Farine, and F. Tondeur, _Thomas-Fermi approach to nuclear mass formula_ , Nucl. Phys. A 528, 1 (1991).
  * (115) National Nuclear Data Center, Nuclear structure and decay data on-line library, Nudat 2.8, 2018. https://www.nndc.bnl.gov/nudat2/.
  * (116) M. Rocchini et al., _Onset of triaxial deformation in 66Zn and properties of its first excited $0^{+}$ state studied by means of Coulomb excition_, Phys. Rev. C 103, 014311 (2021).
  * (117) M. Sugawara et al., _Multiple Coulomb excitation of a 70Ge beam and the interpretation of the $0_{2}^{+}$ state as a deformed intruder_, Eur. Phys. J. A 16, 409 (2003).
  * (118) A. D. Ayangeakaa et al., _Shape coexistence and the role of axial asymmetry in 72Ge_, Phys. Lett. B 754, 254 (2016).
  * (119) Y. Toh et al., _Coulomb excitation of 74Ge beam_, Eur. Phys. J. A 9, 353 (2000).
  * (120) A. D. Ayangeakaa et al., _Evidence for rigid triaxial deformation in 76Ge from a model-independent analysis_, Phys. Rev. Lett. 123, 102501 (2019).
  * (121) E. Clément et al., _Shape coexistence in neutron-deficient krypton isotopes_ , Phys. Rev. C 75, 054313 (2007).
  * (122) A. E. Kavka et al., _Coulomb excitation of 76,80,82Se_, Nucl. Phys. A 593, 177 (1995).
  * (123) E. Clément et al., _Spectroscopic quadrupole moments in 96-98Sr: evidence for shape coexistence in neutron-rich Strontium isotopes at $N=60$_, Phys. Rev. C 94, 054326 (2016).
  * (124) M. Zielińska, _Electromagnetic structure of molybdenum isotopes studied using Coulomb excitation method_ , Ph.D. thesis, Warsaw University, (2005).
  * (125) M. Zielińska et al., _Electromagnetic structure of 98Mo_, Nucl. Phys. A 712, 3 (2002).
  * (126) K. Wrzosek-Lipska et al., _Electromagnetic properties of 100Mo: Experimental results and theoretical description of quadrupole degrees of freedom_, Phys. Rev. C 86, 064305 (2012).
  * (127) J. Srebrny et al., _Experimental and theoretical investigations of quadrupole collective degrees of freedom in 104Ru_, Nucl. Phys. A 766, 25 (2006).
  * (128) K. Wrzosek-Lipska et al., _Quadrupole Deformation of 110Cd studied with coulomb excitation_, Acta Phys. Pol. B 51, 789 (2020).
  * (129) C. Fahlander et al., _Quadrupole collective properties of 114Cd_, Nucl. Phys. A 485, 327 (1988).
  * (130) L. E. Svensson et al., _Multiphonon vibrational states in 106,108Pd_, Nucl. Phys. A 584, 547 (1995).
  * (131) L. Morrison et al., _Quadrupole deformation of 130Xe measured in a Coulomb-excitation experiment_, Phys. Rev. C 102, 054304 (2020).
  * (132) C. Y. Wu et al., _Quadrupole collectivity and shapes of Os-Pt nuclei_ , Nucl. Phys. A 607, 178 (1996).
  * (133) M. Zielinska, private communication.
  * (134) S. Raman, C. W. Nestor, and P. Tikkanen, _Transition probability from the ground to the first-excited $2+$ state of even-even nuclides_, At. Data Nucl. Data Tables, 78, 1 (2001).
  * (135) K. Kumar, _Intrinsic quadrupole moments and shapes of nuclear ground states and excited states_ , Phys. Rev. Lett. 28, 249 (1972).
  * (136) D. Cline, _Nuclear shapes studied by Coulomb excitation_ , Ann. Rev. Nucl. Part. Sci. 36, 683 (1986).
  * (137) S. Sels _et al._ , _Shape staggering of midshell mercury isotopes from in-source laser spectroscopy compared with density-functional-theory and Monte Carlo shell-model calculations_ , Phys. Rev. C 99, 044306 (2019).
  * (138) A. Barzakh, _et al._ , _Large Shape Staggering in Neutron-Deficient Bi Isotopes_ , Phys. Rev. Lett. 127, 192501 (2021).
  * (139) N. Chamel, _Self-interaction errors in nuclear energy density functionals_ , Phys. Rev. C 82, 061307 (2010).
  * (140) S. Goriely, N. Chamel, and J. M. Pearson, _Further explorations of Skyrme-Hartree-Fock-Bogoliubov mass formulas. XII. Stiffness and stability of neutron-star matter_ , Phys. Rev. C 82, 035804 (2010).
  * (141) B. Friedman and V. R. Pandharipande, _Hot and cold, nuclear and neutron matter_ , Nucl. Phys. A 361, 502 (1981).
  * (142) R. B. Wiringa, V. Fiks, and A. Fabrocini, _Equation of state for dense nucleon matter_ , Phys. Rev. C 38, 1010 (1988).
  * (143) A. Akmal, V. R. Pandharipande, and D. G. Ravenhall, _Equation of state of nucleon matter and neutron star structure_ , Phys. Rev. C 58, 1804 (1998).
  * (144) Z. H. Li and H.-J. Schulze, _Neutron star structure with modern nucleonic three-body forces_ , Phys. Rev. C 78, 028801 (2008).
  * (145) J. E. Lynn _et al._ , _Chiral three-nucleon interactions in light nuclei, neutron- $\alpha$ scattering, and neutron matter_, Phys. Rev. Lett. 116, 062501 (2016).
  * (146) C. Drischler, K. Hebeler, and A. Schwenk, _Chiral interactions up to next-to-next-to-next-to-leading order and nuclear saturation_ , Phys. Rev. Lett. 122, 042501 (2019).
  * (147) E. Fonseca _et al._ , _Refined Mass and Geometric Measurements of the High-mass PSR J0740+6620_ , Astr. J. Lett. 915, L12 (2021).
  * (148) J. Erler, P. Klüpfel, and P.-G. Reinhard, _A stabilized pairing functional_ , Eur. Phys. J. A 37, 81-86 (2008).
  * (149) C. Rigollet, P. Bonche, H. Flocard, and P.-H. Heenen, _Microscopic study of the properties of identical bands in the $A\sim 150$ mass region_, Phys. Rev. C 59, 3120-3127 (1999).
  * (150) J. Eberz _et al._ , _Spins, moments and mean square charge radii of 104-127In determined by laser spectroscopy_, Nucl. Phys. A 464, 9-28 (1987).
  * (151) A. R. Vernon _et al._ , _Nuclear moments of indium isotopes reveal abrupt change at magic number 82_ , Nature 607, 260-265 (2022).
  * (152) A. Bohr and B. Mottelson, _Nuclear structure, Vol. I: Single-Particle motion_ , (Benjamin, New York, 1969).
  * (153) A. Bohr and B. Mottelson, _Nuclear structure, Vol II: Nuclear Deformations_ , (Benjamin, New York, 1975).
  * (154) L. Bonneau, N. Minkov, D. D. Duc, P. Quentin and J. Bartel, _Effect of core polarization on magnetic dipole moments in deformed odd-mass nuclei_ , Phys. Rev. C 91, 054307 (2015).
  * (155) S. Péru, S. Hilaire, S. Goriely and M. Martini, _Description of magnetic moments within the Gogny Hartree-Fock-Bogolyubov framework: Application to Hg isotopes_ , Phys. Rev. C 104, 024328 (2021).
  * (156) B. Bally, G. Giacalone, and M. Bender, _Structure of 128,129,130Xe through multi-reference energy density functional calculations_, http://arxiv.org/abs/2207.13576.
  * (157) T. Schmidt, _Über die magnetischen Momente der Atomkerne_ , Z. Phys. 106, 358-361 (1934).
  * (158) T. Ericson and W. Weise, _Pions and Nuclei_ , Clarendon Press (Oxford, 1988).
  * (159) A. Roggero, A. Mukherjee, and F. Pederiva, _Constraining the Skyrme energy density functional with quantum Monte Carlo calculations_ , Phys. Rev. C 92, 054303 (2015).
  * (160) T. Duguet, M. Bender, K. Bennaceur, D. Lacroix, and T. Lesinski, _Particle-number restoration within the energy density functional formalism: Non-vialibility of terms depending on non-integer powers of the density matrices_ , Phys. Rev. C 79, 044320 (2009).
  * (161) W. Ryssens and M. Bender, _Skyrme pseudopotentials at next-to-next-to-leading order: Construction of local densities and first symmetry-breaking calculations_ , Phys. Rev. C 104, 044308 (2021).