# Identifying primes from entanglement dynamics

A. L. M. Southier  Departamento de Física, Universidade Federal do Paraná,
81531-990, Curitiba, PR, Brazil  L. F. Santos  Department of Physics,
University of Connecticut, Storrs, CT, USA  P. H. Souto Ribeiro  Departamento
de Física, Universidade Federal de Santa Catarina, 88040-900, Florianópolis,
SC, Brazil  A. D. Ribeiro  Departamento de Física, Universidade Federal do
Paraná, 81531-990, Curitiba, PR, Brazil

###### Abstract

Factorization is the most fundamental way to determine if a number $n$ is
prime or composite. Yet, this approach becomes impracticable when considering
large values of $n$, a difficulty that is exploited by cryptographic
protocols. We propose an alternative method to decide the primality of a
natural number, that is based on the analysis of the evolution of the linear
entanglement entropy. Specifically, we show that a singular behavior in the
amplitudes of the Fourier series of this entropy is associated with prime
numbers. We also discuss how this idea could be experimentally implemented and
examine possible connections between our results and the zeros of the Riemann
zeta function.

Prime numbers have captured the attention of researchers for centuries. In
pure mathematics, given the prominent role of primes in factorizing a positive
integer $n$, much effort has been made to unravel patterns in their
distribution over the set of natural numbers $\mathbb{N}$ [1]. Remarkable
results in this direction involve the counting function $\pi(n)$, which gives
the number of primes less than or equal to $n$. The complete knowledge of
$\pi(n)$ would imply being able to determine the position of each prime number
$p$, since jumps of this function expose their presence. However, achieving
the counting function with satisfactory accuracy for large values of $n$ has
been shown to be practically impossible. In applied mathematics, this lack of
knowledge is exploited for cryptographic protocols [2], such as the RSA
(Rivest-Shamir-Adleman) algorithm. To break the RSA protocol, one needs to
find the prime factors of a huge $n$, which would require the implementation
of Shor’s algorithm in a quantum computer.

Another fascinating result involving primes is their connection with the zeros
of the Riemann zeta function $\zeta(s)$, where $s$ is a complex variable,
defined as [3]

$\zeta(s)\equiv\sum_{n=1}^{\infty}\frac{1}{n^{s}}=\prod_{p}\frac{1}{1-p^{-s}}=\frac{\Gamma(1-s)}{2\pi
i}\int_{\gamma}\frac{(-x)^{s}}{\mbox{e}^{x}-1}\frac{\mbox{d}x}{x}.$ (1)

The summation and the product above converge when the real part of the
variable $s$ satisfies $\Re[s]>1$. The first equality in Eq. (1) is due to
Euler and makes evident the relation between $p$ and $\zeta(s)$. In 1859 [4],
Riemann achieved the expression shown on the right side of Eq. (1), which
involves a line integral along a particular path $\gamma$ [5] and the Gamma
function $\Gamma(1-s)\equiv\int_{0}^{\infty}x^{-s}\mbox{e}^{-x}\mbox{d}x$. The
expression is analytic for all values of $s$, except for a simple pole at
$s=1$. In the region $\Re[s]>1$, the term with the integral recovers the other
two expressions in Eq. (1), so one considers it as their analytic
continuation. Using these ideas, Riemann [4] demonstrated how the complex
zeros $s_{0}$ of $\zeta(s)$ encode the distribution of $p$. Specifically, he
showed that a certain convergent series, running over all $s_{0}$, recovers
the function $\pi(n)$. But for large values of $n$, the number of zeros
required to accurately get the counting function is so large that the use of
the series to identify $p$ becomes intractable. This sophisticated connection
between $\zeta(s)$ and primes is but one of the impressive results of Ref.
[4]. Arguments in that work also gave origin to the Riemann hypothesis, which
conjectures that $\Re[s_{0}]=\frac{1}{2}$ for all nontrivial zeros of
$\zeta(s)$.

