Integrable active atom interferometry

Active SU(1,1) interferometers are designed to enhance phase sensitivity beyond the standard quantum limit. An atomic version of such an interferometer can be constructed by means of a spinor Bose-Einstein condensate with an $F=1$ groundstate manifold in which spin-changing collisions create entangled pairs of $m=\pm1$ atoms. We use Bethe Ansatz techniques to find exact eigenstates and eigenvalues of the Hamiltonian that models such spin-changing collisions. Using these results, we express the interferometer's phase sensitivity, Fisher information, and Hellinger distance in terms of the Bethe rapidities. By evaluating these expressions we study scaling properties and the interferometer's performance in the presence of depletion, i.e., without the idealising assumption necessary for ideal SU(1,1) interferometry.


I. INTRODUCTION
Atom interferometry uses the wave character of atoms, and in particular the superposition principle, to detect phase differences and perform high-precision measurements in a variety of fields, ranging from measurements of the fine structure constant to gravimetry and atomic clocks [1]. The more commonly used passive interferometers use beam splitters that redistribute a conserved number of atoms among two or more modes. Upon recombining the split beams, the phase difference accrued inside the interferometer is measured in the form of interference fringes. The precision to which this phase difference can be measured is specified by the phase sensitivity ∆φ, which is an important characteristic of any interferometer. The larger the number N of atoms measured at the interferometer's output, the lower can the statistical error be pushed and the higher a phase sensitivity can be reached. Assuming at most classical correlations between the (typically uncorrelated) probed events, passive interferometers are known to have phase sensitivities constrained by the standard quantum limit, ∆φ ≥ 1/ √ N, which is essentially a consequence of the central limit theorem [2]. One way to surpass the standard quantum limit is to feed the interferometer with suitably entangled input states. In this case, the Heisenberg limit ∆φ ≥ 1/N, which is a fundamental constraint resulting from Heisenberg's uncertainty principle, may be approached [3].
Another strategy to surpass the standard quantum limit goes back to Yurke,McCall,and Klauder [4], and consists of exchanging passive beam splitters by active components. These active components generate entanglement within the interferometer, and they may have the advantages of being more robust and their experimental realisation being more practical. Originally such active interferometers had been proposed as optical devices, but more recently active atom interferometers have been built [5][6][7][8][9], and their improved sensitivity, beyond the standard quantum limit, has been confirmed. One of these experimental realisations is based on effective three-level systems in a Bose-Einstein condensate of 87 Rb atoms, and uses spin-changing collisions as the active component of the interferometer; see Refs. [6][7][8] for details. The Hamiltonian describing such spin-changing collisions, given in Eq. (2), describes the nonlinear interactions between three species of bosons (corresponding to the three levels effectively taking part in the dynamics).
In this paper we exploit Bethe-Ansatz integrability to compute exact eigenstates and eigenvalues of the Hamiltonian (2) and for semi-analytic calculations of phase sensitivities and other interferometric characteristics. We use these techniques to analyse phase sensitivities, their scaling properties, and other quantities of interest for several variants of the active atom interferometer, which allow us to identify parameter regimes in which the standard quantum limit can be surpassed.

II. BOSONIC THREE-SPECIES HAMILTONIAN WITH SPIN-CHANGING COLLISIONS
A Hamiltonian describing the scattering between atomic hyperfine states |F, κ with κ ∈ {0, ±1} can be written in second-quantized form as where Ψ † κ creates the hyperfine state |F, κ . Sums in (1) are over {0, ±1}, the coefficients Ω αβγδ specify the two-particle interactions between atoms, and T and V denote single-particle kinetic and potential energy terms. Under a number of assumptions [6,10], including δ-interactions due to dominant s-wave scattering and a single-mode approximation for ultracold and strongly confined gases in the absence of atom loss, one can eliminate the spatial part of the Hamiltonian and derive an effective Hamiltonian for only the spin part of the atomic states. Adding the effect of microwave dressing on the atoms [6,11], one obtains a Hamiltonian of the form where a † 0 and a † ± are creation operators of (the spin part of) the three different hyperfine states, λ is a parameter related to the scattering strength, and q quantifies the microwave dressing. The second-last and last terms in the first line of (2) describe so-called spin-changing collisions, transforming a pair of 0bosons into a + and a − boson, or vice versa. See [6,10] for details of the derivation of (2) from (1) and the assumptions and approximations made along the way. Equation (2) is the starting point for all results reported in this paper.
The way in which the three species of bosonic operators occur in the Hamiltonian (2) suggests to introduce the operators and to write the Hamiltonian as The L κ and K κ operators defined in (3a)-(3c) are one-mode, respectively two-mode, representations of the SU(1,1) algebra [12] satisfying a property that will turn out to be beneficial for treating the Hamiltonian by algebraic techniques.

