Physical entanglement between localized orbitals

The goal of the present work is to guide the development of quantum technologies in the context of fermionic systems. For this, we first elucidate the process of entanglement swapping in electron systems such as atoms, molecules or solid bodies. This demonstrates the significance of the number-parity superselection rule and highlights the relevance of localized few-orbital subsystems for quantum information processing tasks. Then, we explore and quantify the entanglement between localized orbitals in two systems, a tight-binding model of non-interacting electrons and the hydrogen ring. For this, we apply the first closed formula of a faithful entanglement measure, derived in (arXiv:2207.03377) as an extension of the von Neumann entropy to genuinely correlated many-orbital systems. For both systems, long-distance entanglement is found at low and high densities η, whereas for medium densities, η≈12 , practically only neighboring orbitals are entangled. The Coulomb interaction does not change the entanglement pattern qualitatively except for low and high densities where the entanglement increases as function of the distance between both orbitals.


I. INTRODUCTION
Entanglement is an essential resource for quantum information processing (QIP) [1] and it has been intensively studied in the context of distinguishable particles.In the era of the second quantum revolution [2,3], however, the new primary platforms for executing QIP tasks will likely comprise atoms and molecules.It is therefore of utmost importance to identify and quantify the physical entanglement (that which can be physically accessed), for systems of indistinguishable particles, particularly electrons.
That being said, quantifying entanglement in electronic systems is by no means a straightforward task.Even more urgently, there is a misalignment between the common treatment and quantification of electronic entanglement [4][5][6][7][8][9][10][11][12][13] and the goal of utilizing it for QIP tasks [14][15][16][17][18].It is therefore the motivation of our present work to communicate and elucidate three key aspects which are preferable if not indispensable for comprehensively achieving the latter: (i) localized molecular orbitals, on which local operations could be performed to realize QIP tasks, (ii) operationally meaningful entanglement quantification, which evaluates precisely the useful entanglement for QIP, and (iii) computable measures for entanglement.
To recognize the underdevelopment of point (i), we recall that common orbitals used in quantum chemistry are typically not localized.Indeed, by forming delocalized molecular orbitals one can reduce the complexity of the electronic ground state problem since these orbitals tend to correlate less with one another compared to fullylocalized atomic orbitals.In fact, optimizing the orbital basis to lower the correlation and thus the computational cost in solving the ground state problem has become a hot topic [19][20][21][22].However, this objective is in direct * c.schilling@lmu.decontradiction with our ambition to use orbital entanglement for QIP purposes.A large family of QIP protocols uses entanglement as a resource to overcome spatial separation [23][24][25][26][27][28], for instance, teleporting a quantum state.If we are to realize these protocols on a molecule, we must head in the direction opposite to the common approach of delocalizing molecular orbitals.The entanglement needs to be as high as possible and distributed across separated local regions, allowing two local agents to harness it with local operations.To achieve this, one must construct and utilize orbitals localized to as few atomic centers as possible, which inevitably comes with a significant computational cost.
The necessity of point (ii) originates from the superselection rules [29,30] (SSR) for fermions.Performing QIP tasks with electronic orbitals essentially comes down to implementing local operations on individual orbitals [31].Yet, not all local operations can be realized due to nature's restriction ensuring that the number-parity is locally preserved [29].This in turns implies that entanglement arising from breaking the number parity, the socalled "fluffy bunny" entanglement [32], cannot be utilized for our purposes.Ignoring in concrete calculations the SSR would lead to a gross overestimation of the physical entanglement [33][34][35], the entanglement that can be utilized, and thus the working capacity of the electronic system at hand for QIP applications.Therefore the importance of addressing the effect of SSR on the available entanglement in the system can hardly be overemphasized.
As for (iii), although there is no shortage of meaningful definitions of entanglement measures, thanks to the abstract mathematical foundation provided by quantum information theory, they are mostly notoriously difficult to calculate.As a consequence, analytic treatment of orbital entanglement in many-body physics is almost always restricted to the so-called block entropy [36][37][38].This quantity, despite its analytic computability [39], suffers from two critical shortcomings.First, it quantifies the entanglement between the block and its complement only when the overall system is in a pure state.Otherwise, it effectively measures the entanglement of the block and the rest of the universe.Second, it is invariant under unitary transformations within the block.In other words, it alone cannot capture any internal structure such as the entanglement between orbitals, which is the sought-after resource for QIP.To fill this vital gap, Ref. [40] recently provided the first closed formula for the faithful entanglement between electronic orbitals.
It is therefore one of the key achievements of this paper to lay out and accentuate these three essential points (i)-(iii) to the quantum science and technology community.By referring to an entanglement swapping protocol, we demonstrate that the useful orbital entanglement for QIP tasks is quantified with SSR taken into account.To showcase the analytic tools provided in Ref. [40], we then study the exactly solvable tight-binding model.There, we break the usual paradigm of restricting analytic studies to the block entropy by performing a fully analytic analysis of the physical entanglement between any two spatial orbitals.Lastly, we upgrade the tight-binding model to the hydrogen ring and analyze the orbitalorbital entanglement therein.By tackling the numerical complexity from both the orbital localization and the Coulomb interaction, we manage to realize potent entanglement between two far separated fully-localized orbitals.
The paper is structured as follows: In Section II, we first review the relevant notions of fermionic entanglement.After establishing the mathematical framework, we demonstrate the necessity of SSRs by discussing an entanglement swapping protocol.Moreover, we restate the analytic formula for orbital entanglement.In Section III we calculate the orbital-orbital entanglement in the tight-binding model and evaluate by analytical means its asymptotic behavior.Finally, in Section IV we explore the entanglement pattern between localized orbitals in the ground state of the hydrogen ring.