The almost mystical status of prime numbers has also reached the physics
community, especially researchers working with quantum mechanics [6, 7]. A
particularly inspiring idea is the Hilbert-Pólya conjecture that could lead to
the proof of the Riemann hypothesis. The conjecture proposes that the
nontrivial zeros of the Riemann zeta function, which supposedly fall on the
critical line $\Re[s_{0}]=\frac{1}{2}$, correspond to the eigenvalues of a
Hermitian Hamiltonian operator. The first substantial evidence supporting this
conjecture arose in the semiclassical studies carried out by Berry and Keating
[8]. They compared the distribution of $s_{0}$ over the critical line with the
equivalent function for the eigenvalues of a given Hamiltonian achieved
through the Gutzwiller trace formula. In this way, they were able to identify
a Hamiltonian operator that fulfills the conjecture. The operator is obtained
by quantizing the classical Hamiltonian $H_{\mathrm{c}l}=xp$, where $x$ and
$p$ are the particle position and momentum, respectively.

Inspired by the Hilbert-Pólya conjecture and the Berry and Keating $xp$-model,
several works [9, 10, 11, 12, 13, 14, 15] have aimed at interpreting,
extending and circumventing technical difficulties of Ref. [8]. Other
contributions looked for alternative physical systems, where the properties of
$\zeta(s)$ could be identified. This has been done using quantum graph theory
[16], many-body systems [17], wave-packet dynamics [18], statistical mechanics
of random energy landscapes [19], random matrix theory [20, 21], quantum
entanglement [22], and there is also an experiment based on a periodically
driven single qubit [23]. In addition to these works, where the Riemann zeta
function is the protagonist to link prime numbers and quantum physics, other
quantum approaches deal directly with primes, studying, for instance, number
factorization [24, 25] and the properties of prime states [26, 27], which are
superpositions of states corresponding to prime numbers.

In the present paper, we show that the dynamics of the linear entanglement
entropy encodes prime and semiprime numbers. Specifically, we demonstrate that
the Fourier modes $c_{n}$ of this evolving function have amplitudes that
present a singular behavior when $n$ corresponds to a prime, and we find
curves $c_{n}^{f}$ for the location of families $f$ of semiprimes. This is
shown for a system of two coupled harmonic oscillators and a system of two
coupled spins. We discuss a possible experimental realization of this idea and
speculate on how to link $c_{n}$ and the counting function $\pi(n)$.

Entanglement dynamics.— We consider a system consisting of two interacting
parts, $A$ and $B$, to which we assign distinct Hilbert spaces,
$\mathcal{H}_{A}$ and $\mathcal{H}_{B}$, respectively. The system is isolated
and prepared in a pure state, $\rho(0)=|\Psi(0)\rangle\langle\Psi(0)|$, which
evolves according to the total Hamiltonian

$H=H_{A}\otimes\mathbb{1}_{B}+\mathbb{1}_{A}\otimes H_{B}+\lambda H_{A}\otimes
H_{B},$ (2)

where $\lambda$ is the coupling strength between the two subsystems. The
linear entanglement entropy,

$S_{L}(t)=1-\mathrm{Tr}\left[\rho_{A}^{2}(t)\right],$ (3)

where $\rho_{A}(t)\equiv\mathrm{Tr}_{B}[\rho(t)]$ is the reduced density
matrix of system $A$, measures the bipartite entanglement in time between $A$
and $B$. When $\rho(t)$ is separable, $S_{L}(t)\\!=\\!0$, otherwise
$0\\!<\\!S_{L}(t)\\!<\\!1$.

Coupled oscillators.— The first Hamiltonian that we consider describes two
coupled harmonic oscillators, where $H_{A}$ is the Hamiltonian for part $A$
with eigenvalues $\hbar\omega_{0}(n_{A}+\frac{1}{2})$ and eigenstates
$|n_{A}\rangle$, and equivalently for part $B$. The initial state is the
product of two canonical coherent states [28],

$|\Psi(0)\rangle=|\alpha_{A}\,\alpha_{B}\rangle\equiv\mbox{e}^{-u/2}\\!\\!\\!\\!\\!\\!\sum_{n_{A},n_{B}=0}^{\infty}\frac{\alpha_{A}^{n_{A}}}{\sqrt{n_{A}!}}\frac{\alpha_{B}^{n_{B}}}{\sqrt{n_{B}!}}|n_{A}\,n_{B}\rangle,$
(4)

where we chose $|\alpha_{A}|^{2}=|\alpha_{B}|^{2}=\frac{u}{2}$. Defining
$\omega\equiv\hbar\lambda\omega_{0}^{2}$,

