# Nonpropagating ghost in covariant $f(Q)$ gravity

Kun Hu<EMAIL_ADDRESS>Institute of Astrophysics, Central China
Normal University, Wuhan 430079, China  Makishi Yamakoshi
<EMAIL_ADDRESS>Department of Physics, Nagoya
University, Nagoya 464-8602, Japan  Taishi Katsuragawa
<EMAIL_ADDRESS>Institute of Astrophysics, Central China Normal
University, Wuhan 430079, China  Shin’ichi Nojiri
<EMAIL_ADDRESS>Department of Physics, Nagoya University,
Nagoya 464-8602, Japan  Kobayashi-Maskawa Institute for the Origin of
Particles and the Universe, Nagoya University, Nagoya 464-8602, Japan  Taotao
Qiu<EMAIL_ADDRESS>School of Physics, Huazhong University of Science and
Technology, Wuhan, 430074, China

###### Abstract

$f(Q)$ gravity is an extension of the symmetric teleparallel equivalent to
general relativity (STEGR). This work shows that based on the scalar-
nonmetricity formulation, a scalar mode in $f(Q)$ gravity has a negative
kinetic energy. This conclusion holds regardless of the coincident gauge
frequently used in STEGR and $f(Q)$ gravity. To study the scalar mode, we
further consider the covariant $f(Q)$ gravity as a special class in higher-
order scalar tensor (HOST) theory and rewrite the four scalar fields, which
play a role of the Stüeckelberg fields associated with the diffeomorphism, by
vector fields. Applying the standard Arnowitt-Deser-Misner (ADM) formulation
to the new formulation of the $f(Q)$ gravity, we demonstrate that the ghost
scalar mode can be eliminated by the second-class constraints, thus ensuring
that $f(Q)$ gravity is a healthy theory.

##  I Introduction

In recent years, gravitational theories with nonmetricity have been actively
discussed Nester and Yo (1999); Beltrán Jiménez  _et al._ (2018a, 2019);
Rünkla and Vilson (2018); Beltrán Jiménez  _et al._ (2020). These theories are
written by the nonmetricity scalar $Q$ as the fundamental geometrical quantity
and in Palatini formalism, where the connection is an independent variable in
addition to the metric. Imposing that the Riemann tensor and torsion tensor
vanish for the connection, one finds that the connection is written by the
four scalar fields Blixt  _et al._ (2023); Beltrán Jiménez and Koivisto
(2022); Adak (2018); Tomonari and Bahamonde (2023). Moreover, one can choose
specific scalar fields so that the connection vanishes, which is called the
coincident gauge, and then the metric is only a dynamical variable in theory.
It is known that in the coincident gauge, the theory with action linear to $Q$
is equivalent to general relativity (GR), and we call it symmetric
teleparallel equivalent to general relativity (STEGR) Beltrán Jiménez  _et
al._ (2018b); Quiros (2022).

Teleparallelism is another well-known concept in the gravitational theory that
does not rely on curvature, where the torsion scalar $T$ is the fundamental
quantity Maluf (2013); Aldrovandi and Pereira (2013); Cai  _et al._ (2016). By
choosing the Weitzenböck connection, which corresponds to the coincident gauge
in STEGR, the spin connection vanishes, and the tetrad is only a dynamical
variable. The theory with action linear to $T$ is equivalent to GR, and we
call it teleparallel equivalent to general relativity (TEGR) Maluf (2013);
Krššák and Saridakis (2016). Compared with teleparallelism, we utilize
symmetric teleparallelism in nonmetricity gravity to reflect the torsionless
condition, that is, symmetric affine connection. Moreover, as in the way of
extending GR to $f(R)$ gravity characterized by the function of the curvature
scalar, TEGR, and STEGR have been extended to $f(T)$ and $f(Q)$ gravity whose
actions include the arbitrary function of $T$ and $Q$ respectively Bahamonde
_et al._ (2015); Cai  _et al._ (2016); Krššák and Saridakis (2016); Bahamonde
_et al._ (2023); Heisenberg (2019); Beltrán Jiménez  _et al._ (2020); Harko
_et al._ (2018); Järv  _et al._ (2018); Rünkla and Vilson (2018); Capozziello
_et al._ (2022, 2023). Those theories have been intensively examined as one of
the modified gravity theories, and new degrees of freedom (DOF) introduced by
the torsion or nonmetricity show various interesting phenomena; for instance,
cosmological models Hu  _et al._ (2023); Qiu  _et al._ (2019); Nájera and
Fajardo (2022); Li and Zhao (2022); Khyllep  _et al._ (2021); Casalino  _et
al._ (2021); Sharma and Sur (2021, 2022, 2023); Mandal  _et al._ (2020a, b);
Mandal and Sahoo (2021); Arora and Sahoo (2022); Gadbail  _et al._ (2022);
Paliathanasis (2023a, b); Dimakis  _et al._ (2022); Paliathanasis (2023c),
black hole solutions Wang  _et al._ (2022); Lin and Zhai (2021); D’Ambrosio
_et al._ (2022); Bahamonde  _et al._ (2022); Calzá and Sebastiani (2023), and
gravitational waves Hohmann  _et al._ (2018, 2019); Soudi  _et al._ (2019);
Abedi and Capozziello (2018).

Recently, the number of DOF in $f(Q)$ gravity has come into the spotlight Hu
_et al._ (2022); D’Ambrosio  _et al._ (2023); Tomonari and Bahamonde (2023);
Heisenberg (2023); Paliathanasis  _et al._ (2023); Dimakis  _et al._ (2021).
In our previous work Hu  _et al._ (2022), we have shown there are eight DOFs
in coincident $f(Q)$ gravity, and the scalar mode cannot propagate 111The
conclusion in Ref. Hu  _et al._ (2022) have been controversial and under
active discussion recently. For example, see Refs. D’Ambrosio  _et al._
(2023); Tomonari and Bahamonde (2023) for details.. However, the existing
works on the Hamiltonian analysis rely on the specific gauge, the so-called
coincident gauge. It is significant to confirm that there are no ghost modes
with negative kinetic energy in the physical DOF without gauge fixing. This
work mainly discusses the ghost scalar mode in the $f(Q)$ gravity. Moreover,
we investigate the ADM formulation of $f(Q)$ gravity in the arbitrary gauge to
consider the origin of the negative kinetic energy.

This paper is organized as follows: We introduce the underlying geometrical
background of the $f(Q)$ gravity in Sec. II. We apply the conformal rescaling
to the action and discuss the ghost scalar mode in the $f(Q)$ gravity by the
scalar-nonmetricity formulation in Sec. III. In Sec. IV, we propose a new
formulation of the $f(Q)$ gravity as the higher-order scalar-tensor theory.
Finally, we conclude this paper and discuss the implications of ghost mode in
Sec. V. Some mathematical details of the paper are given in the Appendices.

Throughout this paper, we use the leading letters of the Latin alphabet
$\left(a,b,c\right)$ running from $0$ to $3$ to label the tangent space-time
coordinates. The Latin indices $\left(i,j,k\right)$ running from $1$ to $3$
represent the ADM spatial indices, and the Greek indices
$\left(\alpha,\beta,\cdots\right)$ running from $0$ to $3$ does the space-time
indices. For clarity of notation, we define symbols as in Table 1.

Table 1: Conventions and notations $\left\\{{}^{\,\alpha}_{\mu\nu}\right\\}$ | Levi-Civita connection
---|---
$\varGamma^{\alpha}_{\mu\nu}$ | General affine Connection
$\nabla_{\mu}$ | Covariant derivative with respect to Levi-Civita connection
${\hat{\nabla}}_{\alpha}$ | Covariant derivative with respect to general affine connection.
$\mathcal{D}_{i}$ | Three dimensional covariant derivative with respect to $h_{ij}$
$R$ | Curvature scalar with respect to general affine connection
$\mathcal{R}$ | Curvature scalar with respect to Levi-Civita connection
$\mathring{Q}$ | Nonmetricity scalar in coincident gauge

##  II Basics of $f(Q)$ Gravity and Coincident Gauge

We briefly review the basic structure of the $f(Q)$ gravity and symmetries. We
discuss the four scalar fields in the affine connection and their role as the
Stüeckelberg fields.

###  II.1 Geometrical foundations

A general affine connection $\varGamma^{\alpha}_{\,\mu\nu}$ can be decomposed
into three parts:

$\displaystyle\varGamma^{\alpha}_{\,\mu\nu}=\left\\{{}^{\,\alpha}_{\mu\nu}\right\\}+K_{\,\mu\nu}^{\alpha}+L_{\,\mu\nu}^{\alpha}\,.$
(1)

Here, $\left\\{{}^{\,\alpha}_{\mu\nu}\right\\}$ is Levi-Civita connection

$\displaystyle\left\\{{}^{\,\alpha}_{\mu\nu}\right\\}=\frac{1}{2}g^{\alpha\lambda}\left(\partial_{\mu}g_{\lambda\nu}+\partial_{\nu}g_{\mu\lambda}-\partial_{\lambda}g_{\mu\nu}\right)\,,$
(2)

$K_{\,\mu\nu}^{\alpha}$ is contortion

$\displaystyle
K^{\alpha}_{\,\mu\nu}=\frac{1}{2}g^{\alpha\lambda}\left(T_{\mu\lambda\nu}+T_{\nu\lambda\mu}+T_{\lambda\mu\nu}\right)\,,$
(3)

and $L_{\,\mu\nu}^{\alpha}$ is disformation

$\displaystyle
L^{\alpha}_{\,\mu\nu}=\frac{1}{2}g^{\alpha\lambda}\left(Q_{\lambda\mu\nu}-Q_{\mu\lambda\nu}-Q_{\nu\lambda\mu}\right)\,.$
(4)

Torsion tensor $T^{\lambda}_{\,\mu\nu}$ and nonmetricity tensor
$Q_{\alpha\mu\nu}$ are defined as

$\displaystyle T^{\alpha}_{\ \mu\nu}$
$\displaystyle\equiv\varGamma^{\alpha}_{\ \mu\nu}-\varGamma^{\alpha}_{\
\nu\mu}\,,$ (5)
$\displaystyle\begin{split}Q_{\alpha\mu\nu}&\equiv\hat{\nabla}_{\alpha}g_{\mu\nu}\\\
&=\partial_{\alpha}g_{\mu\nu}-g_{\nu\sigma}\varGamma^{\sigma}_{\
\mu\alpha}-g_{\sigma\mu}\varGamma^{\sigma}_{\ \nu\alpha}\,.\end{split}$ (6)

Using two types of trace of nonmetricity tensor

$\displaystyle Q_{\alpha}\equiv
Q_{\alpha}{}_{\mu}{}^{\mu}\,,\quad\tilde{Q}_{\alpha}\equiv Q^{\mu}_{\
\mu\alpha}\,,$ (7)

we define the nonmetricity conjugate

$\displaystyle
P^{\alpha\mu\nu}=-\frac{1}{4}Q^{\alpha\mu\nu}+\frac{1}{2}Q^{(\mu\nu)\alpha}+\frac{1}{4}\left(Q^{\alpha}-\tilde{Q}^{\alpha}\right)g^{\mu\nu}-\frac{1}{4}g^{\alpha(\mu}Q^{\nu)}\,,$
(8)

and nonmetricity scalar

$\displaystyle\begin{split}Q&=Q_{\alpha\mu\nu}P^{\alpha\mu\nu}\\\
&=-\frac{1}{4}Q^{\alpha\mu\nu}Q_{\alpha\mu\nu}+\frac{1}{2}Q^{\nu\mu\alpha}Q_{\alpha\mu\nu}+\frac{1}{4}Q^{\alpha}Q_{\alpha}-\frac{1}{2}\tilde{Q}^{\alpha}Q_{\alpha}\,.\end{split}$
(9)

The action of STEGR is

$\displaystyle
S=\int\mathrm{d}^{4}x\left(\sqrt{-g}Q_{\alpha\mu\nu}P^{\alpha\mu\nu}+\lambda_{\alpha}{}^{\beta\mu\nu}R^{\alpha}{}_{\beta\mu\nu}+\lambda_{a}{}^{\mu\nu}T^{\alpha}_{\
\mu\nu}\right)\,.$ (10)

Variation of the action (10) with respect to Lagrange multipliers
$\lambda_{\alpha}{}^{\beta\mu\nu}$ and $\lambda_{a}{}^{\mu\nu}$ generates the
torsion-free condition

$\displaystyle T^{\alpha}_{\,\mu\nu}\stackrel{{\scriptstyle!}}{{=}}0$ (11)

and vanishing curvature

$\displaystyle{R}_{\
\beta\mu\nu}^{\alpha}(\varGamma)\equiv\partial_{\mu}\varGamma^{\alpha}_{\
\nu\beta}-\partial_{\nu}\varGamma_{\ \mu\beta}^{\alpha}+\varGamma_{\
\mu\lambda}^{\alpha}\varGamma_{\ \nu\beta}^{\lambda}-\varGamma_{\
\nu\lambda}^{\alpha}\varGamma_{\
\mu\beta}^{\lambda}\stackrel{{\scriptstyle!}}{{=}}0\,.$ (12)

The former demands the general affine connection to be symmetric to its
lowered two indices, as the latter implies that the connection must have the
form $\varGamma^{\lambda}_{\,\mu\nu}=(A^{-1})^{\lambda}_{\
\alpha}\partial_{\mu}A^{\alpha}_{\ \nu}$ Zhao (2022).

###  II.2 Coincident gauge

Combining two conditions on the affine connection as in Eqs. (11) and (12), we
find that the matrix $A^{\alpha}{}_{\beta}=\partial_{\beta}\xi^{\alpha}$, and
the affine connection consequently takes the following form:

$\displaystyle\varGamma^{\lambda}_{\,\mu\nu}=\frac{\partial
x^{\lambda}}{\partial\xi^{\alpha}}\frac{\partial^{2}\xi^{\alpha}}{\partial
x^{\mu}\partial x^{\nu}}\,.$ (13)

Here, four scalar fields $\xi^{\alpha}(x)$ are the arbitrary functions of
coordinates. Note that the connection is symmetric with respect to two lower
indices because two partial derivatives of $\xi^{\alpha}$ commute each other.
Moreover, one finds $A^{\alpha}_{\mu}$ corresponds to the tetrad $e^{a}_{\mu}$
for Weitzenböck connection in the teleparallel gravity.

