Prime Factorization of Arbitrary Integers with a Logarithmic Energy Spectrum

TL;DR

This paper introduces a quantum-based iterative method for prime factorization using bosonic atoms with a logarithmic energy spectrum, enabling efficient extraction of prime factors through energy measurements.

Contribution

It presents a novel quantum scheme leveraging a logarithmic energy spectrum of bosonic atoms for prime factorization, a significant departure from classical algorithms.

Findings

Probabilistic factor extraction with 50% success rate per measurement

Utilizes a logarithmic energy spectrum to encode factors

Proposes a protocol for prime-only energy spectra

Abstract

We propose an iterative scheme to factor numbers based on the quantum dynamics of an ensemble of interacting bosonic atoms stored in a trap where the single-particle energy spectrum depends logarithmically on the quantum number. When excited by a time-dependent interaction these atoms perform Rabi oscillations between the ground state and an energy state characteristic of the factors. The number to be factored is encoded into the frequency of the sinusoidally modulated interaction. We show that a measurement of the energy of the atoms at a time chosen at random yields the factors with probability one half. We conclude by discussing a protocol to obtain the desired prime factors employing a logarithmic energy spectrum which consists of prime numbers only.

Equations

$$E_{\ell }(L)\equiv \hslash {\omega }_0\ln \left(\frac{\ell }{L}+1\right),\ell =0,1,2,\dots$$

$$E_{\mathrm{total}}^{(n)}\equiv \hslash {\omega }_0\ln \left(\frac{{\prod }_{i=1}^n{p}_i}{L^n}\right)$$

$$E_{\mathrm{total}}^{(n)}=E_{p_1-L}(L)+E_{p_2-L}(L)+\dots +E_{p_n-L}(L)\equiv E_{\mathbf{p}-\mathbf{L}}(L)$$

$$E_{\mathrm{total}}^{(3)}=E_{p_1\cdot p_2-L}(L)+E_{p_2-L}(L)=E_{p_1-L}(L)+E_{p_2^2-L}(L).$$

$$N\equiv \prod _{i=1}^n{p}_i,$$

$$E_{\mathrm{total}}^{(2)}=\hslash {\omega }_0\ln \left(\frac{N}{L^2}\right).$$

$$E_{\mathrm{total}}^{(2)}=E_{N-L}+\hslash {\omega }_0\ln \left(\frac{1}{L}\right).$$

$$E_{\mathrm{total}}^{(k)}\equiv \hslash {\omega }_0\ln \left(\frac{N}{L^k}\right)\text{mbox}to{E}_{\mathrm{total}}^{(k+1)}\equiv \hslash {\omega }_0\ln \left(\frac{N}{L^{k+1}}\right).$$

$$E_{m_i-L}(L)=\hslash {\omega }_0\ln \left(\frac{m_i}{L}\right),i=1,\dots ,k$$

$$\hslash {\omega }_0\ln \left(\frac{N}{L^{k+1}}\right)=\hslash {\omega }_0[\sum _{i=1}^n\ln \left(\frac{m_i}{L}\right)+\ln \left(\frac{1}{L}\right)],$$

$${\hat{H}}_0^{(1)}(x;L){\phi }_{\ell }(x;L)$$

$${\hat{H}}_0^{(k)}(\mathbf{x};L){\phi }_{\mathbf{m}}(\mathbf{x};L)$$

$${\phi }_{\mathbf{m}}(\mathbf{x};L)\equiv \prod _{i=1}^k{\phi }_{m_i}(x_i;L),$$

$$E_{\mathbf{m}}(L)\equiv \sum _{i=1}^k{E}_{m_i}(L).$$

$$E_{\mathrm{total}}^{(k)}=\hslash {\omega }_0\ln \left(\frac{N}{L^k}\right).$$

$$\delta {\hat{V}}^{(k)}(t)\equiv {\gamma }^{(k)}\sin ({\omega }_{\mathrm{ext}}t)v^{(k)}(\hat{\mathbf{x}})$$

$${\omega }_{\mathrm{ext}}\equiv {\omega }_0\ln \left(\frac{N}{L^k}\right)$$

$$i\hslash \frac{d}{dt}|\Psi (t)⟩=\left({\hat{H}}_0^{(k)}+\delta {\hat{V}}^{(k)}(t)\right)|\Psi (t)⟩$$

$$|\Psi (t)⟩=\sum _{\mathbf{m}}{b}_{\mathbf{m}}(t)e^{-i{E}_{\mathbf{m}}t/\hslash }|\mathbf{m}⟩.$$

$$i\hslash {\dot{b}}_{\mathbf{m}}(t)={\gamma }^{(k)}\sin ({\omega }_{\mathrm{ext}}t)\sum _{\mathbf{n}}{e}^{i(E_{\mathbf{m}}-E_{\mathbf{n}})t/\hslash }{W}_{\mathbf{m},\mathbf{n}}{b}_{\mathbf{n}}(t)$$

$$b_0(0)=1,b_{\mathbf{m}}(0)=0\text{mbox}for\sum _{i=1}^k{m}_i>0.$$

$${\dot{b}}_0(t)=\frac{{\gamma }^{(k)}}{2\hslash }\sum _{j=1}^d{W}_{0,{\mathbf{q}}_j-\mathbf{L}}{b}_{{\mathbf{q}}_j-\mathbf{L}}(t)$$

$${\dot{b}}_{{\mathbf{q}}_j-\mathbf{L}}(t)=\frac{{\gamma }^{(k)}}{2\hslash }{W}_{{\mathbf{q}}_j-\mathbf{L},0}{b}_0(t)$$

$$b_0(t)$$

$$b_{{\mathbf{q}}_j-\mathbf{L}}(t)$$