$S_{L}^{osc}(t)=1-\mbox{e}^{-2u}\\!\\!\sum_{j,k,l,m}\\!\frac{(\frac{u}{2})^{j+k+l+m}}{j!\,k!\,l!\,m!}\mbox{e}^{-i\omega
t(j-k)(l-m)}.$ (5)

The sum above runs from 0 to $\infty$ for all indexes. Since the entropy in
Eq. (5) is periodic in time, with period $T=2\pi/\omega$, we use the Fourier
series and arrive at

$S_{L}^{osc}(t)=c_{0}(u)-\sum_{n=1}^{\infty}\;c_{n}(u)\cos(n\omega t),$ (6)

where the coefficients $c_{n}(u)$ are given by

$c_{n}(u)=4\,\mbox{e}^{-2u}\sum_{k,m}\sum_{j>k}\sum_{l>m}\frac{(\frac{u}{2})^{j+k+l+m}}{j!\,k!\,l!\,m!}\delta_{(j-k)(l-m)}^{n},$
(7)

representing the amplitude of the mode with frequency $n\omega$. In Eq. (7),
the Kronecker delta function, $\delta^{b}_{a}$ for integers $a$ and $b$,
prompts the analysis of how $n$ decomposes as a product of two integers. We
introduce the set $\Lambda_{n}$ composed of all distinct positive divisors of
$n$ and note that the non-null terms in Eq. (7) are those for which $(j-k)$
and $(l-m)$ belong to $\Lambda_{n}$ and their product equals $n$. Therefore,

$c_{n}(u)=4\,\mbox{e}^{-2u}\,\sum_{\mu\in\Lambda_{n}}\;I_{\mu}(u)\,I_{\frac{n}{\mu}}(u),$
(8)

where
$I_{\chi}(w)\equiv\sum_{k=0}^{\infty}[k!(\chi+k)!]^{-1}(\frac{w}{2})^{2k+\chi}$
is the modified Bessel function of the first kind. As we show next, $c_{n}(u)$
is very sensitive to the primality of $n$.

Identifying primes.—To show that the entanglement dynamics quantified by
$S_{L}^{osc}(t)$ is a prime identifier, we first assume that $n$ is prime, so
that $\Lambda_{n}=\\{1,n\\}$ and Eq. (8) yields
$c_{n}(u)=c_{n}^{\scriptscriptstyle{(p)}}(u)$, where

$c_{n}^{\scriptscriptstyle{(p)}}(u)\equiv
8\,\mbox{e}^{-2u}I_{1}(u)\,I_{n}(u).$ (9)

In contrast, when $n$ is as a composite number, $\Lambda_{n}$ necessarily
consists of 1, $n$, and, at least, one more integer. By defining the set
$\Lambda^{\prime}_{n}=\Lambda_{n}-\\{1,n\\}$, Eq. (8) becomes

$c_{n}(u)=c_{n}^{\scriptscriptstyle{(p)}}(u)+4\,\mbox{e}^{-2u}\sum_{\mu\in\Lambda^{\prime}_{n}}I_{\mu}(u)\,I_{\frac{n}{\mu}}(u)\geq
c_{n}^{\scriptscriptstyle{(p)}}(u),$ (10)

which holds only for $n>1$.

Inequality (10) is a main result of this work. It shows that the coefficients
$c_{n}(u)$ coincides with $c_{n}^{\scriptscriptstyle{(p)}}(u)$ in Eq. (9) if,
and only if, $n$ is prime, while for composite numbers, the amplitudes are
strictly lower bounded by $c_{n}^{\scriptscriptstyle{(p)}}(u)$. This means
that if we have a way to determine the values of $c_{n}(u)$, other than
through the construction of $\Lambda_{n}$, the application of Eq. (10) can
reveal the primality of $n$. Indeed, $c_{n}(u)$ can be evaluated theoretically
via Eq. (5) and possibly experimentally.

In Fig. 1, we mark the values of $c_{n}(u)$ with blue dots. For each $n$
identified as a prime, the dots are encircled with a black circle. All values
of $c_{n}(u)$ corresponding to primes coincide with the red dotted line that
represents the lower bound $c_{n}^{\scriptscriptstyle{(p)}}(u)$ and can,
therefore, be clearly distinguished. In Fig. 1 (a), we show results for $u=1$,
but $c_{n}(1)$ becomes too small when $n$ is large, so in Fig. 1 (b), we use
$u=10^{3}$ and larger values of $n$.