III. SU(1,1) INTERFEROMETRY IN THE NON-DEPLETED REGIME
In 1986, Yurke, McCall, and Klauder pointed out that interferometers can be characterised by certain Lie groups [4]. The group SU(2) naturally characterizes passive interferometers like Mach-Zehnder and Fabry-Perot devices. In the same paper the authors introduced a class of active interferometers characterized by the group SU(1,1), and they also proposed an optical realization of such a device by means of active elements such as degenerate-parametric amplifiers and four-wave mixers. Strikingly, in contrast to SU(2) interferometers, SU(1,1) interferometers can achieve a phase sensitivity of ∆φ = 1/N, surpassing the standard quantum limit ∆φ = 1/ √ N, even without the use of entangled input states [4].
In the Lie-group-theoretic language, the effect of an SU(1,1) interferometer on an input state |in is written as with K x = (K + + K − )/2 and real parameters β and φ. The first and third exponentials in (6) describe Lorentz boosts in opposite directions, which may correspond to the effect of four-wave mixers in an optical SU(1,1) interferometer. The second exponential in (6) describes a rotation, corresponding to free propagation in the interferometer; see [4] for a detailed explanation of the Lie-group-theoretic description of interferometers. The experimental realization of an SU(1,1) interferometer then hinges on the ability to sequentially implement time evolution under the Hamiltonians −K x , K z , and K x with evolution times β, φ, and β, respectively. The Hamiltonian (4), and hence (2), can be used to approximately implement K x by requiring the following conditions to hold.
(i) The number N of bosons is large.
(ii) Most of the bosons are in the κ = 0 state, a † 0 a 0 ≈ N. (iii) q and λ are chosen such that q/λ ≈ 1 − 2N.
Conditions (ii) and (iii) imply that 4λL z −2λ +q ≈ 0 and hence the first term on the right-hand side of (4) vanishes. Conditions (i) and (ii) imply that a † 0 a 0 is large, which justifies replacing a 0 operators by complex numbers [13,14], a 0 = √ N = a † 0 , which implies L ∓ K ± = (N/2)K ± . The resulting Hamiltonian where the subscript "nd" is for nondepleted, acts nontrivially on the occupation of the ±-boson modes and is one of the building blocks of SU(1,1) interferometric sequence (6). Switching off all interactions, free phase evolution occurs according to Upon shifting ω 0 to zero, (8) becomes, besides an irrelevant constant term, proportional to K z , which realizes the second building block of the interferometer (6). Based on a gas of ultracold Rubidium atoms in a regime that is described by the Hamiltonian (2), an experimental realization of an atomic SU(1,1) interferometer has been achieved by using samples of N ≈ 500 atoms, preparing an initial state where all atoms are in the 0 hyperfine state, i.e., by tuning q/λ to satisfy condition (iii), and by restricting the experiment to evolution times that are short enough for the system to remain in the nondepleted regime where a † 0 a 0 ≈ N [6][7][8]. The measured output in these experiments is the ±mode occupation at the end of the interferometric sequence.