II. CONCEPTS
In this Section, we first review the relevant quantum information concepts and explain how to operationally transfer these concepts to fermionic systems.Secondly, we demonstrate the effect of superselection rules on the accessible entanglement in fermionic systems, by invoking an entanglement swapping protocol.Lastly, we briefly restate the entanglement formula for two orbitals derived in Ref. [40].

A. Locality and entanglement
The notion of locality is the starting point for the theory of entanglement.We consider two distant parties, Alice and Bob, each possessing a local quantum system with a finite-dimensional Hilbert space denoted by H A and H B , respectively.The two systems are local in the sense that an action of one party on their subsystem has no effect on the other.Mathematically speaking, this separation is described by the tensor product structure of both the total Hilbert space H = H A ⊗ H B and the total algebra of observables With the stage set, a quantum state ρ AB shared by Alice and Bob is said to be separable, if and only if it can be prepared by local operations and classical communication (LOCC) [41].Local operations allow Alice and Bob to prepare any product state ρ A ⊗ ρ B , and classical communication can be used to generate convex mixtures of those states.This then gives rise to the convex set of bipartite separable states D sep : Definition II.1 (Bipartite Separability [41]).A bipartite state ρ AB on H A ⊗ H B is separable if and only if it can be expressed as with p i ≥ 0 and i p i = 1.The set of separable states is denoted as D sep .Any state that is not in D sep is necessarily entangled.
A non-separable, or entangled state contains valuable entanglement resource for QIP tasks [23][24][25][26][27][28].This resourcefulness of entanglement is distinct from the mutual information, which contains both quantum and classical correlation [42].Therefore the detection and quantification of entanglement in a bipartite system is of crucial interest.For a pure state ρ AB = |Ψ AB ⟩⟨Ψ AB |, it is easy to verify whether or not it belongs to D sep .It suffices to check if it is a pure product state.Or equivalently one can check the purity of the reduced state The degree of mixedness of ρ A/B thus accounts for the entanglement in ρ AB .This leads to a well-defined entanglement measure for pure states where S(ρ) = −Tr[ρ log(ρ)] is the von Neumann entropy.For mixed states, however, there is no such straightforward analytic check for separability, except for qubitqubit and qubit-qutrit systems [43,44].In fact, even the problem of determining the separability of a general mixed quantum state is NP-hard [45,46], and in practice can only be approximately solved with significant computational efforts using symmetric extension and semidefinite programming [47,48].On the formal level, the degree of entanglement of a mixed state can be defined geometrically as the distance to the closest separable state (here and in the following we occasionally skip the index AB of the joint state) where the generalized distance function is chosen to be the relative entropy for its rich quantum informational meaning.To be more precise, the quantum relative entropy S(ρ∥σ) exponentially suppresses the "difficulty" of telling ρ and σ apart [49], in the spirit of the so-called quantum Sanov's theorem [49].The measure (3) is called the relative entropy of entanglement [50], and will be used as the standard entanglement measure throughout this paper.When applied to pure states, (3) reduces to the pure state measure via the von Neumann entropy (2) [49].

