# Probing axion mediated fermion–fermion interaction by means of entanglement

A. Capolupo<EMAIL_ADDRESS>Dipartimento di Fisica “E.R. Caianiello”
Università di Salerno, and INFN – Gruppo Collegato di Salerno, Via Giovanni
Paolo II, 132, 84084 Fisciano (SA), Italy  G. Lambiase<EMAIL_ADDRESS>Dipartimento di Fisica “E.R. Caianiello” Università di Salerno, and INFN –
Gruppo Collegato di Salerno, Via Giovanni Paolo II, 132, 84084 Fisciano (SA),
Italy  A. Quaranta<EMAIL_ADDRESS>Dipartimento di Fisica “E.R.
Caianiello” Università di Salerno, and INFN – Gruppo Collegato di Salerno, Via
Giovanni Paolo II, 132, 84084 Fisciano (SA), Italy  S. M. Giampaolo
<EMAIL_ADDRESS>Institut Ruđer Bošković, Bijenička cesta 54, 10000 Zagreb,
Croatia

###### Abstract

We propose a new approach in the investigation and detection of axion and
axion–like particles based on the study of the entanglement for two
interacting fermions. We study a system made of two identical fermions with
spin $-1/2$, and we show that fermion–fermion interaction mediated by axions
leads to a non–zero entanglement between the fermions. An entanglement
measurement can reveal the interaction, providing an indirect evidence of the
existence of axions. We discuss how the other interactions affect the
entanglement, and how to isolate the axion contribution. Particular care is
devoted to the analysis of the magnetic dipole–dipole interaction, which turns
out to be, apart from axions, the most relevant contribution to the
entanglement, and we show that it can be suppressed by setting opportunely the
duration of the observation. We also introduce a two–body correlation
function, which could be directly observed in an experiment, and plays the
role of an entanglement witness.

##  I Introduction

The standard model provides a simple explanation for a wide range of phenomena
involving fundamental particles and their interactions. In spite of its
success, it does not represent the definitive theory of elementary particles.
Several phenomena, ranging from particle mixing intro1 ; intro2 ; intro3 ;
intro4 ; intro5  to the quantum features of gravitation Ellis2009 , are beyond
the standard model. Among the shortcomings of the theory is the so–called
_strong CP problem_ in Quantum Chromodynamics (QCD) PDG2019 ; Peccei1977 ;
Peccei2008 . The QCD Lagrangian features a gluon–gluon interaction term that
in principle allows for an arbitrary violation of CP symmetry, whereas no such
violation is observed in strong processes PDG2019 ; Peccei1977 ; Peccei2008 ;
Wilczek1978 ; Weinberg1978 . As a direct consequence of the CP violating term,
one would expect a relatively large neutron electric dipole moment, which the
recent experiments constrain below $3\times 10^{-26}e\,cm$ (see Pendlebury2015
and, for a more recent result Abel2020 ). Similar bounds on the electric
dipole moment of atoms and molecules Graner2016 ; Flambaum2019  pose a further
constraint on the magnitude of the CP violating term. To remedy this
inconsistency, R. Peccei and H. Quinn introduced a new global symmetry
$U_{PQ}(1)$ (called Peccei–Quinn symmetry) that is spontaneously broken
Peccei1977 ; Peccei2008 . As shown by Frank Wilczek and Steven Weinberg
Wilczek1978 ; Weinberg1978 , this results in a new particle, named the axion,
which is the pseudo–Nambu–Goldstone boson of the broken $U_{PQ}(1)$ symmetry
Raffelt2007 .

The scale at which the symmetry breaking occurs $f_{A}$, known as the axion
decay constant, determines, according to the model considered, both the axion
mass $m_{a}\propto\frac{1}{f_{A}}$ and the effective couplings with the
standard matter. The original proposal considered $f_{A}$ close to the
electroweak scale Peccei1977 ; Wilczek1978 ; Weinberg1978 , a hypotesis that
was later ruled out by the experiments. Soon alternative axion models were
devised, notably the KSVZ model KSVZ , featuring heavy quarks carrying a
Peccei–Quinn (PQ) charge, and the DFSZ model DSFZ , in which the ordinary
quarks and additional Higgs doublets carry the PQ charge. They provide a
reference for two large classes of axion models (_hadronic_ and _GUT_ axions).
Today the axion decay costant is estimated to be very large $f_{A}>10^{9}GeV$,
so that the QCD axions (both hadronic and GUT) have to be very light
$m_{a}\sim 10^{-6}-10^{-2}eV$ and very weakly interacting PDG2019 . These
aspects make axions a natural candidate for dark matter.

Moreover, motivated by the search for dark matter components, a variety of
Axion–like–particles (ALPs) has been introduced. They can deviate
significantly from the original Peccei–Quinn proposal, and are not necessarily
tied to the solution of the strong CP problem, but share the nature of axions
and are weakly interacting with the standard matter. In these models the
relation between the coupling constant and the mass of the ALP can differ from
the direct proportionality that characterizes the PQ axions. They range from
masses $m_{a}\sim 10^{-22}$ eV, characteristic of the ultra–light axions
Kim2016 ; Demartino2017 , up to masses of $1$ TeV for the heavy GUT axions
Rubakov1997 . Axions and ALPs are, to date, one of the most credible
explanations for dark matter Marsch2016 ; Lambiase2018 ; Auriol2019 ;
Houston2018 .

Driven by large theoretical interest, several experiments were designed to
prove the existence of ALPs. Perhaps unsurprisingly, considered their
extremely small interaction rates, the experimental search for ALPs has proven
to be very challenging. Many experiments take advantage of the axion–photon
coupling. Among them, searches for polarization anomalies in the light
propagating through a magnetic field (PVLAS) PVLAS , “light shining through a
wall” experiments (OSQAR, ALPS) OSQAR ; OSQAR2 ; ALPS , detection via the
Primakoff effect (CAST) CAST , and haloscope experiments (ADMX) ADMX . More
recently other approaches based on geometric phases Capolupo2015  and QFT
effects in the axion–photon mixing Capolupo2019  were also suggested.
Astrophysical observations and terrestrial experiments over the decades have
restricted the allowed regions in parameter space, and further constraints
might come from the analysis of the axion–nucleon and axion–lepton
interactions, as suggested for instance in Bezerra2016 . Experiments based on
the axion–fermion interaction have also been proposed (QUAX) QUAX . Despite
this, no evidence for the existence of the ALPs has been found up to now.