Figure 1: Logarithm of the coefficients $c_{n}(u)$ as a function of $n$ for
(a) $u=1$ and (b) $u=10^{3}$. The blue dots, connected with a blue solid line
to guide the eye, represent $c_{n}(u)$, and the red points indicate the lower
bound $c_{n}^{\scriptscriptstyle{(p)}}(u)$. For each prime $n$, the blue dots
are encircled with black circles and touch the red dotted line. Green and cyan
lines are for $c_{n}^{\scriptscriptstyle{(f_{2})}}(u)$ and
$c_{n}^{\scriptscriptstyle{(f_{3})}}(u)$, respectively. Black squares
[triangles] enclose the coefficients of the family $f_{2}$ [$f_{3}$]. Inset
plots show a magnified region of the respective main graph.

Semiprimes.— Figure 1 also reveals the square-free semiprimes, including the
integer 2, which we denote by family $f_{2}$. For these numbers,
$\Lambda_{n}=\\{1,2,P_{n},n\\}$, where $P_{n}=n/2\neq 2$ is a prime, and
$c_{n}(u)=c_{n}^{\scriptscriptstyle{(f_{2})}}(u)$, where

$c_{n}^{\scriptscriptstyle{(f_{2})}}(u)\equiv
8\,\mbox{e}^{-2u}\left[I_{1}(u)\,I_{n}(u)+I_{2}(u)\,I_{\frac{n}{2}}(u)\right].$
(11)

For a composite $n$ that presents the divisors: 1, 2, $n/2$, $n$, and at least
one more integer,

$c_{n}^{\scriptscriptstyle{(2)}}(u)\equiv
c_{n}^{\scriptscriptstyle{(f_{2})}}(u)+4\,\mbox{e}^{-2u}{\sum_{\mu\in\Lambda^{{}^{\prime}\scriptscriptstyle{(2)}}_{n}}}I_{\mu}(u)\,I_{\frac{n}{\mu}}(u)>c_{n}^{\scriptscriptstyle{(f_{2})}}(u),$
(12)

where
$\Lambda^{{}^{\prime}\scriptscriptstyle{(2)}}_{n}=\Lambda_{n}-\\{1,2,n/2,n\\}$.
The semiprimes of family $f_{2}$ are identified with black squares in Fig. 1
and they all fall on the curve $c_{n}^{\scriptscriptstyle{(f_{2})}}(u)$,
indicated with a green line. This line consists of a lower bound for any
integer composed by 2, except for $2$, $2^{2}$, and $2^{3}$.

It is straightforward to extend the analysis above to any other family of
semiprimes. The particular case of the square-free semiprimes containing the
number 3, denoted by $f_{3}$, is shown in Fig. 1. The components of family
$f_{3}$ are marked with black triangles and they are exactly located over the
cyan curve $c_{n}^{\scriptscriptstyle{(f_{3})}}(u)$ [29].

Interacting spins.— Analogously to the system of coupled harmonic oscillators,
we now show that two interacting spins with a large quantum number $S$ is
another physical system that can be used to identify primes. In the
Hamiltonian of Eq. (2), we assume that $H_{A}=\hbar\omega_{0}S_{A}^{z}$, where
$S_{A}^{z}$ is the $z$-component of the spin operator of part $A$ and
$S_{A}^{z}|m_{A}\rangle=m_{A}|m_{A}\rangle$, and equivalently for part $B$.
The initial state is the product of spin coherent states [28],
$|\Psi(0)\rangle=|s_{A}\,s_{B}\rangle$, where

$|s_{A}\,s_{B}\rangle\equiv
N_{s}\\!\\!\\!\\!\\!\\!\sum_{n_{A},n_{B}=0}^{2S}\\!\\!\frac{s_{A}^{n_{A}}s_{B}^{n_{B}}}{n_{A}!n_{B}!}\sqrt{\textstyle\binom{2S}{n_{A}}\binom{2S}{n_{B}}}|n_{A}-S,n_{B}-S\rangle,$
(13)

where $N_{s}=(1+u^{2})^{-2S}$ and we chose $|s_{A}|^{2}=|s_{B}|^{2}=u$.

The evolving linear entanglement entropy (3) becomes