B. Fermionic locality and superselection rules
In the previous section, we explained that the notion of entanglement necessitates two ingredients of locality: (1) a tensor product structure of the Hilbert space, and (2) a tensor product structure of the algebra of observables.In the context of fermionic systems, this raises the following questions.Does there also exist a tensor product structure in a fully-antisymmetrized N -fermion Hilbert space, e.g. for a hydrogen molecule with two electrons?And how can one identify two commuting algebras of observables when fermionic creation/annihilation operators anticommute?
The fermionic analogue of the bipartition into distant laboratories corresponds to a choice of splitting the one particle Hilbert space H 1 (or a subspace of it) into two orthogonal subspaces.For instance, H 1 can be realized by an array of Wannier orbitals, or an arrangement of molecular spin orbitals.From here on we shall refer to the elements of an orthonormal basis of H 1 as orbitals.A partition of an orbital basis {|φ i ⟩} d i=1 into two complementing subsets A and B dissolves the one particle Hilbert space into a direct sum 1 .This for example can correspond to separating a lattice into the left and right parts at some cut.Each part generates a local Fock space and together they form a tensor product of the total Fock space F(H 1 ) ∼ = F(H through the following basis identification where n 1 , . . ., n d ∈ {0, 1}, ( †) i denotes the operator that annihilates (creates) a fermion in state |φ i ⟩.
As for the second necessary ingredient, i.e., commuting local algebras of observables, it seems at first sight to be in direct contradiction with the nature of fermionic operators.Suppose Alice and Bob each possess a fermionic mode at the opposite ends of the universe.Due to the anticommutation relation between the creation operators, Alice's action could in principle still influence Bob's mode, and accordingly realize superluminal signaling [51,52].This manifest violation of special relativity is in reality not possible, due to nature's restriction on the local operations by Alice and Bob, which can be formulated as so-called superselection rules [29,30] (SSR).The relevant superselection rules depend on the conserved observables in the theory or even practical limitations.Nonetheless, the parity superselection rule [29] (P-SSR) is ubiquitous: one cannot superpose quantum states with different local fermionic parity.Equivalently, physical local observables must always commute with the local parity operator [32,33,53,54] Π ≡ exp(iπ l nl ) = P+ − P− where the sum is over all local orbital degrees of freedom, P± is the orthogonal projection onto the even(+)/odd(−) subspace, and nl ≡ f † l f l denotes the particle number operator for mode |φ l ⟩.The local observables Ô fulfilling this requirement are polynomials of f † i f † j , f i f j , and f † i f j , and satisfy where G Π denotes the projective map onto the P-SSR compatible sub-algebra.Then, it follows directly that any two observables Â and B referring to two distinct sets of fermionic modes and satisfying locally (7) necessarily commute with each other.In other words, algebras of parity preserving observables on two distinct fermionic subsystems do exhibit a tensor product structure.More generally, any physical quantum channel E (completely positive and trace-preserving map [55]) on the reduced state must satisfy C. Entanglement swapping under SSR Nature's restriction on the implementable local observables has a profound implication on the entanglement in fermionic states.From the perspective of state tomography, coherent terms between different local parity states can never be observed under P-SSR, for the observable needed for such measurement is simply not available by the laws of physics.Therefore, the entanglement arising from SSR-violating coherence, also known as the "fluffy bunny entanglement" [32], is not accessible.This reduction of entanglement due to SSR is proven rigorously in Ref. [33] on an abstract mathematical level.Here we demonstrate it in concrete terms with a swapping protocol.
In Figure 1 we depict schematically an entanglement swapping protocol between two fermionic orbitals (φ A , φ B ) in the state ρ AB , and two qubits (q A , q B ) in the state |00⟩, both shared by Alice and Bob.We shall denote the orbital Fock spaces as F φ A/B and the qubit Hilbert spaces as H q A/B .For both the Fock spaces and the qubit Hilbert spaces, we denote the local basis as B = {|0⟩, |1⟩}.
FIG. 1. Entanglement swapping protocol under SSR between the two electronic orbitals φA and φB, and the two qubit registers qA and qB.See text for more details.
One should keep in mind that for the qubit space B is the computational basis, while for the Fock space, B refers to the occupational basis.The protocol in Figure 1 extracts the entanglement between the fermionic modes onto the qubit registers, on which this entanglement could be used to perform more flexible quantum information protocols.
The protocol consists of Alice swapping the two states on φ A and q A , and Bob doing the same on φ B and q B .Specifically, the SWAP channel S on a bipartite system is defined element-wise as for all i, j, k, l ∈ {0, 1}.In a world without SSR, after the SWAP operation, the final state of the two qubits would be exactly ρ AB and it would accordingly contain exactly the same amount of entanglement as the original two-orbital state.However, the SWAP operation clearly does not preserve the local parities, and hence cannot be realized with physical fermionic operators.Instead, Alice/Bob can only perform the superselected channel To see the difference between S A/B and SA/B , let us first observe the action of the SWAP channel on a pure state , which clearly violates the P-SSR.In a world free of SSR, |ψ A ⟩ would be invariant under the SWAP channel.Under the P-SSR, however, the action of the physically implementable SA consists of three steps: (1) 2 ⊗|+⟩⟨+|.Namely, the coherent terms between different fermionic parities are erased.(2) S A swaps the two tensor states to 2 , which again contains P-SSR violating terms.(3) G Π A then eliminates said terms and the final state is In this example we see evidently that SWAP does not leave |ψ A ⟩ invariant, despite its apparent symmetry.Instead, P-SSR filters out any coherence between different parity sectors, and only the superselected state ρP ] can be transferred onto the qubit register.Moreover, the final state on the fermionic mode φ A is also the superselected variant of the original qubit state.Now that we understand the action of ( 10), we can proceed to write down the action of the protocol S ≡ SA ⊗ SB on the composite Hilbert space where and the same for σP AB .The final state ρP AB on the qubit registers contains the same information as the original two-orbital state ρ AB , except it is rid of the coherence between the different local parity sectors, which does not contribute to the physical entanglement.The decohering map G Π A ⊗G Π B thus acts as an entanglement filter as depicted in Figure 1, and effectively reduces the available entanglement in the fermionic state.Only the entanglement within each local parity sector can be extracted for operational purposes.