In this paper, we propose a different approach to the detection of axions,
based on the study of the entanglement arising between two fermions
interacting via axion exchange. In recent years, theoretical and experimental
analysis of the entanglement properties have found application in the most
disparate fields, from quantum biology Arndt2009  to statistical physics
Vidal2003 ; Amico2008 , and also as a tool to gain insight on fundamental
interactions, like gravity Bose2017 ; Marletto2017 ; Simonov2019 ;
Marshman2019 . These applications stem from the fact that the emergence of
entanglement between two (or more) physical objects is strictly connected to
the presence of a quantum interaction between them acting as a quantum channel
Giampaolo2019 ; Bennett1999 . Consequently, an analysis of the entanglement
properties provides information about the interaction itself.

We focus on the axion–mediated fermion–fermion interaction, which assumes the
familiar form of a Yukawa coupling between the pseudoscalar field and two
fermionic fields. In the non–relativistic limit, the axion–induced interaction
reduces to an effective two–body potential Moody1984 ; Daido2017 , that acts
as a source of entanglement for the two fermions. An entanglement measurement
might then provide evidence for the pseudoscalar interaction, and thus for the
existence of axions/ALPs. Obviously, the two fermions interact with each other
in many other ways, i.e. gravitationally, magnetically, etc. All these
interactions are potentially a source of entanglement. Hence, one of the main
goals of this work is to show how it is possible to extract the axions–induced
entanglement contribution from the others. We choose to quantify the
entanglement through the 2–Renyi entropy Horodecki2002 ; Plenio2007 ;
Nielsen2011 ; Giampaolo2013  that has the advantage to be, in several systems,
directly connected to experimentally accessible quantities Bovino2006 ;
Brydges2019 ; Lesche04 ; Islam15 ; Abanin12 . Since a direct measurement of
the entanglement entropy is challenging in many cases, here we individuate a
two–body spin correlation function that plays the role of an entanglement
witness, and can be more easily accessed.

An entanglement witness is a physical quantity that is strictly related to the
family of states under analysis, with the property of vanishing simultaneously
with the entanglement. The detection of a non–zero value of the entanglement
witness therefore implies that the entanglement is non–zero. Specialized to
our case, under suitable conditions, a non–vanishing witness would signal the
presence of axions, and then provide an indirect evidence for their existence.
For masses in the range of $(10^{-3}-1)$ eV and coupling constants close to
the actual constraints, that are compatible with some ALPs models, the witness
is significantly different from zero.

The paper is organized as follows. We first recap the axion–fermion
pseudoscalar interaction and the corresponding two–fermion potential in the
non–relativistic limit (Sec II). From the knowledge of the potential we
compute the time–dependent entanglement between the two fermions, and we
individuate a specific time at which, in absence of axions, the entanglement
must vanish. If in correspondence with such time the entanglement is different
from zero, the presence of axions is detected (Sec III). Soon after we
introduce an entanglement witness and we present a numerical analysis (Sec
IV), and finally we draw our conclusions (Sec V).

##  II Fermion–Fermion interaction induced by Axions

Let us start by recalling the main features of the axion–mediated
fermion–fermion interaction. The coupling of axions with fermions is described
by a Yukawa pseudoscalar vertex Moody1984 ; Daido2017 . If $\phi$ is the axion
field and $\psi_{1},\psi_{2}$ are the fermion fields, the interaction term
reads

$\mathcal{L}_{INT}=-\sum_{j=1,2}ig_{pj}\phi\bar{\psi}_{j}\gamma_{5}\psi_{j}$
(1)

where $\gamma_{5}$ is the product of Dirac matrices
$i\gamma^{0}\gamma^{1}\gamma^{2}\gamma^{3}$ and $g_{pj}$ are the effective
axion–fermion coupling constants, which depend critically on the fermions
considered and the underlying axion (or ALP) model. Since the couplings are
expected to be small, i.e. $g_{pj}\ll 1$, the scattering amplitudes can be
well approximated by the leading order in the perturbative expansion. For the
scattering
$\psi_{1}(p_{1})\psi_{2}(p_{2})\rightarrow\psi_{1}(p^{{}^{\prime}}_{1})\psi_{2}(p^{{}^{\prime}}_{2})$
we have

$\imath\mathcal{A}\\!=\\!\bar{u}_{1}^{s^{\prime}_{1}}(p^{{}^{\prime}}_{1})g_{p1}\gamma_{5}u_{1}^{s_{1}}(p_{1})\frac{i}{q^{2}\\!-\\!m^{2}}\bar{u}_{1}^{s^{\prime}_{2}}(p^{{}^{\prime}}_{2})g_{p2}\gamma_{5}u_{2}^{s_{2}}(p_{2})\\!\\!$
(2)

where the pseudoscalar free propagator with momentum
$q=p_{1}^{{}^{\prime}}-p_{1}=p_{2}-p_{2}^{{}^{\prime}}$ appears, and $m$ is
the axion mass. Here $u_{i}^{s_{i}}(p_{i})$ are the solutions of the free
Dirac equation in momentum space, for the i-th fermion with momentum $p_{i}$
and spin projection $s_{i}$:

$u^{s_{i}}_{i}(p_{i})=\sqrt{\frac{\omega_{i}+M_{i}}{2\omega_{i}}}\begin{pmatrix}\phi_{s_{i}}\\\
\frac{\boldsymbol{\sigma}\cdot\
\boldsymbol{p}_{i}}{\omega_{i}+M_{i}}\phi_{s_{i}}\end{pmatrix}$ (3)

with $M_{i}$ mass of the $i$-th fermion,
$\omega_{i}=\sqrt{\boldsymbol{p}_{i}^{2}+M_{i}^{2}}$ energy of the $i$-th
fermion and $\phi_{s_{i}}$ normalized two–component spinors. In the
non–relativistic limit $\omega_{i}\approx M_{i}$, these become

$u^{s_{i}}_{i}(p_{i})\approx\begin{pmatrix}\phi_{s_{i}}\\\
\frac{\boldsymbol{\sigma}\cdot\
\boldsymbol{p}_{i}}{2M_{i}}\phi_{s_{i}}\end{pmatrix}\ .$ (4)

Inserting equations (4) in equation (2), one obtains the scattering amplitude
for non–relativistic fermions

$\mathcal{A}\approx\frac{g_{p_{1}}g_{p_{2}}}{\boldsymbol{q}^{2}+m^{2}}\frac{\left(\phi_{s^{\prime}_{1}}^{\dagger}\left(\boldsymbol{\sigma}\cdot\boldsymbol{q}\right)\phi_{s_{1}}\right)\left(\phi_{s^{\prime}_{2}}^{\dagger}\left(\boldsymbol{\sigma}\cdot\boldsymbol{q}\right)\phi_{s_{2}}\right)}{4M_{1}M_{2}}\
.$ (5)