A special choice $\xi^{\alpha}(x)=x^{\alpha}$ makes the connection vanish
$\varGamma^{\lambda}_{\ \mu\nu}=0$. This choice is called the coincident
gauge. Under the coincident gauge, the covariant derivative is replaced by the
partial derivative $\nabla_{\alpha}\mapsto\partial_{\alpha}$, and the
nonmetricity tensor $Q_{\alpha\mu\nu}$, deformation tensor
$L^{\alpha}_{\,\mu\nu}$, nonmetricity scalar $Q$, and the traces $Q_{\alpha}$,
$\tilde{Q}_{\alpha}$ reduce to

$\displaystyle Q_{\alpha\mu\nu}\mapsto\mathring{Q}_{\alpha\mu\nu}$
$\displaystyle=\partial_{\alpha}g_{\mu\nu}\,,$ (14) $\displaystyle
L^{\alpha}_{\,\mu\nu}\mapsto\mathring{L}_{\,\mu\nu}^{\alpha}$
$\displaystyle=-\left\\{{}^{\,\alpha}_{\mu\nu}\right\\}\,,$ (15)
$\displaystyle\begin{split}Q\mapsto\mathring{Q}&=g^{\mu\nu}\left(\mathring{L}_{\
\sigma\mu}^{\alpha}\mathring{L}_{\ \nu\alpha}^{\sigma}-\mathring{L}_{\
\sigma\alpha}^{\alpha}\mathring{L}_{\ \mu\nu}^{\sigma}\right)\\\
&=g^{\mu\nu}\left(\left\\{{}^{\,\alpha}_{\sigma\mu}\right\\}\left\\{{}^{\,\sigma}_{\nu\alpha}\right\\}-\left\\{{}^{\,\alpha}_{\sigma\alpha}\right\\}\left\\{{}^{\,\sigma}_{\mu\nu}\right\\}\right)\,,\end{split}$
(16)
$\displaystyle\begin{split}Q^{\alpha}-\tilde{Q}^{\alpha}\mapsto\mathring{Q}^{\alpha}-\mathring{\tilde{Q}}^{\alpha}&=g^{\mu\nu}\mathring{L}^{\,\alpha}_{\,\mu\nu}-g^{\mu\alpha}\mathring{L}^{\,\nu}_{\,\mu\nu}\\\
&=-g^{\mu\nu}\left\\{{}^{\,\alpha}_{\mu\nu}\right\\}+g^{\mu\alpha}\left\\{{}^{\,\nu}_{\mu\nu}\right\\}\,.\end{split}$
(17)

Writing the Einstein-Hilbert action as follows

$\displaystyle\begin{split}\mathcal{S}_{{\mathrm{EH}}}&=\int
d^{4}x\,\sqrt{-g}\,g^{\mu\nu}\mathcal{R}_{\mu\nu}\left(\left\\{{}^{\,\alpha}_{\mu\nu}\right\\}\right)\\\
&=\int
d^{4}x\,\sqrt{-g}g^{\mu\nu}\left(\left\\{{}^{\,\alpha}_{\sigma\mu}\right\\}\left\\{{}^{\,\sigma}_{\nu\alpha}\right\\}-\left\\{{}^{\,\alpha}_{\sigma\alpha}\right\\}\left\\{{}^{\,\sigma}_{\mu\nu}\right\\}\right)\\\
&\quad+\int
d^{4}x\,\partial_{\alpha}\left[\sqrt{-g}\,(g^{\mu\nu}\left\\{{}^{\,\alpha}_{\mu\nu}\right\\}-g^{\mu\alpha}\left\\{{}^{\,\nu}_{\mu\nu}\right\\})\right]\,\end{split}$
(18)

we find that the Ricci scalar can be expressed by the nonmetricity scalar and
divergence of two independent traces of the nonmetricity tensor, where the
divergence term is reduced to the boundary term after the integration. This
result indicates the equivalence between GR and STEGR up to the boundary term.
Moreover, for general connection, we find the relation between $\mathcal{R}$:

$\displaystyle\mathcal{R}=Q-{\nabla}_{\alpha}\left(Q^{\alpha}-\tilde{Q}^{\alpha}\right)$
(19)

###  II.3 Symmetry and Stüeckelberg field

Under the general coordinate transformation
$\\{X^{\mu}\\}\rightarrow\\{x^{\mu}\\}$, the connection is transformed as

$\displaystyle\varGamma^{\lambda}_{\ \mu\nu}(x)=\frac{\partial
x^{\lambda}}{\partial X^{\alpha}}\frac{\partial X^{\rho}}{\partial
x^{\mu}}\frac{\partial X^{\sigma}}{\partial x^{\nu}}\varGamma^{\alpha}_{\
\rho\sigma}(X)+\frac{\partial x^{\lambda}}{\partial
X^{\alpha}}\frac{\partial^{2}X^{\alpha}}{\partial x^{\mu}\partial
x^{\nu}}~{}.$ (20)

Provided that the connection vanishes in the coordinate system
$\\{X^{\mu}\\}$, the first term vanishes in Eq. (20), and the second term
reproduces Eq. (13). Thus, we can identify a set of scalar fields
$\\{\xi^{\mu}\\}$ as such a special coordinate system $\\{X^{\mu}\\}$. Note
that the above argument partially includes the definition of the local Lorentz
frame, although the space-time is not necessarily flat in the current setup.

Equation (13) also suggests that the index of the scalar fields is contracted
within the connection, and it does not show up as a free index. Therefore, one
can expect that the index $\alpha$ in $\xi^{\alpha}$ is related to the
inertial symmetry rather than the space-time index. To distinguish it from the
space-time index, hereafter, we denote the scalar field as $\xi^{a}$:

$\displaystyle\varGamma^{\lambda}_{\,\mu\nu}=\frac{\partial
x^{\lambda}}{\partial\xi^{a}}\frac{\partial^{2}\xi^{a}}{\partial
x^{\mu}\partial x^{\nu}}\,.$ (21)

We can find that Eq. (21) is invariant under the following transformation for
$\xi^{a}$:

$\displaystyle\xi^{a}(x)\rightarrow\bar{\xi}^{a}(x)=\mathcal{M}^{a}_{\
b}\xi^{b}(x)+\zeta^{a}\,.$ (22)

Here, $\mathcal{M}^{a}_{\ b}$ is a $4\times 4$ nondegenerate constant matrix,
and $\zeta^{a}$ is a constant vector. In the coincident gauge $\xi^{a}=x^{a}$,
the above transformation leads to the linear transformation between two
coordinate systems,

$\displaystyle\bar{x}^{a}=\mathcal{M}^{a}_{\ b}x^{b}+\zeta^{a}\,.$ (23)

Equation (23) represents the affine transformation, and the connection (21) is
purely inertial, which always can transform to zero under Eq. (23). By
imposing on the flat space-time, the Minkowski metric is invariant under the
transformation, we find that matrix $\mathcal{M}^{a}_{\ b}$ corresponds to the
Lorentz transformation, and Eq. (23) is reduced to the Poincaré transformation
in the coincident gauge.

Thus, any two different frames under the coincident gauge are related by
transformation (23), naturally to introduce the metric $f_{ab}(x)$ in the
coincident gauge. Under the coordinate transformation
$\\{\xi^{a}\\}\rightarrow\\{x^{\mu}\\}$, the space-time metric is transformed
as $f_{ab}(x)\rightarrow g_{\mu\nu}(x)$ written in terms of the Jacobi matrix:

$\displaystyle g_{\mu\nu}=\frac{\partial\xi^{a}}{\partial
x^{\mu}}\frac{\partial\xi^{b}}{\partial x^{\nu}}f_{ab}\,.$ (24)

We present an explicit example in Appendix A to better demonstrate how
$f_{ab}$ transforms under Eq. (24). It is apparent that $\xi^{a}$ plays a role
in restoring the general covariance while fixing $\xi^{a}$ breaks the
diffeomorphism explicitly. A similar structure can be found in the relation
between the tetrad and metric if we read $\partial_{\mu}\xi^{a}=e^{a}_{\mu}$
and $f_{ab}=\eta_{ab}$. Moreover, Eq. (24) shows up in the mass term of the de
Rham-Gabadadze-Tolley (dRGT) massive gravity de Rham  _et al._ (2011, 2012);
Hassan and Rosen (2011); Hinterbichler (2012); Molaee and Shirzad (2018) , and
$\xi^{a}=x^{a}$ is called the unitary gauge therein. The above arguments
suggest that we can treat four scalar fields $\xi^{a}$ as the Stüeckelberg
fields associated with the diffeomorphism as discussed in the dRGT massive
gravity. Hereafter we will handle the generic $\xi^{a}$ to investigate the
$f(Q)$ gravity in a covariant way.

##  III Ghost scalar field in $f(Q)$ gravity

In this section, we recall the scalar-nonmetricity formulation of $f(Q)$
gravity and apply the conformal re-scaling transformation to the action. We
suggest that the scalar field is a ghost mode in the scalar-nonmetricity
formulation.

###  III.1 $f(Q)$ gravity and scalar-nonmetricity formulation

Let us consider the $f(Q)$ gravity by replacing $Q$ with a function of $Q$:

$\displaystyle S_{f(Q)}=\int d^{4}x\sqrt{-g}f(Q)\,.$ (25)

By introducing an auxiliary scalar field $\phi$, one can consider the
following action

$\displaystyle S_{f(Q)}=\int
d^{4}x\sqrt{-g}\left\\{f^{\prime}(\phi)Q-\left[\phi
f^{\prime}(\phi)-f(\phi)\right]\right\\}\;.$ (26)

By varying the above action with respect to $\phi$, one gets
$f^{\prime\prime}(\phi)(Q-\phi)=0$. Provided that $f^{\prime\prime}(\phi)\neq
0$, this equation implies $\phi=Q$ and restores the original action. Moreover,
we can redefine the scalar field $\varphi\equiv f^{\prime}(\phi)$ and
$V(\varphi)=\phi f^{\prime}(\phi)-f(\phi)$. The Lagrangian density of our
interest then takes the following form

$\displaystyle S_{f(Q)}=\int d^{4}x\sqrt{-g}\left[\varphi
Q-V(\varphi)\right]~{}.$ (27)

Using Eq. (19), we obtain

$\displaystyle S_{f(Q)}=\int
d^{4}x\sqrt{-g}\left[\varphi\mathcal{R}-V(\varphi)+\varphi\nabla_{\mu}\left(Q^{\mu}-\tilde{Q}^{\mu}\right)\right]~{}.$
(28)

The third term includes the covariant derivative with respect to the Levi-
Civita connection $\nabla_{\mu}$, and we can utilize the ordinary formula for
the divergence,

$\displaystyle\int d^{4}x\sqrt{-g}~{}\nabla_{\mu}A^{\mu}=\int
d^{4}x~{}\partial_{\mu}\left(\sqrt{-g}A^{\mu}\right)\,,$ (29)

and ignore the surface integration. The integration by parts reduces the
original action (25) to the following form:

$\displaystyle S_{f(Q)}=\int
d^{4}x\sqrt{-g}\left[\varphi\mathcal{R}-V(\varphi)-\partial_{\mu}\varphi\cdot\left(Q^{\mu}-\tilde{Q}^{\mu}\right)\right]~{}.$
(30)

###  III.2 Conformal rescaling of $f(Q)$ gravity

We further deform Eq. (30), following the frame transformation established in
$f(R)$ gravity. We consider the following conformal transformation of the
metric:

$\displaystyle g_{\mu\nu}\rightarrow
g^{\prime}_{\mu\nu}=e^{-\Phi}g_{\mu\nu}\,,\quad\Phi=-\ln\varphi~{}.~{}$ (31)

It is worth mentioning that the affine connection is invariant under the
transformation (31) since the connection is written in Eq. (21) by the
Stüeckelberg field, which is independent of the metric. By the conformal
transformation of metric, the Ricci scalar $\mathcal{R}$ with respect to the
Levi-Civita connection is transformed as

$\displaystyle\mathcal{R}=$
$\displaystyle\mathrm{e}^{-\Phi}\left[\mathcal{R}^{\prime}-3g^{\prime\mu\nu}\nabla^{\prime}_{\mu}\partial_{\nu}\Phi-\frac{3}{2}g^{\prime\mu\nu}\left(\partial_{\mu}\Phi\right)\left(\partial_{\nu}\Phi\right)\right]\,.$
(32)

Here, we denote by $\nabla^{\prime}$ the covariant derivative with respect to
the Levi-Civita connection after the transformation.

Moreover, by the transformation, the nonmetricity tensor
$Q_{\alpha\beta\gamma}$ and its traces $Q^{\mu}$ and $\tilde{Q}^{\mu}$ are
transformed as below:

$\displaystyle\begin{split}Q_{\alpha\beta\gamma}\rightarrow
Q^{\prime}_{\alpha\beta\gamma}&=\hat{\nabla}_{\alpha}g^{\prime}_{\beta\gamma}\\\
&=\hat{\nabla}_{\alpha}\left(e^{-\Phi}g_{\beta\gamma}\right)\\\
&=e^{-\Phi}Q_{\alpha\beta\gamma}-g_{\beta\gamma}e^{-\Phi}\partial_{\alpha}\Phi\;,\end{split}$
(33) $\displaystyle\begin{split}Q^{\mu}\rightarrow
Q^{\prime\mu}&=g^{\prime\mu\alpha}g^{\prime\beta\gamma}Q^{\prime}_{\alpha\beta\gamma}\\\
&=e^{2\Phi}g^{\mu\alpha}g^{\beta\gamma}\left(e^{-\Phi}Q_{\alpha\beta\gamma}-g_{\beta\gamma}e^{-\Phi}\partial_{\alpha}\Phi\right)\\\
&=e^{\Phi}Q^{\mu}-4e^{\Phi}g^{\mu\alpha}\partial_{\alpha}\Phi\;,\end{split}$
(34)
$\displaystyle\begin{split}\tilde{Q}^{\mu}\rightarrow\tilde{Q}^{\prime\mu}&=g^{\prime\mu\alpha}g^{\prime\beta\gamma}Q^{\prime}_{\beta\alpha\gamma}\\\
&=e^{2\Phi}g^{\mu\alpha}g^{\beta\gamma}\left(e^{-\Phi}Q_{\beta\alpha\gamma}-g_{\alpha\gamma}e^{-\Phi}\partial_{\beta}\Phi\right)\\\
&=e^{\Phi}\tilde{Q}^{\mu}-e^{\Phi}g^{\mu\alpha}\partial_{\alpha}\Phi\;.\end{split}$
(35)