D. Analytic formula for orbital-orbital entanglement
In the previous section, we concluded that the physical entanglement in a bipartite fermionic state ρ AB in the presence of P-SSR, is precisely the entanglement in its superselected variant ρP ] quantified in the usual manner without P-SSR.In this Section, we will briefly summarize the analytic formula for physical entanglement between two spatial orbitals in electronic systems, derived in Ref. [40].
Consider two orbitals A and B with internal spin-1 2 degrees of freedom.This system is described by a state ρ AB acting on the Fock space The mixedness of ρ AB accounts for the physically relevant scenario where A and B are in contact with other orbitals.A key step in quantifying the entanglement E(ρ P AB ) is understanding its restriction to several smaller subspaces of the total Fock space, after imposing several common symmetries exhibited by realistic condensed matter and quantum chemical systems.
The relevant systems in our case, the tight-binding model and the hydrogen ring, enjoy many symmetries, including (1) particle number, (2) magnetization and (3) reflection symmetry between orbitals A and B (e.g. two orthogonalized 1s orbitals in a hydrogen ring).Under these symmetries, there are only four entangled pure states compatible with local P-SSR: Naturally, if ρP AB enjoys the symmetries (1), (2), and (3), then |Φ ± ⟩ and |Ψ ± ⟩ are eigenstates of ρP AB , with eigenvalues p ± and q ± , respectively.Even so, the exact P-SSR orbital-orbital entanglement, though entirely analytic, has a heavily involved form, which can be found in the Appendix of Ref. [40].Instead, we highlight another related entanglement formula from Ref. [40] in the case where the particle number superselection rule (N-SSR) applies.In this setting, no superposition between different local particle number states is possible and therefore even |Φ ± ⟩ are forbidden.The N-SSR compatible twoorbital state is obtained with the same reasoning via the superselection map where Pn A ,n B = Pn A ⊗ Pn B and Pn A/B is the projection onto the sector of local particle number n A/B on subsystem A/B.The N-SSR entanglement of ρN AB can then be derived, after considerable efforts, to be where Note that r < t, or explicitly is the exact entanglement criterion.Below we will see that the P-and N-SSR orbital-orbital long-distance entanglement on a localized orbital chain are virtually indistinguishable.The simple form of the latter thus allows us to derive analytically asymptotic properties of the physical entanglement between orbitals for both Pand N-SSR.We remark that ( 15) is only applicable for states with the desired symmetry.Applying our formula to symmetry-broken states could lead to erroneous results.For example, according to (15), the entanglement of In the reminder of our paper we will apply the entanglement measure (15) and the corresponding one for the weaker P-SSR to two model systems in order to investigate their orbital-orbital entanglement as a function of, e.g., the electron density and the distance between both localized orbitals.This will provide new insights which may guide the utilization of the abundant entanglement resource for QIP tasks performed in realistic materials by tuning the particle density and interaction strength.

III. FREE FERMIONS
The full benefit of having an analytic formula for orbital entanglement is best shown for systems where the ground state is also analytically available.Therefore, we turn to the exactly solvable free fermion models.Despite the absence of particle interaction, free fermion systems can exhibit various interesting phenomena such as topological phases [56][57][58] and entanglement phase transition [59,60].In this Section, we will focus on the tight-binding model and explore its entanglement pattern with our analytic tools.But instead of the commonly used block entropy [36][37][38], we will analytically derive the physical entanglement between any two spatial orbitals.
A. Reduced states in free fermions systems A general particle number conserving non-interacting fermionic Hamiltonian can be written as a sum of onebody terms The spectrum of H free is analytically solved by finding a unitary matrix U such that d = UhU † is diagonal.We then rewrite Ĥfree as The N -particle ground states of Ĥfree is simply a Slater determinant state of the form with ground state energy N i=1 d ii , assuming the d ii 's are ordered increasingly.
To study the entanglement within the ground state |Ψ N ⟩, a natural object to consider is the reduced density matrix ρ A of the ground state on a sub-lattice A. Its entropy quantifies the entanglement between the subsystem A and its complement, and serves as a widely used analytic tool [36-38, 61, 62].Note that the lattice basis where such partition is performed is not necessarily the one in which h is diagonalized.For the ground states of Hamiltonians of the form (18), according to Wick's theorem the expectation value of any observable is a function of the two-point correlator ⟨Ψ N |f † i f j |Ψ N ⟩ ≡ γ ji .The latter defines the one-particle reduced density matrix (1RDM) γ of the ground state.As ρ A is uniquely defined via the expectations of all physical observables on the subsystem A [63], one can then in principle express ρ A using γ.While the explicit expression can be cumbersome, spectra of ρ A and γ can be easily related via the remarkable shortcut put forward by Peschel [39].The shortcut makes use of two insights: (1) The reduced state , where K is the normalization factor.(2) The defining matrix C = (C ij ) and the 1RDM restricted to the sub-algebra on subsystem A, γ A , can be simultaneously diagonalized.Combining these two ingredients, Peschel discovered that the two matrices are related via This shortcut greatly simplifies the complexity of computing the spectrum of ρ A : one simply needs to find the restricted 1RDM and its eigenvalues.As the von Neumann entropy is only a function of the spectrum, S(ρ A ) can also be readily calculated, without explicit construction of C or ρ A .In general, however, not all correlation quantities can be determined with the spectrum of the reduced states alone.In particular, suppose A contains two orthonormal orbitals.The entanglement between the two orbitals is not invariant under rotations between them and cannot be retrieved using the same method for S(ρ A ).In the next Section, we will break the current paradigm by calculating explicitly the orbital-orbital entanglement in the tight-binding model, using the analytic formula (15).