By Fourier transforming the amplitude (5) into real space, one finds the
two–body potential due to axion exchange Daido2017

$\displaystyle V(\boldsymbol{r})\\!$ $\displaystyle\\!=\\!$
$\displaystyle\\!-\frac{g_{p_{1}}g_{p_{2}}e^{-mr}}{16\pi
M_{1}M_{2}}\left[{\boldsymbol{\sigma}}_{1}\cdot{\boldsymbol{\sigma}}_{2}\left(\frac{m}{r^{2}}+\frac{1}{r^{3}}+\frac{4}{3}\pi\delta^{3}({\bf{r}})\right)\right.$
(6)
$\displaystyle\;\;\;\;\;\;\;\;\;\left.-\left({\boldsymbol{\sigma}}_{1}\cdot\boldsymbol{\hat{r}}\right)\left({\boldsymbol{\sigma}}_{2}\cdot\boldsymbol{\hat{r}}\right)\left(\frac{m^{2}}{r}+\frac{3m}{r^{2}}+\frac{3}{r^{3}}\right)\right]\
,$

where $r$ ($\boldsymbol{\hat{r}}$) stands for the modulus (the unit vector) of
the relative distance between the fermions and $\delta^{3}({\bf{r}})$ is the
Dirac delta, while $\boldsymbol{\sigma}_{i}$ is the three–dimensional vector
of Pauli operators defined on the $i$–th fermion. For two identical
non–relativistic fermions one has $M_{1}=M_{2}=M$ and
$g_{p_{1}}=g_{p_{2}}=g_{p}$, yielding the interaction Hamiltonian

$\displaystyle H_{p}\\!$ $\displaystyle\\!=\\!$
$\displaystyle\\!-\frac{g_{p}^{2}e^{-mr}}{16\pi
M^{2}}\left[{\boldsymbol{\sigma}}_{1}\cdot{\boldsymbol{\sigma}}_{2}\left(\frac{m}{r^{2}}+\frac{1}{r^{3}}+\frac{4}{3}\pi\delta^{3}({\bf{r}})\right)\right.$
(7)
$\displaystyle\;\;\;\;\;\;\;\;\;\left.-\left({\boldsymbol{\sigma}}_{1}\cdot\boldsymbol{\hat{r}}\right)\left({\boldsymbol{\sigma}}_{2}\cdot\boldsymbol{\hat{r}}\right)\left(\frac{m^{2}}{r}+\frac{3m}{r^{2}}+\frac{3}{r^{3}}\right)\right].$

In the following we will always consider $r$ large enough that the contact
term proportional to $\delta^{3}(\boldsymbol{r})$ can be neglected. Assuming,
in addition, that $\boldsymbol{\hat{r}}$ coincides with the $z$–direction, we
obtain

$\displaystyle H_{p}\\!$ $\displaystyle\\!=-\\!$
$\displaystyle\\!\frac{g_{p}^{2}e^{-mr}}{16\pi
M^{2}r^{3}}\left[m^{2}r^{2}\,\sigma^{z}_{1}\sigma^{z}_{2}+\left(mr+1\right)\mathcal{O}\right]$
(8)

where the operator $\mathcal{O}$ is defined as
$\mathcal{O}=2\sigma^{z}_{1}\sigma^{z}_{2}-\sigma^{x}_{1}\sigma^{x}_{2}-\sigma^{y}_{1}\sigma^{y}_{2}$.
The equation (8) is the interaction, due to axion exchange, between two
identical non–relativistic fermions, and represents the starting point for our
analysis.

##  III Dynamics of entanglement

Since we wish to analyze the entanglement between two fermions due to the
axion–mediated interaction (8), we first need to establish whether the latter
can induce entanglement in a fermionic system. More precisely, considering two
fermions interacting through the Hamiltonian $H_{p}$ of eq. (8) and initially
prepared in a fully separable (non–entangled) state, we need to determine
whether their state develops a non-vanishing entanglement under the action of
$H_{p}$. There exist precise conditions that a Hamiltonian has to fulfill in
order to produce entanglement: 1) a non completely degenerate spectrum; 2) the
impossibility to be reduced to the sum of local terms acting separately on
every single object Giampaolo2019 ; Simonov2019 ; Bennett1999 . As it is easy
to check, the Hamiltonian in eq. (8) fulfills both the requirements, and shall
induce entanglement on the two fermion system. This fact can be used to reveal
axions and ALPs through the analysis of the entanglement properties of the two
fermion state.

Therefore we focus on a system of two identical spin–$\frac{1}{2}$ fermions
(for instance electrons or neutrons) and study the time evolution of its
entanglement properties. We write the state of the system as

$\Psi(\boldsymbol{r}_{1},\boldsymbol{r}_{2},s_{1},s_{2};t)=R(\boldsymbol{R},\boldsymbol{r};t)\psi(s_{1},s_{2};t)$
(9)

where the spatial wave–function $R$ depends on the center of mass position
$\boldsymbol{R}$ and the relative position of the two fermions
$\boldsymbol{r}=\boldsymbol{r}_{1}-\boldsymbol{r}_{2}$, while the spin
wave–function $\psi(s_{1},s_{2};t)$ is a state vector in the product space
$H^{spin}_{1}\otimes H^{spin}_{2}$ of the spin Hilbert spaces associated to
the two particles. In order to simplify the analysis, we assume, as a first
approximation, that the spatial wave–function $R$ is sharply peaked at a given
value of the distance $r=|\boldsymbol{r}_{1}-\boldsymbol{r}_{2}|$, and so
remains during the time interval of interest. We then consider the distance
$r$ in eq. (8) as a parameter, and the Hamiltonian $H_{p}$ as an operator
acting on the spin state alone. Of course, the full wave–function $\Psi$ must
be antisymmetric under particle exchange. Since we shall consider symmetric
spin states $\psi(s_{1},s_{2};t)$, the spatial wave–function
$R(\boldsymbol{R},\boldsymbol{r};t)$ must be antysimmetric.

At $t=0$ we assume that the spin state of the whole system is fully separable,
i.e. that it can be written as the tensor product of two states each defined
on a single fermion. In other words we have that, at $t=0$, the state of the
system can be written as