In Eq. (33), the general affine connection in the covariant derivative
$\hat{\nabla}$ are not transformed. From the above, we can rewrite the
original nonmetricity tensor $Q$ by new ones $Q^{\prime}$ and scalar field
$\Phi$ after the conformal transformation,

$\displaystyle\begin{split}Q^{\mu}&=e^{-\Phi}Q^{\prime\mu}+4g^{\mu\alpha}\partial_{\alpha}\Phi\\\
&=e^{-\Phi}Q^{\prime\mu}+4e^{-\Phi}g^{\prime\mu\alpha}\partial_{\alpha}\Phi\;,\end{split}$
(36)
$\displaystyle\begin{split}\tilde{Q}^{\mu}&=e^{-\Phi}\tilde{Q}^{\prime\mu}+g^{\mu\alpha}\partial_{\alpha}\Phi\\\
&=e^{-\Phi}\tilde{Q}^{\prime\mu}+e^{-\Phi}g^{\prime\mu\alpha}\partial_{\alpha}\Phi\;.\end{split}$
(37)

Thus, the third term in Eq. (30) is given as

$\displaystyle\begin{split}\partial_{\mu}\varphi\cdot\left(Q^{\mu}-\tilde{Q}^{\mu}\right)&=-e^{-\Phi}\partial_{\mu}\Phi\cdot
e^{-\Phi}\left(Q^{\prime\mu}+4g^{\prime\mu\alpha}\partial_{\alpha}\Phi-\tilde{Q}^{\prime\mu}-g^{\prime\mu\alpha}\partial_{\alpha}\Phi\right)\\\
&=-e^{-2\Phi}\left[\partial_{\mu}\Phi\cdot\left(Q^{\prime\mu}-\tilde{Q}^{\prime\mu}\right)+3\partial^{\alpha}\Phi\partial_{\alpha}\Phi\right]\;.\end{split}$
(38)

By substituting Eqs. (32) and (38) into Eq. (30), the action of $f(Q)$ gravity
is written by the Einstein gravity with the minimally coupled scalar field:

$\displaystyle\begin{split}S&=\int
d^{4}x\sqrt{-g^{\prime}}e^{2\Phi}\left\\{\mathrm{e}^{-2\Phi}\left[\mathcal{R}^{\prime}-3g^{\prime\mu\nu}\nabla^{\prime}_{\mu}\partial_{\nu}\Phi-\frac{3}{2}g^{\prime\mu\nu}\left(\partial_{\mu}\Phi\right)\left(\partial_{\nu}\Phi\right)\right]\right.\\\
&\qquad\qquad\qquad\qquad\qquad\left.-V(\Phi)+e^{-2\Phi}\left[\partial_{\mu}\Phi\cdot\left(Q^{\prime\mu}-\tilde{Q}^{\prime\mu}\right)+3\partial^{\alpha}\Phi\partial_{\alpha}\Phi\right]\right\\}\\\
&=\int
d^{4}x\sqrt{-g^{\prime}}\,\left[\mathcal{R}^{\prime}+\frac{3}{2}\partial^{\alpha}\Phi\partial_{\alpha}\Phi-U(\Phi)+\partial_{\mu}\Phi\cdot\left(Q^{\prime\mu}-\tilde{Q}^{\prime\mu}\right)\right]\;.\end{split}$
(39)

Here, we defined $U(\Phi)\equiv e^{2\Phi}V(\Phi)$ and ignored the divergence
term.

In Eq. (39), the sign of the kinetic term of scalar field $\Phi$ is positive,
indicating that it is ghost mode. The conformal transformation of the first
and second terms in Eq. (30) reproduces the minimally coupled scalar field
$\Phi$, which is the well-known result in $f(R)$ gravity. However, since the
third term, the nonmetricity tensor, in Eq. (30) gives rise to a new kinetic
term with the opposite sign, the sign of the kinetic term of the scalar field
is finally reversed.

As a classical theory, the ghost fields have negative kinetic energy and
generate the instability of the system. In a quantum context, the energy is
bounded below by the vacuum, even if there is a ghost. The ghosts, however,
generate negative norm states, which give the negative probability; therefore,
the existence of the ghosts conflicts with the Copenhagen interpretation of
the quantum theory, and the theory is not physically acceptable. In the case
of the gauge theory, if we quantize the system by using the BRS symmetry
Becchi  _et al._ (1976), the ghosts appear in the combinations of zero norm
states in the physical states, which satisfy the constraints coming from the
BRS symmetry or gauge symmetry, and therefore the negative norm states are
eliminated Kugo and Ojima (1978, 1979). The situation in $f(Q)$ gravity could
be similar. Although there appear to be ghosts even in $f(Q)$ gravity, the
propagation of the ghost fields is excluded as the physical degrees of freedom
by the constraints, as we will see in the next section.

##  IV Reformulation of $f(Q)$ gravity

In this section, we reformulate the action of $f(Q)$ gravity. Regarding the
$f(Q)$ gravity as a higher-derivative scalar-tensor theories Belenchia  _et
al._ (2018); Langlois and Noui (2016); Gao  _et al._ (2019); Kimura  _et al._
(2017); Crisostomi  _et al._ (2016); Ben Achour  _et al._ (2016), we develop
ADM formulation of the covariant action and discuss the origin of the ghost
scalar mode.

###  IV.1 Scalar-vector-tensor formulation

The nonmetricity scalar $Q$ includes the metric and its first derivative as in
Eqs. (6) and (9), while the connection includes the first and second
derivatives of Stüeckelberg fields as in Eq. (21). Thus, the Lagrangian
density of $f(Q)$ gravity (27) is symbolically given as

$\displaystyle\sqrt{-g}\mathcal{L}\,\left(\partial_{\nu}\xi^{a},\partial_{\mu}\partial_{\nu}\xi^{a},g_{\mu\nu},\partial_{\beta}g_{\mu\nu}\right)\,.$
(40)

To handle the higher-order derivatives in the Lagrangian density, we can
introduce new variables into the $f(Q)$ theory in a way that does not alter
the theory Motohashi  _et al._ (2016). By introducing four space-time vectors
$A^{a}_{\;\mu}$ as new variables, we can rewrite the Lagrangian density (40)
as 222 Taking the variation of (41) with respect to $\omega_{\;a}^{\nu}$, we
can get $\partial_{\nu}\xi^{a}=A_{\;\nu}^{a}$. By inserting this relation back
into the equation of motion, it is easy to check that the Lagrangian (41) is
equivalent to (40).

$\displaystyle\sqrt{-g}\mathcal{L}^{eq}=\sqrt{-g}\mathcal{L}\,(A^{a}_{\;\nu},\partial_{\mu}A^{a}_{\;\nu},g_{\mu\nu},\partial_{\beta}g_{\mu\nu})+\omega_{\;a}^{\nu}(\partial_{\nu}\xi^{a}-A_{\;\nu}^{a})\,.$
(41)

We note that by rewriting $\partial_{\nu}\xi^{a}=A_{\;\nu}^{a}$, the affine
connection Eq. (21) goes back the original form suggested by Eq. (12).

Using Eq. (41), we can reformulation the action Eq. (30) into following form
333We present the detail calculation in Appendix B,

$\displaystyle S_{f(Q)}^{eq}$ $\displaystyle=\int
d^{4}x\,\left[\sqrt{-g}\,\left(\varphi\mathcal{R}-V\right)-\sqrt{-g}\partial_{\alpha}\varphi\cdot\left(Q^{\alpha}-\tilde{Q}^{\alpha}\right)+\omega_{\;a}^{\nu}(\partial_{\nu}\xi^{a}-A_{\;\nu}^{a})\right]$
$\displaystyle=S_{f(R)}-\int
d^{4}x\,\left\\{\sqrt{-g}\,\partial_{\alpha}\varphi\left[g^{\mu\nu}(A^{-1})^{\alpha}_{\
a}-g^{\mu\alpha}(A^{-1})^{\nu}_{\ a}\right]\cdot\nabla_{\mu}A^{a}_{\
\nu}-\omega^{\nu}_{\ a}\left(\partial_{\nu}\xi^{a}-A^{a}_{\
\nu}\right)\right\\}$ $\displaystyle=S_{f(R)}-\int
d^{4}x\,\left[\sqrt{-g}\,C_{a}^{\alpha\mu\nu}\partial_{\alpha}\varphi\nabla_{\mu}A^{a}_{\
\nu}-\omega^{\nu}_{\ a}\left(\partial_{\nu}\xi^{a}-A^{a}_{\
\nu}\right)\right]\,.$ (42)

To simplify the expression in Eq. (42), we defined a tensor
$C_{a}^{\alpha\mu\nu}$ as

$\displaystyle C_{a}^{\alpha\mu\nu}\equiv g^{\mu\nu}B^{\alpha}_{\
a}-g^{\mu\alpha}B^{\nu}_{\ a}\quad\text{where}\quad B^{\alpha}_{\
a}=(A^{-1})^{\alpha}_{\ a}\,.$ (43)

Where the term $S_{f(R)}$ in action (42) is the scalar-tensor description of
$f(R)$ gravity,

$\displaystyle S_{f(R)}=\int
d^{4}x\,\sqrt{-g}\left[\varphi\mathcal{R}-V(\varphi)\right]\,.$ (44)

We note that the nonmetricity tensor in Eq. (30) induces the coupling between
the scalar field $\varphi$ and vector field $A^{a}_{\ \mu}$ in Eq. (42). This
term corresponds to the origin of the ghost scalar mode after the conformal
transformation of metric, as we observed in the previous section.

###  IV.2 $3+1$ decomposition of the action

In this subsection, we investigate the origin of the ghost mode under the
Lagrangian formulation by performing the $3+1$ decomposition of the action
(42) in terms of ADM variables. We decompose space-time into 3-dimensional
spacelike hypersurfaces $\Sigma_{t}$ and their normal vector $n^{\mu}$, which
satisfies the normalization condition $n^{\mu}n_{\mu}=-1$

$\displaystyle n_{\mu}$ $\displaystyle=\left(-N,0,0,0\right)\,,\
n^{\mu}=\left(1/N,-N^{i}/N\right)\,.$ (45)

$N$ is the lapse function and $N^{i}$ is the shift vector lying on the
$\Sigma_{t}$. We also introduce the three-dimensional induced metric
$h_{\mu\nu}$ defined by

$\displaystyle\begin{split}h_{\mu\nu}&=g_{\mu\nu}+n_{\mu}n_{\nu}=\left(\begin{array}[]{cc}N_{i}N^{i}&N_{i}\\\
N_{i}&h_{ij}\end{array}\right)\,,\\\
h^{\mu\nu}&=g^{\mu\nu}+n^{\mu}n^{\nu}=\left(\begin{array}[]{cc}0&0\\\
0&h^{ij}\end{array}\right)\,.\end{split}$ (46)

Inversely, we can express the metric in terms of $N$, $N^{i}$, and $h_{ij}$,

$\displaystyle\begin{split}g_{\mu\nu}&=\left(\begin{array}[]{cc}-N^{2}+N_{i}N^{i}&N_{i}\\\
N_{i}&h_{ij}\end{array}\right)\,,\\\
g^{\mu\nu}&=\left(\begin{array}[]{cc}-N^{-2}&N^{-2}N^{i}\\\
N^{-2}N^{i}&h^{ij}-\frac{N^{i}N^{j}}{N^{2}}\end{array}\right)\,.\end{split}$
(47)

Moreover, we can decompose the space-time vector $A_{\mu}^{a}$ into tangent
$\bar{A}_{\mu}^{a}$ and normal part $A_{*}^{a}$ separately:

$\displaystyle
A_{\mu}^{a}=\bar{A}_{\mu}^{a}-A_{*}^{a}n_{\mu}\quad\text{with}\quad
A_{*}^{a}=A_{\alpha}^{a}n^{\alpha}\,.$ (48)

For the tangent part, we have $A_{0}^{a}=\bar{A}_{0}^{a}+A_{*}^{a}N$ and
$A_{i}^{a}=\bar{A}_{i}^{a}$. Provided by these two relations, the normal
component of $A_{\mu}^{a}$ can be expressed by

$\displaystyle\begin{split}A_{*}^{a}=&A_{0}^{a}\cdot n^{0}+A_{i}^{a}n^{i}\\\
&=\frac{1}{N}(\bar{A}_{0}^{a}+A_{*}^{a}N-N^{i}\bar{A}_{i}^{a})\,,\end{split}$
(49)

from which we get the relation $\bar{A}_{0}^{a}=N^{i}\bar{A}_{i}^{a}$. After
straightforward calculations presented in Appendix C, one can obtain the
different components of the covariant derivative of four vectors
$A_{\mu}^{a}$:

$\displaystyle\begin{split}\nabla_{0}A_{0}^{a}&=N\dot{A}_{*}^{a}-\mathcal{K}_{ij}\left(A_{*}^{a}N^{i}N^{j}+N\bar{A}^{ia}N^{j}+N\bar{A}^{ja}N^{i}\right)\\\
&-N^{i}A^{a}_{\ast}\mathcal{D}_{i}N-N^{i}\bar{A}^{a}_{j}\mathcal{D}_{i}N^{j}-N\bar{A}^{ia}\mathcal{D}_{i}N+N^{i}\dot{\bar{A}}_{i}^{a}\,,\\\
\nabla_{i}A_{0}^{a}&=-\mathcal{K}_{ij}\left(A_{*}^{a}N^{j}+N\bar{A}^{ja}\right)+N\mathcal{D}_{i}A_{*}^{a}+N^{j}\mathcal{D}_{i}\bar{A}_{j}^{a}\,,\\\
\nabla_{0}A_{i}^{a}&=\dot{\bar{A}}_{i}^{a}-\mathcal{K}_{ij}\left(A_{*}^{a}N^{j}+N\bar{A}^{ja}\right)-A^{a}_{\ast}\mathcal{D}_{i}N-\bar{A}^{a}_{j}\mathcal{D}_{i}N^{j}\,,\\\
\nabla_{i}A_{j}^{a}&=\mathcal{D}_{i}\bar{A}_{j}^{a}-A_{*}^{a}\mathcal{K}_{ij}\,.\end{split}$
(50)