B. Formula for orbital-orbital entanglement in a tight-binding model
We consider the following tight-binding Hamiltonian for electrons on a periodic chain with L lattice sites with periodic boundary condition imposed.Here, f ( †) lσ annihilates(creates) an electron with spin σ at the Wannier orbital/lattice site labeled by l.Each Wannier orbital can host up to two electrons with different spin.The discrete translation invariance of the system allows us to diagonalise the Hamiltonian via a Fourier transform In this basis, the Hamiltonian ( 21) is diagonal, and the one-electron spectrum is given by the dispersion relation where in case of L even only one of the signs ± L−1 2 should be chosen.As the system is non-interacting, the N -electron ground state is the configuration where the N lowest energy levels are filled.For simplicity, we shall assume from now on that the number of electrons is either N = 0, or N = 4k max + 2 where k max ≥ 0 denotes the highest occupied momentum.In the former case, the ground state is simply the vacuum.In the latter, the Nparticle ground state |Ψ N ⟩ is uniquely characterized by the Fermi level where Θ is the Heaviside step function.With ( 22) and ( 24) at hand we are ready to compute the matrix elements of the 1RDM in the spatial basis as When l = l ′ , all summands in (25) become 1, and thus we obtain where η = N/(2L) is the particle filling fraction.This agrees with the physical intuition that the electrons move freely along the chain and appear at every site with uniform probability.When l ̸ = l ′ , by introducing ω ≡ 2π(l−l ′ )

L
, we arrive at In the thermodynamic limit where N, L → ∞ while η = N/(2L) is kept fixed, the off-diagonal elements of the 1RDM become lim where Next, we use Wick's theorem to calculate the quantities needed for the entanglement between the orbitals l and l ′ from the elements of the 1RDM derived above.The relevant operators for calculating r and t in the analytic formula (15) are the orthogonal projection P ll ′ −− onto the local odd parity subspace, and the orthogonal projections P ll ′ Ψ± onto |Ψ ± ⟩, all defined on the orbitals l and l Here, nlσ = f † lσ f lσ denotes the spin-σ particle number operator on orbital l.After some straightforward but lengthy calculation, we obtain the following expressions for the parameters r and t, required in Eq. ( 15), where With this, the N-SSR entanglement can be neatly expressed as provided the entanglement criteria ( 17) (here it translates to A < 2B) is met, and otherwise it is zero.We present in Figure 2 the orbital-orbital entanglement as a function of the filling fraction η, at various orbital separations d, together with the entanglement when only P-SSR is assumed, also obtained analytically (formula can be found in Ref. [40]).

C. Discussion of orbital-orbital entanglement in tight-binding model
In this section we discuss the orbital-orbital entanglement in the tight-binding model based on (32).Figure 2 reveals the rich entanglement structure of the Slater determinant ground state.At different particle filling fraction η, the orbital-orbital entanglement behaves in qualitatively distinct ways, which we will now account for in detail.
First, we notice that the graph of orbital-orbital entanglement is symmetric around the half-filling point η = 1/2.This is due to the unbroken particle-hole symmetry of the Hamiltonian (21) in the ground states.Moreover, using (32), it is easy to prove that the entanglement as a function of η is symmetric about η = 1/2.Because of this, it suffices to consider only 0 ≤ η ≤ 1/2, and treat it as a particle or hole filling fraction.In order to illustrate more details of this η-regime, we use a log-log scale in the right panel of Figure 2.
Second, the analytic properties of the orbital-orbital entanglement highly depend on the combination of orbital separation d and the filling fraction η.The common feature here is the vanishing of orbital-orbital entanglement for all d at both, η = 0 and η = 1.In these two cases the ground state is invariant under any orbital rotation and contains zero correlation in any orbital basis.For two neighboring orbitals (d = 1), their entanglement is maximized exactly when the chain is half filled (η = 1/2), for both P-and N-SSR entanglement.Around half-filling, the P-SSR entanglement is significantly higher than the N-SSR one.This is not surprising, as N-SSR imposes a stronger restriction on the accessible entanglement.Yet, at small particle or hole filling fractions, the two become essentially the same.The reason for this is that in such regimes the weight in the P 1,1 sector dominates over the weight in the P 2,0 +P 0,2 sector, as the particles/holes tend to avoid sitting on the same orbital due to the repulsive interaction.This weight domination effect is stronger when the orbitals are further apart.For d ≥ 2, the Pand N-SSR entanglement is practically indistinguishable.Thus we can restrict the investigation of analytic properties to the latter by referring to (32) for all practical purposes here, and also in the following section.In this case, two maxima occur which appear symmetrically on the η-axis about η = 1/2 (see the left panel of Figure 2) due to the particle-hole symmetry of (32).At lower η, more and farther separated orbital pairs become entangled, though this entanglement grows vanishingly small with the filling fraction.
Third, we find that entanglement across extremely long distances can be engineered by tuning the filling fraction to be infinitesimally small.In the right panel of Figure 2, we see more clearly on the log-log scale that for each separation d, there exists a critical value η(d), below which the two orbitals are always entangled.At the same time, for any filling fraction η, there is a minimal disentangling distance d min (η), beyond which all orbital pairs are in a separable state.This phenomenon of completely vanishing entanglement is called the sudden death of entanglement [64].Geometrically, it means nothing else than that the two-orbital reduced quantum state evolves along a trajectory which at some point enters the convex set of separable states.We shall explore this phenomenon in more details in the next section.

D. Limiting regimes and sudden death of entanglement
In this section we investigate the asymptotic behavior of the orbital-orbital entanglement in various limiting regimes.Previously we observed that interesting phenomena such as long-distance entanglement and sudden death of entanglement occur when the particle filling η or hole filling 1 − η is small.Moreover, the P-and N-SSR entanglement in this region are practically identical.Therefore we can focus on the latter in our analysis.
For a fixed orbital separation d, when η ≪ d −1 , we get from ( 28) W (d, η) = η+O(η 3 ).The N-SSR entanglement between two orbitals separated by d becomes i.e., the N-SSR orbital-orbital entanglement is quadratic in η for small η.In log-log scale this relation becomes linear For the sake of completeness, we also verify that the P-SSR entanglement indeed reduces to the N-SSR one, i.e.

