HHydrogenic entanglement
Sofia Qvarfort , , Sougato Bose , Alessio Serafini QOLS, Blackett Laboratory, Imperial College London, SW7 2AZ London, UnitedKingdom Department of Physics and Astronomy, University College London, Gower Street, WC1E6BT London, United KingdomE-mail: [email protected]
Abstract.
Is there any entanglement in the simplest ubiquitous bound system? We studythe solutions to the time-independent Schr¨odinger equation for a Hydrogenic system anddevise two entanglement tests for free and localised states. For free Hydrogenic systems,we compute the Schmidt basis diagonalisation for general energy eigenstates, and for aHydrogenic system localised to a three-dimensional Gaussian wavepacket, we demonstratethat measuring its second moments is sufficient for detecting entanglement. Our resultsapply to any system that exhibits Hydrogenic structure.
1. Introduction
Are the electron and proton in a Hydrogen atom entangled? In 1926, Erwin Schr¨odingersuccessfully predicted the spectral energies of Hydrogen by solving the wave-equation for anelectron wavefunction in a potential-well created by the positively charged nucleus [1, 2].While the Hydrogenic solutions to Schr¨odinger’s equation have now been known for almosta century, the question of whether the two subsystems are entangled has hardly beeninvestigated. This is most likely due to the prevalence of the Born–Oppenheimer (BO)approximation [3], which explicitly assumes that the motion of two subsystems with vastlydifferent masses can be treated separately. While it sometimes follows that entanglementis explicitly removed in the process of applying the BO approximation, this is in fact onlytrue if the approximation is exact [4]. Indeed, entanglement can be retained as degrees ofthe approximation is applied, and can even be used as a measure of the validity of the BOapproximation [5].In this work, we forego the BO approximation completely in order to study theentanglement of the exact Hydrogenic solutions. The ubiquity of Hydrogen in physical,chemical, and biological systems makes it is one of the most well-studied physical systems,and irrespective of the question as to whether entanglement – if present in Hydrogenic systems– has any applications, it is important to know whether such a basic textbook entity of physicsis intrinsically entangled. a r X i v : . [ qu a n t - ph ] S e p ydrogenic entanglement 𝑟 (a) 𝑟 (b) Figure 1:
Sketch of a bipartite Hydrogenic system. The relative degree-of-freedom wavepackettakes the Hydrogenic solutions and is parametrised by the variable r = | r − r | , where r and r are position vectors for the two subsystems. The centre-of-mass wavepacket is either (a) a freewavepacket, or (b) a localised three-dimensional Gaussian wavepacket. In recent years, the study of quantum information-processing tasks has demonstratedthe importance of entanglement to quantum computing [6], quantum cryptography [7],and quantum sensing [8]. A number of fundamental questions, such as whether Natureis described by collapse theories [9], or whether gravity is a quantum force [10], are aidedby the quantification and detection of entanglement [11, 12, 13, 14, 15]. The ability toexperimentally determine entanglement for Hydrogenic systems could have implications forthe application of Hydrogenic atoms, artificial atoms [16, 17], exiton states [18, 19], or possiblyeven mesoscopic systems, such as levitated nanospheres [20], which could be engineered tointeract via a central potential [21]. In fact, the entanglement of similar systems, such aselectron-electron entanglement in Helium and Helium-like atoms [22, 23, 24], and proton-proton entanglement in Hydrogen molecules [25], has been recently explored.A rudimentary step in certifying the entanglement in Hydrogenic systems has beentaken by Tommasini et al . [26], who showed that entanglement in a free
Hydrogenic groundstate system can be verified by diagonalising the state through the Fourier transform andcomputing its Schmidt coefficients. In the laboratory, however, the assumption that theHydrogenic system is free is no longer accurate. Rather, any prepared state will be localisedto a finite spatial volume, and thus can no longer be diagonalised by the Fourier transform.We must therefore devise an appropriate entanglement test that holds even for localisedstates, which is what we set out to do here.In this work, we first perform the Schmidt basis diagonalisation for arbitrary energyeigenstates of the Hydrogenic system with a free centre-of-mass wavefunction. We find thatthe spread in the Schmidt basis scales with n − , where n is the principal quantum number,which implies that Hydrogenic systems found in highly excited states are less entangled.These results apply to any system that can be approximated as free, which includes ydrogenic entanglement
2. The Hydrogenic solutions
We begin by considering two systems with a joint wavefunction Ψ( r , r ) that is parametrisedby the position vectors r and r . We allow the two systems to interact via a central potentialof the form V ( r ) = α/ | r | , where r = r − r and α > ‡ . The time-independent Schr¨odinger equation that describes thesystem is given by (cid:20) − (cid:126) m ∇ − (cid:126) m ∇ − V ( r ) (cid:21) Ψ( r , r ) = i (cid:126) ∂∂t Ψ( r , r ) . (1)As is well known, this equation can be split into a relative and a centre-of-mass equationresulting from the following change of coordinates: R = m r + m r m + m , r = r − r , (2)where m and m are the masses of the particles. Note that the centre-of-mass and relativemomenta P = ( p + p ) and p = ( m p − m p ) / ( m + m ) are conjugate to R and r .The momentum coordinates must be chosen such that the canonical commutator relationsare preserved at all times. See Appendix A for details on how the separation of variables isperformed. The separated Schr¨odinger equations become − (cid:126) M ∇ R ϕ ( R ) = E CM ϕ ( R ) , (3) (cid:20) − (cid:126) µ ∇ r − V ( r ) (cid:21) ψ ( r ) = E rel ψ ( r ) , (4) ‡ We only consider attractive potentials, which means that V ( r , r ) is positive when the minus sign isexplicitly included in the Schr¨odinger equation in Eq. (A.2) ydrogenic entanglement µ = m m / ( m + m ) is the reduced mass of the two-particle system, and where M = m + m is the total mass. The two decoupled time-independent Schr¨odinger equationsallow for solutions comprised of a relative wavefunction ψ nlm ( r ), where n, l , and m are theprincipal, angular and magnetic quantum numbers, and a centre-of-mass wavefunction ϕ ( R )that takes the form of an infinite plain wave. Later on, we shall also consider the more realisticcase where the centre-of-mass wavepacket is localised, and study how the entanglement isaffected by its finite volume in space.The full state Ψ( r , R ) of the system can now be written in terms of the relative coordinate r and centre-of-mass coordinate R as the following separable state:Ψ( r , R ) = ψ nlm ( r ) ϕ ( R ) . (5)The Hydrogenic wavefunction ψ nlm ( r ) is given by ψ nlm ( r ) = R nl ( r ) Y ml ( θ, φ ) , (6)where R nl ( r ) is the radial wavefunction, which depends on the Laguerre polynomials, and Y ml ( θ, φ ) are the spherical harmonics. These functions are described in full in Appendix A.While the state in Eq. (5) is separable in the { r , R } basis, it cannot be written as aseparable state in the original { r , r } basis. As a result, the two subsystems are entangled,and it remains to determine how such entanglement can be detected.Before we proceed to examine the entanglement of the system, we present the momentumrepresentation of the Hydrogenic wavefunctions, since these will be of great use to us later.They correspond to the Fourier transform of ψ nlm ( r ), and are given by [29]:˜ ψ nlm ( k ) = F nl ( k ) Y ml ( θ, φ ) , (7)where F nl ( k ) is a function of the Gegenbauer polynomials [30], and Y ml ( θ, φ ) are againthe spherical harmonics, and ( k, θ, φ ) are the wavevector coordinates. More details on thissolution and especially the Gegenbauer polynomials can again be found in Appendix A.
3. Entanglement of the free Hydrogen atom
We start by considering the case when the centre-of-mass wavefunction ϕ ( R ) in Eq. (5)corresponds to that of a free wavefunction (see Figure 1a): ϕ ( R ) = 1 √ V e − i K · R , (8)where V is the free space normalisation volume, which we take to infinity at the end of eachcalculation, R is given in Eq. (2) and K = ( k + k ) is the centre-of-mass wavevector of thestate. The full eigenstate, including the Hydrogenic part is then given byΨ( r , r ) = Ψ( r , R ) = ψ nlm ( r ) 1 √ V e − i K · R . (9) ydrogenic entanglement K , one can perform a Galillean transformation into a frame where the centre-of-massis stationary, corresponding to the product of local unitary operations e im K · r / ( m + m ) ⊗ e im K · r / ( m + m ) . As one should expect, the entanglement cannot therefore depend on K ,which can be set to zero without loss of generality.Transforming into a frame with K = 0 leaves us with the following state:Ψ( r , R ) = 1 √ V ψ nlm ( r ) . (10)In the remainder of this section, we compute the Schmidt basis associated with theseeigenstates and generalise the results of Ref. [26] to arbitrary energy eigenstates. We beginby setting out a few preliminaries concerning the Schmidt basis decomposition for continuousvariable systems. A necessary and sufficient test to verify whether a pure state is entangled or not consists incasting it in its Schmidt basis and computing its Schmidt rank. If the Schmidt rank is greaterthan one, the state is entangled. The Schmidt basis and the Schmidt rank are straight-forwardto introduce for discrete systems, but they have also been studied for continuous variablesystems in full generality [31]. We here provide a short review of the Schmidt diagonalisationfor discrete systems and its generalisation to continuous systems, with particular referenceto the wavefunction notation we are adopting.It may be shown that a choice of local bases exist such that an arbitrary bipartite state | Ψ (cid:105) can be expanded in a basis formed by pairs of distinct, orthonormal local vectors, whichwe shall call | w k (cid:105) = | u k (cid:105) ⊗ | v k (cid:105) : | Ψ (cid:105) = (cid:88) k c k | w k (cid:105) = (cid:88) k c k ( | u k (cid:105) ⊗ | v k (cid:105) ) , (11)where the Schmidt coefficients c k may be taken to be positive and real. The local Schmidtbases are nothing but the eigenbases of the local density operators.It is interesting to note that the local state of a subsystem, defined as ˆ (cid:37) = Tr [ | Ψ (cid:105) (cid:104) Ψ | ] = (cid:80) k c k | u k (cid:105) (cid:104) u k | , satisfies the following eigenvalue equation:( ˆ (cid:37) ⊗ ˆ ) | w k (cid:105) = c k | w k (cid:105) . (12)Sometimes it is easier to solve this equation for the basis | w k (cid:105) , and thereby reconstruct thestate in Eq. (11) than to perform the usual Gram–Schmidt diagonalisation procedure.The link between wavefunction notation and bra-ket notation is the inner product, wherethe wavefunctions are scalar functions of the state space, allowing to represent the latter inthe (improper) eigenvectors of the position operator. For example: ψ ( x ) = (cid:104) x | ψ (cid:105) , (13) ydrogenic entanglement | x (cid:105) are the position eigenstates and | ψ (cid:105) is the state of the system.In order to introduce the continuous analogue of the Schmidt basis note that, giventhe wavefunction for a bipartite state Ψ( x, y ), where x and y denote the positions of thesubsystems (assumed here to be one-dimensional, for simplicity), the analogue of the densitymatrix is given by (cid:37) ( x (cid:48) , y (cid:48) , x, y ) = Ψ ∗ ( x (cid:48) , y (cid:48) )Ψ( x, y ) , (14)which is normalised as follows: (cid:90) d x (cid:90) d y (cid:37) ( x, y, x, y ) = 1 . (15)The density matrix in the position representation (cid:37) ( x (cid:48) , x ) of the traced-out subsystem isthen given by (cid:37) ( x (cid:48) , x ) = (cid:90) d y Ψ ∗ ( x (cid:48) , y )Ψ( x, y ) . (16)We concluded before that finding the Schmidt basis for discrete states amounts to solvingthe eigenvalue equation in Eq. (12). The continuous analogue of the eigenvalue problem isthe following integral equation [31], otherwise known as a Fredholm equation of the secondkind [32]: (cid:90) d x (cid:48) (cid:37) ( x (cid:48) , x ) φ i ( x (cid:48) ) = λ i φ i ( x ) , (17)where { λ i } is the set of eigenvalues of (cid:37) ( x (cid:48) , x ) parametrised by the index i .Since we in this work are interested in bipartite systems in three spatial dimensions, wewrite the density matrix of the free Hydrogenic state in Eq. (10) as (cid:37) ( r (cid:48) , r (cid:48) , r , r ) = Ψ ∗ ( r (cid:48) , r (cid:48) ) Ψ( r , r ) . (18)This allows us to rewrite the Fredholm equation in Eq. (17) in terms of a three-dimensionalsubsystem state (cid:37) ( r (cid:48) , r ), given by (cid:37) ( r (cid:48) , r ) = (cid:90) d r Ψ ∗ ( r (cid:48) , r ) Ψ( r , r ) , (19)whence the Fredholm equation becomes (cid:90) d r (cid:48) (cid:37) ( r (cid:48) , r ) φ k ( r (cid:48) ) = λ k φ k ( r ) , (20)where we have indexed φ k with the (possibly) continuous index k . In terms of the full state,the Fredholm equation in Eq. (17) becomes (cid:90) d r (cid:48) d r Ψ ∗ ( r (cid:48) , r ) Ψ( r , r ) φ k ( r (cid:48) ) = λ k φ k ( r ) . (21)To obtain the diagonalisation in three dimensions, Eq. (21) must be solved for φ k . Weanticipate that finding the solution in general will be an extremely challenging task. Inthe next Section, however, we show that Eq. (21) has a simple solution when the system isinvariant under spatial translations. ydrogenic entanglement Before considering the Schmidt diagonalisation of the free Hydrogenic wavefunctions, weproceed with a proof showing that translationally invariant (“ homogeneous ”) states arediagonalised by the Fourier transform, which was first demonstrated in Ref. [26]. We definetranslationally invariant by the fact that displacing the two subsystems by the same distancein space leaves the wavefunction invariant, up to a global phase.A translationally invariant wavefunction must depend on the variable r alone, and noton R , since the former is invariant under an equal translation of r and r whilst the lattermay be modified arbitrarily by such translations. This implies that the local density function (cid:37) ( r (cid:48) , r ), whose eigenbasis determines the Schmidt decomposition, is a function only of thedifference of its arguments, since Eq. (19) in this case reads (cid:37) ( r (cid:48) , r ) = (cid:90) d r Ψ ∗ ( r (cid:48) − r )Ψ( r − r ) (22)= (cid:90) d r Ψ ∗ ( r )Ψ( r − r (cid:48) + r ) = (cid:37) ( r (cid:48) − r ) , where we have made the variable substitution r → r (cid:48) − r , treating r (cid:48) as a constant. TheFourier transform and inverse of (cid:37) ( r − r (cid:48) ) may then be expressed in terms of a single(three-dimensional) variable: ˜ (cid:37) ( k ) = (cid:90) d r e − i k · r (cid:37) ( r ) , (23) (cid:37) ( r ) = 1(2 π ) (cid:90) d k e i k · r ˜ (cid:37) ( k ) , (24)where we have defined r = ( r − r (cid:48) ) and its conjugate variable k = ( k − k (cid:48) ) § .We then recall the Fredholm eigenvalue equation for continuous systems given in Eq. (20).Inserting the translationally invariant state and the momentum eigenstate ansatz φ k = e i k · r (cid:48) ,we find: (cid:90) d r (cid:48) (cid:37) ( r − r (cid:48) ) e i k · r (cid:48) = 1(2 π ) (cid:90) d k (cid:48) ˜ (cid:37) ( k (cid:48) ) e i k (cid:48) · r (cid:90) d r (cid:48) e i ( k − k (cid:48) ) · r (cid:48) = ˜ (cid:37) ( k ) e i k · r , (25)where we have used the Fourier transform of (cid:37) ( r ) shown in Eq. (24). This demonstrates thatthe local density operator is diagonal in the local momentum eigenbasis, which thus form itsassociated Schmidt basis, with Schmidt coefficients given by the Fourier transform ˜ (cid:37) ( k ). Inother words, a translationally invariant state is always diagonalised by its Fourier transform. We have shown that the local density operator of homogeneous state is always diagonalisedby the Fourier transform. A free Hydrogenic system, whose wavefunction is shown in § Note that we have added the bar to the local variables to differentiate them from the relative and centre-of-mass position and momentum variables introduced earlier ydrogenic entanglement r = ( r − r ) between the two systems.In order to compute the Schmidt coefficients of a Hydrogenic eigenfunction, and hencequalify its entanglement, we need first to trace out one of the two particles, as follows (cid:37) ( r − r (cid:48) ) = (cid:90) d r Ψ ∗ ( r (cid:48) , r )Ψ( r , r )= 1 V (cid:90) d r ψ ∗ nlm ( r (cid:48) − r ) ψ nlm ( r − r )= 1 V (cid:90) d y ψ ∗ nlm ( y ) ψ nlm ( r + y ) . (26)Here, through the substitution y = r (cid:48) − r , we have again shown that the traced-out Hydrogenic state is translation invariant, since it only depends on the local variable r = r − r (cid:48) . As we proved in the previous Section, the Fourier transform ˜ (cid:37) ( k ) correspondsto the Schmidt basis of the Hydrogenic system, where k = k − k (cid:48) is the Fourier complementto r .The state is separable if and only if the probability distribution ˜ (cid:37) ( k ) is a delta-function,which provides one with a criterion to test entanglement in this context. This correspondsto the discrete analogue of having one single non-zero Schmidt coefficient.Notice that, even if the continuous analogue of the Schmidt coefficients are known, anactual quantification of the entanglement is difficult in this case. The local von Neumannentropy is in fact not a good quantifier for such continuous systems as, being the Shannonentropy of the local state’s spectrum, it is ill defined if the latter is continuous, as is always thecase for our homogeneous systems, which are diagonal in the (continuous) local momentumbasis. The differential entropy [33], sometimes introduced as the continuous analogue of thevon Neumann entropy, is problematic (it can even be negative!), and lacks the propertiesto qualify as a bona fide entanglement monotone (cid:107) . It may even be shown that, in infinitedimension, states with diverging von Neumann entropy may be found arbitrarily close, inthe trace norm topology, to any quantum state [34].The linear entropy S Lin is sometimes favoured as an entanglement quantifier because itis easier to compute. However, besides being an unwieldy quantity in our case, the locallinear entropy is not endowed with an operational interpretation.It will therefore be convenient to adopt a different quantifier of entanglement, whichallows us to illustrate the behaviour of Hydrogenic entanglement more clearly. We optto evaluate an entanglement quantifier already suggested in [26], namely one given by thestandard deviation of the Schmidt function:∆ k = (cid:113)(cid:10) k (cid:11) − (cid:10) k (cid:11) . (27)Eq. (27) does in a sense quantify the deviation from the delta function that characterisesseparable states: if ∆ k = 0, the state is separable, and if ∆ k >
0, the state is entangled. (cid:107)
The differential entropy should not be confused with the relative entropy , for which the discrete andcontinuous version have identical properties ydrogenic entanglement k . The continuous analogue of computing expectation valueswith the trace operation is given by (cid:10) ¯ k (cid:11) = (cid:90) d r (cid:90) d r (cid:48) δ ( r − r (cid:48) ) (cid:0) − i ∇ ( r − r (cid:48) ) (cid:1) (cid:37) ( r − r (cid:48) ) . (28)We then insert the Fourier transform in Eq. (24) to find (cid:10) ¯ k (cid:11) = − i (2 π ) (cid:90) d r (cid:90) d r (cid:48) (cid:90) d k δ ( r − r (cid:48) ) ˜ (cid:37) ( k ) ∇ ( r − r (cid:48) ) e i ¯ k · ( r − r (cid:48) )= V (2 π ) (cid:90) d k k ˜ (cid:37) ( k ) , (29)where the integration volume V appears because we integrate over all space. This expectationvalue (cid:10) ¯ k (cid:11) is calculated by integrating over all vectors, which averages to zero: (cid:10) k (cid:11) = V (2 π ) (cid:90) d k k ˜ (cid:37) ( k ) = 0 . (30)The variance (cid:10) k (cid:11) , on the other hand, is given by (cid:10) k (cid:11) = V (2 π ) (cid:90) d k k ˜ (cid:37) ( k ) . (31)Using the momentum representation of the Hydrogenic wavefunctions, shown in Eq. (7), andrealising that the free system is proportional to V − , we find (cid:10) k (cid:11) = 1(2 π ) (cid:90) ∞ d k k [ F nl ( k )] (cid:90) π d θ sin θ (cid:90) π d φ | Y ml ( θ, φ ) | . (32)The integrals over θ and φ satisfy the normalisation condition for the spherical harmonics(see Eq. (A.19) in Appendix A), and we are left with (cid:10) ¯ k (cid:11) = 1(2 π ) (cid:90) ∞ d¯ k ¯ k | F nl (¯ k ) | . (33)The integral over ¯ k returns a well-known result from atomic physics. It admits the standardsolution [29]: (cid:10) ¯ k (cid:11) = 1(2 π ) n a . (34)We prove this relation, which extends what was previously known about entanglement forthe Hydrogenic ground state to any energy eigenstate [26], in Appendix B.2.We conclude that the standard deviation of the local momentum, related to the localwavevector by ¯ p = (cid:126) ¯ k , reads ∆¯ p = 1(2 π ) / (cid:126) na = 1(2 π ) / αµn (cid:126) , (35) ydrogenic entanglement a = (cid:126) / ( αµ ) was inserted (recall that µ isthe reduced mass of the system and that α is the potential interaction strength, which hasdimensions of an energy times a length).The standard deviation of the continuous local spectrum only depends on the principalquantum number n , and not on l and m . This is the case since the operations thattransform sectors with the same n and different l and m are local unitary operations(which is particularly clear for the azimuthal number m , whose different values are relatedby rotations), and thus cannot affect properties related to the Schmidt spectrum andentanglement.On the other hand, the distribution of local momentum is more spread out for lowerprincipal quantum numbers, which is non-trivial and suggests that, in practice, lowerquantum states would be more favourable for the detection of quantum correlations inHydrogenic systems. Notice that the fact that we can qualify the presence or absence ofentanglement so liberally (by any positive value of the standard deviation!) is just an artefactof the fact that our theoretical finding applies to ideal (pure) eigenstates.Further, our quantifier grows with the interaction strength, as one should expect, andwith the reduced mass, which indicates that, at a specific given total mass , the entanglementis maximum for two equally distributed masses. However, the eigenstate of a Hydrogen atomwould be about twice more entangled than the corresponding positronium eigenstate, sincethe large mass of the proton nucleus factors out of µ and leaves µ = m e , where m e is theelectron mass. Compare this with the reduced mass of the positronium eigenstate (a boundstate of one electron and one positron), which reads µ = m e m p / ( m e + m p ) = m e /
2. We alsonote that the expression in Eq. (35) is not independent of the choice of Fourier normalisation,however, as this is merely a scaling factor, the choice does not significantly affect the results.Before we proceed to investigate localised Hydrogenic states, let us briefly return tothe linear entropy S Lin of a Hydrogenic system. We find that the linear entropy convergesfor a free Hydrogenic systems, although, as mentioned above, it does not admit a directoperational interpretation as a quantifier of the entanglement. We merely provide this resultfor completeness.In Appendix C, we show that the linear entropy S Lin for a free Hydrogenic system isgiven by S Lin = 1 − V (cid:90) d¯ k (cid:2) F nl (¯ k ) (cid:3) [ Y lm ( θ, φ )] . (36)While we find that the integral in Eq. (36) converges (see the rather lengthy expression inEq. (C.24)), we conclude that when V → ∞ , the linear entropy tends to S Lin = 1 for allenergy eigenstates.
4. Entanglement of a localised Hydrogenic system
The approximation of the system as free breaks down when the overall centre-of-masswavepacket is of a similar width compared with the characteristic length scale a of the ydrogenic entanglement R is given by ϕ ( R ) = 1 π / b / e − R · R / (2 b ) , (37)where b has units of length and encodes the standard deviation of the wavepacket at thespecific moment in time for which we consider the system. The full state of a localisedHydrogenic system therefore becomesΨ( r , R ) = ψ ( r ) ϕ ( R ) . (38)It is clear that the state appears entangled even in the Gaussian wavepacket alone, since ϕ ( R ) in Eq. (37) cannot be written as a product of functions of r and r . We thereforeexpect contributions to entanglement from both the relative and centre-of-mass degrees offreedom.We proceed now to study the entanglement criterion for this state. For simplicity, werestrict the analysis to the case of equal masses, such that m = m . In order to devise an entanglement test for a localised Hydrogenic state, we adopt somewell-established techniques developed in the study of quantum information with continuousvariables [35]. In particular, a necessary and sufficient entanglement criterion based on thesecond order statistical moments of the canonical operator exist for two-mode Gaussianstates [27, 36, 37], which is however sufficient to detect the entanglement of any quantum ydrogenic entanglement ¶ .Such a criterion is based on a general necessary test for the positivity of the partialtranspose (PPT) at the level of second moments, and on the observation that, in turn, PPT isnecessary for the separability of arbitrary states [38, 39]. A necessary criterion for separabilityis equivalent to a sufficient criterion for its converse, namely quantum entanglement.An entanglement test on the second moments can be carried out as follows. Given asuitable canonical basis of self-adjoint operators ˆ r = (ˆ x , ˆ p , ˆ x , ˆ p . . . ˆ x N , ˆ p N ) T for N modes,the covariance matrix σ is given by σ = Tr (cid:2) { ˆ r , ˆ r T } ˆ (cid:37) (cid:3) , (39)where ˆ (cid:37) is a Gaussian (or non-Gaussian) state, and where the transpose T is taken in theouter-product sense, including all possible pairs of operators to form a real, symmetric matrixof correlations.For continuous variable systems, partial transposition is equivalent to changing the signof the second canonical variable. This implies that half of the off-diagonal elements in thecovariance matrix, namely all elements that contain an odd power of the second momentumvariable, gain a minus sign [35] (we demonstrate this explicitly in Appendix D). Any physical(i.e., derived from a trace-class, positive density operator) covariance matrix σ must satisfythe uncertainty principle that can be cast as − (Ω σ Tp ) ≥ ˆ , where Ω is the symplectic formdefined in this basis as Ω = n (cid:77) j =1 ω , with ω = (cid:32) − (cid:33) . (40)The sufficient entanglement test can thus be stated as ∃ j : ˜ ν j < , (41)where ˜ ν j are the N symplectic eigenvalues of the partially transposed covariance matrix σ Tp .The symplectic eigenvalues are defined as the square roots of the eigenvalues of − (Ω σ Tp ) ,which are at least two-fold degenerate (so that there are N independent ones in the case onhand).Hence, in order to demonstrate entanglement for the Hydrogenic wavefunction (which areglobally non-Gaussian due to the relative coordinate contribution, even if the centre-of-masscoordinate is taken in a Gaussian wave-packet, which in turn would amount to a Gaussianstate in the quantum optics nomenclature), it suffices to show that one of the symplecticeigenvalues of its covariance matrix satisfies ˜ ν j <
1. We therefore need to evaluate its secondmoments, which is what we do in the following section. ¶ Indeed, a non-Gaussian state that is entangled might not appear entangled through examination of itssecond moments alone. As already mentioned, the test we will adopt is therefore sufficient but not necessary for determining entanglement. ydrogenic entanglement The bipartite state in Eq. (38) has support in three spatial dimensions. This setting yieldssix independent modes (one for each spatial direction for each of the subsystems). It followsthat the covariance matrix is a 12 ×
12 matrix (with two variables per mode). To constructthe matrix, we must compute the expectation values and variances for the position andmomentum variables in each mode.In the original partition in terms of the r and r coordinates for the original twosubsystems, we construct the following (dimensionfull) basis of coordinates: X = ( x , p x , y , p y , z , p z , x , p x , y , p y , z , p z ) T . (42)It would be cumbersome to compute the first and second moments directly in the basis X of Eq. (42), since the variables r and r are mixed between the Hydrogenic eigenstatesand the Gaussian wavepacket. However, we can simply evaluate the statistical momentsin the { r , p , R , P } basis, separately for the Hydrogenic wavefunction and the Gaussianwavepacket, and then apply the symplectic transformation S that relates the local coordinatesto the decoupled ones on the covariance matrix itself, which transforms, by congruence, as σ (cid:48) = SσS T . With the S shown in Eq.(D.3) in Appendix D, the transformed basis vector inEq. (42) becomes X (cid:48) = ( x, p x , y, p y , z, p z , X, P X , Y, P Y , Z, P Z ) T , (43)where x, y, z and p x , p y , p z are the relative position and momentum coordinates given by (forequal masses) x = x − x , and p x = p x − p x , (44)and so on. Similarly, X, Y, Z and P X , P Y , P Z are the centre-of-mass position and momentumcoordinates given by X = x + x , and P X = p x + p x , (45)and so on. It is possible to find a symplectic transformation that takes unequal masses intoaccount, but we do not do so here.We are now ready to compute the first and second moments of the Hydrogenic states. Toensure that they are dimensionless, we rescale the position coordinates by the characteristiclength scale, the reduced Bohr radius a (see Appendix A), and the momentum coordinatesby the characteristic momentum (cid:126) /a .Using identities that involve the Legendre polynomials and known relations such as theKramer–Pasternak relation [40, 41], we obtain the following first and second moments in the { r , R } basis. See Appendix E for the full, rather lengthy calculations. The dimensionless(where we have chosen to normalise all quantities with respect to the reduced Bohr radius ydrogenic entanglement a ) expectation values and variances for the relative position variables of the Hydrogeniceigenstates are (cid:104) x/a (cid:105) = (cid:104) y/a (cid:105) = (cid:104) z/a (cid:105) = 0 , (46)and, since (cid:104) x /a (cid:105) = (cid:104) y /a (cid:105) : (cid:10) x /a (cid:11) = n (5 n − l ( l + 1) + 1)2 l + l + m − l − l + 3) , (47) (cid:10) z /a (cid:11) = n (5 n − l ( l + 1) + 1)2 1 − l − l + 2 m − l − l . The expectation values of the relative momentum variables are (cid:104) p x a / (cid:126) (cid:105) = (cid:104) p y a / (cid:126) (cid:105) = (cid:104) p z a / (cid:126) (cid:105) = 0 , (48)and the variances are given by, with (cid:104) p x a / (cid:126) (cid:105) = (cid:10) p y a / (cid:126) (cid:11) : (cid:10) p x a / (cid:126) (cid:11) = 1 n l + l + m − l − l + 3) , (cid:10) p z a / (cid:126) (cid:11) = 1 n − l − l + 2 m − l − l . (49)As for the Gaussian wavepacket, the centre-of-mass expectation values and variances read (cid:104) X/a (cid:105) = (cid:104) Y /a (cid:105) = (cid:104) Z/a (cid:105) = 0 , (cid:10) X /a (cid:11) = (cid:10) Y /a (cid:11) = (cid:10) Z /a (cid:11) = b a , (cid:104) P X a / (cid:126) (cid:105) = (cid:104) P Y a / (cid:126) (cid:105) = (cid:104) P Z a / (cid:126) (cid:105) = 0 , (cid:10) P X a / (cid:126) (cid:11) = (cid:10) P Y a / (cid:126) (cid:11) = (cid:10) P Z a / (cid:126) (cid:11) = a b . (50)We now promptly transform back to the original basis X = S T X (cid:48) , and σ = S − σ (cid:48) S − . Wethen apply the PPT criterion and compute the symplectic eigenvalues ˜ ν j of σ Tp . They aregiven by, in terms of the previously computed variances,˜ ν = (cid:113) (cid:104) x /a (cid:105) (cid:104) P X a / (cid:126) (cid:105) , ˜ ν = 4 (cid:113) (cid:104) X /a (cid:105) (cid:104) p x a / (cid:126) (cid:105) , ˜ ν = (cid:113) (cid:104) y /a (cid:105) (cid:104) P Y a / (cid:126) (cid:105) , ˜ ν = 4 (cid:113) (cid:104) Y /a (cid:105) (cid:10) p y a / (cid:126) (cid:11) , ˜ ν = (cid:113) (cid:104) z /a (cid:105) (cid:104) P Z a / (cid:126) (cid:105) , ˜ ν = 4 (cid:113) (cid:104) Z /a (cid:105) (cid:104) p z a / (cid:126) (cid:105) . (51)Due to the rotational symmetry of the system in the x - y plane, the symplectic eigenvalues ydrogenic entanglement ν = ˜ ν and ˜ ν = ˜ ν . The four unique eigenvalues are given by˜ ν (1) nlm = a n b (cid:114) ( l + l + m −
1) (5 n − l ( l + 1) + 1)4 l + 4 l − , ˜ ν (2) nlm = 2 √ ba n (cid:114) l + l + m − l + 4 l − , ˜ ν (5) nlm = a n b (cid:114) (2 l + 2 l − m −
1) (5 n − l ( l + 1) + 1)4 l + 4 l − , ˜ ν (6) nlm = 2 √ ba n (cid:114) l + 2 l − m − l + 4 l − . (52)To detect entanglement, we merely require that one of the eigenvalues becomes smaller thanunity: ν ( j ) nlm < j = 1 , , ,
6. The first symplectic eigenvalue ˜ ν (1) nlm indicates that thelikelihood of the PPT criterion being violated increases with b , in fact, as b → ∞ , ˜ ν j → ν (2) nlm has the opposite dependence of a and b , but scales as n − ,which means that the larger n becomes, the more entangled the state appears. We alsonote that all eigenvalues depend on the ratio a /b (which is independent of the choice ofnormalisation with either a or b ).For the Hydrogenic ground state, with n = 1 and l = m = 0, we find that the eigenvaluesbecome even more degenerate with ˜ ν (5)100 = ˜ ν (1)100 and ˜ ν (6)100 = ˜ ν (2)100 since rotational symmetry isrestored. The remaining unique eigenvalues are˜ ν (1)100 = (cid:114) a b , ˜ ν (2)100 = (cid:114) ba . (53)These values depend purely on the ratio of the Hydrogenic characteristic length-scale and theGaussian wavepacket spread a /b . Entanglement can however be detected for most valuesthrough examination of these two eigenvalues. When a ∼ b , the first eigenvalue alwaysdetects entanglement. For the special case when a /b ∼ √
2, entanglement can briefly notbe detected at all, until a /b > (cid:112) /
3. This implies that it would be difficult to verify theentanglement of systems localised to scales that are of the same order as a .We have plotted the function min(˜ ν (1)100 , ˜ ν (2)100 ) in Figure 2 to demonstrate the regimeswhere entanglement can be detected. The blue area indicates where either of the eigenvaluesdip below 1, and where entanglement can be inferred. The light yellow area indicates whereboth eigenvalues are larger than 1, which means that entanglement cannot be inferred.We conclude that in order to apply this entanglement test, one must measure theexpectation values and variances of the relative and centre-of-mass coordinates in the x and z spatial directions to compute the two symplectic eigenvalues in Eq. (53). It is worthmentioning that this test, in contrast to the the previous one based on the relative momentumvariance, would be sufficient also for noisy, mixed states (albeit for such states it would likelybe more difficult to violate). ydrogenic entanglement Plot of min(˜ ν (1)100 , ˜ ν (2)100 ) shown in Eq. (53) as a function of the reduced Bohr radius a and the Gaussian wavepacket spread b . The plot-range has been chosen merely to demonstratethe validity of the entanglement test. The blue area indicates where either of the two eigenvaluesdips below 1, and where entanglement can be detected. The light yellow area indicates where botheigenvalues are larger than 1, and where entanglement cannot be inferred.
5. Discussion
We derived two expressions that can serve as entanglement tests for a free Hydrogenic system(see Eq. (35)) and a localised Hydrogenic system (see Eq. (52)), respectively. Here, we brieflydiscuss our results and their potential application to detecting entanglement between twosubsystems that interact via a central potential.The first additional system that might exhibit Hydrogenic structure is positronium,which is the bound state of one electron and one positron. With a half-life of 0.1244 ns [42],it is a highly unstable system. We note that detecting entanglement in positronium mightbe extremely challenging, however it may also constitute one of the most straight-forwardroutes towards detecting entanglement between matter and anti-matter systems.We speculate that another candidate family of systems that might exhibit Hydrogenicstructure could be mesoscopic systems that interact through a central potential, such asa Coulombic potential. Recent advances in the ability to control the number of charges onlevitated silica spheres [43] or implant changes through the addition of nitrogen-vacancy (NV)centres [44] provide an excellent means for highly sophisticated operations in the laboratory.It is currently an open question as to whether mesoscopic systems would at all readily exhibitHydrogenic structure, but should this be the case, the methods developed here may be used ydrogenic entanglement
6. Conclusions
In this work, we proposed two entanglement tests for free and localised bipartite systems withHydrogenic structure. We computed the spread of the Schmidt spectrum for an arbitraryenergy eigenstate of a delocalised (free) Hydrogenic system, and showed that it scales with theinverse of the principal quantum number n − . For localised systems, we demonstrated thatmeasuring the variances of the relative and centre-of-mass coordinates suffices to demonstratethat the system is entangled. We believe that these results could potentially aid the study ofentanglement in systems that display Hydrogenic characteristics, such as interacting matter–anti-matter systems or freely-falling mesoscopic systems.A number of open questions remain, such as the extension of these results to dynamicalprocesses and the formation of Hydrogenic states starting from highly localised systems.One may also ask how the entanglement of a Hydrogenic system changes as a functionof an external potential, like, for example, a magnetic field. The inclusion of an externalpotentials has previously been shown to drastically change the entanglement contents ofrelated systems [47]. We leave these questions to future work. Acknowledgments
We warmly thank David Cassidy, Michael R. Vanner, Alfred Harwood, Ryan Marshman,Carlo Sparaciari, Alexander Ferrier, Gavin Morley, Peter F. Barker, A. Douglas K. Plato,Dennis R¨atzel, David Edward Bruschi, Alessio Belenchia, Tania Monteiro, Ivette Fuentes,and Markus Aspelmeyer for helpful discussions and comments. We would also like to conveya special thanks to the referees, whose careful reading of this work was most helpful. SQ isfunded by an EPSRC Doctoral Prize Fellowship. ydrogenic entanglement Bibliography [1] Schr¨odinger E 1926
Annalen der physik https://doi.org/10.1002/andp.19263851302 [2] Schr¨odinger E 2003
Collected papers on wave mechanics vol 302 (American Mathematical Society)[3] Born M and Oppenheimer R 1927
Annalen der physik https://doi.org/10.1002/andp.19273892002 [4] Izmaylov A F and Franco I 2017
Journal of Chemical Theory and Computation https://doi.org/10.1021/acs.jctc.6b00959 [5] Bouvrie P A, Majtey A P, Tichy M C, Dehesa J S and Plastino A R 2014 The European Physical JournalD
346 URL https://doi.org/10.1140/epjd/e2014-50349-2 [6] Jozsa R and Linden N 2003
Proceedings of the Royal Society of London. Series A: Mathematical, Physicaland Engineering Sciences https://doi.org/10.1098/rspa.2002.1097 [7] Ekert A K 1991
Physical Review Letters
661 URL https://doi.org/10.1103/PhysRevLett.67.661 [8] Degen C L, Reinhard F and Cappellaro P 2017
Reviews of Modern Physics https://doi.org/10.1103/RevModPhys.89.03500 [9] Bassi A, Lochan K, Satin S, Singh T P and Ulbricht H 2013 Reviews of Modern Physics
471 URL https://doi.org/10.1103/RevModPhys.85.471 [10] Kiefer C 2007 Why quantum gravity?
Approaches to fundamental physics (Springer) pp 123–130 URL https://doi.org/10.1007/978-3-540-71117-9_7 [11] Belli S, Bonsignori R, D’Auria G, Fant L, Martini M, Peirone S, Donadi S and Bassi A 2016
PhysicalReview A https://doi.org/10.1103/PhysRevA.94.012108 [12] Bose S, Mazumdar A, Morley G W, Ulbricht H, Toroˇs M, Paternostro M, Geraci A A, Barker P F, KimM S and Milburn G 2017 Physical Review Letters https://journals.aps.org/prl/abstract/10.1103/PhysRevLett.119.240401 [13] Marletto C and Vedral V 2017
Physical Review Letters https://journals.aps.org/prl/abstract/10.1103/PhysRevLett.119.240402 [14] Belenchia A, Wald R M, Giacomini F, Castro-Ruiz E, Brukner ˇC and Aspelmeyer M 2018
PhysicalReview D https://doi.org/10.1103/PhysRevD.98.126009 [15] Marshman R J, Mazumdar A and Bose S 2019 arXiv preprint arXiv:1907.01568 URL https://arxiv.org/abs/1907.01568 [16] Ashoori R C 1996
Nature
413 URL https://doi.org/10.1038/379413a0 [17] Fujisawa T, Austing D G, Tokura Y, Hirayama Y and Tarucha S 2002
Nature
278 URL https://doi.org/10.1038/nature00976 [18] Wannier G H 1937
Physical Review
191 URL https://doi.org/10.1103/PhysRev.52.191 [19] Mott N F and Littleton M J 1938
Transactions of the Faraday Society https://pubs.rsc.org/en/content/articlepdf/1938/tf/tf9383400485 [20] Millen J, Fonseca P Z G, Mavrogordatos T, Monteiro T S and Barker P F 2015 Physical Review Letters https://doi.org/10.1103/PhysRevLett.114.123602 [21] Qvarfort S, Bose S and Serafini A 2018 arXiv preprint arXiv:1812.09776
URL https://arxiv.org/abs/1812.09776 [22] Dehesa J S, Koga T, Yanez R J, Plastino A R and Esquivel R O 2011
Journal of Physics B: Atomic,Molecular and Optical Physics https://doi.org/10.1088/0953-4075/45/1/015504 [23] Benenti G, Siccardi S and Strini G 2013 The European Physical Journal D
83 URL https://doi.org/10.1140/epjd/e2013-40080-y [24] Lin Y C, Lin C Y and Ho Y K 2013
Physical Review A (2) 022316 URL https://link.aps.org/doi/10.1103/PhysRevA.87.022316 [25] Ding L and Schilling C 2020 arXiv preprint arXiv:2001.04858 URL https://arxiv.org/abs/2001.04858 [26] Tommasini P, Timmermans E and de Toledo Piza A F R 1998
American Journal of Physics ydrogenic entanglement URL https://doi.org/10.1119/1.18977 [27] Simon R 2000
Physical Review Letters https://doi.org/10.1103/PhysRevLett.84.2726 [28] Shchukin E and Vogel W 2005 Physical Review Letters https://doi.org/10.1103/PhysRevLett.95.230502 [29] Bransden B H, Joachain C J and Plivier T J 2003 Physics of atoms and molecules (Pearson EducationIndia)[30] Erd´elyi A, Magnus W, Oberhettinger F, Tricomi F G and Bateman H 1955
Higher transcendentalfunctions vol 3 (New York)[31] Parker S, Bose S and Plenio M B 2000
Physical Review A https://doi.org/10.1103/PhysRevA.61.032305 [32] Arfken G B, Weber H J and Harris F E 2012 Mathematical methods for physicists (Elsevier) ISBN978-0-12-384654-9 URL https://doi.org/10.1016/C2009-0-30629-7 [33] Cover T M and Thomas J A 2012
Elements of information theory (John Wiley & Sons)[34] Eisert J, Simon C and Plenio M B 2002
Journal of Physics A: Mathematical and General https://doi.org/10.1088/0305-4470/35/17/307 [35] Serafini A 2017 Quantum Continuous Variables: A Primer of Theoretical Methods (CRC Press)[36] Duan L M, Giedke G, Cirac J I and Zoller P 2000
Physical Review Letters https://doi.org/10.1103/PhysRevLett.84.2722 [37] Lami L, Serafini A and Adesso G 2018 New Journal of Physics https://doi.org/10.1088/1367-2630/aaa654 [38] Peres A 1996 Physical Review Letters https://doi.org/10.1103/PhysRevLett.77.1413 [39] Horodecki P 1997 Physics Letters A https://doi.org/10.1016/S0375-9601(97)00416-7 [40] Pasternack S 1937
Proceedings of the National Academy of Sciences of the United States of America
91 URL https://doi.org/10.1073/pnas.23.2.91 [41] Kramers L H 1939
Zeitschrift f¨ur Physik
International Journal of Modern Physics A https://doi.org/10.1142/S0217751X04020142 [43] Moore D C, Rider A D and Gratta G 2014 Physical Review Letters https://doi.org/10.1103/PhysRevLett.113.251801 [44] Neukirch L P, Von Haartman E, Rosenholm J M and Vamivakas A N 2015
Nature Photonics
653 URL [45] Rajantie A 2012
Contemporary Physics https://doi.org/10.1080/00107514.2012.685693 [46] Reis J T and Sauter W K 2017 Physical Review D https://doi.org/10.1103/PhysRevD.96.075031 [47] Nazmitdinov R G, Simonovi´c N S, Plastino A and Chizhov A V 2012 Journal of Physics B:Atomic, Molecular and Optical Physics https://doi.org/10.1088/0953-4075/45/20/205503 [48] Griffiths D J and Schroeter D F 2018 Introduction to quantum mechanics (Cambridge University Press)[49] Avery J S 2012
Hyperspherical harmonics: applications in quantum theory vol 5 (Springer Science &Business Media)[50] Podolsky B and Pauling L 1929
Physical Review
109 URL https://doi.org/10.1103/PhysRev.34.109 [51] Inghoff T, Fritzsche S and Fricke B 2001
Computer physics communications https://doi.org/10.1016/S0010-4655(01)00218-1 [52] San Kim D, Kim T and Rim S H 2012
Advances in Difference Equations
219 URL https://doi.org/10.1186/1687-1847-2012-219 [53] Gram-Hansen B and Ossipov A 2015 An investigation into the creation of entanglement mediated byinteraction URL https://bayesianbrad.github.io/papers/DIS.pdf ydrogenic entanglement [54] Samaddar S N 1974 Mathematics of Computation https://doi.org/10.1090/S0025-5718-1974-0328157-4 ydrogenic entanglement Appendix A. Hydrogenic solutions in position and momentum space
In this Appendix, we introduce the Hydrogenic solutions to the time-independent Schr¨odingerequation in three dimensions for two particles that interact through a central potential. Wedescribe here the general procedure for separating the Schr¨odinger equation in terms of thecentre-of-mass R and relative coordinate r , which will be used to describe the interaction.We begin by presenting the procedure in position space, and then move on to review thecorresponding solutions in momentum space. Appendix A.1. Hydrogenic wavefunctions in position space
We follow the standard derivation presented in Appendix 8 of Ref. [29]. The Schr¨odingerequation for N interacting particles in first quantisation with coordinates r i and masses m i ,is given by − (cid:126) N (cid:88) j =1 m j ∇ j Ψ( r , r , . . . , r N , t ) + V ( r , r , . . . , r N , t )Ψ( r , r , . . . r N , t )= i (cid:126) ∂∂t Ψ( r , r . . . r N , t ) , (A.1)where m j is the mass of particle j , ∇ j is the Laplacian operator, and where V ( r , r , . . . , r N , t )is the sum of central potentials affecting the particles. For two interacting particles, theequation reduces to (cid:20) − (cid:126) m ∇ − (cid:126) m ∇ − V ( r , r ) (cid:21) Ψ( r , r ) = i (cid:126) ∂∂t Ψ( r , r ) , (A.2)where the explicit form of the potential is given by V ( r , r ) = α | r − r | = α | r | = V ( r ) , (A.3)where α > r = r − r .For the Coulomb potential, we find that α = q q / (4 π(cid:15) ), and for a Newtonian gravitationalinteraction, α = Gm m .The two-particle Schr¨odinger equation in Eq. (A.2) is solved by moving to a centre-of-mass coordinate R and a relative coordinate r , defined as R = m r + m r m + m , (A.4) r = r − r , (A.5)(with conjugate momenta P = ( p + p ) and p = ( m p − m p ) / ( m + m ), defining acanonical quadruple). We then introduce the reduced mass µ = m m m + m and the total mass ydrogenic entanglement M = m + m in order to write the original coordinates as r = R + µm r , (A.6) r = R − µm r . (A.7)It can then be shown that the gradient operators can be written as ∇ = ∇ r + µm ∇ R , (A.8) ∇ = − ∇ r + µm ∇ R , (A.9)which allows us to rewrite the two-particle Schr¨odinger equation as (cid:20) − (cid:126) m + m ) ∇ R − (cid:126) µ ∇ r − V ( r ) (cid:21) Ψ( r , R ) = E Ψ( r ; R ) . (A.10)Since the potential only depends on r and not on R , we can separate the variables into theSchr¨odinger equation above, as per Ψ( r , R ) = ψ ( r ) ϕ ( R ), obtaining − (cid:126) M ∇ R ϕ ( R ) = E CM ϕ ( R ) , (A.11) (cid:20) − (cid:126) µ ∇ r + V ( r ) (cid:21) ψ ( r ) = E rel ψ ( r ) , (A.12)where E CM is the centre-of-mass energy and E rel is the relative energy. The total energy ofthe system is E = E CM + E rel . (A.13)Eq. (A.11) describes the free centre-of-mass wavefunction (the overall wavepacket) of thetwo-particle system and admits wave-like solutions. Eq. (A.12), on the other hand, admitsthe Hydrogenic solutions ψ nlm ( r ) with E rel being the usual Hydrogen energy spectrum [48].The solutions to the relative equation are given by ψ nlm ( r ) = (cid:115)(cid:18) na (cid:19) ( n − l − n [( n + l )!] e − r/na (cid:18) rna (cid:19) l (cid:2) L l +1 n − l − (2 r/na ) (cid:3) Y ml ( θ, φ ) , (A.14)where n, l, m are the principal, angular and magnetic quantum numbers, respectively. r, θ and φ are elements of the relative coordinate vector: r = ( r, θ, φ ) T , and where L l +1 n − l − (2 r/na )is the associated Laguerre polynomial, defined by L pq − p ( x ) ≡ ( − p (cid:18) ddx (cid:19) p L q ( x ) , (A.15)and where L q ( x ) denotes the so-called Laguerre polynomial, given by L q ( x ) ≡ e x (cid:18) ddx (cid:19) q (cid:0) e − x x q (cid:1) . (A.16) ydrogenic entanglement a is a length scale set by the interaction coefficient α , correspondingto the reduced Bohr radius: a = (cid:126) µ α . (A.17)Returning to the functions shown in Eq. (A.14), the spherical harmonics Y ml ( θ, φ ) are givenby Y ml ( θ, φ ) = ( − m (cid:115) (2 l + 1)4 π ( l − m )!( l + m )! P lm (cos θ ) e imφ , (A.18)where P lm (cos θ ) are the associated Legendre polynomials without the Condon–Shortleyphase [32]. The spherical harmonics satisfy the following orthogonality relation (cid:90) π d θ sin θ (cid:90) π d φ Y m (cid:48) ∗ l (cid:48) ( θ, φ ) Y ml ( θ, φ ) = δ mm (cid:48) δ ll (cid:48) . (A.19)We make frequent use of this relation in the following Appendices. This concludes oursummary of the Hydrogenic solutions in position space. Appendix A.2. Hydrogenic wavefunctions in momentum space
The Hydrogenic wavefunctions in momentum space ˜ ψ ( k ), where the relation between thewavevector k and the momentum is p = (cid:126) k , are defined as the Fourier transform of theposition space wavefunctions: ˜ ψ ( k ) = (cid:90) d r e − i k · r ψ ( r ) . (A.20)Note that we denote this k with a bar in the main text to differentiate it from another relativequantity, but we omit this here for notational clarity. The explicit form of the Hydrogenicwavefunctions in momentum space can be computed to be (see Appendix 5 in [29]):˜ ψ ( k ) = F nl ( k ) Y ml ( θ, φ ) , (A.21)where F nl ( k ) = (cid:20) π ( n − l − n + l )! (cid:21) / n l +2 l ! (cid:113) a n l a l k l ( n a k + 1) l +2 C l +1 n − l − (cid:18) n a k − n a k + 1 (cid:19) , (A.22)and where C αN ( x ) are the Gegenbauer polynomials [30], which are defined by the relation(1 − xs + s ) − α = ∞ (cid:88) N =0 C αN ( x ) s N , (A.23)with | s | <
1. The Gegenbauer polynomials are special cases of the more general Jacobipolynomials and play a role that is analogous to the role played by Legendre polynomials inthe theory of the three-dimensional spherical harmonics [49]. ydrogenic entanglement (cid:90) − d x (1 − x ) α − / C ( α ) n ( x ) C ( α ) n (cid:48) = 2 − α π Γ( n + 2 α ) n !( n + α )[Γ( α )] δ nn (cid:48) , (A.24)where α is a dummy index, not to be mixed up with the coupling constant defined earlier.This is yet another relation that is used many times in what follows.The normalisation of the momentum wavefunctions reads (cid:90) d k | ˜ ψ ( k ) | = 1 , (A.25)which, since the spherical harmonics are orthogonal (see Eq. (A.19)), implies that (cid:90) ∞ d k k [ F nl ( k )] = 1 . (A.26)This concludes the summary of the Hydrogenic wavefunctions in momentum space. Wereturn to these expressions in Appendix B, where we compute the variance of the relativemomentum variable. Appendix B. Normalisation and expectation values for the momentumHydrogenic wavefunction
In this Appendix, we prove Eq. (34) in the main text, which provides a relation for thevariance of the relative momentum coordinate. This, in turn, corresponds to the continuousanalogue of the Schmidt coefficients of a free Hydrogenic system, and is therefore a key resultin this work. We begin by proving the orthogonality relation and the normalisation conditionof the Gegenbauer polynomials, as this will become useful in proving Eq. (34).We used the notation k in the main text to refer to the local wavevector rather thanthat related to the relative momentum (which is defined as k = ( m k − m k ) / ( m + m )).We here remove the bar for notational convenience. Appendix B.1. Normalisation of the Hydrogenic momentum wavefunction
As discussed in Appendix A, the momentum representation of the Hydrogenic wavefunctionsis given in Eq. (A.21). We begin by proving the normalisation relation in Eq. (A.26), whichwe reprint here for brevity: (cid:90) ∞ d k k | F nl ( k ) | = 1 . (B.1)We reprint the expression for F nl ( k ) here for convenience: F nl ( k ) = (cid:20) π ( n − l − n + l )! (cid:21) / n l +2 l ! (cid:113) a n l a l k l ( n a k + 1) l +2 C l +1 n − l − (cid:18) n a k − n a k + 1 (cid:19) , (B.2) ydrogenic entanglement ψ ( k ) has units of length / due to the appearance of (cid:112) a in Eq. (B.2). This is to ensure that the normalisation is dimensionless.Our goal is to prove the relation in Eq. (B.1). We start by writing out the expression infull: (cid:90) ∞ d k k | F nl ( k ) | = 2 π ( n − l − n + l )! n l +4 ( l !) a × (cid:90) ∞ d k k n l a l k l ( n a k + 1) l +4 (cid:20) C l +1 n − l − (cid:18) n a k − n a k + 1 (cid:19)(cid:21) . (B.3)We now multiply and divide by ( a k ) and 4 l +2 in order to collect the following terms: (cid:90) ∞ d k k | F nl ( k ) | = 2 π ( n − l − n + l )! 2 l ( l !) a × (cid:90) ∞ d k k (cid:18) n a k ( n a k + 1) (cid:19) l +2 (cid:20) C l +1 n − l − (cid:18) n a k − n a k + 1 (cid:19)(cid:21) . (B.4)We now define the new variable x = n a k − n a k + 1 . (B.5)The variable substitution implies that k = 1 a n (cid:114) x − x , and 4 n a k ( n a k + 1) = 1 − x , (B.6)and the derivative becomes dxdk = 4 n a k ( n a k + 1) = 1 − x k . (B.7)With this substitution, the limits become {− , } , which is what we require for thenormalisation condition of the Gegenbauer polynomials. We insert this into Eq. (B.4) tofind (cid:90) ∞ d k k | F n,l ( k ) | = 2 π ( n − l − n + l )! 2 l ( l !) a (cid:90) − d x k − x k (cid:0) − x (cid:1) l +2 (cid:2) C l +1 n − l − ( x ) (cid:3) = 2 π ( n − l − n + l )! 2 l ( l !) n (cid:90) − d x (cid:114) − x x (cid:0) − x (cid:1) l +1 (cid:2) C l +1 n − l − ( x ) (cid:3) = 2 π ( n − l − n + l )! 2 l ( l !) n (cid:90) − d x (1 − x )(1 − x ) l +1 / (cid:2) C l +1 n − l − ( x ) (cid:3) . (B.8)The last integral in Eq. (B.8) can now be written as two separate terms by multiplying outthe first bracket: (cid:90) ∞ d k k | F n,l ( k ) | = 2 π ( n − l − n + l )! 2 l ( l !) n (cid:90) − d x (1 − x ) l +1 / (cid:2) C l +1 n − l − ( x ) (cid:3) − π ( n − l − n + l )! 2 l ( l !) n (cid:90) − d x x (1 − x ) l +1 / (cid:2) C l +1 n − l − ( x ) (cid:3) . (B.9) ydrogenic entanglement α → l + 1 and n → n − l −
1, we find (cid:90) − d x (1 − x ) l +1 / C ( l +1) n − l − ( x ) C ( l +1) n − l − ( x ) = 2 − l − π Γ( n + l + 1) n ( n − l − l + 1)] . (B.10)When we include the prefactor, the full first term of Eq. (B.9) evaluates to2 π ( n − l − n + l )! 2 l ( l !) n (cid:90) − d x (1 − x ) l +1 / (cid:2) C l +1 n − l − ( x ) (cid:3) = 2 π ( n − l − n + l )! 2 l ( l !) n − l − π Γ( n + l + 1) n ( n − l − l + 1)] . (B.11)Since the Gamma function is given by Γ( n ) = ( n − n , we are ableto write 2 π ( n − l − n + l )! 2 l ( l !) n − l − π ( n + l )! n ( n − l − l !) = 1 . (B.12)This is the normalisation we require. It remains to show that the term on the left-hand sideof Eq. (B.11) vanishes. As noted in Ref. [50], to complete the proof, we can use the followingrecursion formula for the Gegenbauer polynomials: x C αn ( x ) = 12( n + α ) (cid:2) ( n + 1) C αn +1 ( x ) + ( n + 2 α − C αn − ( x ) (cid:3) . (B.13)With this relation, the second integral in Eq. (B.11) (without its prefactor and with theappropriate index substitutions) can be written as (cid:90) − d x x (1 − x ) l +1 / (cid:2) C l +1 n − l − ( x ) (cid:3) = 12 n (cid:90) − d x (1 − x ) l +1 / (cid:2) ( n − l ) C l +1 n − l ( x ) + ( n + l ) C l +1 n − l − ( x ) (cid:3) C l +1 n − l − ( x ) . (B.14)Since the indices of the Gegenbauer polynomials inside and outside the bracket differ, theintegral vanishes due to the orthogonality of the Gegenbauer polynomials in Eq. (A.24). Thisconcludes our proof of Eq. (B.1). Appendix B.2. Proof of Hydrogenic wavefunction momentum variance
In Section 3.3, we showed that the variance of the relative momentum variable k (referred toas ¯ k in the text) corresponds to the continuous analogue of the Schmidt coefficients of thesystem, which is given by the concise expression in Eq. (35). Here, we prove the relation inEq. (35).We wish to show that (cid:90) ∞ d k k [ F nl ( k )] = 1 a n , (B.15) ydrogenic entanglement F nl ( k ) is given in Eq. (B.2). This allows us to write the integral as (cid:90) ∞ d k k [ F nl ( k )] = 2 π ( n − l − n + l )! n l +4 ( l !) a × (cid:90) ∞ d k k n l a l k l ( n a k + 1) l +4 (cid:20) C l +1 n − l − (cid:18) n a k − n a k + 1 (cid:19)(cid:21) . (B.16)By now multiplying and dividing by a and 4 l +2 , we can write Eq. (B.16) as (cid:90) ∞ d k k F nl ( k ) F ∗ nl ( k ) = 2 π ( n − l − n + l )! 2 l ( l !) a (B.17) × (cid:90) ∞ d k (cid:18) n a k ( n a k + 1) (cid:19) l +2 (cid:20) C l +1 n − l − (cid:18) n a k − n a k + 1 (cid:19)(cid:21) . To evaluate the integral, we perform the same variable substitution as in Eq. (B.5). Giventhe relationships in Eqs. (B.6) and (B.7), and the new limits x ∈ {− , } , we find (cid:90) ∞ d k k F nl ( k ) F ∗ nl ( k ) = 2 π ( n − l − n + l )! 2 l ( l !) a (cid:90) − d x k − x (cid:0) − x (cid:1) l +2 (cid:2) C l +1 n − l − ( x ) (cid:3) = 2 π ( n − l − n + l )! 2 l ( l !) a (cid:90) − d x n (cid:114) x − x (1 − x ) l +1 (cid:2) C l +1 n − l − ( x ) (cid:3) = 2 π ( n − l − n + l )! 2 l ( l !) a n (cid:90) − d x (1 + x ) (cid:0) − x (cid:1) l +1 / (cid:2) C l +1 n − l − ( x ) (cid:3) . (B.18)The last integral in Eq. (B.18) can again be written in terms of two separate terms: (cid:90) ∞ d k k F nl ( k ) F ∗ nl ( k ) = 2 π ( n − l − n + l )! 2 l ( l !) a n (cid:90) − d x (cid:0) − x (cid:1) l +1 / (cid:2) C l +1 n − l − ( x ) (cid:3) (B.19)+ 2 π ( n − l − n + l )! 2 l ( l !) a n (cid:90) − d x x (cid:0) − x (cid:1) l +1 / (cid:2) C l +1 n − l − ( x ) (cid:3) . As in the previous section, the second term vanishes due to the mismatched indices. Then,by using the orthogonality condition for the Gegenbauer polynomials in Eq. (B.10), we findthat the first term evaluates to2 π ( n − l − n + l )! 2 l ( l !) a n (cid:90) − d x (cid:0) − x (cid:1) l +1 / (cid:2) C l +1 n − l − ( x ) (cid:3) = 2 π ( n − l − n + l )! 2 l ( l !) a n − l − π Γ( n + l + 1) n ( n − l − l + 1)] . (B.20)Noting that Γ( n ) = ( n − (cid:90) ∞ d k k F nl ( k ) F ∗ nl ( k ) = 2 π ( n − l − n + l )! 2 l ( l !) a n − l − π ( n + l )! n ( n − l − l !) = 1 n a , (B.21)which is the expression we wanted. This concludes the proof of Eq. (35). ydrogenic entanglement Appendix C. Entropy measures applied to free Hydrogenic systems
The goal of this Appendix is to compute the linear entropy of a generic Hydrogenic energyeigenstate. We provide this result for completeness and we emphasise that it is not possible toinvoke an operational interpretation of the linear entropy as a quantifier of the entanglement.The continuous variable analogue of the linear entropy S Lin is S Lin = 1 − (cid:90) d r (cid:90) d r (cid:48) (cid:37) ( r , r (cid:48) ) (cid:37) ( r (cid:48) , r ) , (C.1)where (cid:37) ( r (cid:48) , r ) is the continuous analogue of the subsystem density matrix, given by tracingout one of the variables: (cid:37) ( r (cid:48) , r ) = (cid:90) d r Ψ ∗ ( r (cid:48) , r )Ψ( r , r ) . (C.2)The form of the linear entropy in Eq. (C.1) is the continuous analogue of the discrete matrixmultiplication. Inserting the Fourier transform Eq. (24) of the state (cid:37) ( r (cid:48) , r ) into Eq. (C.1)gives S Lin ( (cid:37) ) = 1 − (cid:90) d r (cid:90) d r (cid:48) (cid:37) ( r , r (cid:48) ) (cid:37) ( r (cid:48) , r )= 1 − (cid:90) d r (cid:90) d r (cid:48) π ) (cid:90) d k ˜ (cid:37) ( k ) e i k · ( r − r (cid:48) ) π ) (cid:90) d k (cid:48) ˜ (cid:37) ( k (cid:48) ) e i k (cid:48) ( r (cid:48) − r ) = 1 − π ) (cid:90) d k ˜ (cid:37) ( k ) (cid:90) d k (cid:48) ˜ (cid:37) ( k (cid:48) ) (cid:90) d r e i r ( k − k (cid:48) ) (cid:90) d r (cid:48) e i r (cid:48) · ( k (cid:48) − k ) = 1 − (cid:90) d k ˜ (cid:37) ( k ) ˜ (cid:37) ( k ) δ (0)= 1 − V (cid:90) d k [ ˜ (cid:37) ( k )] , (C.3)where we have defined the Dirac delta function as δ ( k (cid:48) − k ) = 1(2 π ) (cid:90) d x e i x · ( k (cid:48) − k ) , and δ (0) = V. (C.4)With the definition of the momentum space Hydrogenic wavefunctions in Eq. (A.21), we find S Lin = 1 − V (cid:90) d k [ F nl ( k )] | Y ml ( θ, φ ) | . (C.5)To the authors’ knowledge, this integral has no known standard solutions. Instead, we aimto derive a closed-form expression that can be evaluated for specific choices of n , l and m .To simplify the evaluation of the integral, we write it in terms of a radial integral I k and anangular integral I θ , such that S Lin = 1 − V I rad I ang . (C.6)We now proceed to evaluate each integral separately. ydrogenic entanglement Appendix C.1. The angular integral
Let us start by evaluating the angular integral, which involves four spherical harmonics: I ang = (cid:90) π d θ sin θ (cid:90) π d φ | Y ml ( θ, φ ) | . (C.7)To evaluate this integral, we first seek to reduce the order of the spherical harmonics fromfour to two. We use the following Clebsch–Gordon expansion [51] to write Y l ,m ( θ, φ ) Y l ,m ( θ, φ ) (C.8)= (cid:114) (2 l + 1)(2 l + 1)4 π (cid:88) l ,m ( − m (cid:112) l + 1 (cid:32) l l l m m − m (cid:33) (cid:32) l l l (cid:33) Y l ,m ( θ, φ ) , where the 3 × j symbols. Theyare defined in terms of the Clebsch–Gordon coefficients and read: (cid:32) j j j m m m (cid:33) ≡ ( − j − j − m √ j + 1 (cid:104) j m j m | j ( − m ) (cid:105) , (C.9)where (cid:104) j m j m | j ( − m ) (cid:105) denotes the inner product between two eigenstates. The Wigner3- j are zero unless the following conditions are fulfilled: m i ∈ {− j i , − j i + 1 , − j i + 2 , . . . , j i } , for i = 1 , , ,m + m + m = 0 , | j − j | ≤ j ≤ j + j ,j + j + j is an integer . (C.10)When l = l = l and m = m = m , as in our case, we find Y lm ( θ, ϕ ) Y lm ( θ, ϕ ) = (2 l + 1)2 √ π l (cid:88) l (cid:48) ,m (cid:48) ( − m (cid:48) √ l (cid:48) + 1 (cid:32) l l l (cid:48) m m − m (cid:48) (cid:33) (cid:32) l l l (cid:48) (cid:33) Y l (cid:48) ,m (cid:48) ( θ, ϕ ) ≡ (cid:88) l (cid:48) ,m (cid:48) f l,m,l (cid:48) ,m (cid:48) Y l (cid:48) m (cid:48) ( θ, ϕ ) , (C.11)where f l,m,l (cid:48) ,m (cid:48) is defined as f l,m,l (cid:48) ,m (cid:48) = ( − m (cid:48) (2 l + 1) √ l (cid:48) + 12 √ π (cid:32) l l l (cid:48) m m − m (cid:48) (cid:33) (cid:32) l l l (cid:48) (cid:33) . (C.12) ydrogenic entanglement I ang = (cid:90) π d θ sin θ (cid:90) π d φ Y ml ( θ, φ ) Y ml ( θ, φ ) Y m ∗ l ( θ, φ ) Y m ∗ l ( θ, φ )= (cid:90) π d θ sin θ (cid:90) π d φ (cid:88) l (cid:48) ,m (cid:48) f l,m,l (cid:48) ,m (cid:48) Y m (cid:48) l (cid:48) ( θ, φ ) (cid:88) l (cid:48)(cid:48) ,m (cid:48)(cid:48) f ∗ l,m,l (cid:48)(cid:48) ,m (cid:48)(cid:48) Y m (cid:48)(cid:48) ∗ l (cid:48)(cid:48) ( θ, φ )= (cid:88) l (cid:48) ,m (cid:48) ,l (cid:48)(cid:48) ,m (cid:48)(cid:48) f l,m,l (cid:48) ,m (cid:48) f ∗ l,m,l (cid:48)(cid:48) ,m (cid:48)(cid:48) δ l (cid:48) ,l (cid:48)(cid:48) δ m (cid:48) ,m (cid:48)(cid:48) = (cid:88) l (cid:48) ,m (cid:48) f l,m,l (cid:48) ,m (cid:48) f ∗ l,m,l (cid:48) ,m (cid:48) . (C.13)Using the definition of f l,m,l (cid:48) ,m (cid:48) in Eq. (C.12), we expand the expression to find I ang = (cid:88) l (cid:48) ,m (cid:48) f l,m,l (cid:48) ,m (cid:48) f ∗ l,m,l (cid:48) ,m (cid:48) (C.14)= (cid:88) l (cid:48) ,m (cid:48) (2 l + 1) (2 l (cid:48) + 1)4 π (cid:32) l l l (cid:48) m m − m (cid:48) (cid:33) (cid:32) l l l (cid:48) (cid:33) (cid:32) l l l (cid:48) m m − m (cid:48) (cid:33) ∗ (cid:32) l l l (cid:48) (cid:33) ∗ . One of the requirements to ensure that the Wigner 3- j symbols are non-zero states that thesecond row satisfies m + m + m = 0. In our case, this corresponds to 2 m = m (cid:48) . Thisrequirement is only fulfilled for a single value of the sum, namely m (cid:48) = 2 m , which cancelsthe sum over m (cid:48) . Hence we obtain I ang = (cid:88) l (cid:48) (2 l + 1) (2 l (cid:48) + 1)4 π (cid:32) l l l (cid:48) m m − m (cid:33) (cid:32) l l l (cid:48) (cid:33) . (C.15)We proceed to compute the radial integral. Appendix C.2. The radial integral
The radial integral is given by I rad = (cid:90) ∞ d k k | F nl ( k ) | . (C.16)Now, F nl ( k ) is given in Eq. (B.2). Just like we did in Appendix B, we rewrite the integral inEq. (C.16) by dividing and multiplying with a k and 4 l +4 to find: I rad = (cid:20) π ( n − l − n + l )! (cid:21) n l +8 ( l !) a (cid:90) ∞ d k k n l a l k l ( n a k + 1) l +8 (cid:20) C l +1 n − l − (cid:18) n a k − n a k + 1 (cid:19)(cid:21) = (cid:20) π ( n − l − n + l )! (cid:21) l ( l !) a (cid:90) ∞ d k k (cid:18) n a k ( n a k + 1) (cid:19) l +4 (cid:20) C l +1 n − l − (cid:18) n a k − n a k + 1 (cid:19)(cid:21) . (C.17) ydrogenic entanglement I rad = (cid:20) π ( n − l − n + l )! (cid:21) l ( l !) a (cid:90) ∞ d k k (cid:18) n a k ( n a k + 1) (cid:19) l +4 (cid:20) C l +1 n − l − (cid:18) n a k − n a k + 1 (cid:19)(cid:21) = (cid:20) π ( n − l − n + l )! (cid:21) l ( l !) a (cid:90) − d x k k − x (cid:0) − x (cid:1) l +4 (cid:2) C l +1 n − l − ( x ) (cid:3) = (cid:20) π ( n − l − n + l )! (cid:21) l ( l !) n a (cid:90) − d x (1 − x ) (cid:0) − x (cid:1) l +1 / (cid:2) C l +1 n − l − ( x ) (cid:3) . (C.18)This integral does not satisfy the Gegenbauer normalisation condition in Eq. (A.24). We mustinstead find another way to evaluate the integral. To this end, there exists an alternativeexpansion of the Gegenbauer polynomials in terms of polynomials, given by [52] C ( λ ) n ( x ) = (cid:32) n + 2 λ − n (cid:33) n (cid:88) k =0 (cid:32) nk (cid:33) (2 λ + n ) k ( λ + 1 / k (cid:18) x − (cid:19) k = n (cid:88) k =0 f k (1 − x ) k , (C.19)where we have defined f k = (cid:32) n + 2 λ − n (cid:33) (cid:32) nk (cid:33) (2 λ + n ) k ( λ + 1 / k ( − k k , (C.20)and where a k is the Pochhammer symbol, which is defined such that a k = a ( a + 1)( a +2) . . . ( a + k − I r = (cid:20) π ( n − l − n + l )! (cid:21) l ( l !) n a (cid:90) − d x (1 − x ) (cid:0) − x (cid:1) l +1 / (cid:2) C l +1 n − l − ( x ) (cid:3) (C.21)= (cid:20) π ( n − l − n + l )! (cid:21) l ( l !) n a (cid:88) a,b,c,d f a f b f c f d (cid:90) − d x (1 − x ) a + b + c + d (cid:0) − x (cid:1) l +1 / . Defining γ = a + b + c + d , we find (cid:90) − d x (1 − x ) γ (cid:0) − x (cid:1) l +1 / = ( γ + 5) F (cid:0) , − γ − , − γ − ; , l + ; 1 (cid:1) l + 3+ F (cid:0) , , − l − ; γ + , γ + 4; 1 (cid:1) γ + 6+ 2 l + γ +6 Γ (cid:0) l + (cid:1) Γ (cid:0) l + γ + (cid:1) Γ(4 l + γ + 8) , (C.22) ydrogenic entanglement F is the generalised hypergeometric function. Together with its prefactor, the radialintegral becomes I rad = (cid:20) π ( n − l − n + l )! (cid:21) ( l !) l n a (cid:34)(cid:32) n + ln − l − (cid:33)(cid:35) × (cid:88) a,b,c,d ( − a + b + c + d a + b + c + d (cid:32) n − l − a (cid:33) ( n + l + 1) a ( l + 3 / a (cid:32) n − l − b (cid:33) ( n + l + 1) b ( l + 3 / b × (cid:32) n − l − c (cid:33) ( n + l + 1) c ( l + 3 / c (cid:32) n − l − d (cid:33) ( n + l + 1) d ( l + 3 / d × (cid:20) l + a + b + c + d +6 Γ (cid:0) l + (cid:1) Γ (cid:0) l + a + b + c + d + (cid:1) Γ(4 l + a + b + c + d + 8)+ F (cid:0) , , − l − ; a + b + c + d + , a + b + c + d + 4; 1 (cid:1) a + b + c + d + 6+ ( a + b + c + d + 5) F (cid:0) , − a + b + c + d − , − a + b + c + d − ; , l + ; 1 (cid:1) l + 3 (cid:21) . (C.23)We are now ready to put everything together. Appendix C.3. Final expression for the linear entropy
Combining the angular integral in Eq. (C.15) and the radial integral in Eq. (C.23), we find S Lin = 1 − V I rad I ang = 1 − n a V (cid:88) l (cid:48) (2 l + 1) (2 l (cid:48) + 1)4 π (cid:32) l l l (cid:48) m m − m (cid:33) (cid:32) l l l (cid:48) (cid:33) × (cid:20) π ( n − l − n + l )! (cid:21) ( l !) l (cid:34)(cid:32) n + ln − l − (cid:33)(cid:35) × (cid:88) a,b,c,d ( − a + b + c + d a + b + c + d (cid:32) n − l − a (cid:33) ( n + l + 1) a ( l + 3 / a (cid:32) n − l − b (cid:33) ( n + l + 1) b ( l + 3 / b × (cid:32) n − l − c (cid:33) ( n + l + 1) c ( l + 3 / c (cid:32) n − l − d (cid:33) ( n + l + 1) d ( l + 3 / d × (cid:20) l + a + b + c + d +6 Γ (cid:0) l + (cid:1) Γ (cid:0) l + a + b + c + d + (cid:1) Γ(4 l + a + b + c + d + 8)+ F (cid:0) , , − l − ; a + b + c + d + , a + b + c + d + 4; 1 (cid:1) a + b + c + d + 6+ ( a + b + c + d + 5) F (cid:0) , − a + b + c + d − , − a + b + c + d − ; , l + ; 1 (cid:1) l + 3 (cid:21) . (C.24) ydrogenic entanglement n , l and m . For the Hydrogenic ground state,for example, the linear entropy becomes S Lin = 1 − V a π , (C.25)where we note that the fraction is dimensionless due to the appearances of the integrationvolume V and the reduced Bohr radius a . Since we take V → ∞ at the end of eachcalculation, we are left with S Lin = 1, which conventionally implies that the entropy isinfinite. However, as mentioned, the linear entropy cannot be interpreted as a measure ofentanglement in this context. A similar result for the ground state linear entropy of theHydrogen atom was reported in a Masters thesis, see Ref. [53].
Appendix D. Covariance matrix and symplectic eigenvalues
In the main text, we noted that successfully detecting entanglement in the first and secondmoments of a Gaussian state through the application of the PPT criterion immediatelyimplies that a non-Gaussian state with the same moments is also entangled [35]. In thisAppendix, we demonstrate in detail how we compute the PPT criterion for Hydrogenicsystems, by first performing a symplectic transformation from the { r , r } -basis into therelative { r , R } basis. We provide the explicit form of the symplectic transformation andcompute the symplectic eigenvalues of the final covariance matrix.Our task is to compute the entanglement between the two subsystems in the { r , r } basis. We define σ in this basis with elements given by σ = (cid:104) X i X j + X j X j (cid:105) − (cid:104) X i (cid:105) (cid:104) X j (cid:105) . (D.1)The diagonal elements of σ are given bydiag[ σ ] = 2 (cid:0)(cid:10) x (cid:11) , (cid:10) p x, (cid:11) , (cid:10) y (cid:11) , (cid:10) p y, (cid:11) , (cid:10) z (cid:11) , (cid:10) p z, (cid:11) , (cid:10) x (cid:11) , (cid:10) p x, (cid:11) , (cid:10) y (cid:11) , (cid:10) p y, (cid:11) , (cid:10) z (cid:11) , (cid:10) p z, (cid:11)(cid:1) . (D.2)Where x j , y j , z j and p x,j , p y,j , p z,j are the position and momentum coordinates of the tworespective systems.To transform into the basis of the subsystems parametrised by r and r , we perform asymplectic transformation S that maps the covariance matrix from the { r , p , r , p } basisto the { r , p , R , P } basis.Assuming that the two masses of the systems are the same, the symplectic transformation S that maps the covariance matrix from the relative and centre-of-mass basis to the subsystem ydrogenic entanglement S = − − − − − −
00 0 0 0 0 1 0 0 0 0 0 1 . (D.3)In language commonly used in the context of quantum optics, this transformation isequivalent to the combination of a 50:50 beam-splitter and a diagonal squeezer.The new covariance matrix is given by σ (cid:48) = SσS T . We prove in the following Appendixthat all the first moments in the relative basis are zero, and that σ (cid:48) is a diagonal matrix withthe following entries: σ (cid:48) = 2 diag (cid:0)(cid:10) x (cid:11) , (cid:10) p x (cid:11) , (cid:10) y (cid:11) , (cid:10) p y (cid:11) , (cid:10) z (cid:11) , (cid:10) p z (cid:11) , (cid:10) X (cid:11) , (cid:10) P X (cid:11) , (cid:10) Y (cid:11) , (cid:10) P Y (cid:11) , (cid:10) Z (cid:11) , (cid:10) P Z (cid:11)(cid:1) , (D.4)where we have defined the relative and centre-of-mass position coordinates as follows: x = x − x , X = x + x ,y = y − y , Y = y + y ,z = z − z , Z = z + z , (D.5)and, similarly, the momentum variables: p x = p x, − p x, , P X = p x, + p x, ,p y = p y, − p y, , P Y = p y, + p y, ,p z = p z, − p z, , P Z = p z, + p z, , (D.6)Once the expectation values of the coordinates in Eqs. (D.5) and (D.6) have been computed,we transform back into the original basis with σ = S − σ (cid:48) S − . To better see the explicitaction of this transformation and the return of coherence to the off-diagonal elements, wewrite σ in terms of block matrices: σ = (cid:32) σ σ σ σ (cid:33) . (D.7) ydrogenic entanglement (cid:37) ( r (cid:48) , r ) = (cid:37) ( r (cid:48) , r ), the covariancematrices for both systems are the same with σ = σ , where σ is given by σ = 2 (cid:10) x (cid:11) + (cid:10) X (cid:11) (cid:10) p x (cid:11) + (cid:10) P X (cid:11) (cid:10) y (cid:11) + (cid:10) Y (cid:11) (cid:10) p y (cid:11) + (cid:10) P Y (cid:11) (cid:10) z (cid:11) + (cid:10) Z (cid:11)
00 0 0 0 0 (cid:10) p z (cid:11) + (cid:10) P Z (cid:11) . (D.8) The off-diagonal elements, which contain the correlations between subsystem 1 and 2, aregiven by σ = 2 (cid:10) X (cid:11) − (cid:10) x (cid:11) (cid:10) P X (cid:11) − (cid:10) p x (cid:11) (cid:10) Y (cid:11) − (cid:10) y (cid:11) (cid:10) P Y (cid:11) − (cid:10) p y (cid:11) (cid:10) Z (cid:11) − (cid:10) z (cid:11)
00 0 0 0 0 (cid:10) P Z (cid:11) − (cid:10) p z (cid:11) . (D.9) We are now ready to perform the partial transposition of σ and compute its symplecticeigenvalues to determine whether the r and r subsystems are entangled. We call ν j thesymplectic eigenvalues of σ and ˜ ν j the symplectic eigenvalues of the partially transposedcovariance matrix σ Tp .In a continuous variable system, partial transposition is equivalent to introducing aminus sign to all the momenta (or all the positions) of a subsystem [35]. This implies thatevery element corresponding to the relative momentum expectation values (cid:104) p x, (cid:105) , (cid:104) p y, (cid:105) , (cid:104) p z, (cid:105) becomes inverted. The diagonal elements of σ (cid:48) Tp remain invariant, and the off-diagonalelements ˜ σ are now given by ˜ σ = 2 (cid:10) X (cid:11) − (cid:10) x (cid:11) (cid:10) p x (cid:11) − (cid:10) P X (cid:11) (cid:10) Y (cid:11) − (cid:10) y (cid:11) (cid:10) p y (cid:11) − (cid:10) P Y (cid:11) (cid:10) Z (cid:11) − (cid:10) z (cid:11)
00 0 0 0 0 (cid:10) p z (cid:11) − (cid:10) P Z (cid:11) . (D.10) We can now compute the symplectic eigenvalues of the partially transposed system. They ydrogenic entanglement i Ω σ (cid:48) Tp , and we find:˜ ν = (cid:113) (cid:104) x (cid:105) (cid:104) P X (cid:105) , ˜ ν = 4 (cid:112) (cid:104) X (cid:105) (cid:104) p x (cid:105) , ˜ ν = (cid:113) (cid:104) y (cid:105) (cid:104) P Y (cid:105) , ˜ ν = 4 (cid:113) (cid:104) Y (cid:105) (cid:10) p y (cid:11) , ˜ ν = (cid:113) (cid:104) z (cid:105) (cid:104) P Z (cid:105) , ˜ ν = 4 (cid:112) (cid:104) Z (cid:105) (cid:104) p z (cid:105) . (D.11)We note that the symplectic eigenvalues do not mix between the modes. Entanglement canthen be established once the first and second moments of the Hydrogenic system have beencomputed, which we do in the next Appendix. Due to the rotational symmetry of the system,these three pairs of quantities will take the same along each spatial direction, so that ourentanglement test will effectively reduce to a two-mode problem. Appendix E. First and second moments of a localised Hydrogenic system
In this Appendix, we compute the first and second moments of the Hydrogenic energyeigenstates for the case when the centre-of-mass wavefunction is a three-dimensional localisedGaussian wavepacket. This derivation is valid for the case when the two masses of thebipartite system are the same with m = m = m . Appendix E.1. Preliminaries
We choose to set the centre-of-mass wavepacket ϕ ( R ) to be a Gaussian wavepacket centredat the origin. It is given by ϕ ( R ) = 1 π / b / e − R · R / (2 b ) . (E.1)Then the full state is given byΨ( r , R ) = ψ nlm ( r ) 1 π / b / e − R · R / (2 b ) , (E.2)where ψ nlm ( r ) are the Hydrogenic eigenstates shown in Eq. (A.14). We proceed to computeall the elements of the covariance matrix. Given the definitions of the relative and centre-of-mass coordinates in Eqs. (D.5) and (D.6), full two-mode covariance matrix is a diagonal12 ×
12 matrix that is explicitly given by σ = 2 (cid:10) x (cid:11) − (cid:104) x (cid:105) , σ = 2 (cid:10) X (cid:11) − (cid:104) X (cid:105) ,σ = 2 (cid:10) p x (cid:11) − (cid:104) p x (cid:105) , σ = 2 (cid:10) P X (cid:11) − (cid:104) P X (cid:105) ,σ = 2 (cid:10) y (cid:11) − (cid:104) y (cid:105) , σ = 2 (cid:10) Y (cid:11) − (cid:104) Y (cid:105) ,σ = 2 (cid:10) p y (cid:11) − (cid:104) p y (cid:105) , σ , = 2 (cid:10) P Y (cid:11) − (cid:104) P Y (cid:105) ,σ = 2 (cid:10) z (cid:11) − (cid:104) z (cid:105) , σ , = 2 (cid:10) Z (cid:11) − (cid:104) Z (cid:105) ,σ = 2 (cid:10) p z (cid:11) − (cid:104) p z (cid:105) , σ , = 2 (cid:10) P Z (cid:11) − (cid:104) P Z (cid:105) . (E.3) ydrogenic entanglement σ i,j = (cid:104) X i X j + X j X i (cid:105) − (cid:104) X i (cid:105) (cid:104) X j (cid:105) , (E.4)for i (cid:54) = j . However, since the wavefunction is a product of terms that depend separatelyon each such variable, we only have to compute the diagonal elements. The advantage fromworking in this basis now becomes evident. Depending on whether the variable of interestbelongs to the relative or the centre-of-mass coordinates, we find that the other part of thestate does not influence the result. In other words, if we are computing a moment of one ofthe centre-of-mas variables R j , then we find (cid:104) R j (cid:105) = (cid:90) R d R Ψ ∗ ( r , R ) R j Ψ( r , R ) , (E.5)where R denotes the integration over the full spatial domain of the bipartite state. Whenthe state separates in terms of the relative and centre-of-mass coordinates as Ψ( r , R ) = ψ nlm ( r ) ϕ ( R ), the part of the state that does not relate to the expectation value in questionsimply satisfies the normalisation relation. For example, if we wish to compute theexpectation value (cid:104) R j (cid:105) = (cid:90) d x d y d z | ψ nlm ( r ) | (cid:90) d X d Y d Zϕ ∗ ( R ) R j ϕ ( R ) , (E.6)where the first integral is just the normalisation condition for ψ nlm ( r ), we are left with thetask to evaluate the Gaussian integral. Finally, we remark that since the transformation S , shown in Eq. (D.3), between the { r , r } basis and the { r , R } basis is symplectic, theJacobian for the substitution of variables in the integrals is equal to unity. Appendix E.2. Expectation values and variances for Gaussian centre-of-mass wavepacket
If the Gaussian wavepacket is centred at the origin, its expectation values are zero. This isthe case for the state we consider, which is shown in Eq. (E.1), and thus we find (cid:104) X (cid:105) = (cid:104) Y (cid:105) = (cid:104) Z (cid:105) = 0 . (E.7)The variances, on the other hand, are given by (cid:10) X (cid:11) = 1 π / b (cid:90) ∞−∞ d X (cid:90) ∞−∞ d Y (cid:90) ∞−∞ d Z X e − ( X + Y + Z ) /b = 1 π / b (cid:90) ∞−∞ d Y (cid:90) ∞−∞ d Ze − ( Y + Z ) /b (cid:90) ∞−∞ d X X e − X /b . (E.8)The first integral is equal to (cid:90) ∞−∞ d Y (cid:90) ∞−∞ d Ze − ( Y + Z ) /b = b π, (E.9) ydrogenic entanglement (cid:90) ∞−∞ d X X e − X /b = √ πb . (E.10)So we are left with (cid:10) X (cid:11) = 1 π / b b π √ π b = b . (E.11)This element will be the same for all other quantities. In summary, we have (cid:10) X (cid:11) − (cid:104) X (cid:105) = b , (cid:10) Y (cid:11) − (cid:104) Y (cid:105) = b , (cid:10) Z (cid:11) − (cid:104) Z (cid:105) = b . (E.12)Next, we wish to compute the centre-of-mass momentum expectation values. The momentumvariable in the X direction is given by P X = − i (cid:126) ddX , (E.13)and similarly for P Y and P Z . Taking the expectation value, we find (cid:104) P X (cid:105) = 1 π / b (cid:90) ∞−∞ d X e − ( X + Y + Z ) / (2 b ) (cid:18) − i (cid:126) ddX (cid:19) e − ( X + Y + Z ) / (2 b ) = − i (cid:126) π / b (cid:90) ∞−∞ d X e − ( X + Y + Z ) / (2 b ) (cid:18) − Xb (cid:19) e − ( X + Y + Z ) / (2 b ) = i (cid:126) π / b (cid:90) ∞−∞ d X X e − ( X + Y + Z ) /b = 0 , (E.14)where the integral vanishes because the function is odd. Due to the symmetry of thewavepacket, it follows all the expectation values of the momentum coordinates are zero.We proceed to compute the variance of the momentum variable. We find (cid:68) ˆ P X (cid:69) = 1 π / b (cid:90) ∞−∞ d X (cid:90) ∞−∞ d Y (cid:90) ∞−∞ d Z e − ( X + Y + Z ) / (2 b ) (cid:18) − i (cid:126) ddX (cid:19) e − ( X + Y + Z ) / (2 b ) = − (cid:126) π / b (cid:90) ∞−∞ d Y (cid:90) ∞−∞ d Z e − ( Y + Z ) /b (cid:90) ∞−∞ d X (cid:18) X b − b (cid:19) e − X /b . (E.15)The first two integrals over Y and Z again evaluate to that in Eq. (E.9). The second integralbecomes (cid:90) ∞−∞ d X (cid:18) X b − b (cid:19) e − X /b = − √ π (cid:114) b . (E.16) ydrogenic entanglement (cid:10) P X (cid:11) = (cid:126) π / b √ π (cid:114) b b π = (cid:126) b . (E.17)Again, due to symmetry, we find (cid:104) P X (cid:105) = (cid:104) P Y (cid:105) = (cid:104) P Z (cid:105) . So in summary, the momentumvariances become (cid:10) P X (cid:11) − (cid:104) P X (cid:105) = (cid:126) b , (cid:10) P Y (cid:11) − (cid:104) P Y (cid:105) = (cid:126) b , (cid:10) P Z (cid:11) − (cid:104) P Z (cid:105) = (cid:126) b . (E.18)We proceed to consider the Hydrogenic wavefunctions. Appendix E.3. Position expectation values of the Hydrogenic subsystems
We wish to compute the expectation values (cid:104) x (cid:105) , (cid:104) y (cid:105) , and (cid:104) z (cid:105) . We recall that the Hydrogenicwavefunctions are given by ψ nlm ( r ) = (cid:115)(cid:18) na (cid:19) ( n − l − n [( n + l )!] e − r/na (cid:18) rna (cid:19) l (cid:2) L l +1 n − l − (2 r/na ) (cid:3) Y ml ( θ, φ ) . (E.19)We can rewrite ψ nlm ( r ) in terms of a radial and an angular part ψ nlm ( r ) = R nl ( r ) Y ml ( θ, φ ) , (E.20)where we have defined R nl ( r ) = (cid:115)(cid:18) na (cid:19) ( n − l − n [( n + l )!] e − r/na (cid:18) rna (cid:19) l (cid:2) L l +1 n − l − (2 r/na ) (cid:3) . (E.21)In what follows, we will often make use of the absolute value | Y ml ( θ, φ ) | , which is independentof the variable φ . To emphasise this fact, we write this as | Y ml ( θ, φ ) | = Y ml ( θ ) . (E.22)We also note that many of the expressions that we wish to compute are easier to evaluatein spherical polar coordinates. The transformation from Cartesian coordinates to sphericalpolar coordinates is given by x = r sin θ cos φ,y = r sin θ sin φ,z = r cos θ, (E.23) ydrogenic entanglement x d y d z = r sin θ d r d θ d φ .We proceed to compute the expectation values of the position variables x , y , and z .Starting with (cid:104) x (cid:105) , the expression that we must calculate is given by (cid:104) x (cid:105) = (cid:90) ∞−∞ d X (cid:90) ∞−∞ d Y (cid:90) ∞−∞ d Z | ϕ ( R ) | (cid:90) ∞−∞ d x (cid:90) ∞−∞ d y (cid:90) ∞−∞ d z ψ ∗ nlm ( r ) x ψ nlm ( r )= (cid:90) ∞−∞ d x (cid:90) ∞−∞ d y (cid:90) ∞−∞ d z ψ ∗ nlm ( r ) x ψ nlm ( r ) . (E.24)Then, using the fact that x = r sin θ cos φ , we write (cid:104) x (cid:105) = (cid:90) ∞ d r r (cid:90) π d θ sin θ (cid:90) π d φ ψ ∗ nlm ( r ) r sin θ cos φ ψ nlm ( r )= (cid:90) ∞ d r r R nl ( r ) (cid:90) π d θ sin θ (cid:90) π d φ cos φ | Y ml ( θ, φ ) | = (cid:90) ∞ d r r R nl ( r ) (cid:90) π d θ sin θ Y ml ( θ ) (cid:90) π d φ cos φ = 0 , (E.25)which follows because the final integral over cos φ is zero. Next, we evaluate (cid:104) y (cid:105) , where wecan use the same trick and write y = r sin θ sin φ . We have (cid:104) y (cid:105) = (cid:90) ∞−∞ d X (cid:90) ∞−∞ d Y (cid:90) ∞−∞ d Z | ϕ ( R ) | (cid:90) ∞−∞ d x (cid:90) ∞−∞ d y (cid:90) ∞−∞ d z ψ ∗ nlm ( r ) y ψ nlm ( r )= (cid:90) ∞−∞ d x (cid:90) ∞−∞ d y (cid:90) ∞−∞ d z ψ ∗ nlm ( r ) y ψ nlm ( r )= (cid:90) ∞ d r r (cid:90) π d θ sin θ (cid:90) π d φ ψ ∗ nlm ( r ) r sin θ sin φ ψ nlm ( r )= (cid:90) ∞ d r r R nl ( r ) (cid:90) π d θ sin θ Y ml ( θ ) (cid:90) π d φ sin φ = 0 , (E.26)again, because the last integral over φ is zero. Finally, the expectation value (cid:104) z (cid:105) , with z = r cos θ , becomes (cid:104) z (cid:105) = (cid:90) ∞−∞ d X (cid:90) ∞−∞ d Y (cid:90) ∞−∞ d Z | ϕ ( R ) | (cid:90) ∞−∞ d x (cid:90) ∞−∞ d y (cid:90) ∞−∞ d z ψ ∗ nlm ( r ) z ψ nlm ( r )= (cid:90) ∞−∞ d x (cid:90) ∞−∞ d y (cid:90) ∞−∞ d z ψ ∗ nlm ( r ) z ψ nlm ( r )= (cid:90) ∞ d r r (cid:90) π d θ sin θ (cid:90) π d φ ψ ∗ nlm ( r ) r cos θ ψ nlm ( r )= (cid:90) ∞ d r r R nl ( r ) (cid:90) π d θ sin θ cos θ Y ml ( θ ) (cid:90) π d φ = 2 π (cid:90) ∞ d r r R nl ( r ) (cid:90) π d θ sin θ cos θ Y ml ( θ ) . (E.27) ydrogenic entanglement θ . Given theexpression for the spherical harmonics in Eq. (A.18), we write (cid:90) π sin θ cos θ d θ Y ml ( θ ) = 2 l + 14 π ( l − m )!( l + m )! (cid:90) π d θ sin θ cos θP ml (cos θ ) P ml (cos θ ) . (E.28)We then perform the substitution cos θ = u . Thus d u = − sin θ d θ , leading to d θ = − d u/ sin θ .The limits (0 , π ) become (1 , − (cid:90) π d θ sin θ cos θP ml (cos θ ) P ml (cos θ ) = − (cid:90) − d u u P ml ( u ) P ml ( u ) = (cid:90) − d u u P ml ( u ) P ml ( u ) . (E.29)We can now use the following recurrence relation for the associated Legendre polynomials: xP ml ( x ) = 1(2 l + 1) (cid:0) ( l − m + 1) P ml +1 ( x ) + ( l + m ) P ml − ( x ) (cid:1) . (E.30)Inserting this into Eq. (E.29) yields (cid:90) π d θ sin θ cos θ Y ml ( θ ) (E.31)= 2 l + 14 π ( l − m )!( l + m )! 1(2 l + 1) (cid:90) − d u (cid:0) ( l − m + 1) P ml +1 ( u ) + ( l + m ) P ml − ( u ) (cid:1) P ml ( u ) . However, the orthogonality relation for the associated Legendre polynomials states that (cid:90) − d u P ml ( u ) P ml (cid:48) ( u ) = 2( l + 1)!(2 l + 1)( l − m )! δ ll (cid:48) , (E.32)which means that the indices must match in order for the integral in Eq. (E.31) to be non-zero. In our case, we have differing indices, and thus the integral is zero. This same result isalso noted in Ref. [54]. Therefore, in summary, (cid:104) x (cid:105) = (cid:104) y (cid:105) = (cid:104) z (cid:105) = 0 , (E.33)which is to be expected since the Hydrogenic subsystem is not expected to display anynet movement in any spatial direction. We proceed with calculating the variances of theHydrogenic subsystem. Appendix E.4. Variances of the Hydrogenic subsystems
Our goal is to compute the following three variances: (cid:10) x (cid:11) , (cid:10) y (cid:11) , and (cid:10) z (cid:11) . (E.34)To do so, we first present some preliminary results that will aid our calculations. ydrogenic entanglement Appendix E.4.1. Preliminaries
We begin by introducing the Kramer–Pasternack relation,which provides a closed-form expression for any powers of the expectation value of (cid:104) r q (cid:105) . Ingeneral, the expectation value of the variable r q is given by (cid:104) r q (cid:105) = (cid:90) ∞ d r r ( r ) q R nl ( r ) , (E.35)where q is an integer number, and where R nl ( r ) is the radial wavefunction. The Kramer–Pasternack relation reads [40, 41]4( q + 1) (cid:104) r q (cid:105) − n (2 q + 1) (cid:10) r q − (cid:11) + n q [(2 l + 1) − q ] (cid:10) r q − (cid:11) = 0 . (E.36)As will become clear below, we are interested in the second order expression obtained throughthe Kramer–Pasternack relation. It is given by (cid:104) ψ nlm | r | ψ nlm (cid:105) = a n (5 n − l ( l + 1) + 1)2 . (E.37)We now proceed to compute the variances of x , y , and z one by one. In doing so, we willagain make frequent use of the transformation from Cartesian coordinates to spherical polarcoordinates shown in Eq. (E.23).Before we proceed, we also wish to evaluate the following integral, which will appear anumber of times in the calculations below: (cid:10) sin θ Y ml ( θ ) (cid:11) = (cid:90) π d θ sin θ Y ml ( θ ) . (E.38)Note the extra factor of sin θ in the integral, which arises due to the inclusion of the Jacobian.Using the explicit expressions for the spherical harmonics in Eq. (A.18), we write: (cid:10) sin θ Y ml ( θ ) (cid:11) = 2 l + 14 π ( l − m )!( l + m )! (cid:90) π d θ sin θ P ml (cos θ ) P ml (cos θ ) . (E.39)As we did above, we now let u = cos θ , so that d θ = − d u/ sin θ . We obtain (cid:10) sin θ Y ml ( θ ) (cid:11) = − l + 14 π ( l − m )!( l + m )! (cid:90) − d u (1 − u ) P ml ( u ) P ml ( u )= 2 l + 14 π ( l − m )!( l + m )! (cid:90) − d u (1 − u ) P ml ( u ) P ml ( u ) . (E.40)This integral can be divided into two parts: (cid:10) sin θ Y ml ( θ ) (cid:11) = 2 l + 14 π ( l − m )!( l + m )! (cid:18)(cid:90) − d u P ml ( u ) P ml ( u ) − (cid:90) − d u u P ml ( u ) P ml ( u ) (cid:19) . (E.41)The first integral satisfies the orthogonality relation for the associated Legendre polynomials,which reads (cid:90) − d u P ml ( u ) P ml ( u ) = 2( l + m )!(2 l + 1)( l − m )! , (E.42) ydrogenic entanglement (cid:10) sin θ Y ml ( θ ) (cid:11) = 2 l + 14 π ( l − m )!( l + m )! (cid:18) l + m )!(2 l + 1)( l − m )! − (cid:90) − d u u P ml ( u ) P ml ( u ) (cid:19) = 12 π − l + 14 π ( l − m )!( l + m )! (cid:90) − d u u P ml ( u ) P ml ( u ) . (E.43)We now use the recurrence relation in Eq. (E.30) to write the remaining integral in Eq. (E.43)as (cid:90) − d u u P ml ( u ) P ml ( u ) = (cid:90) − d u (cid:18) l − m + 12 l + 1 P ml +1 ( u ) + l + m l + 1 P ml − ( u ) (cid:19) (E.44) × (cid:18) l − m + 12 l + 1 P ml +1 ( u ) + l + m l + 1 P ml − ( u ) (cid:19) = (cid:90) − d u (cid:34)(cid:18) l − m + 12 l + 1 (cid:19) P ml +1 ( u ) P ml +1 ( u ) + (cid:18) l + m l + 1 (cid:19) P ml − ( u ) P ml − ( u ) (cid:35) , where the cross-terms vanish because of the orthogonality relation in Eq. (E.42). Then, weuse the same relation to find that (cid:90) − d u (cid:34)(cid:18) l − m + 12 l + 1 (cid:19) P ml +1 ( u ) P ml +1 ( u ) + (cid:18) l + m l + 1 (cid:19) P ml − ( u ) P ml − ( u ) (cid:35) (E.45)= (cid:18) l − m + 12 l + 1 (cid:19) l + 1 + m )!(2( l + 1) + 1)( l + 1 − m )! + (cid:18) l + m l + 1 (cid:19) l − m )!(2( l −
1) + 1)( l − − m )! . Simplifying this expression and multiplying it by the prefactor of the integral, we find (cid:10) sin θ Y ml ( θ ) (cid:11) = 2 l + 14 π ( l − m )!( l + m )! (cid:20)(cid:18) l + 1 − m l + 1 (cid:19) l + 1 + m )!(2 l + 3)( l + 1 − m )! + (cid:18) l + m l + 1 (cid:19) l − m )!(2 l − l − − m )! (cid:21) = 14 π (2 l + 1) (cid:20) ( l − m )!( l + m )! ( l + 1 − m ) l + 1 + m )!(2 l + 3)( l + 1 − m )!+ ( l − m )!( l + m )! ( l + m ) l − m )!(2 l − l − − m )! (cid:21) . (E.46)Now rename p = l + m and q = l − m to find (cid:10) sin θ Y ml ( θ ) (cid:11) = 14 π (2 l + 1) (cid:20) q ! p ! ( q + 1) p + 1)!(2 l + 3)( q + 1)! + q ! p ! ( p ) p − l − q − (cid:21) . (E.47)We then note the following simplifications:( p + 1)! p ! = p + 1 , p !( p + 1)! = 1 p + 1 , ( p − p ! = 1 p , p !( p − p . (E.48) ydrogenic entanglement (cid:10) sin θ Y ml ( θ ) (cid:11) = 14 π (2 l + 1) (cid:20) q ! p ! ( q + 1) p + 1)!(2 l + 3)( q + 1)! + q ! p ! ( p ) p − l − q − (cid:21) = 14 π (2 l + 1) (cid:20) q + 1 ( p + 1) 22 l + 3 ( q + 1) + p q p l − (cid:21) = 12 π (2 l + 1) (cid:20) ( p + 1)( q + 1) 12 l + 3 + pq l − (cid:21) . (E.49)Inserting the original expressions p = l + m and q = l − m , we simplify and find that (cid:10) sin θ Y ml ( θ ) (cid:11) = l + l + m − π (2 l − l + 3) . (E.50)With this relation and the Kramer–Pasternack relation for (cid:104) r (cid:105) in Eq. (E.37), we are readyto compute the remaining variances for the Hydrogenic wavefunctions. Appendix E.4.2. Calculating (cid:104) x (cid:105) The integral for this expectation value is given by (cid:10) x (cid:11) = (cid:90) ∞−∞ d x (cid:90) ∞−∞ d y (cid:90) ∞−∞ d z x | ψ nlm ( r ) | . (E.51)Again splitting the wavefunction into a radial part R nl ( r ) and an angular part | Y ml ( θ, φ ) | = Y ml ( θ ), and using the fact that x = r sin θ cos φ , we find (cid:10) x (cid:11) = (cid:90) ∞ d r r R nl ( r ) (cid:90) π d θ sin θ Y ml ( θ ) (cid:90) π d φ cos φ = π (cid:90) ∞ d r r R nl ( r ) (cid:90) π d θ sin θ Y ml ( θ ) . (E.52)We then proceed to use the Kramer–Pasternack relation in Eq. (E.37) to evaluate the radialintegral, and the result for the angular integral listed in Eq (E.50) to find (cid:10) x (cid:11) = a n (5 n − l ( l + 1) + 1)2 l + l + m − l − l + 3) . (E.53) Appendix E.4.3. Calculating (cid:104) y (cid:105) We now proceed with (cid:104) y (cid:105) , where we recall that y = r sin θ sin φ . We find (cid:10) y (cid:11) = (cid:90) ∞ d r r (cid:90) π d θ sin θ (cid:90) π d φ y | ψ nlm ( r ) | = (cid:90) ∞ d r r R nl ( r ) (cid:90) π d θ sin θ (cid:90) π d φ sin φ | Y ml ( θ, φ ) | . (E.54)We start with the last integral, which, since | Y ml ( θ, φ ) | = Y ml ( θ ) is independent of φ , becomes (cid:90) π d φ sin φ = π. (E.55) ydrogenic entanglement (cid:10) y (cid:11) = a n (5 n − l ( l + 1) + 1)2 l + l + m − l − l + 3) . (E.56) Appendix E.4.4. Computing (cid:104) z (cid:105) This is a bit different, since z = r cos θ . We find (cid:10) z (cid:11) = (cid:90) ∞ d r r R nl ( r ) (cid:90) π d θ sin θ r cos θ Y ml ( θ ) (cid:90) π d φ = 2 π (cid:90) ∞ d r r R nl ( r ) (cid:90) π d θ sin θ cos θ Y ml ( θ ) . (E.57)Starting with the angular integral, and using the expression for the spherical harmonics inEq. (A.18), we write (cid:90) π d θ sin θ cos θ Y ml ( θ ) = 2 l + 14 π ( l − m )!( l + m )! (cid:90) π d θ sin θ cos θ P ml (cos θ ) P ml (cos θ ) . (E.58)We again perform the substitution u = cos( θ ), and we find d u = − d θ sin( θ ), so thatd θ = − d u/ sin( θ ), which gives (cid:90) π d θ sin θ cos θ P ml (cos θ ) P ml (cos θ ) = − (cid:90) − d u u P ml ( u ) P ml ( u )= (cid:90) − d u u P ml ( u ) P ml ( u ) . (E.59)We already obtained the answer to this quantity in Eq. (E.45). The result is2 l + 14 π ( l − m )!( l + m )! (cid:90) − d u u P ml ( u ) P ml ( u )= 2 l + 14 π ( l − m )!( l + m )! (cid:20)(cid:18) l − m + 12 l + 1 (cid:19) l + 1 + m )!(2( l + 1) + 1)( l + 1 − m )!+ (cid:18) l + m l + 1 (cid:19) l − m )!(2( l −
1) + 1)( l − − m )! (cid:21) , (E.60)which can be simplified, so that we ultimately find (cid:90) π d θ sin θ cos θ Y ml ( θ ) = 12 π − l − l + 2 m − l − l . (E.61)Then again using the Kramer–Pasternack relation in Eq. (E.37), we find that (cid:10) z (cid:11) = a n (5 n − l ( l + 1) + 1)2 1 − l − l + 2 m − l − l . (E.62)We now have all the covariance matrix elements, and we are ready to compute theentanglement. ydrogenic entanglement Appendix E.4.5. Summary
In summary, we have that (cid:10) x (cid:11) = a n (5 n − l ( l + 1) + 1)2 l + l + m − l − l + 3) , (cid:10) y (cid:11) = a n (5 n − l ( l + 1) + 1)2 l + l + m − l − l + 3) , (cid:10) z (cid:11) = a n (5 n − l ( l + 1) + 1)2 1 − l − l + 2 m − l − l . (E.63)The fact that (cid:104) z (cid:105) differs from the other two variances is reasonable given the identificationof the rotation axis. Appendix E.5. Expectation values of p x , p y and p z The expectation values can be quickly computed by realising that we obtain the same integralsover the spherical harmonics that we already evaluated in Section Appendix E.3.We start with the first expression (cid:104) p x (cid:105) , however we will first compute (cid:104) k x (cid:105) and thenmultiply by (cid:126) to obtain the momentum. This quantity is easier to compute in the momentumbasis. We therefore use the Fourier transform in Eq. (23) and Eq. (24), to write (cid:104) k x (cid:105) = (cid:90) ∞ d x (cid:90) ∞ d y (cid:90) ∞ d z ψ ∗ nlm ( r ) (cid:18) − i ∂∂x (cid:19) ψ nlm ( r )= − i (2 π ) / (cid:90) d k (cid:48) ˜ ψ ∗ ( k (cid:48) ) (cid:90) d k ik x ˜ ψ ( k ) (cid:90) ∞ d x (cid:90) ∞ d y (cid:90) ∞ d z e i ( k − k (cid:48) ) · r = − i (cid:90) d k (cid:48) ˜ ψ ∗ ( k (cid:48) ) (cid:90) d k ik x ˜ ψ ( k ) δ ( k − k (cid:48) )= (cid:90) d k ˜ ψ ∗ ( k ) k x ˜ ψ ( k ) . (E.64)The Hydrogen wavefunction in the momentum representation is given in Eq. (A.21) in termsof the wavevector k . Using this expression, and moving to spherical coordinates ( k, θ, φ ),where k x = k cos θ sin φ , we find (cid:104) k x (cid:105) = (cid:90) ∞−∞ d k x (cid:90) ∞−∞ d k y (cid:90) ∞−∞ d k z ˜ ψ ∗ ( k ) k x ˜ ψ ( k )= (cid:90) ∞ d k k (cid:90) π d θ sin θ (cid:90) π d φ k cos θ sin φ | F nl ( k ) | | Y ml ( θ, φ ) | = (cid:90) ∞ d k k | F nl ( k ) | (cid:90) π d θ sin θ cos θ Y ml ( θ ) (cid:90) π d φ sin φ, (E.65)where we again used the fact that | Y ml ( θ, φ ) | = Y ml ( θ ), which is independent of φ . However,we again note that the last integral is again zero, and thus (cid:104) p x (cid:105) = 0. In fact, since the angularintegral is the same for the position and momentum expectation values, we conclude thatwe also have (cid:104) p y (cid:105) = (cid:104) p z (cid:105) = 0. Intuitively, this is reasonable, since a non-zero momentumexpectation value would mean that the system has a net non-zero motion in one of thedirections. ydrogenic entanglement Appendix E.6. Variance of p x , p y and p z We proceed to compute the variances of the relative momentum variables. In the wavevectorbasis, we find that (cid:10) k x (cid:11) = (cid:90) ∞−∞ d kp x (cid:90) ∞−∞ d kp y (cid:90) ∞−∞ d k z ˜ ψ ∗ nlm ( k ) k x ˜ ψ nlm ( k )= (cid:90) ∞ d p p | F nl ( p ) | (cid:90) π d θ sin θ (cid:90) π d φ cos φ | Y ml ( θ, φ ) | = π (cid:90) d p p | F nl ( p ) | (cid:90) d θ sin θ Y ml ( θ, φ ) . (E.66)By using the result for the weighted Gegenbauer polynomials in Eq. (B.15), and theexpectation value for the angular integral in Eq. (E.50), we find that (cid:10) k x (cid:11) = (cid:18) a n (cid:19) l + l + m − l − l + 3) . (E.67)By symmetry, and using the results we derived in the previous section, the other valuesbecome: (cid:10) k y (cid:11) = (cid:18) a n (cid:19) l + l + m − l − l + 3) , (E.68)and, using the result in Eq. (E.61), we have (cid:10) k z (cid:11) = 2 π (cid:18) a n (cid:19) (cid:10) sin θ cos θ Y ml ( θ ) (cid:11) = (cid:18) a n (cid:19) − l − l + 2 m − l − l . (E.69) Appendix E.7. Summary of expectation values and variances
The expectation values and variances for the centre-of-mass variables are (cid:104) X (cid:105) = (cid:104) Y (cid:105) = (cid:104) Z (cid:105) = 0 , (cid:10) X (cid:11) = (cid:10) Y (cid:11) = (cid:10) Z (cid:11) = b , (cid:104) P X (cid:105) = (cid:104) P Y (cid:105) = (cid:104) P Z (cid:105) = 0 , (cid:10) P X (cid:11) = (cid:10) P Y (cid:11) = (cid:10) P Z (cid:11) = (cid:126) b . (E.70) ydrogenic entanglement (cid:104) x (cid:105) = (cid:104) y (cid:105) = (cid:104) z (cid:105) = 0 , (cid:10) x (cid:11) = (cid:10) y (cid:11) = a n (5 n − l ( l + 1) + 1)2 l + l + m − l − l + 3) , (cid:10) z (cid:11) = a n (5 n − l ( l + 1) + 1)2 1 − l − l + 2 m − l − l , (cid:104) p x (cid:105) = (cid:104) p y (cid:105) = (cid:104) p z (cid:105) = 0 , (cid:10) p x (cid:11) = (cid:10) p y (cid:11) = (cid:18) (cid:126) a n (cid:19) l + l + m − l − l + 3) , (cid:10) p z (cid:11) = (cid:18) (cid:126) a n (cid:19) − l − l + 2 m − l − l .l .