$S_{L}^{spin}(t)=1-\\!\\!\\!\\!\\!\sum_{j,k,l,m=0}^{2S}\\!\\!\\!\\!\xi_{j,k,l,m}\,u^{j+k+l+m}\mbox{e}^{-i\omega
t(j-k)(l-m)},$ (14)

where
$\xi_{j,k,l,m}\equiv\,{\textstyle\binom{2S}{j}\binom{2S}{k}\binom{2S}{l}\binom{2S}{m}}(1+u)^{-8S}$.
The Fourier series of Eq. (14) has the same structure as Eq. (6), but with the
coefficients

$\bar{c}_{n}(u)=4\sum_{k,m=0}^{2S-1}\sum_{j>k}^{2S}\sum_{l>m}^{2S}\xi_{j,k,l,m}\,u^{j+k+l+m}\delta_{(j-k)(l-m)}^{n}.$
(15)

The nonvanishing terms in the equation above are again those for which
$\mu\equiv(j-k)$ and $\nu\equiv(l-m)$ belong to the set $\Lambda_{n}$ of
divisors of $n$. However, this system is bounded, so the upper limits of the
above summations add new constraints to $\mu$ and $\nu$. Given a value for $k$
and $m$, we have that $\mu\in\mathcal{I}_{k}\equiv[1,2S-k]$ and
$\nu\in\mathcal{I}_{m}\equiv[1,2S-m]$. In addition, since
$\delta_{\mu\nu}^{n}$ implies that $\nu=n/\mu$, we also get that
$\mu\in\mathcal{I}_{m,n}\equiv[\frac{n}{2S-m},n]$. Therefore,

$\bar{c}_{n}(u)=4\sum_{k,m=0}^{2S-1}\,\sum_{\mu\in\bar{\Lambda}^{\scriptscriptstyle{(k,m)}}_{n}}\\!\\!\\!\bar{\xi}\,u^{2k+2m+\mu+n/\mu},$
(16)

where $\bar{\xi}\equiv\xi_{k+\mu,k,m+\frac{n}{\mu},m}$ and
$\bar{\Lambda}_{n}^{\scriptscriptstyle{(k,m)}}$ is the set of divisors of $n$
that satisfies all constraints commented above,

$\bar{\Lambda}_{n}^{\scriptscriptstyle{(k,m)}}\equiv\Lambda_{n}\cap\mathcal{I}_{k}\cap\mathcal{I}_{m,n}.$
(17)

Contrary to $S_{L}^{osc}(t)$ in Eq. (6), there is a finite number of Fourier
modes in $S_{L}^{spin}(t)$ due to the finite Hilbert space for the spin
system. For $n>4S^{2}$, we have that
$\mathcal{I}_{k}\cap\mathcal{I}_{m,n}=\emptyset$ for any value of $k$ and $m$,
so $\bar{\Lambda}_{n}^{\scriptscriptstyle{(k,m)}}=\emptyset$ and
$\bar{c}_{n}(u)=0$. In what follows, to assess the primality of $n$, we
conveniently divide the interval $1<n\leq 4S^{2}$ in two regions of interest,
region I ($1<n\leq 2S$) and region II ($2S<n\leq 4S$). In the complementary
range $4S<n\leq 4S^{2}$, one cannot find a behavior that distinguishes primes
from composite numbers. There, even though all amplitudes of the Fourier modes
for a prime $n$ vanish, the same can also happen to some composite numbers.

Region I ($1<n\leq 2S$): A careful inspection of Eq. (17) reveals that
$1\in\bar{\Lambda}^{(k,m)}_{n}$ for all values of $k$ and $0\leq m\leq 2S-n$,
while $n\in\bar{\Lambda}^{(k,m)}_{n}$ for all values of $m$ and $0\leq k\leq
2S-n$. Applying these results to Eq. (16), we find that, for $n$ prime,
$\bar{c}_{n}(u)=\bar{c}^{{\scriptscriptstyle(p)}}_{n,I}(u)$, where

$\bar{c}_{n,I}^{{\scriptscriptstyle(p)}}(u)\equiv
8\,(1+u)^{-8S}G_{1}(u)G_{n}(u),$ (18)