IV. PHASE SENSITIVITY AND FISHER INFORMATION
Interferometers detect relative phase shifts, which can be determined by measuring interferences of the output beams. A key figure of merit for any interferometer is the phase sensitivity quantified by (∆φ) 2 . For an interferometer where phase shifts of the observable η are measured, a direct way to calculate the phase fluctuations is via the error propagation formula [4] where expectation values η as well as variances (∆η) 2 = η 2 − η 2 are calculated with respect to the state |out . Both η and (∆η) 2 acquire φ-dependences through the interferometric sequence (6) that generates |out . The phase sensitivity of the ideal SU(1,1) interferometer (6) can be calculated analytically [4]. For an initial state satisfying (9a) and (9b), this analytical result simplifies to [15] (∆φ where is the number of ± atoms "inside" the interferometer, i.e., after evolution under only the rightmost exponential in the interferometric sequence (6). For φ = 2πn with n ∈ N, the phase sensitivity in (12) scales like 1/η 2 1 asymptotically for large η 1 , which surpasses the standard quantum limit 1/η 1 . In Refs. [7,8], the phase sensitivity of an SU(1,1) atom interferometer has been determined experimentally via the error propagation formula (11) by measuring η and (∆η) 2 over a range of phases φ, confirming a precision beyond the standard quantum limit.
As an alternative to the error propagation formula (11), one can use the Cramer-Rao bound [16] (∆φ) 2 1 to estimate the phase sensitivity in terms of the Fisher information [2] where is a probability distribution over the occupation numbers η. |η/2 denotes the Fock state made up of an equal number η/2 of + and − atoms (which is one of the elements of the Fock basis (30) defined later).
To determine F I directly from its definition (15), the derivative of P η with respect to φ needs to be computed, which requires knowledge of the functional dependence of P η on φ. Experimentally determining this functional dependence is challenging, but can be circumvented by resorting to the Hellinger distance [6,17] which does not contain derivatives of P η . Taylor-expanding dH 2 φ,φ+∆ for small ∆, one finds and this leading-order proportionality permits to infer F I from measurements of dH 2 φ,φ+∆ [6,17].

V. SCC INTERFEROMETRY
The theoretical analysis of an SU(1,1) atom interferometer realized by means of the three-species bosonic Hamiltonian (2) reviewed in Sec. III is mostly based on the simplifying assumptions (i)-(iii) made in that section, which give rise to the nondepleted Hamiltonian (7). Effects beyond this idealized Hamiltonian are inevitably present in real experiments but are harder to deal with theoretically. More realistically, the time evolution under K x operators in the idealized interferometric sequence (6) needs to be replaced by evolution under the full three-species Hamiltonian, where the notation H(q, λ) highlights a certain choice of the parameters q and λ in (2). Even the free phase evolution exp(−iuH φ ) under the Hamiltonian H φ in Eq. (8) may be seen as an idealization, as it may be difficult to switch off the interaction strength λ sufficiently fast. Instead, in the experimental realization of the interferometer reported in Refs. [7,8], an approximately free evolution is achieved by switching the interacting Hamiltonian (2) to a large value of q, such that the interferometric sequence is given by with q q, λ. In this limit, the Hamiltonian (2) is given by to leading order in the small parameter λ/q , justifying the claim of a quasifree phase evolution [18]. Unlike in the case of the ideal SU(1,1)-interferometer (6), the phase φ does not feature as a parameter in the interferometric sequence (20). Instead, φ is expected to be approximately proportional to the quasifree evolution time u, at least for those parameter values for which the device indeed functions as an interferometer. We will come back to this issue, and determine the proportionality constant between φ and u, in Sec. VIII.
Effects beyond the ideal SU(1,1) interferometric sequence (6) are certainly harder to deal with theoretically. However, we show in the following that the full Hamiltonian (4) is amenable, for arbitrary parameter values, to an exact analytic treatment by means of Bethe Ansatz techniques. These methods allow us to calculate phase sensitivities and the Hellinger distance essentially analytically for either of the interferometric protocols (19) or (20). Since these protocols are based on the Hamiltonian (4) that models spin-changing collisions (SCC), we refer to both sequences (19) and (20) as SCC interferometry.

VI. BETHE ANSATZ SOLUTION
The Hamiltonian (4) satisfies the conditions of a Richardson-Gaudin model [19][20][21] and its exact eigenstates and eigenvalues can be determined by techniques that fall into the broader class of the algebraic Bethe Ansatz. More specifically, our model is an example of a bosonic pairing model of the type analysed in Ref. [19]. We will make use of results from that study, suitably adapted, in what follows. The appendix contains more information on the derivation of these results, and on how the SCC Hamiltonian fits into the general pairing model formalism.
The starting point for obtaining the solution is the observation that the L + and K + operators in (3b) can be used to span subspaces of the Fock space of the bosonic system, where and each subspace is labelled by the seniorities ν 0 ∈ {0, 1} and ν 1 ∈ N 0 in the multi-index ν = (ν 0 , ν 1 ). The seniorities can be interpreted as the numbers of unpaired 0-bosons and ±bosons, respectively. The subspaces defined in (22) are closed under the application of the L κ and K κ operators defined in (3a)-(3c), We use Richardson's Ansatz [22] |ψ s := for states with a specified seniority ν and number of pairs n.
Our aim is to determine the rapidities e sα in such a way that the states (25) are eigenstates of the Hamiltonian (4). The index s in (25) labels different eigenstates of H, each of which is specified by a different set {e sα } α=1,...,n of rapidities. In Ref. [19] it is shown that if, for all α = 1, . . . , n, the rapidities satisfy the Richardson equation holds with In this way the task of determining the eigenstates and eigenvalues of the Hamiltonian (4) is reduced to finding the roots of a set of n coupled nonlinear algebraic equations (26). These equations have n + 1 (in general different) sets of roots, each of which corresponds to a different eigenstate of H. For the present bosonic case these rapidities are known to be real [19,23]. The numerical calculation of the roots of (26) is done by a mapping to the zeros of special polynomials [24] and computation of these zeros by standard numeric libraries.

VII. OUTPUT STATES IN TERMS OF RAPIDITIES
Our aim is to use the Bethe Ansatz solution (25)-(28c) to compute the state |out at the end of either of the interferometric sequences (19) or (20). This output state, in turn, can then be used to calculate expectation values like out|η|out , or the probabilities P η in (16) that are required for computing the Fisher information (15) or the Hellinger distance (17). The three exponentials occurring in (19) or (20) are most easily evaluated in the respective eigenbases of the Hamiltonians. To this aim, we introduce transformations between the relevant bases, which then allow us to evaluate the sequence of three time evolutions successively and write |out in terms of transformation matrix elements and phase factors.
Guided by the experimental realizations [7,8], we assume the initial state to consist of 2n bosons in the 0-hyperfine state, This state is in the seniority sector ν 0 = 0 = ν 1 of the Hilbert space and, since all evolution operators in (19) and (20) conserve seniority as well as boson number, the system will remain in that sector throughout the interferometric sequence. One distinguished basis is therefore the orthonormalized 2n-boson Fock basis {|k } n k=0 with where N k = 2 2(k−n) [2(n − k)]!(k!) 2 , and it follows that in this basis |in = |0 . The eigenstates (25) of the Hamiltonian (4), which form a second distinguished basis {|ψ s } n s=0 for describing the interferometer, have been derived in Sec. VI. In our selected seniority sector these states can be written by means of a binomial-type expansion as where N s = n k=0 N k |α sk | 2 , with where the elements (34) define an (n + 1) × (n + 1) unitary transformation matrix. By means of this transformation, the initial state (29) can be written in the energy eigenbasis as |in = n s=0 c * sn |ψ s .
As a first application, we calculate the expectation value of the number of ±-bosons produced by the first (rightmost) of the three evolution factors in the interferometric sequences (19) and (20). Using the basis transformation (33) and the fact that the Fock states (30) are eigenstates of η, one obtains the expression where E r and E s are eigenenergies (28a) of H. Analytical expressions for other expectation values can be calculated in a similar manner. The nontrivial ingredients on the right-hand side of Eq. (36) are the c-coefficients (34), which, via α sk defined in (32), depend on the rapidities e sα . The rapidities are determined numerically from the Richardson equations (26) as outlined at the end of Sec. VI. Figure 1 shows the dependence of the expectation value η 1 (36) on the duration t of the first "active" phase of the interferometric sequences (19) or (20). The larger η 1 , the more ±-bosons are available "inside" the interferometer as a resource of entanglement, which is essential for the surpassing of the standard quantum limit. We observe that, in the regime where λ and q are of similar magnitude and for sufficiently long seeding times t, a substantial number of ±-bosons is produced, fluctuating roughly around N/2. For large values of q, the production of ±-bosons is strongly suppressed and η oscillates in an approximately sinusoidal fashion around a value much smaller than N/2. This confirms, as discussed towards the end of Sec. IV, that time evolution under the Hamiltonian H(q, λ) with q λ can be used to approximate free phase evolution, as proposed in the interferometric sequence (20).
Similar to the derivation of Eq. (36), the output state at the end of the full interferometric sequence (19) can be evaluated by repeated transformations between the Fock basis (30) and the energy basis (31), yielding For the interferometric sequence (20), in which free phase evolution is replaced by evolution under H(q , λ) with q large, we transform, in addition to the eigenbasis (31) of H(q, λ), also to the eigenbasis {|ψ s } n s=0 of H(q , λ). Expansion coefficients α and basis transformation coefficients c are defined analogous to their non-primed counterparts (32) and (34). With this notation, by repeated transformations between the Fock basis (30), the energy basis (31), and the primed energy basis, the output state at the end of the interferometric sequence (20) can be written as By expanding |ψ q in (37) or (39) in the Fock basis (30), which is an eigenbasis of η, we can evaluate the expectation values η 0 ≡ out| η |out and η ≡ out | η |out at the end of the respective interferometric sequences. Figure 2 (left) shows, for two choices of the seeding time t, η 0 as a function of the dwell time u. In the case of a short seeding time t = 0.006 (blue) we observe clear interference fringes with an oscillation period of 0.003. As in the case of the ideal SU(1,1) interferometer (Eq. (9.28) in Ref. [4]) the fringes are approximately sinusoidal, which makes this regime particularly suitable for interferometry. For longer seeding time t = 0.03 (orange) the same fundamental period of 0.003 is observed, but with higher frequency contributions superimposed. Figure 2 (right) shows similar data, but for the expectation value η calculated for the interferometric The blue curve uses a seeding time t = 0.006, which, according to Fig. 1, corresponds to a production of approximately 1.5 pairs of ±-bosons at the end of the first (rightmost) exponential in the sequence (19). The orange curve is for t = 0.03, corresponding to the creation of approximately 21 pairs. Right: As in the left plot, but showing the expectation value η ≡ out | η |out calculated for the interferometric sequence (20) with quasifree phase evolution and q = 1000.
sequence (20) with quasifree phase evolution. In this case, the strict periodicity is spoiled, which is particularly evident for the example with the larger seeding time t (orange).

VIII. PHASE SENSITIVITY IN TERMS OF RAPIDITIES
To calculate the phase sensitivities of the interferometric sequences (19) or (20), we relate the dwell time u to the phase φ. This is achieved by numerically determining η 0 or η as a function of u. If, as for the parameter values λ, q, and q in Fig. (2) (right), a roughly periodic dependence on u is observed, then the angular frequency Ω of the oscillatory behaviour can be read off, and we identify φ = Ωu. The output states |out (37) and |out (39) depend on the dwell time u, and hence on the interferometric phase φ, only through the exponentials in the coefficients (38) and (40), respectively. Taking derivatives with respect to φ, as required for the calculation of the phase sensitivity (11), can therefore be done analytically. Using the output state (19) with coefficients (38) to calculate the expectation value out| η |out and then taking its derivative with respect to φ, the phase sensitivity (∆φ) 2 (11) can be expressed in terms of the rapidities e sα , For the interferometer with quasifree time evolution (20) the same formula holds with coefficients x replaced by x . Figure 3 (left) shows the sensitivity (∆φ) 2 of the quasifree interferometer (20) as a function of the phase φ. The sensitivity exhibits pronounced minima when the phase is around multiples of 2π, affirming that this is where the interferometer, like its ideal SU(1,1) counterpart, performs at its most precise. We will choose φ = 0 for all numerical explorations from here on. As illustrated in the plot, the minimum value of (∆φ) 2 is well below the standard quantum limit 1/η 1 and approaches fairly closely, but does not quite reach, the sensitivity 1/[η 1 (η 1 + 2)] of the ideal SU(1,1) interferometer (12). Similar behavior is found for the interferometric sequence (19) with free phase evolution (not shown).

IX. HELLINGER DISTANCE IN TERMS OF RAPIDITIES
As an alternative method for estimating the phase sensitivity, we show in the following how to express the Hellinger distance (17) in terms of the rapidities e sα , which in turn yields an estimate of (∆φ) 2 via Eqs. (14) and (18). Similarly, the Fisher information (15) could be expressed in terms of the rapidities and then used to estimate the phase sensitivity, but we do not pursue this approach further.
To compute the probabilities P η (φ) = | η/2|out | 2 in the definition (17) of the Hellinger distance, we write x q |ψ q = n p,q=0 where we have used Eqs. (37) and (33). Calculating the modulus squared of this result and plugging it into (17), one obtains an expression for the Hellinger distance in terms of the rapidities, which, while being lengthy, is fairly straightforward to evaluate numerically. Based on this computation and making use of Eqs. (14) and (18), we define as a proxy for the phase sensitivity (∆φ) 2 . In the following we investigate the dependence of ( ∆φ) 2 on parameters in the Hamiltonian (2) and in the interferometric sequences (19) or (20), with the aim of singling out the parameter regime of optimal performance of the interferometer. We checked that the numerical results reported in this section are insensitive to the choice of the parameter ∆ in the definition of the Hellinger distance (17), as long as it is much smaller than 2π. Figure 3 (right) shows the sensitivity ( ∆φ) 2 as a function of the phase φ.
The plot uses the same parameter values as for (∆φ) 2 in Fig. 3 (left). While the plots show minor differences, the results are generally in good agreement, confirming that ( ∆φ) 2 is a valid, and numerically more stable, proxy for the phase sensitivity (∆φ) 2 .
The main interest in active interferometers, like the ideal SU(1,1) interferometer (6) or the SCC interferometers (19) and (20), lies in their phase fluctuations having the potential to surpass the standard quantum limit, and potentially approach the Heisenberg limit, without the need for entangled input states |in . Both the standard quantum limit ∼ 1/η 1 and the Heisenberg limit ∼ 1/η 2 1 are expressed in terms of the number η 1 of seeded ±-bosons after the first (rightmost) exponential in the interferometric sequences (19) or (20). To assess the influence of the seeding on the performance of the interferometer, we show in Fig. 4 the phase sensitivity ( ∆φ) 2 as a function of η 1 for the free interferometric sequence (19) (red) and the quasifree interferometric sequence (20) (blue). In both cases the sensitivity decays monotonically with η 1 , showing phase fluctuations that behave qualitatively similar to, but are slightly larger than, those of the ideal SU(1,1) interferometer (green line in Fig. 4). When operating the interferometer at short seeding times t, which results in small values of η 1 , we find, as expected, a phase sensitivity very close to that of the ideal SU(1,1) interferometer. For larger values of η 1 , deviations from the ideal case become visible (see insert of Fig. 4), but The red line shows ( ∆φ) 2 for the interferometric sequence (19) with free phase evolution for parameter values N = 100, λ = 1, and q = 4/3. The blue lines are obtained with the same parameter values, but using the interferometric sequence (20) with quasifree phase evolution and parameter values q = 125, 250, 500, and 1000 (top to bottom). All those curves are well below the standard quantum limit 1/η 1 (orange). The phase sensitivity of the free interferometer (red) is very close to the Heisenberg limit 1/[η 1 (η 1 + 2)] (green) of the ideal SU(1,1) interferometer (12), and the quasifree interferometer (blue) approaches that limit with increasing q . The inset shows the same data on a logarithmic scale.

X. CONCLUSIONS
tion spoils the periodicity of the interferometric fringes and operates at slightly inferior phase sensitivity compared to the case of free phase evolution, our results clearly indicate that the SCC interferometer with quasifree phase evolution can successfully function with a phase sensitivity well below the standard quantum limit and, for suitable parameter values, close to the Heisenberg limit accessible by the ideal SU(1,1) interferometer proposed by Yurke, McCall, and Klauder [4]. While we exploited integrability of the SCC Hamiltonian in order to elegantly and efficiently calculate quantities of interest by expressing them in terms of the Bethe rapidities, integrability does, to our understanding, not affect the performance characteristics of the SCC interferometer. However, the techniques developed in the present paper are general and may potentially be of use for a broad range of applications in quantum dynamics as well as equilibrium physics. A sudden quench protocol and the subsequent equilibration dynamics under the SCC Hamiltonian is an example where integrability is expected to play an important role, and where the Bethe-Ansatz techniques developed in this paper may be brought to use in future work.
The numerical evaluations of the Bethe rapidities, or quantities derived from them, reported in this paper are for moderate boson numbers of N = 100. This particle number can be reached, and exceeded, on a regular desktop computer at the time of writing. Besides the polynomial mapping we used for the calculation of the Bethe rapidities, other numerical approaches have been reported in the literature, and also more efficient methods for the computation of overlaps of Bethe eigenstates are known [25]. Here, we did not make a concerted effort to reach larger sizes by following any of these and instead opted to focus on conceptual aspects, but we expect that at least an order of magnitude in system size can be gained with a bit of effort. To reach even larger system sizes, and possibly even perform analytical calculations in the large-N limit, the SCC Hamiltonian (4) expressed in terms of the generators of the group SU(1,1)⊗SU(1,1) constitutes a promising starting point for Holstein-Primakoff expansions or other analytic approaches.
Note added: When adding the finishing touches to the paper we became aware of the recent preprint Ref. [26] that uses integrability for the analysis of a passive atom interferometer.

ACKNOWLEDGMENTS
In an early phase of the project the authors have benefited from helpful discussions with Markus Oberthaler, Frederik Scholtz, and Helmut Strobel.

Appendix: The SCC Hamiltonian as a solvable pairing model
Here we provide a brief account of the origin of the Richardson equation in (26), and of the expressions for the eigenvalues of the SCC Hamiltonian (4) in Eqs. (28a)-(28c). The key observation is that the SCC Hamiltonian falls within the class of solvable models studied in Ref. [19]. In the language of that paper we are dealing with a particular two-level bosonic pairing model. We will label the levels using l ∈ {0, 1}, where l = 0 corresponds to the 0-boson mode, and l = 1 to the ±boson modes. The degeneracies of the levels are denoted by Ω 0 = 1 and Ω 1 = 2. The L ± and K ± operators in (3a) and (3b) create and destroy pairs of bosons in the l = 0 and l = 1 levels respectively. Following Ref. [19] we introduce the operators where g, X ll , and Y ll are scalar parameters. The choice of these parameters is dictated by two requirements, namely that [R 0 , R 1 ] = 0 and that the SCC Hamiltonian (4) can be expressed as a function of these two commuting operators. The former condition is met by setting with η 0 and η 1 arbitrary unequal real numbers. This yields the so-called rational model of Ref. [19]. The second requirement is satisfied by choosing g = 2λ/q and η 0 = −η 1 = 1/2, (A.3) which allows the Hamiltonian (4) to be expressed as This implies that R 0 and R 1 are conserved charges of H, and so the eigenstates of H are the simultaneous eigenstates of R 0 and R 1 , R 0 |ψ s = r 0s |ψ s , R 1 |ψ s = r 1s |ψ s . (A.5) Here s labels the various eigenstates. Richardson's Ansatz [22] provides an explicit form for these states as Inserting η 0 = −η 1 = 1/2 into these expressions yields Eqs. (26) and (28a)-(28c).