Subsequently, by using the same approach, the $3+1$ decomposition of
$B_{a}^{\alpha}$ is given by

$\displaystyle
B_{a}^{0}=-B_{a}^{*}\frac{1}{N}~{},~{}B_{a}^{i}=\bar{B}_{a}^{i}+B_{a}^{*}\frac{N^{i}}{N}~{},~{}$
(51)

the calculation details of Eq. (51) and the $3+1$ decomposition of the
interaction term $C_{a}^{\alpha\mu\nu}$ have also shown explicitly in Appendix
C. It is worth mentioning that since the $A_{a}^{\mu}$ and $B_{a}^{\mu}$ are
not independent space-time vectors, they are related by the
$B_{a}^{\mu}A_{\nu}^{a}=\delta^{\mu}_{\;\nu}$, which lead to the following
important relations:

$\displaystyle
A_{*}^{a}B_{a}^{*}=A_{\alpha}^{a}n^{\alpha}B_{a}^{\beta}n_{\beta}=\delta_{\;\alpha}^{\beta}n^{\alpha}n_{\beta}=-1~{},~{}$
(52)

and

$\displaystyle
B_{a}^{\mu}A_{\mu}^{a}=B_{a}^{0}A_{0}^{a}+B_{a}^{i}A_{i}^{a}=\bar{B}_{a}^{i}\bar{A}_{i}^{a}-B_{a}^{*}A_{*}^{a}=4~{},~{}$
(53)

from which, one concludes that the tangent parts of those vectors also enjoy
the relation $\bar{B}_{a}^{i}\bar{A}_{j}^{a}=\delta^{i}_{j}$. In particular,
we have

$\displaystyle\bar{A}_{j}^{a}B_{a}^{*}=A_{j}^{a}(B_{a}^{\alpha}n_{\alpha})=-\delta_{j}^{0}N=0~{}.~{}$
(54)

It shows the tangent part of $A_{a}^{\mu}$ and the normal part of
$B_{a}^{\mu}$ are orthogonal to each other and vice versa, i.e.,
$A_{*}^{a}\bar{B}_{a}^{i}=0$. Putting together all our results and utilizing
relation (52)-(54), we finally obtain the full $3+1$ decomposition of the
action (42),

$\displaystyle\begin{split}\mathcal{L}^{eq}_{f(Q)}&=\,\mathcal{L}_{f(R)}-\mathcal{L}_{\text{int}}+\omega_{a}^{\nu}(\partial_{\nu}\xi^{a}-A_{\nu}^{a})\\\
&=\,N\sqrt{h}\varphi\left({}^{(3)}\mathcal{R}+\mathcal{K}^{ij}\mathcal{K}_{ij}-\mathcal{K}^{2}-\frac{V(\varphi)}{\varphi}\right)\\\
&\quad+\sqrt{h}\dot{\varphi}\left(-\frac{1}{N}\bar{B}^{i}_{a}\dot{\bar{A}}^{a}_{i}+B^{\ast}_{a}\mathcal{D}_{i}\bar{A}^{ia}+\frac{N^{i}}{N}\bar{B}^{j}_{a}\mathcal{D}_{i}\bar{A}^{a}_{j}+\frac{1}{N}\mathcal{D}_{i}N^{i}\right)+\sqrt{h}\dot{A}^{a}_{\ast}\bar{B}^{i}_{a}\mathcal{D}_{i}\varphi\\\
&\quad+\sqrt{h}\mathcal{D}_{i}\varphi\mathcal{D}^{i}N+\sqrt{h}\frac{N^{i}}{N}\left(\bar{B}^{j}_{a}\dot{\bar{A}}^{a}_{j}-\mathcal{D}_{j}N^{j}-N^{j}\bar{B}^{k}_{a}\mathcal{D}_{j}\bar{A}^{a}_{k}\right)\mathcal{D}_{i}\varphi\\\
&\quad-\sqrt{h}N\left(B^{\ast}_{a}\mathcal{D}^{i}A^{a}_{\ast}+\bar{B}^{i}_{a}\mathcal{D}^{j}\bar{A}^{a}_{j}-\bar{B}_{ja}\mathcal{D}^{i}\bar{A}^{ja}\right)\mathcal{D}_{i}\varphi\\\
&\quad-\sqrt{h}N^{i}\left(B^{\ast}_{a}\mathcal{D}_{j}\bar{A}^{ja}\mathcal{D}_{i}\varphi+\bar{B}^{j}_{a}\mathcal{D}_{i}A^{a}_{\ast}\mathcal{D}_{j}\varphi\right)\\\
&\quad+\omega_{a}^{0}\left(\dot{\xi}^{a}-NA_{*}^{a}-N^{i}\bar{A}^{a}_{i}\right)+\omega_{a}^{i}(D_{i}\xi^{a}-\bar{A}_{i}^{a})\,,\end{split}$
(55)

where $\mathcal{K}_{ij}$ is extrinsic curvature of the hypersurface,

$\displaystyle\mathcal{K}_{ij}=\frac{1}{2N}\left(\dot{h}_{ij}-D_{i}N_{j}-D_{j}N_{i}\right)\,.$
(56)

In the Lagrangian formulation, the field equations obtained from (55) involve
only first-order time derivatives of $\varphi$ and ${A}_{*}^{a}$. However,
wave functions require second-order derivatives. In other words, the action
(55) strongly indicates these scalar fields cannot propagate as physical DOF.
Furthermore, we can glimpse the prospect of the Hamiltonian formulation. The
canonical momentum variable with respect to the scalar field $\varphi$ is
given as

$\displaystyle
p=\sqrt{h}\left[-\frac{1}{N}\bar{B}^{i}_{a}\dot{\bar{A}}^{a}_{i}+B^{\ast}_{a}\mathcal{D}_{i}\bar{A}^{ia}+\frac{N^{i}}{N}\bar{B}^{j}_{a}\mathcal{D}_{i}\bar{A}^{a}_{j}+\frac{1}{N}\mathcal{D}_{i}N^{i}\right]\,.$
(57)

Moreover, we consider the constraint from the Lagrange multiplier
$\omega^{\nu}_{\ a}$. We obtain $\nabla_{\mu}A^{a}_{\
\nu}=\nabla_{\nu}A^{a}_{\ \mu}$ after we use the constraint $A^{a}_{\
\mu}=\partial_{\mu}\xi^{a}$ Langlois and Noui (2016). Since the dummy indices
$\mu,\nu$ are symmetric, we find that $\nabla_{i}A_{0}^{a}$ equals to
$\nabla_{0}A_{i}^{a}$, which induces a relation
$\dot{\bar{A}}_{i}^{a}=\mathcal{D}_{i}A_{0}^{a}$. If we use
$\dot{\bar{A}}^{a}_{i}=\mathcal{D}_{i}A^{a}_{0}$ and
$\mathcal{D}_{i}\bar{A}^{a}_{j}=\mathcal{D}_{j}\bar{A}^{a}_{i}$, Eq. (57) is
reduced to the following form:

$\displaystyle
p\approx\sqrt{h}(B_{a}^{\ast}D_{i}\bar{A}^{ia}-\bar{B}_{a}^{i}D_{i}A_{*}^{a})~{}\,,$
(58)

which generates constraints between phase space variables on the hypersurfaces
$\Sigma_{t}$. In other words, we can eliminate the phase space variables
$\varphi$ and its corresponding conjugate momentum $p$ by using Eq. (58) and
other constraints. It indicates that the ghost scalar mode discussed in Sec.
III is nonphysical; thus, we can conclude that the covariant $f(Q)$ gravity
theory is free from the ghosts.

##  V Conclusions and discussions

In this work, we have revisited the theoretical structure of the $f(Q)$
gravity and discussed the ghost scalar mode in the physical DOF. Four scalar
fields $\xi^{a}$ in the affine connection play a role of the Stüeckelberg to
restore the diffeomorphism, similar to the dRGT massive gravity, and the
coincident gauge in the nonmetricity gravity corresponds to the unitary gauge
in the dRGT massive gravity de Rham (2014); Hassan  _et al._ (2012). In the
covariant $f(Q)$ gravity, we have considered the scalar-nonmetricity
formulation and conformal re-scaling of the metric, where this method was well
established in $f(R)$ gravity to study the scalar field in the theory. We have
shown that the scalar field $\varphi$, which stems from the $f(Q)$ functional
DOF, has the negative kinetic energy and is the ghost mode. This result is
consistent with earlier work Beltrán Jiménez and Koivisto (2021), and the
opposite sign in front of the kinetic term of the scalar field originates from
the nonmetricity tensor when we rewrite the nonmetricity scalar by the Ricci
scalar.

To further investigate the ghost mode, we have developed a new formulation of
the $f(Q)$ gravity, inspired by the technique to address the higher derivative
of the scalar field of the HOST theory. Introducing the vector field
$A^{a}_{\mu}$ to rewrite the derivative of the four Stüeckelberg field
$\partial_{\mu}\xi^{a}$, we have shown that the covariant $f(Q)$ gravity
written by the scalar field $\varphi$, vector fields $A^{a}_{\mu}$, and the
tensor field $g_{\mu\nu}$ can be equivalent to the HOST theory. The structure
of the corresponding Lagrangian is relatively simple and explicit, equals to
$f(R)$ term plus an interaction term, which might be reminiscent of a special
class of Horndeski theories in four dimensions, i.e., $L^{\text{H}}_{4}$
Horndeski (1974); Kobayashi (2019). However, it is worth mentioning that the
nonquadratic interacting term in Lagrangian (42) might indicate that the
vectors $A^{a}_{\mu}$ are nonpropagate modes. The coupling between the scalar
and vector fields in the new formulation rephrases the opposite sign of the
scalar field from the nonmetricity tensor. Moreover, we have also applied
$3+1$ decomposition to the new formulation of the $f(Q)$ gravity. By using ADM
variables, we are able to write the full Lagrangian (42) in a relatively
compact form, which significantly simplifies the computation of the
Hamiltonian for further study. At last, we found that in the context of
canonical structure, the ghost scalar mode can be removed by using
constraints; thus, it ceases to propagate, and the theory is free from the
ghost.

We make several remarks on the ghost scalar field in $f(Q)$ gravity. The ghost
is a dynamical DOF with a negative norm and often cancels the physical DOF in
quantum theory. Even in classical theory, if we consider the trace or
determinant of the coefficients of the kinetic term, ghosts may lower the rank
of the kinetic matrix. Therefore, if we incorrectly assume all the DOFs have
the positive norm, the number of the DOFs looks different from what it really
is. It is mandatory to handle the Hamiltonian analysis and to evaluate the
physical DOF, assuming that ghost mode can exist in the theory. The new
formulation we developed will be helpful to understand the theoretical
structure of the covariant $f(Q)$ gravity. In the future, it would be
necessary to confirm the conjecture that whether it suffers from Ostrogradsky
instability brought by higher-order derivatives Kobayashi (2019); Gleyzes  _et
al._ (2015), and we need to check the DOF of covariant $f(Q)$ gravity by
performing complete Hamiltonian analysis. Finally, we comment on the ghost
mode at the quantum level. Although we have confirmed that the ghost scalar is
not propagating at the classical level, such a ghost might propagate at the
quantum level. Nonpropagation of the ghost scalar at the classical level
relies on the symmetry $\nabla_{\mu}A^{a}_{\ \nu}=\nabla_{\nu}A^{a}_{\ \mu}$
after we use the constraint $A^{a}_{\ \mu}=\partial_{\mu}\xi^{a}$. It is
intriguing to investigate the relation between the propagation of the ghost
scalar and the symmetry induced by the constraint in terms of the quantum
field theory.

###### Acknowledgements.

T.K. is supported by the National Key R&D Program of China (2021YFA0718500)
and by Grant-in-Aid of Hubei Province Natural Science Foundation (2022CFB817).
T.Q. and K.H. are supported by the National Key Research and Development
Program of China under Grant No. 2021YFC2203100, and the National Science
Foundation of China under Grant No. 11875141. S.N. was partially supported by
the Maria de Maeztu Visiting Professorship at the Institute of Space Sciences,
Barcelona.

##  Appendix A THE GENERAL COORDINATE TRANSFORMATION: AN EXAMPLE

One can use Eq. (24) to freely transform from tangent space-time coordinates
to arbitrary space-time coordinates. In order to better illustrate this point,
we consider the nonmetricity tensor.

$\displaystyle\begin{split}Q_{\alpha\mu\nu}&=\hat{\nabla}_{\alpha}g_{\mu\nu}=\partial_{\alpha}g_{\mu\nu}-g_{\nu\sigma}\varGamma_{\;\mu\alpha}^{\sigma}-g_{\sigma\mu}\varGamma_{\;\nu\alpha}^{\sigma}\\\
&=\partial_{\alpha}g_{\mu\nu}-g_{\nu\sigma}\left(\frac{\partial
x^{\sigma}}{\partial\xi^{c}}\frac{\partial^{2}\xi^{c}}{\partial
x^{\mu}\partial x^{\alpha}}\right)-g_{\sigma\mu}\left(\frac{\partial
x^{\sigma}}{\partial\xi^{c}}\frac{\partial^{2}\xi^{c}}{\partial
x^{\nu}\partial x^{\alpha}}\right)\\\
&=\partial_{\alpha}g_{\mu\nu}-\left(\frac{\partial\xi^{a}}{\partial
x^{\nu}}\frac{\partial\xi^{b}}{\partial x^{\sigma}}f_{ab}\right)\frac{\partial
x^{\sigma}}{\partial\xi^{c}}\cdot\partial_{\alpha}\frac{\partial\xi^{c}}{\partial
x^{\mu}}-\left(\frac{\partial\xi^{a}}{\partial
x^{\sigma}}\frac{\partial\xi^{b}}{\partial
x^{\mu}}f_{ab}\right)\left(\frac{\partial
x^{\sigma}}{\partial\xi^{c}}\cdot\partial_{\alpha}\frac{\partial\xi^{c}}{\partial
x^{\nu}}\right)\\\
&=\partial_{\alpha}g_{\mu\nu}-\frac{\partial\xi^{c}}{\partial
x^{\nu}}f_{ac}\cdot\partial_{\alpha}\frac{\partial\xi^{a}}{\partial
x^{\mu}}-\frac{\partial\xi^{a}}{\partial
x^{\mu}}f_{ac}\cdot\partial_{\alpha}\frac{\partial\xi^{c}}{\partial
x^{\nu}}\\\
&=\partial_{\alpha}g_{\mu\nu}-\left(\frac{\partial\xi^{c}}{\partial
x^{\nu}}\cdot\partial_{\alpha}\frac{\partial\xi^{a}}{\partial
x^{\mu}}+\frac{\partial\xi^{a}}{\partial
x^{\mu}}\cdot\partial_{\alpha}\frac{\partial\xi^{c}}{\partial
x^{\nu}}\right)f_{ac}\\\
&=\partial_{\alpha}g_{\mu\nu}-\partial_{\alpha}\left(\frac{\partial\xi^{c}}{\partial
x^{\nu}}\frac{\partial\xi^{a}}{\partial
x^{\mu}}f_{ac}\right)+\left(\partial_{\alpha}f_{ac}\right)\cdot\frac{\partial\xi^{c}}{\partial
x^{\nu}}\frac{\partial\xi^{a}}{\partial x^{\mu}}\\\
&=\frac{\partial\xi^{c}}{\partial x^{\nu}}\frac{\partial\xi^{a}}{\partial
x^{\mu}}\frac{\partial\xi^{d}}{\partial
x^{\alpha}}\mathring{Q}_{dac}\,.\end{split}$ (59)