with
$G_{\chi}(w)\equiv\sum_{k=0}^{2S-\chi}\binom{2S}{k}\binom{2S}{k+\chi}w^{2k+\chi}$.
On the other hand, for a composite number $n>1$ inside region I, we have that
$\bar{c}_{n}(u)=\bar{c}_{n}^{{\scriptscriptstyle(p)}}(u)+\bar{d}_{n,I}(u)$,
where

$\bar{d}_{n,I}(u)\equiv
4{\sum_{k,m=0}^{2S-1}\,\sum_{\mu\in{\bar{\Lambda}^{\prime\scriptscriptstyle{(k,m)}}_{n}}}}\,\bar{\xi}\,u^{2k+2m+\mu+n/\mu}$
(19)

and
${\bar{\Lambda}^{\prime\scriptscriptstyle{(k,m)}}_{n}}\equiv\Lambda^{\prime}_{n}\cap\mathcal{I}_{k}\cap\mathcal{I}_{m,n}$,
where $\Lambda^{\prime}_{n}$ was defined above Eq. (10). One can show that at
least for $k=m=0$, there is a positive non-null term in Eq. (19). Therefore,
for a composite number, we necessarily have
$\bar{c}_{n}(u)>\bar{c}_{n,I}^{{\scriptscriptstyle(p)}}(u)$ and this
inequality can be employed in the search for prime numbers in region I.

Region II ($2S<n\leq 4S$): In this case, prime modes are not included in the
Fourier series of Eq. (14). Indeed, for a prime number $n$, we now have that
$1\not\in\mathcal{I}_{m,n}$ and $n\not\in\mathcal{I}_{k}$, which implies that
$\bar{\Lambda}_{n}^{\scriptscriptstyle{(k,m)}}=\emptyset$ and
$\bar{c}_{n}(u)=0$. For a composite $n$ in the same interval, similarly to
region I, there is at least one integer in $\Lambda_{n}$, in addition to 1 and
$n$, so that $\bar{c}_{n}(u)>0$. So we can distinguish primes from composite
numbers by checking if $\bar{c}_{n}(u)=0$ or $\bar{c}_{n}(u)\neq 0$.

Figure 2: Setup for the identification of prime numbers. The symbols BS, HWP,
PBS, M, KM, and PD refer to beam splitter, half-wave plate, polarizing beam
splitter, mirror, Kerr medium, and photodiode, respectively. See text for
details.

Experimental proposal.— The Hamiltonian for the two coupled harmonic
oscillators studied above also describes the interaction between two optical
fields via a Kerr nonlinear medium [30, 31]. In Fig. 2, we show the sketch of
an experimental setup for implementing this optical system. A laser beam of
frequency $\omega_{0}$ and linear polarization at $45^{\circ}$ is sent to a
polarizing beam splitter (PBS). The vertically polarized ($\updownarrow$)
component goes directly to the homodyne detection scheme [32], shown in the
lower right corner of the figure, to act as the local field. The horizontally
polarized ($\leftrightarrow$) component passes through a half waveplate (HWP),
which rotates the polarization to $45^{\circ}$, and enters an unbalanced Mach-
Zehnder interferometer, identified in Fig. 2 with dashed lines. In the
interferometer, after the PBS, the $\updownarrow$ component beam propagates
through the short path, and the $\leftrightarrow$ one goes through the long
course. The path difference is longer than the coherence length, so that the
recombined $\leftrightarrow$ and $\updownarrow$ beams at the interferometer
output PBS are separable and no longer result in a pure mode with linear
diagonal polarization. Next, these two beams are injected in the nonlinear
Kerr medium (KM) of length $L$ and Kerr optical nonlinearity $\chi^{(3)}$.
During the propagation time inside the KM, each beam will experience a
modified index of refraction due to the action of the other beam. This is how
the coupling between the oscillators is physically implemented. By varying the
length $L$, one changes the interaction time $t$. After passing KM, the
$\leftrightarrow$ and $\updownarrow$ beams are split again. To identify the
prime numbers, we measure one of the beams and ignore the other, which is
equivalent to performing the trace over one of the interacting systems. We can
analyze any of the two beams because of the interaction symmetry. In Fig. 2,
we choose the $\updownarrow$ beam, which goes to the homodyne detector, where
quantum state tomography is performed to reconstruct the density matrix
$\rho_{A}$. With the reduced density matrix, we calculate the linear
entanglement entropy. After several experimental realizations spanning the
needed values of $t$, curves like those in Fig. 1 can be generated.