log(E(ρ
Eqs. ( 34) and ( 35) explain the asymptotic linear behavior in Figure 2 (right panel) of both the P-SSR (red) and N-SSR (blue) entanglement on the log-log scale, with slope 2 valid for any separation d.This relation is more robust for small separation d, as for large separation the orbital-orbital entanglement quickly deviates from the linear asymptote and plunges into sudden death as η increases.On the other hand, for small η and d < d min (η), the orbital-orbital entanglement is only a function of η, but not of d.All orbital pairs below the critical separation are equally entangled with one another.Due to the particle-hole symmetry, the orbital-orbital entanglement at small hole fraction 1 − η ≪ d −1 is also quadratic in 1−η.Above the critical separation d min (η), all orbital pairs are disentangled.Based on Figure 3, d min (η) seems to diverge as η approaches 0. To understand the asymptotic behavior of d min (η) at small η, we analyze the entanglement criteria (17).Recall from the previous section that the two orbitals are in a separable state if A and B in (31) satisfy the relation A ≥ 2B, which in explicit terms reads Although ( 36) cannot be solved analytically, in the limit of large d using the expansion W (d, η) ∼ (πd) −1 we can derive an analytic estimation for the minimal disentangling separation d min (η) In log-log scale, this relation becomes linear In Figure 3, we plot the exact minimal disentangling distance obtained by numerically locating the sudden death of entanglement using (36), together with the analytic estimation (37).We can see the agreement is excellent, even up to the order 10 −1 in η.
To summarize Section III, building on Ref. [40], we derived an analytic expression for the entanglement between any two spatial orbitals in the tight-binding model (21), as a function of both inter-orbital separation d and particle filling fraction η.In particular, Eq. ( 15) revealed the existence of long-distance entanglement when the filling fraction is low or close to 1.Moreover, our asymptotic analysis for small η revealed that this entanglement between orbitals is almost independent of the orbital separation d, and quadratic in the particle filling η for small η (and quadratic in hole filling 1 − η for small 1 − η).For a fixed filling η, the entangling range is still finite, as we also observed the sudden death of entanglement when the two orbitals are separated beyond a critical distance.Utilizing our analytic expression, the leading order of the minimal disentangling distance d min is extracted to be diverging as η −1 when η is small, and (1−η) −1 when η is close to 1.

IV. HYDROGEN RING
In the previous section, we thoroughly explored the orbital-orbital entanglement of the tight-binding model (21) for free electrons, using the analytic tools provided in Ref. [40].This rises the question of how interactions between electrons influence the orbital-orbital entanglement.Therefore we study in this section the hydrogen ring H 16 .In order to motivate more the choice of that system, we recall that the hydrogen ring (or open chain) is a common system for benchmarking ground state methods [65][66][67][68][69].This is due to the fact that it contains a significant amount of entanglement and thus aims at simulating electron correlation effects in realistic materials.Yet, the strong correlation therein presents a significant numerical challenge.Moreover, this and related ring models can be experimentally realized in the context of ultracold atom experiments [70][71][72][73][74].
The hydrogen ring shares many similarities with the tight-binding chain: a localized orbital basis can be established, and there is a local hopping generated by the overlapping atomic orbitals.As a molecular system, however, it additionally contains a periodic nuclear potential, and the Coulomb interaction between the electrons which is the source of computational complexity.It is thus of considerable interest to ask to what extend the entanglement effects in the hydrogen ring can be explained by the insights we gained from the tight-binding model (21) of non-interacting electrons in the previous section.Although the orbital-orbital correlation structure of similar systems has been analyzed [75], the orbital-orbital entanglement structure is yet to be determined.

A. Model and Hamiltonian
We consider a finite periodic hydrogen ring, defined by the uniform nearest internuclear distance R, with the minimal STO-3G basis set.This basis contains one 1s orbital at each atomic center and it is the minimal setting for describing the dissociation of hydrogen molecules.
The electronic Hamiltonian contains three terms: the nuclear potential, the kinetic energy, and the Coulomb interaction.The corresponding Hamiltonian contains oneand two-electron terms ĥ and V , respectively, and takes the form where Here, ϕ i 's are orthonormal molecular orbitals which form a complete basis of the one-particle Hilbert space spanned by the non-orthogonal atomic orbitals.M is the total number of nuclear centers with charges Z m .In the case of H 16 , Z m = 1 for all m = 1, 2, . . ., 16 and X m are defined as so that the 16 hydrogen atoms are arranged periodically on a ring, with the straight-line distance between two nearest neighbor atoms equal to R. From here on, we shall use atomic units (denoted by a.u.) for both distance and energy, which correspond to Bohr and Hartree, respectively.The matrix elements of the one-and two-electron Hamiltonians are explicitly orbital basis dependent.For the purpose of QIP, we choose a set of localized orthonormal orbitals, obtained by symmetrically orthogonalizing [76] the local 1s atomic orbitals.In Figure 4, we plot these localized orbitals for different nearest neighbor distances R. Locally they resemble the atomic orbitals, but they necessarily have finite contributions from atomic orbitals on other nuclear centers to ensure orthogonality.At R = 1, each localized orbital at the corresponding atomic center still contains significant weight of the atomic orbitals on the nearest neighbors, while for R ≥ 3, the localized orbitals are virtually single-centered.Do these localized orbitals exhibit a local interaction structure similar to the tight-binding chain?To answer this question, we plot the hopping terms h ij for i ̸ = j in R = 1 (a.u.)

B. Numerical methods
The calculation of the ground state of the hydrogen ring is done in three steps: (1) using the atomic wavefunction overlap matrix provided by Molpro [77][78][79], we symmetrically orthogonalize [76] the atomic orbitals into the localized orbitals, for nearest neighbor separation R = 1, 2, 3, 5 (a.u.).( 2) With a preceding Hartree-Fock calculation, we then transform the Hartree-Fock canonical orbitals into the localized orbitals as a post Hartree-Fock step, and compute the matrix elements of the oneand two-electron Hamiltonians followed by outputting them in FCIDUMP [80] file format, all within the framework of Molpro.(3) With the FCIDUMP file, quantum chemistry density matrix renormalization group (QC-DMRG) calculations are performed to solve for the ground states, using the SyTen [81,82] package, originally created by Claudius Hubig.The calculations are carried out for various electron numbers, initialized with random configuration states.The DMRG convergence is particularly challenging, because of both the strong intrinsic correlation in the system and the additional entanglement due to the localized orbital basis representation.For the data we use, the maximal bond dimension is set at m = 1000.The final energies are chemically accurate, verified by the energy improvement below the order 10 −3 (a.u.) in the DMRG calculations performed with bond dimension 2000.

C. Plots and results
In this section we present our main findings on the entanglement structure of the hydrogen ring H 16 .We tabulate in Figure 6 the P-SSR entanglement between two localized orbitals for various nearest neighbor separations R and different numbers of electrons N in the ground state.For readability we exploit the discrete rotational symmetry of the system and only plot the entanglement associated to the localized orbital labeled 1 (represented as filled red circle).Entanglement involving other localized orbitals can be understood by relabeling them.Exact values of the P-and N-SSR entanglement can be found in Table I in Appendix A. With these data, we observe three consequential features of the orbital-orbital entanglement between the localized orbitals in H 16 .
First, we observe that at half-filling (N = 16), the orbital-orbital entanglement is in complete agreement with the tight-binding model (21), namely only neighboring localized orbitals are entangled.In fact, this agreement extends beyond half-filling, as there is a region of filling fraction around 1/2 where all orbital pairs separated by more than one lattice constant are disentangled, exactly as in Figure 2.This is particularly surprising for systems with small nearest neighbor distance R, since there the hopping strength beyond nearest neighbors is rather considerable.
Second, at low or high filling, the hydrogen ring also displays long-distance entanglement similar to the tightbinding model.This phenomenon here differs from the latter, however, by preferring to entangle localized orbitals that are farther apart.Recall that in the tightbinding model, the orbital-orbital entanglement is uniform when long-distance entanglement is realized.This departure from the non-interacting scenario is best exemplified by the N = 2 column in Figure 6.In the ground state, due to the strong Coulomb repulsion, the two electrons would like to situate themselves as far apart as possible, forming singlet pairs at distant localized orbitals.Singlet states of such type, usually observed in a dissociated molecule, contain significant amount of entanglement.Singlet states realized on orbitals located closer to each other are of course possible, but they contribute to the total ground state with smaller weights.These two observations put together qualitatively explain the system's preference for long-distance entanglement at low particle filling.Furthermore, we notice that despite the lack of particle-hole symmetry in the system, the orbital entanglement in the N = 30 column in Figure 6 displays a similar behavior as in N = 2.This can be qualitatively understood with the observation that at very low (or very high) filling, the particle-particle (or hole-hole) interaction becomes negligible compared to the one-particle terms.In the extreme case where the system hosts only one particle or one hole, the corresponding Hamiltonians are related to each other by a minus sign.Consequently, the two ground states are identical up to a global phase, and therefore share the same orbital entanglement pattern.For orbitals separated by large distances, this near particle-hole symmetry is even more apparent.
Third, the long-distance entanglement in the hydrogen ring is far more potent than that in the tight-binding model.In Eq. ( 33), we see that the long-distance entanglement in the tight-binding chain is quadratic in η for small η, which is obviously not the case here.This is best demonstrated by the entanglement between orbitals 1 and 9, which form the farthest separated pair, at N = 2.The entanglement between them is of the order 10 −2 , which is the same order as the nearest neighbor entanglement.In the tight-binding model, the available entanglement at this filling fraction would yield at the order 10 −3 .This yet again validates that interacting electronic systems contain significant entanglement, and the entanglement is shown to be naturally distributed across distant regions, making it ideal for realizing quantum information protocols.
To summarize the present section, we took the interaction between electrons into account by studying the hydrogen ring, and analyzed the orbital-orbital entanglement therein.We discovered that the entanglement between the localized orbitals on the hydrogen ring can to a large degree be rationalized using the much simpler tight-binding model: close to half-filling only nearest neighbor orbitals can be entangled, but away from half-filling, long-distance entanglement appears.Yet, in contrast to the vanishingly weak long-distance entanglement in the tight-binding model, the entanglement between two localized orbitals sitting on far-separated nuclei is much more potent, at an order comparable to the nearest neighbor entanglement.This surprising feature, conjectured to be a result of Coulomb interaction, suggests that chain-like molecules can indeed be considered as reservoirs of entanglement.

V. SUMMARY AND CONCLUSIONS
It was the goal of our work to guide the development of quantum technologies in the era of the ongoing second quantum revolution.For this, we identified three key aspects which are beneficial if not indispensable for the successful utilization of orbital entanglement in atoms and molecules for quantum information processing (QIP): (i) orbital localization, (ii) operationally meaningful entanglement quantification, and (iii) computable entanglement measures.To elucidate those three aspects, we first highlighted the necessity of superselection rules (SSR) with an entanglement swapping protocol.This confirmed that only the SSR-compatible entanglement in fermionic systems can be extracted from individual molecules and used for QIP tasks.We then forwent the commonly used block entropy and calculated by analytical means the physical entanglement between spatial orbitals for the exactly solvable tight-binding model.In that sense, we also overcame the restriction of previous entanglement analyses to unitarily invariant settings, i.e. bipartitions of the system into macroscopically large and thus impractical subsystems.Finally, a numerically exact study of the hydrogen ring confirmed that the same entanglement analysis can be realized in continuous molecular systems, provided a scheme for localizing the molecular orbitals is implemented.
In the tight-binding model, we found a rich entanglement structure despite its simplicity.Around half-filling, the entanglement is strictly local, i.e., orbitals are entangled only if they are located at neighboring atomic centers.At very low and very high particle filling, however, the system displays long-distance entanglement between orbitals, even though the hopping is strictly between nearest neighbors.To investigate this phenomenon, we extracted the asymptotic behavior of this long-distance entanglement, and found it to be quadratic in the filling fraction η at small values of η.This low entanglement is shared uniformly within the entangling distance.The size of the entangling region, characterized by the minimal disentangling distance d min , diverges as the inverse of η.
For the interacting hydrogen ring H 16 we solved for the ground state in the localized orbital basis.This process was by itself a numerical challenge, for both the intrinsic complexity of the interacting ground state, and the extrinsic complexity from the localized orbital representation.After obtaining the ground state, we found the entanglement pattern to be similar to the tight-binding model.Only nearest neighbors are entangled near halffilling and long-distance entanglement appears as the filling approaches 0 or 1.What distinguishes the hydrogen ring from the tight-binding model is that in the former the long-distance entanglement is much more potent.In the hydrogen ring, the ground state actually prefers to place more orbital-orbital entanglement in far separated orbital pairs, due to the Coulomb repulsion between the electrons.This effect becomes evident as the molecule starts to dissociate and the entanglement between two localized orbitals sitting at the opposite sides of the ring grows stronger, making the hydrogen ring model a more suited one for preparing and storing long-distance entanglement.
Our results confirm that molecular systems offer indeed excellent prospects for providing valuable entanglement resources required in QIP tasks.At the same time, they invite a systematic exploration by both communities of quantum science and technology, and of quantum chemistry, under the framework put forward in this work.Open challenges include tailoring an orbital localization scheme that maximizes the orbital entanglement and in that sense maximizes the resourcefulness of molecular systems [83].Another crucial point is to take into account the practical restrictions on operations and measurements on localized orbitals [31], resulting in more inaccessible coherent superposition terms in the two-orbital reduced state, similar to the effect of the parity or particle number SSRs.Such technical limitations enforce an effective SSR.Based on this, the actual entanglement resource available in atoms and molecules will need to be further reassessed, following the same treatment as in this work.Finally, we anticipate the challenge of extending the scope of our entanglement analysis to the case of multiple copies of the molecular system.There, the local operations do not need to respect the SSR on the individual orbitals anymore but just on the collection of all local orbitals, thus offering more entanglement resource.Yet, both quantifying the entanglement in the multi-copy quantum state, and experimentally manipulating many orbitals at the same time, are immensely difficult tasks to accomplish.

FIG. 2 .
FIG.2.Entanglement between Wannier orbitals/lattice sites separated by inter-orbital distance 1 (solid), 2 (dashed), 10 (dashed-dotted), and 100 (dotted), under P-SSR (red) and N-SSR (blue) in the ground state of the tight-binding Hamiltonian (21) as a function of the particle/hole density η.In the left panel, the dashed-dotted and dotted curves look essentially flat.They can be seen more clearly in the log-log plot in the right panel.

FIG. 4 .
FIG. 4. Values of the localized orbitals' wavefunctions in H16 in the STO-3G basis, for various nearest neighbor distances R (a.u.).Nuclear centers are represented as black dots.The horizontal range is 16 a.u., and the vertical range 6 a.u..

Figure 5 .
Figure 5. Recall that the hydrogen atoms are arranged on a ring, and therefore orbital 1 and 16 are in fact nearest neighbors, despite appearing on the opposite ends of the axes in Figure5.At small separation, e.g.R = 1 (a.u.), hopping between distant localized orbitals is possible, but the hopping strength decays exponentially as internuclear distance increases.As we stretch the hydrogen ring, the hopping strength decays faster with the orbital separation.Finally at R = 5 (a.u.), the hopping is strictly between nearest neighbors, just as in the tightbinding model(21).