We have used the condition $\mathring{Q}_{abc}=\partial_{a}f_{bc}$ in the last
step. One may notice Eq. (59) is nothing but the transformation of
nonmetricity tensor between arbitrary gauge and coincident gauge.

##  Appendix B SCALAR-VECTOR-TENSOR REPRESENTATION OF $f(Q)$ GRAVITY

We show the calculation detail in Eq. (42) up to $S_{f(R)}$ and Lagrange
multiplier $\omega_{\;a}^{\nu}(\partial_{\nu}\xi^{a}-A_{\;\nu}^{a})$:

$\displaystyle\int
d^{4}x\,\left[\sqrt{-g}\partial_{\alpha}\varphi\cdot(Q^{\alpha}-\tilde{Q}^{\alpha})\right]\,.$
(60)

Computing the integrand, we find the traces of the nonmetricity tensor are
written by metric $g_{\mu\nu}$ and vector field $A^{a}_{\ \mu}$ as follows:

$\displaystyle Q^{\alpha}-\tilde{Q}^{\alpha}$
$\displaystyle=g^{\beta\alpha}g^{\mu\nu}(\partial_{\beta}g_{\mu\nu}-\partial_{\mu}g_{\beta\nu})+\left[g^{\mu\nu}(A^{-1})_{\;a}^{\alpha}-g^{\mu\alpha}(A^{-1})_{\;a}^{\nu}\right]\partial_{\mu}A_{\;\nu}^{a}\,.$
(61)

Using the following two relations

$\displaystyle\begin{split}\partial_{\beta}g_{\mu\nu}&=g_{\lambda\nu}\left\\{{}^{\,\lambda}_{\mu\beta}\right\\}+g_{\lambda\mu}\left\\{{}^{\,\lambda}_{\nu\beta}\right\\}\,,\\\
\nabla_{\mu}A_{\nu}^{a}&=\partial_{\mu}A_{\nu}^{a}-\left\\{{}^{\,\alpha}_{\mu\nu}\right\\}A_{\alpha}^{a}\,,\end{split}$
(62)

we can compute the first term in Eq. (61) as

$\displaystyle\begin{split}g^{\beta\alpha}g^{\mu\nu}\left(\partial_{\beta}g_{\mu\nu}-\partial_{\mu}g_{\beta\nu}\right)&=\left(g^{\beta\alpha}g^{\mu\nu}-g^{\mu\alpha}g^{\beta\nu}\right)\partial_{\beta}g_{\mu\nu}\\\
&=\left(g^{\beta\alpha}g^{\mu\nu}-g^{\mu\alpha}g^{\beta\nu}\right)\left(g_{\lambda\nu}\left\\{{}^{\,\lambda}_{\mu\beta}\right\\}+g_{\lambda\mu}\left\\{{}^{\,\lambda}_{\nu\beta}\right\\}\right)\\\
&=g^{\beta\alpha}\left\\{{}^{\,\lambda}_{\lambda\beta}\right\\}-g^{\beta\nu}\left\\{{}^{\,\alpha}_{\nu\beta}\right\\}~{},\end{split}$
(63)

and the second term in Eq. (61) as

$\displaystyle\begin{split}\left[g^{\mu\nu}(A^{-1})_{\,a}^{\alpha}-g^{\mu\alpha}(A^{-1})_{\,a}^{\nu}\right]\cdot\partial_{\mu}A_{\;\nu}^{a}&=\left[g^{\mu\nu}(A^{-1})_{\,a}^{\alpha}-g^{\mu\alpha}(A^{-1})_{\,a}^{\nu}\right]\cdot\left(\nabla_{\mu}A_{\;\nu}^{a}+\left\\{{}^{\,\lambda}_{\mu\nu}\right\\}A_{\,\lambda}^{a}\right)\\\
&=\left[g^{\mu\nu}(A^{-1})_{\;a}^{\alpha}\cdot\left\\{{}^{\,\lambda}_{\mu\nu}\right\\}A_{\,\lambda}^{a}-g^{\mu\alpha}(A^{-1})_{\;a}^{\nu}\cdot\left\\{{}^{\,\lambda}_{\mu\nu}\right\\}A_{\,\lambda}^{a}\right]\\\
&\qquad+\left[g^{\mu\nu}(A^{-1})_{\;a}^{\alpha}-g^{\mu\alpha}(A^{-1})_{\;a}^{\nu}\right]\cdot\nabla_{\mu}A_{\;\nu}^{a}\\\
&=g^{\mu\nu}\cdot\left\\{{}^{\,\alpha}_{\mu\nu}\right\\}-g^{\mu\alpha}\cdot\left\\{{}^{\,\nu}_{\mu\nu}\right\\}\\\
&\qquad+\left[g^{\mu\nu}(A^{-1})_{\;a}^{\alpha}-g^{\mu\alpha}(A^{-1})_{\;a}^{\nu}\right]\cdot\nabla_{\mu}A_{\;\nu}^{a}~{}.\end{split}$
(64)

Finally, we obtain

$\displaystyle Q^{\alpha}-\tilde{Q}^{\alpha}$
$\displaystyle=\left[g^{\mu\nu}(A^{-1})_{\;a}^{\alpha}-g^{\mu\alpha}(A^{-1})_{\;a}^{\nu}\right]\cdot\nabla_{\mu}A_{\;\nu}^{a}~{},$
(65)

and thus

$\displaystyle\begin{split}\int
d^{4}x\,\left[\sqrt{-g}\partial_{\alpha}\varphi\cdot(Q^{\alpha}-\tilde{Q}^{\alpha})\right]=\int
d^{4}x\,\left[\sqrt{-g}\partial_{\alpha}\varphi\left[g^{\mu\nu}(A^{-1})_{\;a}^{\alpha}-g^{\mu\alpha}(A^{-1})_{\;a}^{\nu}\right]\cdot\nabla_{\mu}A_{\;\nu}^{a}\right]\,.\end{split}$
(66)

##  Appendix C ADM DECOMPOSITION OF THE ACTION

This Appendix aims to apply ADM decomposition to the action (42). In order to
facilitate the calculation, we will divide our goal into two parts
$\nabla_{\mu}A_{\nu}^{a}$, $C_{a}^{\alpha\mu\nu}$, and compute them
separately.

####  C.0.1 $3+1$ decomposition of $\nabla_{\mu}A_{\nu}^{a}$

Utilizing the expressions of the Christoffel symbols (2) in terms of ADM
quantities (47), we have

$\displaystyle\begin{split}\left\\{{}^{\,0}_{00}\right\\}&=\frac{1}{N}\left(N^{i}N^{j}\mathcal{K}_{ij}+\dot{N}+N^{i}\mathcal{D}_{i}N\right)\,,\\\
\left\\{{}^{\,k}_{00}\right\\}&=NN^{i}\left(2h^{jk}-\frac{N^{j}N^{k}}{N^{2}}\right)\mathcal{K}_{ij}+\dot{N^{k}}\\\
&\qquad-\frac{N^{k}}{N}\dot{N}+N^{i}\mathcal{D}_{i}N^{k}+N\left(h^{ki}-\frac{N^{k}N^{i}}{N^{2}}\right)\mathcal{D}_{i}N\,,\\\
\left\\{{}^{\,0}_{0i}\right\\}&=\frac{1}{N}\left(N^{j}\mathcal{K}_{ij}+\mathcal{D}_{i}N\right)\,,\\\
\left\\{{}^{\,j}_{0i}\right\\}&=N\left(h^{jk}-\frac{N^{j}N^{k}}{N^{2}}\right)\mathcal{K}_{ik}+\mathcal{D}_{i}N^{j}-\frac{N^{j}}{N}\mathcal{D}_{i}N\,,\\\
\left\\{{}^{\,0}_{ij}\right\\}&=\frac{1}{N}\mathcal{K}_{ij}\,,\\\
\left\\{{}^{\,k}_{ij}\right\\}&=-\frac{N^{k}}{N}\mathcal{K}_{ij}+^{3}\left\\{{}^{\,k}_{ij}\right\\}\,.\end{split}$
(67)

where $\mathcal{K}_{ij}$ is extrinsic curvature of the hypersurface defined by
Eq. (56). From these expressions, we compute each component of the covariant
derivative of vector $\nabla_{\mu}A_{\nu}^{a}$:

$\displaystyle\begin{split}\nabla_{0}A_{0}^{a}&=N\dot{A}_{*}^{a}-\mathcal{K}_{ij}\left(A_{*}^{a}N^{i}N^{j}+N\bar{A}^{ia}N^{j}+N\bar{A}^{ja}N^{i}\right)\\\
&-N^{i}A^{a}_{\ast}\mathcal{D}_{i}N-N^{i}\bar{A}^{a}_{j}\mathcal{D}_{i}N^{j}-N\bar{A}^{ia}\mathcal{D}_{i}N+N^{i}\dot{\bar{A}}_{i}^{a}\,,\\\
\nabla_{i}A_{0}^{a}&=-\mathcal{K}_{ij}\left(A_{*}^{a}N^{j}+N\bar{A}^{ja}\right)+N\mathcal{D}_{i}A_{*}^{a}+N^{j}\mathcal{D}_{i}\bar{A}_{j}^{a}\,,\\\
\nabla_{0}A_{i}^{a}&=\dot{\bar{A}}_{i}^{a}-\mathcal{K}_{ij}\left(A_{*}^{a}N^{j}+N\bar{A}^{ja}\right)-A^{a}_{\ast}\mathcal{D}_{i}N-\bar{A}^{a}_{j}\mathcal{D}_{i}N^{j}\,,\\\
\nabla_{i}A_{j}^{a}&=\mathcal{D}_{i}\bar{A}_{j}^{a}-A_{*}^{a}\mathcal{K}_{ij}\,.\end{split}$
(68)

####  C.0.2 $3+1$ decomposition of $B_{a}^{\alpha}$ and $C_{a}^{\alpha\mu\nu}$

First, we apply $3+1$ decomposition to vector field $B_{a}^{\alpha}$. From the
property $B_{a}^{\alpha}=\bar{B}_{a}^{\alpha}-B_{a}^{*}n^{\alpha}$, we find

$\displaystyle\begin{split}B_{a}^{0}=&\bar{B}_{a}^{0}-B_{a}^{*}\frac{1}{N}\,,\\\
B_{a}^{i}=&\bar{B}_{a}^{i}+B_{a}^{*}\frac{N^{i}}{N}\,.\\\ \end{split}$ (69)

The norm part of the vector $B_{a}^{\alpha}$ is given as

$\displaystyle\begin{split}B_{a}^{*}&=B_{a}^{\alpha}n_{\alpha}=B_{a}^{0}n_{0}\\\
&=-N\left(\bar{B}_{a}^{0}-B_{*}\frac{1}{N}\right)\\\
&=-N\bar{B}_{a}^{0}+B_{a}^{*}\,,\end{split}$ (70)

from which we can conclude $\bar{B}_{a}^{0}=0$. Thus, we obtain

$\displaystyle B_{a}^{0}=-B_{a}^{*}\frac{1}{N}\,,\
B_{a}^{i}=\bar{B}_{a}^{i}+B_{a}^{*}\frac{N^{i}}{N}\,.$ (71)

Inserting Eq. (71) into Eq. (43), we can compute each component of
$C_{a}^{\alpha\mu\nu}$:

$\displaystyle\begin{split}C^{000}_{a}&=g^{00}B^{0}_{a}-g^{00}B^{0}_{a}=0\,,\\\
C^{00i}_{a}&=g^{0i}B^{0}_{a}-g^{00}B^{i}_{a}=\frac{N^{i}}{N^{2}}B^{0}_{a}+\frac{1}{N^{2}}B^{i}_{a}=\frac{1}{N^{2}}\bar{B}^{i}_{a}\,,\\\
C^{0i0}_{a}&=g^{i0}B^{0}_{a}-g^{i0}B^{0}_{a}=0\,,\\\
C^{i00}_{a}&=g^{00}B^{i}_{a}-g^{0i}B^{0}_{a}=-\frac{N^{i}}{N^{2}}B^{0}_{a}-\frac{1}{N^{2}}B^{i}_{a}=-\frac{1}{N^{2}}\bar{B}^{i}_{a}\,,\\\
C^{ij0}_{a}&=g^{j0}B^{i}_{a}-g^{ji}B^{0}_{a}=-\left(h^{ij}-\frac{N^{i}N^{j}}{N^{2}}\right)B^{0}_{a}+\frac{N^{j}}{N^{2}}B^{i}_{a}=\frac{h^{ij}}{N}B^{\ast}_{a}+\frac{N^{j}}{N^{2}}\bar{B}^{i}_{a}\,,\\\
C^{i0j}_{a}&=g^{0j}B^{i}_{a}-g^{0i}B^{j}_{a}=\frac{N^{j}}{N^{2}}B^{i}_{a}-\frac{N^{i}}{N^{2}}B^{j}_{a}=\frac{1}{N^{2}}\left(N^{j}\bar{B}^{i}_{a}-N^{i}\bar{B}^{j}_{a}\right)\,,\\\
C^{0ij}_{a}&=g^{ij}B^{0}_{a}-g^{i0}B^{j}_{a}=\left(h^{ij}-\frac{N^{i}N^{j}}{N^{2}}\right)B^{0}_{a}-\frac{N^{i}}{N^{2}}B^{j}_{a}=-\frac{h^{ij}}{N}B^{\ast}_{a}-\frac{N^{i}}{N^{2}}\bar{B}^{j}_{a}\,,\\\
C^{ijk}_{a}&=g^{jk}B^{i}_{a}-g^{ji}B^{k}_{a}=\left(h^{jk}-\frac{N^{j}N^{k}}{N^{2}}\right)B^{i}_{a}-\left(h^{ij}-\frac{N^{i}N^{j}}{N^{2}}\right)B^{k}_{a}\\\
&=\frac{1}{N}\left(h^{jk}N^{i}-h^{ij}N^{k}\right)B^{\ast}_{a}+\left(h^{jk}-\frac{N^{j}N^{k}}{N^{2}}\right)\bar{B}^{i}_{a}-\left(h^{ij}-\frac{N^{i}N^{j}}{N^{2}}\right)\bar{B}^{k}_{a}\,.\end{split}$
(72)