Conclusion.— This paper shows that by analyzing the bipartite entanglement in
time, one can identify prime and semiprime numbers in $\mathbb{N}$. The main
ingredient is a Hamiltonian composed of two parts, $A$ and $B$, with
equidistant energy levels for $H_{A}$ and $H_{B}$. We discussed how this idea
could be implemented taking advantage of an existing experimental setup.

Our work resonates with the Hilbert-Pólya conjecture in the sense of proposing
a physically measurable quantity to determine prime numbers. We speculate that
there may be a way to connect our results with $\zeta(s)$. In fact, by
defining $A(n,u)\equiv c_{n}(u)-c^{\scriptscriptstyle(p)}_{n}(u)$ for a fixed
$u$, the task of counting primes along $\mathbb{N}$ is equivalent to counting
the zeros of $A(n,u)$. This alternative method to build $\pi(n)$ could be used
to study the zeros of $\zeta(s)$ by means, for example, of a hypothetical
inversion of the Riemann series that connects $s_{0}$ with $\pi(n)$, shedding
new light on the Riemann hypothesis.

Our last remark concerns the semiclassical approaches originated from the
ideas of Berry and Keating [8]. The semiclassical version of the linear
entanglement entropy was addressed in [33, 34] for the same two physical
systems treated here. In those papers, entanglement is reproduced by summing
over sets of classical trajectories determined by the solutions of a given
transcendental equation. The purely quantum formalism presented here
encourages the re-examination of those works aiming at connecting those
solutions with the distribution of primes.

## Acknowledgements

This study was financed by the Brazilian agencies: Fundação de Amparo à
Pesquisa do Estado de Santa Catarina (FAPESC – DOI 501100005667), Coordenação
de Aperfeiçoamento de Pessoal de Nível Superior (CAPES – DOI 501100002322),
Conselho Nacional de Desenvolvimento Científico e Tecnológico (CNPq – DOI
501100003593), and Instituto Nacional de Ciência e Tecnologia de Informação
Quântica (INCT-IQ 465469/2014-0). LFS was supported by the NSF CCI grant
(Award Number 2124511). The authors also thank Gisele T. Paula and João V. P.
Poletto for their valuable discussion at the very beginning of the work.