$${\Omega }_k\equiv \frac{{\gamma }^{(k)}}{2\hslash }\sqrt{\sum _{j=1}^d{W}_{0,{\mathbf{q}}_j-\mathbf{L}}^2}.$$

$$\sum _{i=0}^d{b}_{{\mathbf{q}}_i-\mathbf{L}}(t_{\mathrm{m}}{)}^2={\sin }^2({\Omega }_k{t}_{\mathrm{m}}).$$

$$P_T\equiv ⟨{\sin }^2({\Omega }_{k}t){⟩}_T\equiv \frac{1}{T}{\int }_0^T{\sin }^2({\Omega }_{k}t)dt=\frac{1}{2}-\frac{\sin (2{\Omega }_{k}T)}{4{\Omega }_{k}T}$$

$$P_{\mathrm{Shor}}=1-\frac{1}{2^{m-1}}.$$

$${\Omega }_n=\frac{{\gamma }^{(n)}}{2\hslash }{W}_{0,\mathbf{p}-\mathbf{L}}.$$

$$\Omega_{n}\propto\left\{\begin{array}[]{r@{\quad\quad}l}N^{-1/2}&\mbox{for $n$ even},\\ N^{-(1/2+{(2L-1)/(2n^{2})})}&\mbox{for $n$ odd}\end{array}\right.$$

Full text

This Accepted Manuscript is available for reuse under a CC BY-NC-ND 3.0 licence after the 12 month embargo period provided that all the terms of the licence are adhered to.

This is a peer-reviewed, un-copyedited version of an article published in Journal of Physics B: Atomic, Molecular and Optical Physics. IOP Publishing Ltd is not responsible for any errors or omissions in this version of the manuscript or any version derived from it. The Version of Record is available online at 10.1088/1361-6455/aa9957.

1 Introduction

The decomposition of a positive integer into a product of prime factors is a difficult problem in number theory. The hardest one is to factor semiprimes, i.e. numbers composed of two primes, provided they are different and quite large. No wonder that this property of semiprimes finds applications in cryptography.

For example, the famous RSA algorithm [1] uses a public key in order to code a message. However, only the owner of the private key can quickly decode it because among the steps needed to achieve this secret key from the public one is the factorization of a very large semiprime. All known factoring algorithms run in non-polynomial time on classical computers thus prohibiting unauthorized decoding in a reasonable time. But on a large ideal quantum computer Shor’s celebrated factorizing algorithm [2] takes only polynomial time and is therefore expected to break the RSA scheme in the future.

We have studied [3, 4] the factorization of semiprimes using one- and two-dimensional traps111Such traps may be realized in the future using methods from the emerging field of atomtronics [5] as e.g. the technique of ”painted potentials” [6] or ”adiabatic potentials” [7], both with a logarithmic energy spectrum [8]. In the one-dimensional case two interacting bosonic atoms are excited into a state with an energy determined by the number to be factored and the measured single-particle energies yield the factors. For the two-dimensional case a single atom suffices with the factors being given by the two energies contained in the two motional degrees of freedom.

In the present article we extend our method to the factorization of an arbitrary natural number into an a priori unknown number of prime factors222 The distribution of the number of factors of an arbitrary integer was studied in [9]. Here we are using again a one-dimensional trap with a logarithmic energy spectrum [10]. The lack of knowledge, however, leads us to a completely new protocol of factorization.

Our article is organized as follows. In Section 2 we first introduce the logarithmic energy spectrum and discuss the distribution of a given energy which depends on the product of a given number of prime factors onto a collection of single-particle states called factor states. We then employ these results to develop an iterative scheme which finally leads to the prime factors. Section 3 discusses the realization of our idea by bosonic atoms moving in a trap with an appropriate one-dimensional potential. Moreover, we solve the Schrödinger equation that governs the dynamics of the bosons under the influence of a time-dependent perturbation within the rotating-wave approximation. In Section 4 we first show that after a measurement of the single-particle energies at randomly chosen times one of the factor states is found with a probability of approximately one half. If, however, additional information like e.g. the number of prime factors is available the factor state can be found with a probability of about unity. Limitations of our method caused by the experimental conditions are discussed in Section 5.

Finally, we append a study of a different factorization scheme in Section 6 provided by a rather exotic potential which, nevertheless, has many advantages over the protocol discussed before. The realization of the potential, however, remains an open question. We close our article in Section 7 with a short summary while background material and calculational details are collected in four appendices.

2 Essential idea

In the present section we first introduce the logarithmic energy spectrum and discuss its special role in finding the factors of an integer. We then turn to the distribution of a given energy onto several subsystems. This discussion constitutes the foundation for our iterative factorization protocol.

2.1 Distribution of logarithmic energy onto subsystems

Our factorization scheme is based on a logarithmic energy spectrum of the type

[TABLE]

with $E_{0}(L)=0$. Here, the constant $L$ plays the role of a scaling parameter and $\hbar\omega_{0}$ is the unit of energy.

Given $n$ prime numbers $p_{i}$ which form the vector ${\bf p}\equiv(p_{1},p_{2},\ldots,p_{n})$ we now consider the following problem: Can an energy

[TABLE]

determined by the product of the $n$ prime factors $p_{i}$ be distributed onto $m$ non-vanishing energies of type (1), and if so under which condition is this energy decomposition unique?

In order to answer this question we first note that we can cast $E_{\rm total}^{(n)}$ into the sum

[TABLE]

of logarithmic energies (1). The $n$ energies on the right-hand side belong to the spectrum (1) if the scaling parameter $L$ is positive and integer, while all energies $E_{p_{i}-L}$ are positive if all primes $p_{i}$ obey $p_{i}>L$.

Next we discuss the number $m$ of energies defined by (1) onto which we want to distribute the total energy $E_{\rm total}^{(n)}$. If $m<n$ the desired distribution is possible but not in a unique way.

For example, the energy $E_{\rm total}^{(3)}(L)$ defined by (2) and determined by $n=3$ prime factors with $N=p_{1}\cdot p_{2}^{2}$ can be distributed in two ways onto $m=2$ energies (1) because here

[TABLE]

If, however, the number $n$ of prime factors of $N$ agrees with the number $m$ of logarithmic energies (1) the distribution of $E_{\rm total}^{(n)}$ is unique because the fundamental theorem of arithmetics [11] guarantees that the prime factorization of a natural number is unique and so is the decomposition on the right-hand side of (3)333Recall that unity is defined not to be prime. To exclude a factor unity from the product in (2) we have to choose a scaling parameter $L\geq 2$..

The case $m>n$ is trivial since the product $\prod_{1}^{n}\,p_{i}$ cannot be decomposed into more than $n$ factors. Obviously the same is true for the distribution of the energy (2).

In A.3 it will turn out that for our purposes we need stronger conditions: All factors have to obey $p_{i}>L$. Moreover, $L$ has to be odd which means $L\geq 3$. Together with these additional assumptions the results of the present section are the central tools for our factorization protocol, which we shall outline in the next section.

2.2 Factorization protocol

As we shall see below for every step in our iterative protocol we need an ensemble of $k$ identical systems each with an energy spectrum (1) characterized by a single quantum number $\ell\geq 0$. The integer $k$ is stepwise increased from $2$ to $n+1$. The subsystems together form a larger system with $k$ degrees of freedom which constitutes the basis for our factorization procedure.

We want to perform a complete factorization of a given integer

[TABLE]

where neither the number $n$ nor the values of the prime factors $p_{i}$ are known. However, to be able to apply the results of the last section $N$ must not contain factors $p_{i}\leq L$. In the case of a scaling parameter $L=3$, for example, this requirement is unimportant because factors $2$ or $3$ are easy to recognize and to remove from the number $N$ to be factored.

To find the prime factors $p_{i}$ we perform a sequence of measurements. In the first step we compose a system from two subsystems with spectrum (1) and try to prepare it in a state with total energy

[TABLE]

Then we measure the energies of the two subsystems. This step may serve as a primality test. Indeed, if $N$ happens to be prime the energy (5) can be written as

[TABLE]

The second term on the right-hand side evidently is not part of the single-particle spectrum (1) and therefore $E_{\rm{total}}^{(2)}$ cannot be distributed onto the two subsystems. So if we fail to prepare the system with energy (5) we have proven that $N$ is prime.

However, if $N$ is not prime the subsystems were found with energies $E_{m_{1}-L}$ and $E_{m_{2}-L}$, respectively. Clearly, $m_{1}$ and $m_{2}$ are factors of $N$ but we still do not know if they are prime or composite.

Since we want to find the prime factors of $N$, we iterate the procedure while changing the number of subsystems from $k$ to $k+1$ and the total energy from

[TABLE]

In every step the energies

[TABLE]

of the $k$ subsystems were measured and the factors $m_{i}$ are found by inverting (7).

We continue this process until $k$ is equal to the number $n$ of prime factors. Since $n$ is unknown we need to find it. For this purpose we increase $k$ once more. However, now the new energy cannot be distributed on the $k+1$ subsystems as is shown by the expression

[TABLE]

where the right-hand side is not a sum of energies of (1). Hence, the $k$ factors $m_{i}$ found in the iteration before are the $n$ prime factors $p_{i}$ of $N$.

We emphasize that our protocol comprises $n$ steps until the $n$ prime factors of $N$ are found. One could think of a different protocol which iterates the first step onto the two factors $m_{1}$ and $m_{2}$. Each one of them either produces two new factors or is identified as prime. Generally more steps are needed than in our protocol e.g. in the case of mutually different prime factors the number of steps is equal to $2n-1$.

In the remainder of our article we describe how the system and the subsystems used in the protocol above are realized as well as how the system is transferred into a state with a definite energy. Moreover, we show how to find out if this transfer is possible or not.

3 Model system for factorization protocol

In the preceding section we have outlined our idea of factoring an arbitrary integer consisting of an a priori unknown number of primes. We now describe in detail a model system achieving this task. For this purpose we first introduce the relevant time-independent Schrödinger equation of $n$ particles each moving in a one-dimensional potential leading to a logarithmic energy spectrum. We then couple these motions by a time- and position-dependent interaction. The number to be factored is encoded in the frequency of this perturbation. The resulting infinite system of coupled equations reduces to a finite system which we can solve. For a brief introduction into bosonic quantum states and the details of this reduction we refer to A and B, respectively.

3.1 One-dimensional motion of $n$ particles

The subsystem with energy spectrum (1) can be realized by a particle with mass $\mu$ which moves along the x-axis in a one-dimensional potential $V^{(1)}=V^{(1)}(x;L)$. Its real-valued energy wave functions $\varphi_{\ell}\equiv\varphi_{\ell}(x;L)$ obey the time-independent Schrödinger equation

[TABLE]

For a given scaling parameter $L$ we can construct numerically the potential $V^{(1)}=V^{(1)}(x;L)$ such that the solutions of (3.1) just reproduce the spectrum (1) for the energies $E_{\ell}(L)$. Here we have indicated by the argument $L$ after the semicolon that the potential $V^{(1)}$, as well as the wave functions $\varphi_{\ell}$ depend on the scaling parameter $L$. Note that the wave functions $\varphi_{\ell}=\varphi_{\ell}(x;L)$ are even (odd) for even (odd) integer indices $\ell$.

Our iteration algorithm to obtain $V^{(1)}$ is based on the Hellmann-Feynman theorem and is described in a previous article [8]. In Fig. 1 we show $V^{(1)}$ together with the eigenfunctions $\varphi_{\ell}$ for $0\leq\ell\leq 6$ for the case $L=3$. No degeneracy is present in this one-dimensional problem.

In order to realize the composite system described in Section 2.2 we consider an ensemble of $k$ non-interacting particles each one moving along the $x$-axis under the influence of the potential $V^{(1)}(x;L)$. Using a shorthand notation ${\bf x}\equiv(x_{1},x_{2},\ldots,x_{k})$ for the $k$ coordinates $x_{i}$ and ${\bf m}\equiv(m_{1},m_{2},\ldots,m_{k})$ for the $k$ quantum numbers $m_{i}$ of the $k$ subsystems the separable solution of the $k$-dimensional Schrödinger equation

[TABLE]

reads

[TABLE]

with

[TABLE]

In the next section we shall show that the factorization protocol described above can be executed with the help of the $k$-particle system (3.1) provided we are able to prepare it in a state with energy (12). In the remainder of this article we suppress for the sake of simplicity in notation the scaling parameter $L$ in the argument of the wave functions and the energies.

3.2 Time-dependent perturbation

According to the factorization protocol described in Section 2.2 we need to prepare the system with Hamiltonian $\hat{H}^{(k)}_{0}$ defined by (3.1) in a state of total energy

[TABLE]

We assume that at time $t=0$ our system is in the $k$-particle ground state $|{{\bf{0}}}\rangle$ and the perturbation

[TABLE]

is switched on. The external frequency

[TABLE]

is chosen such that the energy $\hbar\omega_{\rm ext}$ agrees with $E_{\rm total}^{(k)}$ and is determined by the number $N$ we want to decompose into $k$ factors $q_{j}$ which may or may not be prime.

The spatial part $v^{(k)}\equiv v^{(k)}({\hat{\bf x}})$ of the interaction specified in A.2 couples all $k$ degrees of freedom. The strength of $\delta{\hat{V}}^{(k)}$ given by (14) is characterized by the constant $\gamma^{(k)}$.

We derive the equations of motion for the amplitudes $b_{\bf m}(t)$ for a transition from the ground state $|{\bf 0}\rangle$ to the eigenstate $|{\bf m}\rangle$ by expanding the solution $|{\Psi(t)}\rangle$ of the Schrödinger equation

[TABLE]

into the eigenstates $|{{\bf{m}}}\rangle$ of the unperturbed Hamiltonian $\hat{H}^{(k)}_{0}$, that is

[TABLE]

Since the energies $E_{\bf m}$ given by (12) are the eigenvalues of ${\hat{H}}^{(k)}_{0}$ we arrive at the coupled system

[TABLE]

which has to be solved for the initial conditions

[TABLE]

The matrix elements $W_{{\bf m},{\bf n}}\equiv\langle{\bf m}|v^{(k)}(\hat{\bf x})|{\bf n}\rangle$ are built from the spatial part $v^{(k)}({\hat{\bf x}})$ of the interaction (14) and the unperturbed eigenstates $|{\bf m}\rangle$.

In the present article we consider identical bosonic atoms. Therefore, the state vectors $|{\bf m}\rangle$, the probability amplitudes $b_{\bf m}$, and the matrix elements $W_{{\bf m},{\bf n}}$ are understood as ”bosonic” ones in the sense defined by A. Moreover, the summation in (16) and (3.2) is performed according to (41). In the next section we show that the system under the influence of the perturbation (14) indeed can be found in a state with energy (13).

3.3 Solution of coupled equations

In B we use the so-called secular or rotating wave approximation (RWA) [12] to reduce the infinite system (3.2) to a finite number of equations with constant coefficients and derive the $d+1$ coupled equations

[TABLE]

and

[TABLE]

with $j=1,\ldots d$. The degeneracy $d$ is discussed in C.

Equations (18) and (19) capture the essential features of the dynamics of the $k$ boson system. Together with the initial values (17) and with the help of the symmetry relation (51) they lead to the solutions

[TABLE]

with the frequency

[TABLE]

Equations (20) and (21) represent our main result: The system of $k$ bosons performs a Rabi oscillation, i.e. it oscillates between the ground state $|{\bf 0}\rangle$ and a superposition of the $d$ factor states $|{{\bf q}_{i}-{\bf L}}\rangle$ with one single frequency $\Omega_{k}$ we henceforth call Rabi frequency444The problem how to exclude a vanishing of $\Omega_{k}$ is postponed to A.3.. We emphasize that the solutions (20) and (21) of the system of differential equations (3.2) are exact within the rotating wave approximation.

The system can be found with probability $|b_{{\bf q}_{i}-{\bf L}}(t)|^{2}$ in the factor state $|{{\bf q}_{i}-{\bf L}}\rangle$ and, at times equal to an odd multiple of $\pi/(2\Omega_{k})$ with certainty. Consequently, a measurement could be performed at these times provided the Rabi frequency $\Omega_{k}$ would be known.

4 Measurement strategies

Unfortunately, an estimate of the Rabi frequency $\Omega_{k}$ requires additional information, namely the knowledge of the degeneracy $d$ which in turn depends on the number $n$ of prime factors and on their multiplicities $\nu_{i}$. In this article we assume that in general these parameters are unknown.

Therefore, we now present an alternate approach to determine the single particle energies (1) of the $k$-boson system. Moreover, we also sketch how a Rabi frequency can be determined if additional information is available.

Indeed, two measurement strategies offer themselves: (i) We make a measurement at a random time in the interval from $t=0$ where the perturbation $\delta V^{(k)}$ is switched on, and the time $T$. The probability to find any one of the states $|{{\bf q}_{i}-{\bf L}}\rangle$ approaches one half for $\Omega_{k}^{-1}\ll T$, and (ii) a measurement at a well-defined time $t_{n}\equiv\pi/(2\Omega_{n})$ which requires knowledge of the Rabi frequency $\Omega_{n}$ of the $n$-boson system. Here, the probability to find the unique factor state $|{{\bf p}-{\bf L}}\rangle$ is unity. In D we derive an asymptotic expression for $\Omega_{n}$.

4.1 Probabilistic approach

The measurement is performed at time $t$ chosen at random from a time interval $[\,0,T]$. According to (21) the time-dependent probability to find any of the $d$ factor states $|{{\bf q}_{i}-\bf{L}}\rangle$ reads

[TABLE]

In Fig. 2 we show the average probability

[TABLE]

to find a factor state as a function of the duration $T$ of the time interval. Evidently, the best choice for $T$ is the decoherence time where we tacitly assume $T\gg\Omega_{k}^{-1}$. Then $P_{T}\approx 1/2$ and the ground state $|{\bf 0}\rangle$ and any of the factor states emerge with probability $1/2$, respectively. If a factor state $|{{\bf q}_{j}-{\bf L}}\rangle$ is found we have succeeded, but if the ground state is obtained the procedure has to be repeated until either a factor state is found, or a transition into this state is excluded with a sufficiently large confidence interval.

We recall that the famous Shor algorithm [2, 13] can find a factor of an odd integer $N$ composed from $m$ different prime factors of arbitrary multiplicity with a probability

[TABLE]

We, however, find with probability of approximately $1/2$ all $k$ factors of $N$ in a single run as long as the number of bosons $k$ does not exceed the number $n$ of prime factors of $N$. It is amazing and amusing but seems to be sheer coincidence that for integers composed from two different primes both protocols succeed with about the same probability $1/2$.

4.2 Additional information

If the number $n$ of prime factors $p_{i}$ is known we modify the factorization protocol of Section 2.2 which now comprises one single step only. We use $n$ bosons in the trap and an external frequency $\omega_{\rm ext}\equiv\omega_{0}\ln(N/L^{n})$ to excite a Rabi oscillation between the ground state $|{\bf 0}\rangle$ and the non-degenerate factor state $|{\bf p-L}\rangle$. Again a random measurement results with a probability $P_{T}\approx 1/2$ in the single-particle energies (1) and hence in the desired $n$ prime factors $p_{i}$.

If, in addition, the multiplicities $\nu_{j}$ of the $n$ prime factors with $\sum_{j}\nu_{j}=n$ are known we may estimate the Rabi frequency

[TABLE]

The system is excited with $\omega_{\rm ext}$ as above and the measurement now is made not at a random time but at time $t_{n}\equiv\pi/(2\Omega_{n})$. Hence, the system is found in the factor state $|{\bf p-L}\rangle$ with high probability $|b_{\bf p-L}(t_{n})|^{2}$ given by (21), and the prime factors are determined as before.

In D we give an example how to calculate the Rabi frequency (22). Here we only give the scaling properties

[TABLE]

of the asymptotic expression of $\Omega_{n}$ for large $N$ derived for a very simple pattern of factorization and refer for the prefactors and calculational details to D. We emphasize that according to (27) for an odd number $n$ of factors the Rabi frequency decreases faster with increasing $N$.

5 Limitations

According to Ref. [12] it is a shortcoming of the RWA that transitions into off-resonant states are more probable when the difference between these states and the excitation energy $\hbar\omega_{0}\ln(N/L^{n})$ becomes comparable to or even smaller than the energy $\hbar\Omega_{n}$ of the Rabi oscillation555 We tacitly assume that all transitions are allowed..

Therefore, in order to neglect transitions into off-resonant states we have the condition

[TABLE]

which apparently is not difficult because the Rabi frequency $\Omega_{n}$ defined by (22) can be made arbitrarily small by decreasing the strength $\gamma^{(n)}$ of the perturbation (14).

Unfortunately, a second condition arises from Section 4.1 where a measurement of the energies of the $n$ bosons was performed at times $t_{m}$ chosen randomly from a time interval of length $T$. To find a factor state with probability $\approx 1/2$ requires the condition $\Omega_{n}T\gg 1$ following from (24).

Clearly, the system has to be free of decoherence during the time interval $[0,T]$ hence, the Rabi frequency $\Omega_{n}$ has to obey the inequality

[TABLE]

determined by the decoherence time $T_{\rm dec}$, indicating that in contrast to (28) resulting from the suppression of the off-resonant transitions $\Omega_{n}$ must have a lower limit as well.

To demonstrate that these conflicting conditions can be fulfilled simultaneously we refer to D where the Rabi frequency $\Omega_{n}$ was calculated for the simple example $N=p^{n}$. Having in mind that $\Omega_{n}$ depends on both, the integer to be factored $N$ and the interaction strength $\gamma^{(n)}$, we determined numerically values of $N$ and $\gamma^{(n)}$ which cause the Rabi frequency $\Omega_{n}$ to obey the inequalities

[TABLE]

as well as

[TABLE]

The condition (31) was motivated by Fig. 2 and the requirement that the average probability $P_{T}$ to find a factor state has reached the asymptotic value $1/2$.

At points in the plane spanned by $\gamma^{(4)}/\hbar\omega_{0}$ and $N$ below the curve marked with ”RWA” displayed in Fig. 3 the Rabi frequency $\Omega_{4}$ obeys the inequality (30) while at those below the curve indicated by ”dec” $\Omega_{n}$ satisfies $(\ref{RabiDec})$. The shaded region consists of the points where the Rabi frequency satisfies both conditions, (30) and (31), respectively.

Our calculation was performed with the help of (71) and (73) and with parameters $L=3$, $n=4$, $T_{\rm dec}=2\,\,\,\rm sec$, $\nu_{0}\equiv\omega_{0}/(2\pi)=5\,\,{\rm kHz}$.

Figure 3 clearly demonstrates that here the maximal $N$ that can be factored with our method is about $1.2\cdot 10^{4}$. Hence, we conjecture that the upper limit of integers with a different prime decomposition which we factored with the trap parameters given above will be of the same order of magnitude.

6 Factoring with a logarithmic spectrum of primes

In this section we propose a different factorization scheme which uses a quite exotic potential and has advantages compared to the protocol discussed in extenso above. Because there are many analogies to the protocol of Section 2.2 we keep this section short.

6.1 Potential for logarithmic spectrum of primes

We consider an ensemble of non-interacting bosons with a single-particle energy spectrum111 Tran and Bhaduri [14] have determined the ground state fluctuations of ensembles of bosons in a trap with a spectrum given by (32).

[TABLE]

that is the energy levels depend on the prime numbers $p_{\ell}$ with $p_{1}=2,p_{2}=3,p_{3}=5,\ldots$. The ground state with energy $E_{0}=\hbar\omega_{0}\ln p_{0}=0$, for $p_{0}=1$ was added because we need it in our factorization protocol222Recall that unity is by definition not a prime.. We emphasize that in contrast to the energy spectrum (1) which covers all integers $\ell$ we now restrict it solely to primes.

We have employed methods described in [8] to find [15] the potential $V_{p}=V_{p}(x)$ displayed in Figure 4 together with the $14$ lowest energy levels of the spectrum (32). Due to the intricate structure of $V_{p}$ it may be difficult, however, to realize a trap with such a potential.

6.2 Factorization algorithm

Next we develop a factorization protocol analogous to the one of Section 2.2 but now based on the spectrum (32). As before we consider an integer

[TABLE]

containing the product of $n$ unknown prime factors $p_{i}$. The number $n$ of factors is unknown as well.

From arguments already used in Section 2.1 it follows that an ensemble of $n$ non-interacting bosons with total energy

[TABLE]

can be hosted by the spectrum (32) in a unique way. But it is not possible that $m$ bosons with $m<n$ can be in a state with total energy (34) because in this case at least one of them has to be in a single particle state characterized by a product of primes. It is easy to see that, for example, an energy (34) where the integer $N$ is determined by $n=3$ prime factors such as e.g. $N=p^{3}$ cannot be distributed onto $m=2$ bosons with single-particle spectrum (32) because the state $E_{N=p^{3}}$ is not present in the two-particle spectrum.

In contrast, an ensemble of $m$ bosons with $m>n$ can be in a state with energy (34) because $n$ of them are in the appropriate excited states with $E_{\ell>0}$ while the $m-n$ excess particles remain in the ground state of zero energy333 Note that this behavior is opposite to the one described in Section 2.2 where $m$ bosons with $m<n$ can be in a state with total energy (5) while $m$ bosons with $m>n$ cannot..

Our protocol to factor the integer $N$ defined by (33) is as follows. A system of $m$ bosons is prepared in the $m$-particle ground state which is even under reflection of all position coordinates of the bosons. Like in Section 3.2 an adequate perturbation causes a transition into a state from which the factors of $N$ can be extracted.

Two conditions have to be fulfilled for a successful transition: (i) The factor state has to be even as well, and (ii) we need to have enough bosons to cover the primes contained in $N$. We now address both conditions.

First we recall that the first excited state $\varphi_{1}$ with energy $E_{1}=\hbar\omega_{0}\ln 2$ is odd. Therefore, the state with energy

[TABLE]

has opposite parity to the state with energy $E_{N}$.

Using $m$ bosons we choose the periodic perturbation

[TABLE]

where

[TABLE]

with a spatial part $v^{(m)}({\bf x})$ given by (48) and two periodic terms. Only one of them is able to cause Rabi oscillations between the ground state and a state which has an energy $E_{N}$ or $E_{2N}$, respectively. It does not matter which one of the states is even because from any of them the factors of $N$ can be deduced as described in Section 2.2.

Next we turn to the second condition stating that enough bosons have to be provided to factor $N$. While the number 2 is the smallest prime we choose the number

[TABLE]

of bosons because $2N$ can never have more than $m$ factors.

With these conditions fulfilled we know with certainty that a Rabi oscillation is present in the system. Unfortunately, we cannot estimate the Rabi frequency and we have to content ourselves with a probability of $1/2$ to find a factor state when we measure at an arbitrary time the single particle energies and determine all prime factors of $N$ in a single run. However, if there are not enough bosons in the trap no Rabi oscillation can be excited with the perturbation (36) and the $m$-boson ground state is found with a probability of 1.

7 Summary

In the present article we have proposed two methods to find all prime factors of an arbitrary integer $N$ based on the quantum dynamics of a suitable number of identical bosonic atoms moving in a one-dimensional trap. Both methods rely on exciting a collective vibratory motion by a periodic time-dependent perturbation whose frequency is determined by the number to be factored. In both cases the atoms get transfered to a state where the energies of the individual atoms represent the factors.

Our first method employs the trap potential displayed in Fig. 1 which is designed to obtain the logarithmic energy spectrum (1). Moreover, we assume that the number $n$ of prime factors of $N$ is unknown. Using an iterative scheme the state characterized by the $n$ single-particle energies is obtained after $n$ steps where the probability to be successful in a single run is one half. If, however, the number $n$ of prime factors is known then this procedure is modified such that the factor state can be obtained after one step only.

Our second method is based on the potential of Fig. 4 where according to (32) the energies of the individual atoms depend logarithmically on the prime numbers. It is superior to the first one because with enough atoms in the trap the transition into the factor state can take place in a single step.

Since this potential may be difficult - if ever - to be constructed we have focused on a thorough discussion of the first method.

Appendix A Bosonic state vectors and matrix elements

In the present appendix we briefly summarize the essential ingredients of the quantum mechanical description of an ensemble of bosonic atoms. Here we concentrate on the state vectors and the matrix elements.

A.1 State vectors and probability amplitudes

State vectors of $n$ identical bosonic atoms must be symmetric under permutations of the atoms, and in the one-dimensional case this requirement is equivalent to permutations of their quantum numbers. The state vectors are built from product states by the well-known formula

[TABLE]

Here the $n$ quantum numbers $k_{i}$ are divided into $m$ groups containing $\nu_{i}$ identical numbers each with $\sum\nu_{i}=n$, and the sum is over all permutations $P\{{\bf k}\}$ of the indices $k_{i}$ in the product states.

The normalization factor $\cal N$ is given by [16]

[TABLE]

since the sum in (38) contains $n!/\nu_{1}!,\ldots,\nu_{m}!\equiv{\cal N}(\{\nu_{\bf k}\})^{-2}$ terms.

Next we represent the $n$-particle state

[TABLE]

given by (16) in a basis of bosonic state vectors $|{\bf k}\rangle_{\rm B}$ defined by (38).

When we sum over bosonic states we have to avoid counting the same state several times. Therefore, a ”bosonic” sum is defined as

[TABLE]

that is, only summands with indices $k_{1}\leq k_{2},k_{2}\leq k_{3},\ldots k_{n-1}\leq k_{n}$ are taken into account, and the ”ordinary” sum $\sum_{\bf k}$ can be expressed by the identity

[TABLE]

The probability amplitudes $b_{\bf k}(t)$ are solutions of the system (3.2) with coefficients $W_{\bf i,k}$ which are symmetric under the exchange of indices (51) and, as a consequence, the amplitudes $b_{\bf k}(t)$ enjoy the same property. With the help of (38), (40), and (42) it is now easy to derive the expansion

[TABLE]

where the bosonic probability amplitudes are given by

[TABLE]

The square of the absolute value of the bosonic probability amplitude $|b_{\bf k}^{\rm B}(t)|^{2}$ denotes the probability $P_{\bf k}(t)$ that after a simultaneous measurement of the energies of the bosons performed at the time $t$ one of the bosons is found in a state with energy $E_{k_{1}}$, another boson with energy $E_{k_{2}}$ and so on until a last one with $E_{k_{n}}$.

A.2 Symmetric interaction

If the matrix element $A_{\bf jk}\equiv\langle{\bf j\,}|{\hat{A}}|{\,\bf k}\rangle$ of an arbitrary operator $\hat{A}$ is symmetric under permutations of the indices $\{j_{i}\}$ and $\{k_{i}\}$, respectively, then the expression for the corresponding bosonic matrix element is simply

[TABLE]

and reads

[TABLE]

As an example of this relation we consider the matrix element

[TABLE]

for a short-range collective $n$-particle interaction

[TABLE]

which couples all degrees of freedom of the $n$ bosons and is symmetric, that is, $v^{(n)}({\bf x})=v^{(n)}({\bf-x})$. Here we have introduced powers of

[TABLE]

as to make the matrix element $W_{\rm{\bf 0,k}}$ dimensionless.

Using the real-valued product states $\varphi_{\bf j}({\bf x})$ defined by (11) and the interaction (48) it is easy to derive the expression

[TABLE]

Evidently, the matrix element $W^{\bf B}_{{\bf 0},{\bf q}_{j}-{\bf L}}$ is symmetric under any permutation of the indices $\{0_{\ell},({\bf q}_{j})_{i}-L\}$, and hence

[TABLE]

The sub- or superscript B used here is omitted in the main body of our article because we only deal with bosons.

A.3 Odd scaling parameter

We note that the matrix element vanishes if the product in (50) contains an odd number $n$ of odd wave functions $\varphi_{({\bf q}_{j})_{i}-L}(x)$. Since this feature emerges for every ${\bf q}_{j}$ with $j=1\ldots d$ the Rabi frequency $\Omega_{n}$ defined by (22) also vanishes: a transition from the even ground state into an odd factor state cannot take place through the interaction (48).

We recall that in Section 2.1 we have postulated that the scaling factor $L\geq 3$ has to be odd. Moreover, all prime factors obey $p_{j}>L$. Then any factor $q_{i}$, prime or composite, is odd and consequently the difference $({\bf q}_{j})_{i}-L$ is even and the integrand in (50) is even as well: $W_{{\bf 0},{\bf q}_{j}-{\bf L}}$ never vanishes whatever factors $({\bf q}_{j})_{i}>L$ appear in $N$ if we have chosen an odd scaling parameter $L$ in the single-particle spectrum (1).

Appendix B Reduction of the coupled system of equations

In this appendix we reduce the infinite coupled system of equations (3.2) to a finite one with constant coefficients. This simplification is possible due to the special choice of the logarithmic energy spectrum and the application of the rotating wave approximation.

We consider $k$ bosons with the single-particle spectrum (1) and assume first that $k$ is less or equal to the number $n$ of prime factors $p_{j}$ of $N$ defined by (4). Here we represent $N$ as a product

[TABLE]

of $k$ factors $q_{j}$ which may or may not be prime.

We substitute the resonance condition (15) into the term with all $k$ integers $m_{i}=0$ of (3.2), and after little algebra we arrive at the relation

[TABLE]

The first phase factor in the curly brackets of (53) is unity for $n_{j}=q_{j}-L$, and for any permutation of the factors $\{q_{j}\}$. However, according to (41) only the permutation with $q_{1}\leq q_{2}\leq q_{3}\ldots\leq q_{k}$ needs to be taken into account. In the second phase factor in the curly brackets no terms survive because for $N>L^{k}$ the argument of the logarithm can never be unity.

As long as $k<n$ the representation (52) of $N$ may have more than one solution $\{q_{1},q_{2},\ldots q_{k}\}$. If, for example, $N$ has $n=3$ prime factors which have to be distributed onto $k=2$ bosons, as $N=p_{1}\cdot p_{2}^{\,2}$, two solutions ${\bf q}_{1}=(p_{1},p_{2}^{\,2})$ and ${\bf q}_{2}=(p_{1}\cdot p_{2},p_{2})$ exist.

However, if $k=n$ the representation (52) has a unique solution and all the factors $q_{j}=p_{j}$ are prime as follows from the fundamental theorem discussed in Section 2.1. Thus, we are left with the equation

[TABLE]

where $d$ denotes the number of solutions ${\bf q}_{j}$ of (52), that is, the degeneracy discussed in C.

Finally, for $k>n$ all phase factors in the curly brackets of the right-hand side of (53) oscillate and average out according to the rotating wave approximation. Consequently, the equation

[TABLE]

together with the initial conditions (17) immediately yields $b_{\bf 0}(t)\equiv b_{\bf 0}(0)$, i.e. the $k$ bosons remain in their ground state $|{\bf 0}\rangle$.

In the remainder of this article we concentrate on the case $k\leq n$.

Next we address the system (3.2) with ${\bf m}={\bf q}_{j}-{\bf L}$ and derive for every solution ${\bf q}_{j}$ of (52) the equation

[TABLE]

The first phase factor in curly brackets always oscillates since the integer $N^{2}$ in the numerator of the argument of the logarithm is built from $2n$ prime factors where each one obeys444The first phase factor would survive for some integers $k_{j}$ if at least half of the factors $p_{j}$ were equal to the scaling parameter $L$. the condition $p_{j}>L$ as required in Section 2.1.

The second phase factor is unity for all indices $n_{j}=0$ leading to $d$ differential equations

[TABLE]

Equations (54) and (57) describe the main features of quantum dynamics induced by the time- and position-dependent perturbation $\delta V^{(k)}$. In Section 3 we solve these equations.

Appendix C Degeneracy

The integer $d$ introduced in Section 3.3 denotes the number of ways the equation

[TABLE]

can be solved for $k$ positive integers $q_{j}$ regardless of their order. In (58) the number $N$ is built from $n$ prime factor $p_{i}$. Of course, for $k=n$ the factors $q_{j}$ are identical to the prime factors $p_{i}$ due to the fundamental theorem of arithmetics, hence $d=1$ as discussed in section 2.1. Unfortunately, for $k\neq n$ a closed-form expression for $d$ exists only in two situations.

In the first case all $p_{j}$ are mutually different and we have [17]

[TABLE]

where $S(n,k)$ denotes the Stirling number of the second kind, that is the number of partitions of a set of $n$ distinct elements into exactly $k$ nonempty subsets.

In the second case, if all $p_{i}$ are equal i.e. $N=p^{n}$, we have [17]

[TABLE]

where the restricted partition $p_{k}(n)$ denotes the number of partitions of $n$ into at most $k$ parts.

In all other cases $d$ has to be constructed ”by hand”.

Appendix D Matrix elements

In this appendix we provide an elementary example for estimating a Rabi frequency by considering the case where all $n$ non-interacting bosons are in the same single-particle state $|{p-L}\rangle$ with $p^{n}=N$ and ${\cal N}(\{\nu_{\bf p}\})=1$. Now the matrix element (50) is given by

[TABLE]

For the wave functions appearing in the integrand we use the following approximations:

(i) We replace the ground state wave function $\varphi_{0}$ by that of a harmonic oscillator

[TABLE]

with an inverse characteristic length $\alpha_{\rm eff}\equiv(\mu\omega_{\rm eff}/\hbar)^{1/2}$. As detailed in [8] close to the origin the potential $V^{(1)}(x)$ is approximately that of an harmonic oscillator with an effective frequency

[TABLE]

Here we have changed the notation of [8] to conform to Eq. (1) of the present article.

(ii) For the even wave function $\varphi_{\ell}$ we use the simplified WKB approximation [3, 18]

[TABLE]

valid near the center of the trap where

[TABLE]

Under these assumptions, and aside from a prefactor, the right-hand side of the matrix element $W_{\bf 0,p-L}$ given by (61) is transformed into the integral

[TABLE]

where we have introduced the variable

[TABLE]

and the constant

[TABLE]

For even $n$ we use the expansion

[TABLE]

to convert (65) into a sum of standard integrals

[TABLE]

while the coefficients $d_{j}(n)$ can be found e.g. in [19].

Equation (69) shows that asymptotically for large $N$ i.e. for large $\bar{\beta}$ only the term with $j=0$ has to be taken into account.

With the help of the relation

[TABLE]

it is a straightforward but lengthy algebra to derive the result

[TABLE]

with

[TABLE]

For an odd number of factors $n$ a similar procedure leads to

[TABLE]

Because here the factor state $|{\bf p-L}\rangle$ is non-degenerate the Rabi frequency

[TABLE]

is readily determined for both cases, even and odd, respectively.

We thank M. G. Boshier, M. Freyberger and K. Vogel for stimulating discussions on this topic. WPS is grateful to the Hagler Institute for Advanced Study at Texas A&M University for a Faculty Fellowship and to Texas A&M University AgriLife Research for its support. The research of the IQST is financially supported by the Ministry of Science, Research and Arts Baden-Württemberg.