## References

  * Nester and Yo (1999) J. M. Nester and H.-J. Yo, Chin. J. Phys. 37, 113 (1999), arXiv:gr-qc/9809049 .
  * Beltrán Jiménez  _et al._ (2018a) J. Beltrán Jiménez, L. Heisenberg, and T. S. Koivisto, JCAP 08, 039 (2018a), arXiv:1803.10185 [gr-qc] .
  * Beltrán Jiménez  _et al._ (2019) J. Beltrán Jiménez, L. Heisenberg, and T. S. Koivisto, Universe 5, 173 (2019), arXiv:1903.06830 [hep-th] .
  * Rünkla and Vilson (2018) M. Rünkla and O. Vilson, Phys. Rev. D 98, 084034 (2018), arXiv:1805.12197 [gr-qc] .
  * Beltrán Jiménez  _et al._ (2020) J. Beltrán Jiménez, L. Heisenberg, T. S. Koivisto, and S. Pekar, Phys. Rev. D 101, 103507 (2020), arXiv:1906.10027 [gr-qc] .
  * Blixt  _et al._ (2023) D. Blixt, A. Golovnev, M.-J. Guzman, and R. Maksyutov, (2023), arXiv:2306.09289 [gr-qc] .
  * Beltrán Jiménez and Koivisto (2022) J. Beltrán Jiménez and T. S. Koivisto, Int. J. Geom. Meth. Mod. Phys. 19, 2250108 (2022), arXiv:2202.01701 [gr-qc] .
  * Adak (2018) M. Adak, Int. J. Geom. Meth. Mod. Phys. 15, 1850198 (2018), arXiv:1809.01385 [gr-qc] .
  * Tomonari and Bahamonde (2023) K. Tomonari and S. Bahamonde, (2023), arXiv:2308.06469 [gr-qc] .
  * Beltrán Jiménez  _et al._ (2018b) J. Beltrán Jiménez, L. Heisenberg, and T. Koivisto, Phys. Rev. D 98, 044048 (2018b), arXiv:1710.03116 [gr-qc] .
  * Quiros (2022) I. Quiros, Phys. Rev. D 105, 104060 (2022), arXiv:2111.05490 [gr-qc] .
  * Maluf (2013) J. W. Maluf, Annalen Phys. 525, 339 (2013), arXiv:1303.3897 [gr-qc] .
  * Aldrovandi and Pereira (2013) R. Aldrovandi and J. G. Pereira,  _Teleparallel Gravity: An Introduction_ (Springer, 2013).
  * Cai  _et al._ (2016) Y.-F. Cai, S. Capozziello, M. De Laurentis, and E. N. Saridakis, Rept. Prog. Phys. 79, 106901 (2016), arXiv:1511.07586 [gr-qc] .
  * Krššák and Saridakis (2016) M. Krššák and E. N. Saridakis, Class. Quant. Grav. 33, 115009 (2016), arXiv:1510.08432 [gr-qc] .
  * Bahamonde  _et al._ (2015) S. Bahamonde, C. G. Böhmer, and M. Wright, Phys. Rev. D 92, 104042 (2015), arXiv:1508.05120 [gr-qc] .
  * Bahamonde  _et al._ (2023) S. Bahamonde, K. F. Dialektopoulos, C. Escamilla-Rivera, G. Farrugia, V. Gakis, M. Hendry, M. Hohmann, J. Levi Said, J. Mifsud, and E. Di Valentino, Rept. Prog. Phys. 86, 026901 (2023), arXiv:2106.13793 [gr-qc] .
  * Heisenberg (2019) L. Heisenberg, Phys. Rept. 796, 1 (2019), arXiv:1807.01725 [gr-qc] .
  * Harko  _et al._ (2018) T. Harko, T. S. Koivisto, F. S. N. Lobo, G. J. Olmo, and D. Rubiera-Garcia, Phys. Rev. D 98, 084043 (2018), arXiv:1806.10437 [gr-qc] .
  * Järv  _et al._ (2018) L. Järv, M. Rünkla, M. Saal, and O. Vilson, Phys. Rev. D 97, 124025 (2018), arXiv:1802.00492 [gr-qc] .
  * Capozziello  _et al._ (2022) S. Capozziello, V. De Falco, and C. Ferrara, Eur. Phys. J. C 82, 865 (2022), arXiv:2208.03011 [gr-qc] .
  * Capozziello  _et al._ (2023) S. Capozziello, V. De Falco, and C. Ferrara, Eur. Phys. J. C 83, 915 (2023), arXiv:2307.13280 [gr-qc] .
  * Hu  _et al._ (2023) K. Hu, T. Paul, and T. Qiu, (2023), arXiv:2308.00647 [hep-th] .
  * Qiu  _et al._ (2019) T. Qiu, K. Tian, and S. Bu, Eur. Phys. J. C 79, 261 (2019), arXiv:1810.04436 [gr-qc] .
  * Nájera and Fajardo (2022) A. Nájera and A. Fajardo, JCAP 03, 020 (2022), arXiv:2111.04205 [gr-qc] .
  * Li and Zhao (2022) M. Li and D. Zhao, Phys. Lett. B 827, 136968 (2022), arXiv:2108.01337 [gr-qc] .
  * Khyllep  _et al._ (2021) W. Khyllep, A. Paliathanasis, and J. Dutta, Phys. Rev. D 103, 103521 (2021), arXiv:2103.08372 [gr-qc] .
  * Casalino  _et al._ (2021) A. Casalino, B. Sanna, L. Sebastiani, and S. Zerbini, Phys. Rev. D 103, 023514 (2021), arXiv:2010.07609 [gr-qc] .
  * Sharma and Sur (2021) M. K. Sharma and S. Sur, (2021), arXiv:2102.01525 [gr-qc] .
  * Sharma and Sur (2022) M. K. Sharma and S. Sur, Int. J. Mod. Phys. D 31, 2250017 (2022), arXiv:2112.08477 [astro-ph.CO] .
  * Sharma and Sur (2023) M. K. Sharma and S. Sur, Phys. Dark Univ. 40, 101192 (2023), arXiv:2112.14017 [astro-ph.CO] .
  * Mandal  _et al._ (2020a) S. Mandal, D. Wang, and P. K. Sahoo, Phys. Rev. D 102, 124029 (2020a), arXiv:2011.00420 [gr-qc] .
  * Mandal  _et al._ (2020b) S. Mandal, P. K. Sahoo, and J. R. L. Santos, Phys. Rev. D 102, 024057 (2020b), arXiv:2008.01563 [gr-qc] .
  * Mandal and Sahoo (2021) S. Mandal and P. K. Sahoo, Phys. Lett. B 823, 136786 (2021), arXiv:2111.10511 [gr-qc] .
  * Arora and Sahoo (2022) S. Arora and P. K. Sahoo, Annalen Phys. 534, 2200233 (2022), arXiv:2206.05110 [gr-qc] .
  * Gadbail  _et al._ (2022) G. N. Gadbail, S. Mandal, and P. K. Sahoo, Phys. Lett. B 835, 137509 (2022), arXiv:2210.09237 [gr-qc] .
  * Paliathanasis (2023a) A. Paliathanasis, (2023a), arXiv:2309.14669 [gr-qc] .
  * Paliathanasis (2023b) A. Paliathanasis, Phys. Dark Univ. 41, 101255 (2023b), arXiv:2304.04219 [gr-qc] .
  * Dimakis  _et al._ (2022) N. Dimakis, A. Paliathanasis, M. Roumeliotis, and T. Christodoulakis, Phys. Rev. D 106, 043509 (2022), arXiv:2205.04680 [gr-qc] .
  * Paliathanasis (2023c) A. Paliathanasis, (2023c), arXiv:2310.16357 [gr-qc] .
  * Wang  _et al._ (2022) W. Wang, H. Chen, and T. Katsuragawa, Phys. Rev. D 105, 024060 (2022), arXiv:2110.13565 [gr-qc] .
  * Lin and Zhai (2021) R.-H. Lin and X.-H. Zhai, Phys. Rev. D 103, 124001 (2021), arXiv:2105.01484 [gr-qc] .
  * D’Ambrosio  _et al._ (2022) F. D’Ambrosio, S. D. B. Fell, L. Heisenberg, and S. Kuhn, Phys. Rev. D 105, 024042 (2022), arXiv:2109.03174 [gr-qc] .
  * Bahamonde  _et al._ (2022) S. Bahamonde, J. Gigante Valcarcel, L. Järv, and J. Lember, JCAP 08, 082 (2022), arXiv:2206.02725 [gr-qc] .
  * Calzá and Sebastiani (2023) M. Calzá and L. Sebastiani, Eur. Phys. J. C 83, 247 (2023), arXiv:2208.13033 [gr-qc] .
  * Hohmann  _et al._ (2018) M. Hohmann, M. Krššák, C. Pfeifer, and U. Ualikhanova, Phys. Rev. D 98, 124004 (2018), arXiv:1807.04580 [gr-qc] .
  * Hohmann  _et al._ (2019) M. Hohmann, C. Pfeifer, J. Levi Said, and U. Ualikhanova, Phys. Rev. D 99, 024009 (2019), arXiv:1808.02894 [gr-qc] .
  * Soudi  _et al._ (2019) I. Soudi, G. Farrugia, V. Gakis, J. Levi Said, and E. N. Saridakis, Phys. Rev. D 100, 044008 (2019), arXiv:1810.08220 [gr-qc] .
  * Abedi and Capozziello (2018) H. Abedi and S. Capozziello, Eur. Phys. J. C 78, 474 (2018), arXiv:1712.05933 [gr-qc] .
  * Hu  _et al._ (2022) K. Hu, T. Katsuragawa, and T. Qiu, Phys. Rev. D 106, 044025 (2022), arXiv:2204.12826 [gr-qc] .
  * D’Ambrosio  _et al._ (2023) F. D’Ambrosio, L. Heisenberg, and S. Zentarra, (2023), arXiv:2308.02250 [gr-qc] .
  * Heisenberg (2023) L. Heisenberg, (2023), arXiv:2309.15958 [gr-qc] .
  * Paliathanasis  _et al._ (2023) A. Paliathanasis, N. Dimakis, and T. Christodoulakis, (2023), arXiv:2308.15207 [gr-qc] .
  * Dimakis  _et al._ (2021) N. Dimakis, A. Paliathanasis, and T. Christodoulakis, Class. Quant. Grav. 38, 225003 (2021), arXiv:2108.01970 [gr-qc] .
  * Zhao (2022) D. Zhao, Eur. Phys. J. C 82, 303 (2022), arXiv:2104.02483 [gr-qc] .
  * de Rham  _et al._ (2011) C. de Rham, G. Gabadadze, and A. J. Tolley, Phys. Rev. Lett. 106, 231101 (2011), arXiv:1011.1232 [hep-th] .
  * de Rham  _et al._ (2012) C. de Rham, G. Gabadadze, and A. J. Tolley, Phys. Lett. B 711, 190 (2012), arXiv:1107.3820 [hep-th] .
  * Hassan and Rosen (2011) S. F. Hassan and R. A. Rosen, JHEP 07, 009 (2011), arXiv:1103.6055 [hep-th] .
  * Hinterbichler (2012) K. Hinterbichler, Rev. Mod. Phys. 84, 671 (2012), arXiv:1105.3735 [hep-th] .
  * Molaee and Shirzad (2018) Z. Molaee and A. Shirzad, Nucl. Phys. B 933, 248 (2018), arXiv:1709.04847 [hep-th] .
  * Becchi  _et al._ (1976) C. Becchi, A. Rouet, and R. Stora, Annals Phys. 98, 287 (1976).
  * Kugo and Ojima (1978) T. Kugo and I. Ojima, Phys. Lett. B 73, 459 (1978).
  * Kugo and Ojima (1979) T. Kugo and I. Ojima, Prog. Theor. Phys. Suppl. 66, 1 (1979).
  * Belenchia  _et al._ (2018) A. Belenchia, M. Letizia, S. Liberati, and E. D. Casola, Rept. Prog. Phys. 81, 036001 (2018), arXiv:1612.07749 [gr-qc] .
  * Langlois and Noui (2016) D. Langlois and K. Noui, JCAP 07, 016 (2016), arXiv:1512.06820 [gr-qc] .
  * Gao  _et al._ (2019) X. Gao, M. Yamaguchi, and D. Yoshida, JCAP 03, 006 (2019), arXiv:1810.07434 [hep-th] .
  * Kimura  _et al._ (2017) R. Kimura, A. Naruko, and D. Yoshida, JCAP 01, 002 (2017), arXiv:1608.07066 [gr-qc] .
  * Crisostomi  _et al._ (2016) M. Crisostomi, K. Koyama, and G. Tasinato, JCAP 04, 044 (2016), arXiv:1602.03119 [hep-th] .
  * Ben Achour  _et al._ (2016) J. Ben Achour, D. Langlois, and K. Noui, Phys. Rev. D 93, 124005 (2016), arXiv:1602.08398 [gr-qc] .
  * Motohashi  _et al._ (2016) H. Motohashi, K. Noui, T. Suyama, M. Yamaguchi, and D. Langlois, JCAP 07, 033 (2016), arXiv:1603.09355 [hep-th] .
  * de Rham (2014) C. de Rham, Living Rev. Rel. 17, 7 (2014), arXiv:1401.4173 [hep-th] .
  * Hassan  _et al._ (2012) S. F. Hassan, R. A. Rosen, and A. Schmidt-May, JHEP 02, 026 (2012), arXiv:1109.3230 [hep-th] .
  * Beltrán Jiménez and Koivisto (2021) J. Beltrán Jiménez and T. S. Koivisto, Universe 7, 143 (2021), arXiv:2104.05566 [gr-qc] .
  * Horndeski (1974) G. W. Horndeski, Int. J. Theor. Phys. 10, 363 (1974).
  * Kobayashi (2019) T. Kobayashi, Rept. Prog. Phys. 82, 086901 (2019), arXiv:1901.07183 [gr-qc] .
  * Gleyzes  _et al._ (2015) J. Gleyzes, D. Langlois, F. Piazza, and F. Vernizzi, Phys. Rev. Lett. 114, 211101 (2015), arXiv:1404.6495 [hep-th] .