$\displaystyle\ket{\psi(0)}$ $\displaystyle=$
$\displaystyle\ket{\varphi(0)}_{1}\otimes\ket{\varphi(0)}_{2}$ (10)
$\displaystyle\ket{\varphi(0)}_{i}$ $\displaystyle=$
$\displaystyle\cos(\theta)\ket{\uparrow}_{i}+e^{\imath\phi}\sin(\theta)\ket{\downarrow}_{i}$

where we have dropped the spin arguments $s_{1},s_{2}$ and have switched to a
more convenient Dirac notation. Here $\ket{\uparrow}_{i}$ and
$\ket{\downarrow}_{i}$ denote the eigenstates of the magnetic moment along the
direction joining the two fermions ($z$-direction). Being the state of Eq.
(10) separable, the entanglement vanishes for $t=0$.

If the two fermions would interact with each other only through axions, a
non–vanishing entanglement would directly signal their presence, but this is
not the case. The two fermions generally interact with each other through
several channels, and any of these is potentially a source of entanglement,
depending on the distance $r$ and on the states that one considers.

As we are interested in highlighting the entanglement due to axions, these
interactions produce an unwanted contribution that has to be minimized in
order to successfully reveal the axion–induced entanglement. Then we must
devise a setting in which these additional sources of entanglement are
suppressed or significantly reduced.

At first, let us deal with the weak and the strong nuclear interactions. They
are relevant only for very small distances and, hence, assuming, in our setup,
a relative distance large enough ($r>10^{-12}m$), we can neglect them
altogether. Secondly, let us consider the gravitational and the electrostatic
(if the fermions have a non–vanishing electric charge) interactions. Their
action on the evolution of the spin state cannot induce entanglement and
amounts to a global phase factor. Indeed, as we have previously said, a
fundamental requirement for the interaction to induce entanglement in a
system, is that the spectrum of the associated Hamiltonian is not fully
degenerate. Particles such as electrons, protons, or neutrons are
characterized by precise values of charge and mass that do not depend on their
spin states. Therefore, any state depending solely on the spin of the
particles, like the one we are considering, will react to gravitational and
electrostatic interactions in the same way and to not generate entanglement.

Another unwanted entanglement source is the dipole–dipole magnetic
interaction, whose Hamiltonian reads

$\displaystyle H_{\mu}=-\frac{1}{4\pi
r^{3}}\frac{g^{2}{q_{e}}^{2}}{16\,M^{2}}\mathcal{O}\,,$ (11)

where $g$ is the g–factor and $q_{e}$ is the charge of the electron. As done
in Eq. (8), we have dropped a contact term proportional to
$\delta^{3}(\boldsymbol{r})$ in Eq. (11), assuming that $r$ is large enough.
Being associated with the exchange of massless photons, this interaction is
not confined at short–range and it is sensible to the different spin states in
a way that cannot be reduced to the action of local operators. Hence its
contribution to the entanglement is different from zero. However, in the
following we will show that, by properly setting the time interval of the
entanglement measurement, the contribution to entanglement due to the
dipole–dipole magnetic interaction can be neglected and the entanglement
reduces to the axion contribution alone.

Starting from the initial state $\ket{\psi(0)}$, the state at $t>0$ can be
obtained as $\ket{\psi(t)}=U(t)\ket{\psi(0)}$, where the time evolution
operator $U(t)$ is unitary since we consider our system to be closed (we
neglect any other interaction with the surrounding world). $U(t)$ can be
written as $U(t)=\exp[-\imath t(H_{T})]$, where the total Hamiltonian $H_{T}$
is the sum of the magnenetic $(H_{\mu})$ and the axion term $(H_{p})$, i.e.
$H_{T}=H_{p}+H_{\mu}$, and reads

$\displaystyle\\!H_{T}\\!$ $\displaystyle\\!=\\!$
$\displaystyle\\!-\frac{A}{r^{3}}\left[\mathcal{O}+Be^{-mr}\left(m^{2}r^{2}\,\sigma^{z}_{1}\sigma^{z}_{2}+\mathcal{O}\left(mr+1\right)\right)\right].\,\,\,\,$
(12)

In eq. (12) the parameter $A=\frac{g^{2}{q_{e}}^{2}}{64\pi M^{2}}$ is the
strength of the magnetic interaction while
$B=\frac{4g_{p}^{2}}{g^{2}q_{e}^{2}}=\frac{4g_{p}^{2}}{\alpha g^{2}}$
quantifies the relative weight of the axion interaction and $\alpha$ denotes
the fine structure constant.

Since the operator $U(t)$ is unitary for any time $t\geq 0$, the state
$\ket{\psi(t)}$ remains a pure state, although in general, differently from
$\ket{\psi(0)}$, it is entangled. The amount of entanglement between the two
fermions in $\ket{\psi(t)}$ can be quantified using different measures. In the
present work we consider the $2$–Renyi entropy Nielsen2011 ; Horodecki2002 ;
Plenio2007 ; Giampaolo2013 , that has the advantage, with respect to the other
entropy–based entanglement measures to be associated, at least in some
experimental devices, to experimentally accessible quantities Lesche04 ;
Islam15 ; Abanin12 ; Bovino2006 ; Brydges2019 . The $2$–Renyi entropy is
defined as $S_{2}=-\ln(\mathcal{P}(\rho_{i}(t)))$, where
$\rho_{i}(t)=\mathrm{Tr_{j\neq i}}(\ket{\psi(t)}\bra{\psi(t)})$ is the reduced
density matrix obtained projecting $\ket{\psi(t)}$ on the Hilbert space
defined on one of the two femions and
$\mathcal{P}(\rho_{i}(t))=\mathrm{Tr}\rho^{2}_{i}(t)$ is the purity of
$\rho_{i}(t)$. In our system the $2$-Renyi entropy reads

$\displaystyle S_{2}(t)\\!$ $\displaystyle\\!=\\!$
$\displaystyle\\!-\\!\ln\left[\\!1\\!-\frac{\sin^{4}(2\theta)}{2}\sin^{2}(\Gamma
t)\right]$ (13)

where
$\Gamma=\frac{6A}{r^{3}}\left[1+\frac{B}{3}e^{-rm}\left(3+3rm+r^{2}m^{2}\right)\right]$.
As we can see from the expression of $\Gamma$, the entanglement derives from
both the dipole–dipole magnetic interaction (due to $\frac{6A}{r^{3}}$) and
from the presence of the axions (due to the term
$\frac{2A}{r^{3}}\left[Be^{-rm}\left(3+3rm+r^{2}m^{2}\right)\right]$ ). Given
the form of the $2$–Renyi entropy (13), there exist certain times at which the
entanglement is only due to the axion–mediated interaction. Indeed, by setting
$t=nt^{*}$, with $n$ a positive integer and $t^{*}=\frac{\pi r^{3}}{6A}$, we
obtain:

$\displaystyle S_{2}(nt^{*})$ $\displaystyle=$
$\displaystyle-\ln\Bigg{\\{}1-\frac{\sin^{4}(2\theta)}{2}\times$
$\displaystyle\times$
$\displaystyle\sin^{2}\left[n\pi\left(1+\frac{B}{3}e^{-mr}(3+3mr+m^{2}r^{2})\right)\right]\Bigg{\\}}\,$
$\displaystyle\simeq$
$\displaystyle\frac{\sin^{4}(2\theta)}{2}n^{2}\pi^{2}B^{2}e^{-2mr}\\!\left(1+mr+\frac{m^{2}r^{2}}{3}\right)^{2}.$

In equation (III), the dependence on the dipole–dipole magnetic interaction
has disappeared. Recalling that $B\propto g_{p}^{2}$, we can see at once that
if there is no axion–mediated interaction ($g_{p}\rightarrow 0$), the
$2$-Renyi entropy, and thus the entanglement, vanishes. More precisely, as
$g_{p}\ll 1$, the last line of eq. (III) shows that the entropy is
proportional to $B^{2}$ and then to $g_{p}^{4}$. In the case in which a
non–zero entanglement is detected, in correspondence with the times
$t=nt^{*}$, one can conclude that the former is a consequence of the
axion–induced interaction alone. This would constitute an indirect proof of
the existence of axions.

For a numerical analysis of the $2$-Renyi entropy (III), we focus on ALPs in
the mass range $10^{-3}-1$ eV. These have been considered in the refs.
Bezerra2014 ; Bimonte2016 ; Chen2016 ; Klimchitskaya2015 ; Klimchitskaya2017 ,
where experimental constraints have been obtained on the coupling constant
$g_{p}$ as a function of the ALP mass, both for protons and neutrons. In Fig.
1 we consider two neutrons at a distance $r=1\,nm$ apart, and initial state
given by Eq. (10) with $\theta=\pi/4$ and $\phi=0$. We set, for each value of
the mass in the range $10^{-3}-1$ eV, the coupling constant $g_{p}$ equal to
the threshold value $g_{threshold}$ obtained from the experimental analyses in
the refs. Bezerra2014 ; Bimonte2016 ; Chen2016 ; Klimchitskaya2015 ;
Klimchitskaya2017 . In particular: for the black dot–dashed line, we set
$g_{p}=g_{CP}$, where $g_{CP}$ is the threshold from effective Casimir
pressure measurements Bezerra2014 , and sample values are $g_{CP}=0.0327$ for
$m=10^{-3}eV$, $g_{CP}=0.0348$ for $m=0.05eV$, $g_{CP}=0.0674$ for $m=1eV$;
for the red solid line we set $g_{p}=g_{CF}$, where $g_{CF}$ is the threshold
from measurements of the difference of Casimir forces Klimchitskaya2017 , and
sample values are $g_{CF}=0.007$ for $m=10^{-3}eV$, $g_{CF}=0.012$ for
$m=0.05eV$, $g_{CF}=0.066$ for $m=1eV$; for the blue dashed line we set
$g_{p}=g_{IE}$, where $g_{IE}$ is the threshold from isoelectronic experiments
Klimchitskaya2015 , and sample values are $g_{IE}=0.0036$ for $m=10^{-3}eV$,
$g_{IE}=0.006$ for $m=0.05eV$, $g_{IE}=0.07$ for $m=1eV$. For the three cases,
the $2$-Renyi entropy $S_{2}(t^{*})$ at $t=t^{*}$ is shown in the upper panel
of Fig. 1. All the lines depicted in the figure show a very similar behavior.
At constant distance, larger axion masses imply a stronger damping by the
Yukawa factor $e^{-mr}$. Increasing the mass one eventually arrives at a value
for which the Yukawa damping starts to be relevant. From this point on, a
further increment of the mass of the axions will suppress the entanglement
exponentially.

Figure 1:  (color online) Plot of the $2$–Renyi entropies (upper panel) and of
the Witnesses $C_{yz}$ (lower panel) at $t=t^{*}$, i.e. for $n=1$, between two
neutrons as function of the axions mass $m$ at distance $r=1\,nm$. In the
plots we have considered two neutrons at a distance $r=1\,nm$ apart, and
initial state given by Eq. (10) with $\theta=\pi/4$ and $\phi=0$. We set, for
each value of the mass in the range $10^{-3}-1$ eV, the coupling constant
$g_{p}$ equal to the threshold value $g_{threshold}$ obtained from the
experimental analyses in the refs. Bezerra2014 ; Bimonte2016 ; Chen2016 ;
Klimchitskaya2015 ; Klimchitskaya2017 . In particular: for the black
dot–dashed line we set $g_{p}=g_{CP}$, where $g_{CP}$ is the threshold from
effective Casimir pressure measurements Bezerra2014 , and sample values are
$g_{CP}=0.0327$ for $m=10^{-3}eV$, $g_{CP}=0.0348$ for $m=0.05eV$,
$g_{CP}=0.0674$ for $m=1eV$; for the red solid line we set $g_{p}=g_{CF}$,
where $g_{CF}$ is the threshold from measurements of the difference of Casimir
forces Klimchitskaya2017 , and sample values are $g_{CF}=0.007$ for
$m=10^{-3}eV$, $g_{CF}=0.012$ for $m=0.05eV$, $g_{CF}=0.066$ for $m=1eV$; for
the blue dashed line we set $g_{p}=g_{IE}$, where $g_{IE}$ is the threshold
from isoelectronic experiments Klimchitskaya2015 , and sample values are
$g_{IE}=0.0036$ for $m=10^{-3}eV$, $g_{IE}=0.006$ for $m=0.05eV$,
$g_{IE}=0.07$ for $m=1eV$.