## References

  * [1] P. Ribenboim, The book of prime number records (Springer-Verlag, New York, 1988)
  * [2] N. Koblitz, A course in number theory and cryptography (Springer, New York, 1994)
  * [3] H. M. Edwards, Riemann’s zeta function (Dover, New York, 2001).
  * [4] B. Riemann, Uber die Anzahl der Primzahlen unter einer gegebenen Grösse, Monatsb. der Berliner Akad., 671 (1859).
  * [5] The path of integration $\gamma$, defined in the $x$-complex plane, starts at $x\to\infty+i\epsilon$, with $\epsilon\to 0_{+}$, goes to the origin parallelly to the real axis, encircles it in the counter clock sense, and returns to the final point $x\to\infty-i\epsilon$, also parallelly to the real axis.
  * [6] D. Schumayer and D. A. W. Hutchinson, Physics of the Riemann hypothesis, Rev. Mod. Phys. 83, 307 (2011).
  * [7] M. Wolf, Will a physicist prove the Riemann hypothesis?, Rep. Prog. Phys. 83, 036001 (2020).
  * [8] M. V. Berry and J. P. Keating, The Riemann zeros and eigenvalues asymptotics, SIAM 41, 236 (1999).
  * [9] G. Sierra and P. K. Townsend, Landau levels and Riemann zeros, Phys. Rev. Lett. 101, 110201 (2008).
  * [10] D. Schumayer, B. P. van Zyl, and D. A. W. Hutchinson, Quantum mechanical potentials related to the prime numbers and Riemann zeros, Phys. Rev. E 78, 056215 (2008).
  * [11] G. Sierra and J. Rodríguez-Laguna, $H=xp$ model revisited and the Riemann zeros, Phys. Rev. Lett. 106, 200201 (2011).
  * [12] M. Srednicki, The Berry–Keating Hamiltonian and the local Riemann hypothesis, J. Phys. A: Math. Theor. 44, 305202 (2011).
  * [13] G. Sierra, General covariant $xp$ models and the Riemann zeros, J. Phys. A: Math. Theor. 45, 055209 (2012).
  * [14] M. Wolf, Nearest-neighbor-spacing distribution of prime numbers and quantum chaos, Phys. Rev. E 2014, 022922 (2014).
  * [15] C. M. Bender, D. J. Brody, and M. P. Müller, Hamiltonian for the Zeros of the Riemann Zeta Function, Phys. Rev. Lett. 118, 130201 (2017).
  * [16] J. Kuipers, Q. Hummel, and K. Richter, Quantum graphs whose spectra mimic the zeros of the Riemann zeta function, Phys. Rev. Lett. 112, 070406 (2014).
  * [17] G. Mussardo, A. Trombettoni, and Z. Zhang, Prime suspects in a quantum ladder, Phys. Rev. Lett. 125, 240603 (2020).
  * [18] R. Mack, J. P. Dahl, H. Moya-Cessa, W. T. Strunz, R. Walser, and W. P. Schleich, Riemann $\zeta$ function from wave-packet dynamics, Phys. Rev. A 82, 032119 (2010).
  * [19] Y. V. Fyodorov, G. A. Hiary, and J. P. Keating, Freezing Transition, Characteristic Polynomials of Random Matrices, and the Riemann Zeta Function, Phys. Rev. Lett. 108, 170601 (2012).
  * [20] E. B. Bogomolny and J. P. Keating, Random matrix theory and the Riemann zeros. I. Three- and four-point correlations, Nonlinearity 8, 1115 (1995).
  * [21] P. Leboeuf and A. G. Monastra, Quantum thermodynamic fluctuations of a chaotic Fermi-gas model, Nucl. Phys. A 724, 69 (2003).
  * [22] C. Feiler and W. P. Schleich, Entanglement and analytical continuation: an intimate relation told by the Riemann zeta function, N. J. Phys. 15, 063009 (2013).
  * [23] Ran He, Ming-Zhong Ai, Jin-Ming Cui, Yun-Feng Huang, Yong-Jian Han, Chuan-Feng Li, Tao Tu, C. E. Creffield, G. Sierra, and Guang-Can Guo, Identifying the Riemann zeros by periodically driving a single qubit, Phys. Rev. A 101, 043402 (2020).
  * [24] F. Gleisberg, F. Di Pumpo, G. Wolff, and W. P. Schleich, Prime factorization of arbitrary integers with a logarithmic energy spectrum, J. Phys. B: Atom, Mol. and Opt. Phys. 51, 035009 (2018)
  * [25] J. L. Rosales and V. Martin, Quantum simulation of the integer factorization problem: Bell states in a Penning trap, Phys. Rev. A 97, 032325 (2018).
  * [26] J. I. Latorre and G. Sierra, There is entanglement in the primes, Quant. Info. Comput. 15, 622 (2015).
  * [27] D. García-Martín, E. Ribas, S. Carrazza, J. I. Latorre, and G. Sierra, The Prime state and its quantum relatives, Quantum 4, 371 (2020).
  * [28] J. P. Gazeau, Coherent states in quantum physics (Wiley, New York, 2009).
  * [29] All amplitudes corresponding to family $f_{3}$ satisfy $c_{n}(u)=c_{n}^{\scriptscriptstyle{(f_{3})}}(u)\equiv 8\mbox{e}^{-2u}[I_{1}(u)\,I_{n}(u)+I_{3}(u)\,I_{\frac{n}{3}}(u)]$.
  * [30] M. O. Scully and M. S. Zubairy, Quantum optics (Cambridge University Press, U.K., 1997).
  * [31] N. Imoto, H. A. Haus, and Y. Yamamoto, Quantum nondemolition measurement of the photon number via the optical Kerr effect, Phys. Rev. A 32, 2287 (1985).
  * [32] H. P. Yuen and V. W. S. Chan, Noise in homodyne and heterodyne detection, Opt. Lett. 83, 177 (1983).
  * [33] A. J. S. Lara and A. D. Ribeiro, Entanglement via entangled-boundary-condition trajectories: Long-time accuracy, Phys. Rev. A 100, 042123 (2019)
  * [34] M. V. Scherer and A. D. Ribeiro, Entanglement dynamics of spins using a few complex trajectories, Phys. Rev. A 104, 042222 (2021).