# Sieve of Eratosthenes for Bose-Einstein Condensates in Optical Moiré
Lattices

Dmitry Kouznetsov  KU Leuven, Dept. of Physics and Astronomy, Research unit
Quantum Solid-State Physics, B-3001 Leuven, Belgium  imec, Kapeldreef 75,
B-3001 Leuven, Belgium  Pol Van Dorpe  imec, Kapeldreef 75, B-3001 Leuven,
Belgium  KU Leuven, Dept. of Physics and Astronomy, Research unit Quantum
Solid-State Physics, B-3001 Leuven, Belgium  Niels Verellen  imec, Kapeldreef
75, B-3001 Leuven, Belgium

###### Abstract

We catalog known optical moiré lattices and uncover exotic lattice
configurations following a geometric analog of the ancient sieve of
Eratosthenes algorithm for finding prime numbers. Rich dynamics of Bose-
Einstein condensates loaded into these optical lattices is revealed from
numerical simulations of time-of-flight interference patterns. What sets this
method apart is the ability to tune the periodicity of the optical lattices
without changing the wavelength of the laser, yet maintaining the local
potential at the individual lattices sites. In addition, we discuss the
ability to spatially translate the optical lattice through applying a
structured phase only.

Introduction – Optical lattices (OLs) generated by standing wave laser fields
form a potential landscape into which ultracold atoms can be loaded. The
formation of a Bose-Einstein Condensate (BEC) of bosonic atoms under low
temperatures subject to the symmetries of an OL enables the study of
fundamental models of many-body systems Jaksch  _et al._ (1998); Greiner  _et
al._ (2002); Lewenstein  _et al._ (2007); Greiner and Fölling (2008),
facilitates observation of coherent phenomena such as Bloch oscillations Ben
Dahan  _et al._ (1996) or Wannier-Stark ladders Wilkinson  _et al._ (1996),
and serves as a possible implementation of quantum computation Mandel  _et
al._ (2003).

Over the years, a zoo of OLs has been accumulated due to the myriad ways in
which laser beam configurations can be arranged to form the trapping
potentials. Aside from the familiar square or triangular OLs Grynberg and
Robilliard (2001), quasiperiodic potentials recently piqued the interest of
researchers due to the expected hybrid crystalline/amorphous features Spurrier
and Cooper (2018); Corcovilos and Mittal (2019); Viebahn  _et al._ (2019).
Another example is the Kagome lattice Isakov  _et al._ (2006); Jo  _et al._
(2012) which adds geometric frustration to the system and can lead to the
appearance of a flat band. However, these lattices and the dynamics of the
atomic system are typically treated on an individual basis.

The common parameter of periodic OLs is the periodicity $d$, and it is a key
parameter that dictates the dynamics of the ultracold atom system. From moiré
theory Amidror (2009) the definition of periodicity can be extended to quasi-
periodic OLs in terms of dominant frequency components. The periodicity is
usually tuned by changing the wavelength of the laser. However, the atom-
lattice interaction depends on the detuning between the atomic transition
frequency and the laser field frequency. The effect of changing the wavelength
of the laser can be as drastic as changing the behavior at the high intensity
foci of the OL from trapping atoms to repelling them Bloch  _et al._ (2008).
For one-dimensional systems the lattice periodicity can be modified without
changing the wavelength by controlling the relative angle of the two
interfering laser beams Fallani  _et al._ (2005), but in a two-dimensional
system this approach has the adverse effect of also changing the lattice
symmetry. Moreover, these approaches stretch the individual lattice sites,
modifying the local gradient of the potential. This, in turn, changes the
number of atoms that can be loaded per site and thus distorts the engineered
dynamics of the system. Therefore, a framework for designing OLs with varying
periodicity, yet diffraction limited foci is very valuable. As an additional
benefit, such a framework would catalog known OLs and pave the way for
designing more exotic potential landscapes.

Our recently developed integer lattice method Kouznetsov  _et al._ (2020)
describes an algorithmic approach to computing the orientations of laser beams
for generation of OLs with variable periodicity and symmetry. This technique
combines prime number factorization in the complex plane with moiré theory to
compute wavevector components that can be used to generate desired
interference patterns. In this letter, we consider writing the OL potential in
terms of integer factorization over a number field as defined by the integer
lattice method. Our simulations of BECs in such OLs reproduce matter-wave
interaction patterns reported in literature and naturally lead to a
classification of OLs according to the key parameter in the integer lattice
method: the field norm. We demonstrate the utility of the developed
classification by calculating the set of OLs that would allow for
monochromatic tuning of the OL periodicity.

Since the integer lattice method is strongly linked to the moiré effect, the
proposed classification scheme shows promise in other systems with
controllable superlattice ordering, such as twisted multi-layer graphene Cao
_et al._ (2018); Phinney  _et al._ (2021); Gadelha  _et al._ (2021), twisted
van der Waals materials Wu  _et al._ (2018); Wang  _et al._ (2020a); Andersen
_et al._ (2021) and Anderson localization in photonic moiré lattices Wang  _et
al._ (2020b).

Optical lattice potential in the integer lattice method formalism – Before
writing down the governing wavefunction of the system of cold atoms in an OL,
we will integrate the integer lattice method description of coherent lattices
into vector notation of the optical potential. In Ref. Kouznetsov  _et al._
(2020), the set of orientations of the wavevectors of each light beam
component was introduced as

$P(n)=\\{N(\alpha)=n\,\mid\,\alpha\in\mathbb{Z}[\zeta_{m}]\\}\,,$ (1)

where the ring $\mathbb{Z}[\zeta_{m}]=\\{a+b\zeta_{m}\mid a,b\in\mathbb{Z}\\}$
determines the symmetry of the system by the choice of integer $m$ in
$\zeta_{m}=e^{2\pi i/m}$. The corresponding field norm
$N(\alpha)=\alpha\bar{\alpha}=n$ then selects only specific wavevectors with
an integer magnitude $n$. However, the set $P(n)$ in Eq. (1) contains only
complex numbers and therefore, to give physical meaning to these elements we
introduce the set

$K_{n}\equiv\frac{2\pi}{\lambda}\text{vec}\left(\frac{P(n)}{\sqrt{n}}\right)\,,$
(2)

where $\text{vec}:x+yi\mapsto(x,y)$ simply converts the complex number to
vector notation. A planar arrangement of light beams is the simplest
configuration to form two-dimensional patterns and is shown in Fig. 1a. The
way in which the orientations of these beams is determined from the points
$P(n)$ in the complex plane is illustrated in Fig. 1b. Note that the points
$P(n)$ are normalized to have unit length in Eq. (2) and are subsequently
scaled by the wavelength $\lambda$ of the light beams.

Figure 1: Generation of OLs using the integer lattice method. A planar laser
beam configuration (a) is arranged such that the input beam orientations
correspond to the concyclic points in the triangular integer lattice
$\mathbb{Z}[\zeta_{6}]$ (b). The resulting interference pattern in the
overlapping region of the input beams (c) shows a triangular symmetry.

By introducing $K_{n}$, we can interface the integer lattice method with the
description of the general optical potential for interfering laser beams given
by Popov (1987)

$V_{latt}(\mathbf{r})=V_{0}\left|\sum_{j}\mathcal{E}_{j}\boldsymbol{\epsilon}_{j}e^{-i(\mathbf{k}_{j}\cdot\mathbf{r}+\varphi_{j})}\right|^{2}\,,$
(3)

where $\mathbf{r}$ is the position vector in two dimensions, $V_{0}$ is the
overall strength of the potential, $\mathcal{E}_{j}\in[0,1]$ is the relative
intensity of the laser beam, $\boldsymbol{\epsilon}_{j}$ is the polarization
and $\varphi_{j}$ the phase. The wavevectors $\mathbf{k}_{j}$ can now be
replaced by the calculated wavectors in $K_{n}$ as follows:

$V_{latt}^{\prime}(\mathbf{r})=V_{0}\left|\sum_{\mathbf{k}_{j}\in
K_{n}}e^{-i(\mathbf{k}_{j}\cdot\mathbf{r}+\varphi_{j})}\right|^{2}\,,$ (4)

assuming all beams are linearly polarized in the same direction such that
$\boldsymbol{\epsilon}_{j}$ disappears, and have equal intensity, i.e.
$\mathcal{E}_{j}=1$. The phase $\varphi_{j}$ will come into play when
discussing the translation of the OL. For now, we set all the laser beams to
be in phase, i.e. $\varphi_{j}=0$. An example of a triangular OL generated
from $P(7)$ is shown in Fig. 1c. This final form of the optical potential
facilitates the analysis of BECs loaded in OLs in the context of the integer
lattice method.

Matter-wave interference pattern simulation – The momentum distribution of the
BEC holds key information of the system and is most often experimentally
obtained by observing matter-wave diffraction in time-of-flight imaging. We
therefore target simulating matter-wave interference patterns in the following
discussion.

Aside from the OL potential, we add a harmonic confinement to the BEC, such
that the total potential becomes

$V_{tot}(\mathbf{r})=V_{latt}^{\prime}(\mathbf{r})+\frac{1}{2}\omega|\mathbf{r}|^{2}\,,$
(5)

with trapping frequency $\omega$ in units of the recoil frequency
$\omega_{R}=\hbar|\mathbf{k}|^{2}/2m$ with $m$ the atomic mass. Therefore,
$V_{0}$ in Eq. (4) is in units of the recoil energy $E_{R}=\omega_{R}\hbar$.
For such a system, the weakly interacting bosonic gas is known to be well
described by the time dependent Gross-Pitaevskii equation (GPE) Dalfovo  _et
al._ (1999), which is written in dimensionless form as

$i\partialderivative{\Psi(\mathbf{r},t)}{t}=\left(-\frac{1}{2}\mathbf{\nabla}^{2}+V_{\text{tot}}(\mathbf{r})+g|\Psi(\mathbf{r},t)|^{2}\right)\Psi(\mathbf{r},t)\,,$
(6)

with $g$ being the variable interaction strength parameter. We find the ground
state of the system described by Eq. (6) using an imaginary time evolution
with the Fourier split-step operator method Gaunt (2015) with absorbing
boundary conditions implemented using the QuantumOptics.jl framework Krämer
_et al._ (2018).

Figure 2: Integer lattice method classification of OLs analogous to the sieve
of Eratosthenes algorithm. Each point in the diagram corresponds to a
triangular $(m=6)$ OL generated from an integer $n$. The distinction between
lattices is based on lattice periodicity, with first occuring lattices
$\boldsymbol{\mathbb{\Lambda}}_{6}$ marked red and identical lattices marked
gray. The trend lines of the periodicity are added as a visual aid (dashed
lines). Non-inert prime numbers $\mathbb{P}_{6}$ coincide with $n$ of
$\boldsymbol{\mathbb{\Lambda}}_{6}$ (black squares). Higher order moiré OLs
(arrows) are generated from products of prime numbers, e.g. $91=7\cdot 13$.

Sieve of Eratosthenes – Having all components for simulating bosonic gasses in
OLs in place, we turn to a number theory technique to algorithmically analyze
two-dimensional OL configurations. A good place to start is the result from
the integer lattice method linking coherent OLs to moiré superlattices with
periodicity determined by the prime number factorization of an integer $n$ via
Eq. (1) in the complex plane Kouznetsov  _et al._ (2020). The major
distinguising factor between OLs is their symmetry, which is fixed by choosing
$m$. Furthermore, from moiré theory it is known that the dominant spatial
features are determined by the smallest components in momentum space Amidror
(2009), such that for $\mathbf{k}_{i}\in K_{n}$ the pattern periodicity $d$
can be written as

$d=\min|\mathbf{k}_{i}-\mathbf{k}_{j}|^{-1}\,.$ (7)

This quantity will act as the secondary distinguishing factor between the
generated patterns. To show that this distinction is sufficient, we consider
the distribution of the set of prime numbers, which we will denote by
$\mathbb{P}=\\{p\in\mathbb{N}\mid p\text{ is prime}\\}$. An important
observation is that not all primes remain prime in the complex plane. For
example, $5=(2+i)(2-i)$ is no longer prime in $\mathbb{Z}[\zeta_{4}]$. This
has important consequences for the generated coherent lattices.

Suppose a field norm $n\in\mathbb{P}$ remains prime (known as inert primes
Rotman and Cuoco (2013)) in $\mathbb{Z}[\zeta_{m}]$, then $P(n)=\emptyset$,
since these numbers cannot be split in the complex plane. Let $\mathbb{I}_{m}$
be the set of all such inert primes in $\mathbb{Z}[\zeta_{m}]$. We can then
filter out all these prime numbers that have no valid field norm assossiated
with them, and be left with the set

$\mathbb{P}_{m}=\mathbb{P}\backslash\mathbb{I}_{m}\,.$ (8)