To our knowledge, for ALPs in the mass range $10^{-3}-1eV$, the values
reported in refs. Bezerra2014 ; Bimonte2016 ; Chen2016 ; Klimchitskaya2015 ;
Klimchitskaya2017  represent the strongest model–independent constraints on
axion–nucleon interactions from laboratory experiments. For this class of
ALPs, as the figure 1 shows, the $2$-Renyi entropy is significantly different
from zero, so that the laboratory constraints might be strengthened by several
orders of magnitude from entanglement measurements. Of course, the so obtained
constraints would be model–independent, since no specific axion or ALP model
has been assumed. On the other hand, for QCD axion models there exist several
constraints from astrophysical sources, primarily from supernovae Chang2018 ,
neutron star cooling Hamaguchi2018  and Black Hole superradiance
Arvanitaki2015 . The indicative bound set by supernovae $g_{aNN}^{2}\sim
10^{-19}$ renders our method unviable for QCD axions, at least for present day
technologies. For general ALPs, since mass and coupling constants are
essentially unrelated, and the latter can in principle assume any value, our
approach can strengthen the current laboratory constraints, with the only
limitations presented by the experimental sensitivities and coherence time. It
is worth to note that the idea, here presented, to use entanglement to test
theories of fundamental physics is in line with several recent works, see for
example Ref. Marletto2017 ; Bose2017 ; Hsu2016  that suggest exploiting the
entanglement as a probe for the quantum nature of the gravity. In fact, in
these papers, the spatial wave function has a non trivial evolution, and
consequently, gravitational interactions can induce entanglement.

##  IV The entanglement witness

In the upper panel of Fig. (1), we plot the entanglement entropy between two
neutrons as a function of the axions mass, at distance $r=1nm$ and for ALPs
with coupling constants $g_{p}$ constrained by the analysis of Refs.
Bezerra2014 , Klimchitskaya2017  and Klimchitskaya2015 . Notice that a direct
measurement of the 2–Renyi entropy is in general not easy to accomplish. To
overcome this difficulty we can make use of an entanglement witness, which is
a quantity strictly related to the family of states and to the dynamics of the
system under analysis, whose value is able to signal the presence of
entanglement. Let us set to zero the phase $\phi$ of the initial state in eq.
(10). In this case, the entanglement witness can be identified with the
two–spins correlation functions

$\displaystyle
C_{yz}=\bra{\psi(t)}\sigma_{1}^{y}\sigma_{2}^{z}\ket{\psi(t)}\equiv\bra{\psi(t)}\sigma_{1}^{z}\sigma_{2}^{y}\ket{\psi(t)}.$
(15)

For any choice of $\phi$, one can find a different correlation function
playing the same role. It is straightforward to show that, for any time
greater than zero, the correlation function $C_{yz}$ is given by
$C_{yz}(t)=-\sin(2\theta)\sin\left(\Gamma t\right)$. This function, for any
$\theta\neq\frac{k\pi}{2}$ with $k$ integer, vanishes only when the
entanglement is zero. For $t=nt^{*}=\frac{n\pi r^{3}}{6A}$, the correlation
function $C_{yz}$ depends only on the interaction between axions and fermions
and reduces to

$\displaystyle\\!C_{yz}(nt^{*})\\!\\!$ $\displaystyle\\!=\\!$
$\displaystyle\\!\\!-\\!\sin(2\theta)\\!\sin\\!\left[\\!n\pi\\!\left(\\!1+\frac{B}{3}e^{-\\!mr}(3\\!+\\!3mr\\!+\\!m^{2}r^{2})\\!\right)\\!\right]\;$
(16) $\displaystyle\\!\simeq\\!$
$\displaystyle\\!\\!(-1)^{n+1}\\!\sin(2\theta)n\pi
Be^{-rm}\\!\left(1\\!+\\!rm\\!+\\!\frac{r^{2}m^{2}}{3}\right)\,.$

Therefore, the detection of such a quantity could demonstrate the existence of
axions.

From eq. (16), we can see that the entanglement witness, for $n\pi B\ll 1$,
has a dependence on the interaction strength proportional to $g_{p}^{2}$
rather than to $g_{p}^{4}$ making the signal larger (indeed
$B=\frac{4g_{p}^{2}}{\alpha g^{2}}$). In fact, the entanglement witness of our
system assumes values larger than the entropy, making the detection of the
axion and of the ALPs much more viable. This is shown in the lower panel of
Fig. (1), where plots of the entanglement witness for $t=t^{*}$ are reported,
and compared with those of the $2$-Renyi entropy in correspondence with the
same parameters.

##  V Discussions and Conclusions

In our analysis on the dynamics induced, in a system of two spin-$\frac{1}{2}$
fermions, by the axion–mediated fermion–fermion interaction, we have taken
into account two important constraints: the finiteness of the coherence time
Buchleitner  and the limitation coming from the Yukawa damping factor
$e^{-mr}$, which strongly suppresses the axion–mediated interaction, outside a
limited spatial region $r<\frac{1}{m}$.

Relatively to the first constraint, it is well known that one of the main
problems, when we are interested in the experimental analysis of the
entanglement of a system, is the finiteness of the coherence time. In
realistic systems, the coherent superposition characterizing the pure quantum
states is destroyed by interactions that, unavoidably, every quantum system
shares with the surrounding world. In our computations, we have neglected such
kind of interaction and we assumed that our system is completely isolated from
the rest of the universe. This assumption is realistic since we have
considered a finite time interval which is lesser than the coherence time,
which in modern experimental setups can reach values of order of the second
Abobeih2018  and is continuously extending due to the progress in the
experiments. The characteristic time interval of our system, which has to be
compared with the coherence time, is the minimum time interval needed to
isolate the axion contribution to the entanglement ($t^{*}$). For a system of
two electrons a value of $t^{*}=7\,s$, which is at the limit of today’s
technology, is obtained by considering the relative distance $r=0.1\,\mu m$. A
similar result is achieved for a system of two neutrons by considering
relative distances of order of the nanometer. On the other hand, the Yukawa
damping factor has the only effect of reducing the range of distances $r$ for
which the entanglement is significantly different from zero. This is not a
serious limitation, since the smaller the axion mass, the larger is the
spatial region (with $r<\frac{1}{m}$) where the model is efficient.

The knowledge of $t^{*}$ with a high precision is crucial for the approach
described in the paper. The largest source of error that can affect this
quantity is represented by the uncertainty $\Delta r$ on the distance between
the two fermions. A simple analysis based on the definition of the
entanglement witness in eq. (15) shows that an uncertainty $\Delta t^{*}$ on
$t^{*}$ implies an uncertainty $\Delta C_{yz}(t^{*})$ on the entanglement
that, to lowest order in $B$, is equal to $\Delta C_{yz}(t^{*})=\frac{\Delta
t^{*}}{t^{*}}\sin(2\theta)$. The higher order terms $\mathcal{O}(B^{2})$ can
be safely neglected in virtue of the smallness of the coupling constant
$g_{p}$. Any error $\Delta r$ on the inter–fermion distance then leaads to an
uncertainty $\Delta C_{yz}(t^{*})=3\frac{\Delta r}{r}\sin(2\theta)$. From the
definition of $t^{*}$ we can see that another source of uncertainty is
represented by the fermion magnetic moments, which are known only with finite
precision. Nevertheless, we expect that the uncertainty on the distance
$\frac{\Delta r}{r}$ dominates over the uncertainty on the magnetic moments,
which, in comparison, are known to a high degree of precision Tanabashi2019 .
Given these considerations, a promising framework for the realization of our
approach are optical lattices Bloch2005 ; McGuyer2015 . In this context it is
indeed possible to directly access spin correlation functions Parsons2016 ,
with single site resolution imaging of fermions Haller2015  and a precise
control over the particle spacing.

In conclusion, taking into account the above constraints, we have analyzed the
dynamics induced, in a system of two spin-$\frac{1}{2}$ fermions, by the
axion–mediated fermion–fermion interaction in the non–relativistic regime. We
have shown that it is characterized by the rising of entanglement between the
two fermions. Moreover, we have shown that, by suitably tuning the observation
time $t=nt^{*}$ and the distance between the two fermions, one can get rid of
the contribution given to entanglement by the dipole–dipole interaction of
magnetic origin. In this way, any residual entanglement–entropy can be seen as
a direct consequence of the presence of axions and hence constitutes a proof
of their existence. On the other hand if such entanglement is not detected
this observation can be used to strengthen the current constraints.

In addition, to overcome possible difficulties in direct entropy measurement,
we have introduced a spin–spin correlation function which we have proved to be
a suitable entanglement witness, that vanishes if and only if the entanglement
goes to zero. Such witness has also the advantage to be proportional, for
$t=t^{*}$, to $g_{p}^{2}$ rather than to $g_{p}^{4}$ as the 2–Renyi entropy.
This fact allows extending the range of applicability of our experiment of
several orders of magnitude. The method we propose can likely probe a coupling
constant range $g_{p}=10^{-3}-10^{-1}$ for any axion mass up to $m\simeq
1\,eV$, being particularly efficient for ALPs with low masses and large
coupling constants Marques2018 . For coupling constants below $10^{-3}$ and
masses beyond $m\simeq 1\,eV$ measurements are limited by the current
experimental precision. Improvements in this respect may render wider regions
of parameter space accessible in the next future.

## Acknowledgments

– A.C. G.L. and A.Q. thank partial financial support from MIUR and INFN. A.C.
and G.L. thank also the COST Action CA1511 Cosmology and Astrophysics Network
for Theoretical Advances and Training Actions (CANTATA). SMG acknowledges the
European Regional Development Fund the Competitiveness and Cohesion
Operational Programme (KK.01.1.1.06–RBI TWIN SIN), the Croatian Science Fund
Projects No. IP-2016–6–3347 and IP-2019–4–3321, and the QuantiXLie Center of
Excellence, a project co–financed by the Croatian Government and European
Union through the European Regional Development Fund–the Competitiveness and
Cohesion Operational Programme (Grant KK.01.1.1.01.0004).