The unique factorization theorem states that positive integers greater than
zero can be represented in exactly one way as a product of prime numbers –
essentially describing primes as building blocks of the natural numbers.
Similarly, the factorization of $n\in\mathbb{P}_{m}$ in
$\mathbb{Z}[\zeta_{m}]$ cannot be a composite decomposition of other field
norms due to its primeness. We therefore can leverage the principle of the
sieve of Eratosthenes – an algorithm in which all multiples of a number are
marked iteratively such that only all primes remain Greaves  _et al._ (1997) –
to identify all OLs which appear for the first time, i.e. have the lowest $n$
for a given periodicity $d$. These lattices will be labeled
$\mathbb{\Lambda}_{m}$.

In Fig. 2, each red dot corresponds to one of the first occuring lattices
$\Lambda_{n}\in\mathbb{\Lambda}_{6}$. Lattices generated from different field
norm values $n$ are said to be degenerate if they have the same periodicity
$d$ (see Eq. (7)). These duplicate lattices are iteratively marked grey. For
example, one can easily verify that solving Eq. (1) for $n=1$ and $n=4$
results in the same set of wavevectors due to the normalization factor
$\sqrt{n}$ in Eq. (2), e.g. $K_{1}=K_{4}$. Also, the trend lines of the
periodicity are proportional to $\sqrt{n}$, and are plotted as a visual aid.
From construction, any lattice $\Lambda_{n}$ generated from
$n\in\mathbb{P}_{6}$ will be in $\mathbb{\Lambda}_{6}$. However, some lattices
are constructed from a product of prime numbers (highlighted with arrows in
Fig. 2), and will appear for the first time. For example $\Lambda_{91}$, where
$n=91=13\cdot 7$. In Ref. Kouznetsov  _et al._ (2020) these were identified as
higher order moiré superlattices, since their construction is a superposition
of past lattices. Of course, the exception is $\Lambda_{1}$, since $1$ is not
a prime number. It is important to note that even though $\Lambda_{n}$ have
varying periodicity, the foci are diffraction limited.

Figure 3: Calculated matter-wave interference patterns corresponding to the
momentum distribution of the superfluid in the OLs obtained using the integer
lattice method (see Fig. 2). Each of the momentum distributions (a-f) shows
different dynamics of the system, solely dependent on the choice of $n$. The
dashed triangle (a-c) denotes the decrease of the first Brillouin zone,
corresponding to more dense lattice sites. The parameters used in the
calculations are $\omega=0.08,V_{0}=-6$ and $g=10$.

The analysis of the infinitely many possible OLs is beyond the scope of this
paper. However, we will briefly discuss several key examples for
$\mathbb{\Lambda}_{6}$. In Fig. 3, the matter-wave interference patterns are
plotted that correspond to the momentum distribution of the BEC. First and
foremost, we achieve increased lattice periodicity (closer Bragg peaks in
momentum space), by choosing larger value of $n$. The density of the lattice
sites ranges from dense Fig. 3(a) ($n=1$), to intermediate Fig. 3(b) ($n=7$),
to sparse Fig. 3(c) ($n=13$). Second, contrasting to the regular triangular
OLs, Fig. 3(d) shows additional interference peaks inside the first Brillouin
zone for $n=67$, revealing auxiliary dynamics of the system. These stem from
emerging secondary lattice sites of the OL. Moreover, more exotic OLs which
result in a moiré superlattice with lattice mismatch are shown for $n=109$ in
Fig. 3(e). Finally, in Fig. 3(f), the distribution of the peaks displays a
clear twelvefold rotational symmetry for $n=181$, with the distinctive
structure of a quasicrystal.

Note that the OLs in Fig. 3 all rely on triangular base patterns ($m=6$).
Switching symmetry by fixing a different $m$ will open rich families of OLs
for exploration. For example, setting $m=5$ in Eq. (1) will result in a wide
range of quasi-periodic lattices with tenfold rotational symmetry.

Phase synchronization – The main advantage of relying only on wavevector
orientations to determine the OL symmetry and periodicity is that switching
between different lattices amounts to activating the desired laser beams
($\mathcal{E}_{j}=1$ in Eq. (3)) and deactivating others. This has the benefit
of not relying on changing individual beam orientations or the operating
wavelength to tune the lattice periodicity. Similarly, spatially moving the OL
can be achieved without reorienting the wavevector components, but by applying
structured phase shifts $\varphi_{j}$ (see Eq. (4)) to the input beams.

Figure 4: Phase tuning scheme for spatial translation of the OL. A diagram of
the phase shifts (a) for refocusing of the constructive interference of a
system with three input beams (red dots) from point $A$ to point $B$. The
phase shift can be thought of as a displacement of the plane wavefront;
initial (gray lines) and displaced (black lines). Each phase shift has
magnitude $\varphi_{j}$. Positive (orange bars) or negative (blue bars) phase
shifts (b), when synchronized to follow a cardioid curve (dashed line), enable
the spatial translation of the OL (c).

Displacing the superpattern by shifting each wave component along the
wavevector is a known result from moiré theory Amidror (2009), and can be
readily applied to the beam components. The phase tuning scheme relies on
synchronizing the displacement of the plane wave components such that the
wavefronts maintain the interference pattern at each point in space. This is
illustrated in Fig. 4(a). The structured phase shifts have a magnitude that is
determined by the projection of the wavevectors onto the direction of the
desired displacement. For example, horizontal motion of the moiré superpattern
requires the following phase shifts:

$\varphi_{j}=s\cdot\cos(\arg(k_{j}))\,.$ (9)

The total displacement of the OL can be tuned with the scaling factor $s$.
Visualizing $\varphi_{j}$ in a radial plot unveils the cardioid envelope
(dashed line in Fig. 4(b)). Displacing the OL in arbitrary directions is thus
achieved by orienting the cardioid curve along the corresponding axis.

Tuning the phase of individual beams is overall challenging in an experimental
setting. However, advances in the spatial light modulator technologies have
already shown that tuning a large parameter space is feasible, for example in
arrays of optical tweezers for cold-atom experiments Barredo  _et al._ (2016).
In comparison, the integer lattice method greatly reduces the number of
parameters that need to be tuned to generate and move complex OLs. Moreover,
tuning the phase can be avoided altogether for static OLs, since the
prerequisite in the integer lattice method is that all laser beams are in
phase, which can be readily achieved with a binary amplitude mask.

Conclusion – We describe a design scheme for OLs for ultracold atom research
that is predicted to give rise to rich distributions of particle momenta.
These distributions, characterized by the localized Bragg peaks, are found by
numerically solving the Gross–Pitaevskii equation.

By recognizing that the OL symmetries are intimately linked to prime number
distributions according to the integer lattice method, the possible ground
states of the system are identified by extending the sieve of Eratosthenes
algorithm to the norm of the wavevector orientations. This approach covers the
known (quasi-)periodic OLs and uncovers a wide range of possible OL
configurations, previously unexplored to our knowledge in the scope of cold-
atom physics.

The wavevector orientations of the input laser beams calculated using the
integer lattice method are linked to moiré theory such that tuning the lattice
periodicity can be achieved by simply switching the laser beams on or off,
without the need to change their operating wavelength. Therefore, the method
lends itself to the dynamic study of many-body systems under varying
periodicity of the carrier OL. In addition, the OL can be continuously
displaced by introducing a structured phase across the input beams, further
highlighting the utility of the integer lattice method as the go-to tool for
designing OLs for future cold-atom experiments.

Although finding the most suitable experimental realization is a key future
challenge, the discussed beam arrangements lend themselves to be generated
using standard experimental techniques. As such, it offers an exciting
framework for controlled studies of (quasi-)periodic systems of ultracold
atoms, a major topic of current research. In addition, the two-dimensional
nature of the theoretical framework enables embedding the generated OLs in
planar on-chip next-generation quantum simulation devices.

###### Acknowledgements.

We would like to thank Prof. Jacques Tempere for the helpful discussions. We
acknowledge grant support from FWO Vlaanderen (No. 1SC0321N) to D.K. This work
is part of a project that has received funding from the European Research
Council (ERC) under the European Union’s Horizon 2020 research and innovation
programme (Grant Agreement No. 805222).

## References

  * Jaksch  _et al._ (1998) D. Jaksch, C. Bruder, J. I. Cirac, C. W. Gardiner, and P. Zoller, Physical Review Letters 81, 3108 (1998).
  * Greiner  _et al._ (2002) M. Greiner, O. Mandel, T. Esslinger, T. W. Hänsch, and I. Bloch, Nature 415, 39 (2002).
  * Lewenstein  _et al._ (2007) M. Lewenstein, A. Sanpera, V. Ahufinger, B. Damski, A. Sen(De), and U. Sen, Advances in Physics 56, 243 (2007).
  * Greiner and Fölling (2008) M. Greiner and S. Fölling, Nature 453, 736 (2008).
  * Ben Dahan  _et al._ (1996) M. Ben Dahan, E. Peik, J. Reichel, Y. Castin, and C. Salomon, Physical Review Letters 76, 4508 (1996).
  * Wilkinson  _et al._ (1996) S. R. Wilkinson, C. F. Bharucha, K. W. Madison, Q. Niu, and M. G. Raizen, Physical Review Letters 76, 4512 (1996).
  * Mandel  _et al._ (2003) O. Mandel, M. Greiner, A. Widera, T. Rom, T. W. Hänsch, and I. Bloch, Nature 425, 937 (2003).
  * Grynberg and Robilliard (2001) G. Grynberg and C. Robilliard, Physics Reports 355, 335 (2001).
  * Spurrier and Cooper (2018) S. Spurrier and N. R. Cooper, Physical Review A 97, 043603 (2018).
  * Corcovilos and Mittal (2019) T. A. Corcovilos and J. Mittal, Applied Optics 58, 2256 (2019).
  * Viebahn  _et al._ (2019) K. Viebahn, M. Sbroscia, E. Carter, J.-C. Yu, and U. Schneider, Physical Review Letters 122, 110404 (2019).
  * Isakov  _et al._ (2006) S. V. Isakov, S. Wessel, R. G. Melko, K. Sengupta, and Y. B. Kim, Physical Review Letters 97, 147202 (2006).
  * Jo  _et al._ (2012) G.-B. Jo, J. Guzman, C. K. Thomas, P. Hosur, A. Vishwanath, and D. M. Stamper-Kurn, Physical Review Letters 108, 045305 (2012).
  * Amidror (2009) I. Amidror,  _The Theory of the Moiré Phenomenon: Volume I: Periodic Layers_ , 2nd ed., Computational Imaging and Vision, 38 (Springer, 2009).
  * Bloch  _et al._ (2008) I. Bloch, J. Dalibard, and W. Zwerger, Reviews of Modern Physics 80, 885 (2008).
  * Fallani  _et al._ (2005) L. Fallani, C. Fort, J. E. Lye, and M. Inguscio, Optics Express 13, 4303 (2005).
  * Kouznetsov  _et al._ (2020) D. Kouznetsov, Q. Deng, P. Van Dorpe, and N. Verellen, Physical Review Letters 125, 184101 (2020).
  * Cao  _et al._ (2018) Y. Cao, V. Fatemi, A. Demir, S. Fang, S. L. Tomarken, J. Y. Luo, J. D. Sanchez-Yamagishi, K. Watanabe, T. Taniguchi, E. Kaxiras, R. C. Ashoori, and P. Jarillo-Herrero, Nature 556, 80 (2018).
  * Phinney  _et al._ (2021) I. Y. Phinney, D. A. Bandurin, C. Collignon, I. A. Dmitriev, T. Taniguchi, K. Watanabe, and P. Jarillo-Herrero, Physical Review Letters 127, 056802 (2021).
  * Gadelha  _et al._ (2021) A. C. Gadelha, D. A. A. Ohlberg, C. Rabelo, E. G. S. Neto, T. L. Vasconcelos, J. L. Campos, J. S. Lemos, V. Ornelas, D. Miranda, R. Nadas, F. C. Santana, K. Watanabe, T. Taniguchi, B. van Troeye, M. Lamparski, V. Meunier, V.-H. Nguyen, D. Paszko, J.-C. Charlier, L. C. Campos, L. G. Cançado, G. Medeiros-Ribeiro, and A. Jorio, Nature 590, 405 (2021).
  * Wu  _et al._ (2018) F. Wu, T. Lovorn, E. Tutuc, and A. H. MacDonald, Physical Review Letters 121, 026402 (2018).
  * Wang  _et al._ (2020a) L. Wang, E.-M. Shih, A. Ghiotto, L. Xian, D. A. Rhodes, C. Tan, M. Claassen, D. M. Kennes, Y. Bai, B. Kim, K. Watanabe, T. Taniguchi, X. Zhu, J. Hone, A. Rubio, A. N. Pasupathy, and C. R. Dean, Nature Materials 19, 861 (2020a).
  * Andersen  _et al._ (2021) T. I. Andersen, G. Scuri, A. Sushko, K. De Greve, J. Sung, Y. Zhou, D. S. Wild, R. J. Gelly, H. Heo, D. Bérubé, A. Y. Joe, L. A. Jauregui, K. Watanabe, T. Taniguchi, P. Kim, H. Park, and M. D. Lukin, Nature Materials 20, 480 (2021).
  * Wang  _et al._ (2020b) P. Wang, Y. Zheng, X. Chen, C. Huang, Y. V. Kartashov, L. Torner, V. V. Konotop, and F. Ye, Nature 577, 42 (2020b).
  * Popov (1987) V. N. Popov,  _Functional integrals and collective excitations_ (Cambridge University Press, 1987).
  * Dalfovo  _et al._ (1999) F. Dalfovo, S. Giorgini, L. P. Pitaevskii, and S. Stringari, Reviews of Modern Physics 71, 463 (1999).
  * Gaunt (2015) A. L. Gaunt,  _Degenerate Bose Gases: Tuning Interactions & Geometry_ (University of Cambridge, 2015).
  * Krämer  _et al._ (2018) S. Krämer, D. Plankensteiner, L. Ostermann, and H. Ritsch, Computer Physics Communications 227, 109 (2018).
  * Rotman and Cuoco (2013) J. J. Rotman and A. Cuoco,  _Learning Modern Algebra_ , 1st ed. (American Mathematical Society, 2013).
  * Greaves  _et al._ (1997) G. R. H. Greaves, G. Harman, and M. N. Huxley,  _Sieve Methods, Exponential Sums, and their Applications in Number Theory_ , 1st ed. (Cambridge University Press, 1997).
  * Barredo  _et al._ (2016) D. Barredo, S. de Léséleuc, V. Lienhard, T. Lahaye, and A. Browaeys, Science 354, 1021 (2016).