## References

  * (1) S. M. Bilenky and B. Pontecorvo, Phys. Lett. B 61, 248 (1976).
  * (2) S. M. Bilenky and B. Pontecorvo, Yad. Fiz. 3, 603 (1976).
  * (3) P. K. Kabir, “The CP Puzzle”, Academic Press, London (1968).
  * (4) O. Nachtmann, “Elementary Particle Physics: Concepts and Phenomena”, Springer, Berlin (1990).
  * (5) R. D. Peccei and H. Quinn, Phys. Rev. Lett. 38, 1440 (1977).
  * (6) J. Ellis, Nucl. Phys. A827, 187-198 (2009).
  * (7) Particle Data Group review (Axions) (2019).
  * (8) R.D. Peccei, Lect. Notes Phys. 741, 3 (2008).
  * (9) R. D. Peccei and H. Quinn, Phys. Rev. D 16 1791 (1977).
  * (10) F. Wilczek, Phys. Rev. Lett. 40, 279 (1978).
  * (11) S. Weinberg, Phys. Rev. Lett. 40, 223 (1978).
  * (12) J. M. Pendlebury et al., Phys. Rev. D 92, 092003 (2015).
  * (13) C. Abel et al., Phys. Rev. Lett. 124, 081803 (2020).
  * (14) B. Graner, Y. Chen, E. G. Lindahl, B. R. Heckel, Phys. Rev. Lett. 116, 161601 (2016).
  * (15) V. V. Flambaum, M. Pospelov, A. Ritz, Y. V. Stadnik, arXiv:1912.13129 (2019).
  * (16) G. G. Raffelt, J. Phys. A 40, 6607 (2007).
  * (17) J. E. Kim, Phys. Rev. Lett 43, 103 (1979); M. A. Shifman, A. I. Vainstein, V. I. Zakharov, Nucl. Phys. B 166, 493 (1980).
  * (18) A. P. Zhinitskii, Sov. J. Nucl. phys 31, 260 (1980); M. Dine, F. Fischler, M. Srednicki, Phys. Lett. B 104, 199 (1981).
  * (19) J. E. Kim and D. J. E. Marsch, Phys. Rev. D 93, 025027 (2016).
  * (20) I. De Martino, T. Broadhurst, S. H. H. Tye, T. Chiueh, H. Y. Schive and R. Lazkoz, Phys. Rev. Lett. 119, 221103 (2017)
  * (21) V. A. Rubakov, JETP Lett. 65, 621-624 (1997).
  * (22) D. J. E. Marsch Phys. Rep., 643, 1 (2016).
  * (23) G. Lambiase and S. Mohanty, “The 21-cm axion” arXiv:1804.05318 (2018).
  * (24) A. Auirol, S. Davidson and G. Raffelt, Phys. Rev. D 99, 023013 (2019).
  * (25) N. Houston, C. Li, T. Li, Q. Yang, X. Zhang, Phys. Rev. Lett. 121, 111301 (2018).
  * (26) E. Zavattini et al. (PVLAS Collaboration), Phys. Rev. D 77, 032006 (2008).
  * (27) P. Pugnat et al. (OSQAR Collaboration), Phys. Rev. D 78, 092003 (2008).
  * (28) R. Ballou et al. (OSQAR Collaboration), Phys. Rev. D 92, 092002 (2015).
  * (29) K. Ehret et al. (ALPS Collaboration), Phys. Lett. B 689, issues 4-5, pages 149-155 (2010).
  * (30) S. Aune et al. (CAST Collaboration), Phys. Rev. Lett. 107, 261302 (2011).
  * (31) N. Du et al. (ADMX Collaboration), Phys. Rev. Lett. 120, 151301 (2018).
  * (32) A. Capolupo, G. Lambiase, G. Vitiello, Adv. In High En. Phys. 826051 (2015).
  * (33) A. Capolupo, I. De Martino, G. Lambiase, A. Stabile, Phys. Lett. B 790, 427-435 (2019).
  * (34) V.B. Bezerra, G.L. Klimchitskaya, V.M. Mostepanenko and C. Romero, Phys. Rev. D 94, 035011 (2016).
  * (35) R. Barbieri et al., Phys. Dark Univ. 15, 135-141 (2017).
  * (36) M. Arndt, T. Juffmann and V. Vedral, HFSP Journal 3, 386 (2009).
  * (37) G. Vidal, J. I. Latorre, E. Rico and A. Kitaev, Phys. Rev. Lett. 90, 227902 (2003).
  * (38) L. Amico, R. Fazio, A. Osterloh and V. Vedral Rev. Mod. Phys. 80, 517 (2008).
  * (39) S. Bose, A. Mazumdar, G. W. Morley, H. Ulbricht, M. Toroš, M. Paternostro, A. A. Geraci, P. F. Barker, M. S. Kim, and G. Milburn, Phys. Rev. Lett. 119, 240401 (2017).
  * (40) C. Marletto and V. Vedral, Phys. Rev. Lett. 119, 240402 (2017).
  * (41) R. J. Marshman, A. Mazumdar, S. Bose, “Locality & Entanglement in Table-Top Testing of the Quantum Nature of Linearized Gravity”, arXiv:1907.01568 (2019).
  * (42) K. Simonov, A. Capolupo and S. M. Giampaolo, Eur. Phys. J. C 79, 902 (2019).
  * (43) S. M. Giampaolo and T. Macri, Scientific Reports 9 , 11362 (2019).
  * (44) C. H. Bennett, P. W. Shor, J. A. Smolin and V. A. Thapliyal, Phys. Rev. Lett. 83, 3081 (1999).
  * (45) J.E. Moody and F. Wilczek, Phys. Rev. D 30 130 (1984).
  * (46) R. Daido and F. Takahashi, Phys. Lett. B 772, 127 (2017).
  * (47) P. Horodecki and A. Ekert, Phys. Rev. Lett. 89, 127902 (2002).
  * (48) M. B. Plenio and S. Virmani, Quant. Inf. Comput. 7, 1 (2007).
  * (49) M. A. Nielsen and I. L. Chuang, Quantum Computation and Quantum Information, Cambridge University Press (2010).
  * (50) S. M. Giampaolo, S. Montangero, F. Dell’Anno, S. De Siena, and F. Illuminati, Phys. Rev. B 88, 125142 (2013).
  * (51) F. A. Bovino, G. Castagnoli, A. Ekert, P. Horodecki, C. Moura Alves and A. V. Sergienko, Phys. Rev. Lett. 95, 240407 (2006)
  * (52) T. Brydges, A. Elben, P. Jurcevic, B. Vermersch, C. Maier, B. P. Lanyon, P. Zoller, R. Blatt and C. F. Roos, Science 364, 260 (2019).
  * (53) B. Lesche, Phys. Rev. E 70, 017102 (2004).
  * (54) R. Islam, R. Ma, P. M. Preiss, M. E. Tai, A. Lukin, M. Rispoli and Markus Greiner, Nature 528, 77 (2015).
  * (55) D. A. Abanin and E. Demler Phys. Rev. Lett. 109, 020504 (2012).
  * (56) V. B. Bezerra, G. L. Klimchitskaya, V. M. Mostepanenko and C. Romero, Eur. Phys. J. C 74, 2859 (2014).
  * (57) G. Bimonte, D. Lopez and R. S. Decca, Phys. Rev. B 93, 184434 (2016).
  * (58) G. L. Klimchitskaya and V. M. Mostepanenko, Eur. Phys. J. C 75, 164 (2015).
  * (59) Y. J. Chen, W. K. Tham, D. E. Krause, D. Lopez, E. Fischbach and R. S. Decca, Phys. Rev. Lett. 116, 221102 (2016).
  * (60) G. L. Klimchitskaya and V. M. Mostepanenko, Phys. Rev. D 95 123013 (2017).
  * (61) J. H. Chang, R. Essig, S. D. McDermott, J. High Energ. Phys. 2018, 51 (2018).
  * (62) K. Hamaguchi, N. Nagata, K. Yanagi, J. Zheng, Phys. Rev. D 98, 103015 (2018).
  * (63) A. Arvanitaki, M. Baryakhtar, X. Huang, Phys. Rev. D 91, 084011 (2015).
  * (64) J.-F. Hsu, P. Ji, C. W. Lewandowski, and B. D’Urso, Sci. Rep 6, 30125 (2016)
  * (65) A. Buchleitner, C. Viviescas and M. Tiersch, ”Entanglement and Decoherence - Foundations and Modern Trends”, Lecture Notes in Physics 768 (Springer-Verlag Berlin Heidelberg).
  * (66) M. H. Abobeih, J. Cramer, M. A. Bakker, N. Kalb, M. Markham, D. J. Twitchen, T. H. Taminiau, Nature comm. 9 2552 (2018).
  * (67) M. TanabashiÂ et al.Â (Particle Data Group), Phys. Rev. DÂ 98, 030001 (2018) and 2019 update.
  * (68) I Bloch, ”Ultracold quantum gases in optical lattices”, Nature Phys. 1, pp. 23-30 (2005).
  * (69) B. H. McGuyer, M. McDonald, G. Z. Iwata, M. G. Tarallo, A. T. Grier, F. Apfelbeck, T. Zelevinsky, ”High-precision spectroscopy of ultracold molecules in an optical lattice”, New Journal of Physics 17, 055004 (2015).
  * (70) M. F. Parsons, A.Mazurenko, C. S. Chiu, G. Ji, D. Greif, M. Greiner, Science 16 Sep 2016: Vol. 353, Issue 6305, pp. 1253-1256 (2016).
  * (71) E. Haller, J. Hudson, A. Kelly et al., ”Single-atom imaging of fermions in a quantum-gas microscope”, Nature Phys. 11, pp. 738-742 (2015).
  * (72) G. Marques-Tavares and M. J. Teo, High Energ. Phys. 2018, 180 (2018).