Entanglement phase transition due to reciprocity breaking without measurement or post-selection

Despite its fully unitary dynamics, the bosonic Kitaev chain (BKC) displays key hallmarks of non-Hermitian physics including non-reciprocal transport and the non-Hermitian skin effect. Here we demonstrate another remarkable phenomena: the existence of an entanglement phase transition (EPT) in a variant of the BKC that occurs as a function of a Hamiltonian parameter g, and which coincides with a transition from a reciprocal to a non-reciprocal phase. As g is reduced below a critical value, the post-quench entanglement entropy of a subsystem of size l goes from a volume-law phase where it scales as l to a super-volume law phase where it scales like lN with N the total system size. This EPT occurs for a system undergoing purely unitary evolution and does not involve measurements, post-selection, disorder or dissipation. We derive analytically the entanglement entropy out of and at the critical point for the $l=1$ and $l/N \ll 1$ case.

While MiPT is typically studied in systems where entanglement must be averaged over different random trajectories corresponding to distinct measurement outcomes, recent work has shown that entanglement phase transitions (EPT) can also occur without any stochasticity, in systems evolving under a non-Hermitian Hamiltonian [26][27][28].There is a direct connection to MiPT, as non-Hermitian dynamics is naturally interpreted as arising from measurement dynamics where one post-selects on a specific set of null measurement outcomes.Of particular interest here are studies of EPT in non-Hermitian models exhibiting non-reciprocity (e.g.directional systems where hopping to the right is much stronger than to the left).Kawabata et al. [27] studied an EPT in such a system (two coupled fermionic Hatano-Nelson [29,30] chains) from volume law to area law entanglement scaling.They found that the transition coincided with the transition between a reciprocal and a non-reciprocal phase.The latter could be directly witnessed by the non-Hermitian skin effect (NHSE), a phenomenon occurring in the non-reciprocal phase where all modes localize under open boundary conditions [29][30][31][32][33][34][35][36][37][38].Despite this strik-ing correspondence, one could still view the EPT transition here as being measurement driven, as the strength of non-reciprocity is directly tied to the strength of a post-selected measurement.
Taking inspiration from these previous studies, in this work we ask whether an EPT can occur without any stochasticity and without any need for measurements (post-selected or not).Similar to Ref. [27], we consider a translationally-invariant model that exhibits a transition between reciprocal and non-reciprocal phases (with the latter exhibiting the NHSE).In contrast to that work, our model has a fully Hermitian Hamiltonian and unitary evolution, and there is no need for any kind of measurement, post-selection or dissipation.Our setting is Hermitian, quadratic many-body bosonic Hamiltonians that do not conserve particle number.Such models can exhibit non-reciprocity: despite being fully unitary, the dynamics can nonetheless exhibit directionality at the level of the equations of motion for quadratures [34,[39][40][41].Given that such models can exhibit non-reciprocity transitions, can they also exhibit entanglement transitions despite the lack of any connection to measurements?
We find that the answer to this question is, surprisingly, yes.Our main result is to present the first instance of an EPT in a non-disordered bosonic system under purely unitary dynamics.We stress that this EPT requires no post-selection whatsoever.Furthermore, as this model contains neither randomness nor measurements, one can unambiguously attribute this EPT to reciprocity breaking.Our model of interest is a variant of the the bosonic Kitaev chain (BKC) model introduced in Ref. [34].By tuning a Hamiltonian parameter g, the BKC undergoes a phase transition from a phase where the dynamics of the q and p quadratures are nonreciprocal to a phase where they are reciprocal.At long times after a quench, the reciprocal phase is characterized by a volume law for the entanglement entropy (EE) of a subsystem of size l, i.e. it scales linearly with l.On the other hand, the non-reciprocal phase has even stronger  1)).b.Schematic phase diagram depicting super-volume law and volume law phases that arise as a function of g/∆ (with g, ∆ < w).For g < ∆, the BKC exhibits the NHSE, and has nonreciprocal dynamics leading to super-volume law entanglement scaling (see main text), whereas for g > ∆, the BKC is reciprocal, and exhibits the usual volume law of free systems.c.Log-log plot of long-time-averaged EE for a subsystem formed by the leftmost N/4 sites, divided by system size.Circles are values from numerical simulation, dashed lines are simply a guide to the eye.Since the y axis is divided by N , a slope of 0 (greater than 0) indicates volume law (super-volume law).A clear transition is seen as g is increased above ∆.We fix w = 1, ∆ = 0.25, and the values of g for the lines from blue to purple are 0, 0.2, 0.24, 0.245, 0.25, 0.255, 0.26.d.Scaling collapse of the long-time averaged EE for subsystem of size N/4, with ν = 0.5 in Eq. ( 9).entanglement growth.It is characterized by what we call a super-volume law: the EE of a subsystem of size l scales as lN where N is the total system size.Hence, if take a symmetric bipartition, the EE grows as N 2 .These behaviours are in stark contrast with [27], which also tied entanglement and reciprocity transitions, but found that non-reciprocity is detrimental to entanglement generation, leading to area law behaviour (see Sec. III).We also note that our results (for an unmeasured system) are distinct from the behaviour of explicitly monitored quadratic bosonic systems, which do not exhibit an EPT [42,43].In addition to being of fundamental interest, our setup is also attractive for experiments.Not requiring any measurements nor post-selection greatly simplifies implementation, while the the Hermitian bosonic pairing terms we require can be implemented in a variety of different platforms (they correspond to parametric drives or parametric down-conversion).
The remainder of this paper is organized as follows.In Sec.II, we recall the basic phenomenology of the BKC in both non-reciprocal and reciprocal phases.In Sec.III, we present a numerical demonstration of the EPT.In Sec.IV, we study analytically in depth the special minimal bipartition case l = 1 and show that it captures already the essential features of the EPT.In Sec.V, we extend the l = 1 results to the case l/N ≪ 1 by relying on a local thermalization hypothesis towards a generalized Gibbs ensemble (GGE).Finally, in Sec.VI, we conclude and discuss future directions.

II. MODEL
The BKC describes bosonic modes on a 1D lattice that are coupled by hopping and pairing terms on each nearest-neighbor bond.The Hamiltonian is Ĥ = 1 2 where N is the total number of sites, âj are bosonic operators, [â i , â † j ] = δ ij , and g, w, ∆ are real parameters of the model (see Fig. 1).We will call r := (q 1 , p1 , ..., qN , pN ) T the vector of quadrature operators, qj : Since Ĥ is quadratic, Gaussian states remain Gaussian under time evolution, and are fully specified by their 1point function ⟨r⟩ and their 2N × 2N covariance matrix σ ij = ⟨{r i − ⟨r i ⟩, rj − ⟨r j ⟩}⟩.In the remaining, we will fix the initial state to be the vacuum so that ⟨r⟩ = 0 at all times.The equations of motion (EOMs) for σ close on themselves and are given by where h is the bosonic Bogoliubov-de Gennes (BdG) 2N × 2N matrix defined through Ĥ = rT hr and Ω is the symplectic matrix Ω := N j=1 0 1 −1 0 .Note that the dynamics is completely linear in σ.
The qualitative properties of the BKC are best understood by inspecting the Heisenberg EOMs of r: (3) For g = 0 (i.e.purely imaginary hopping), these EOMs would describe independent, non-reciprocal propagation of the q and p quadratures, with each having opposite directionality.This mimicks the physics of two independent Hatano-Nelson chains [29,30].
We focus throughout in this work on the case w > ∆ [44].In this regime, the model is always dynamically stable for open boundary conditions (OBC), while for periodic boundary conditions (PBC) the system has a transition from stable to unstable as one goes from ∆ < g to g < ∆.The stability of this bosonic system for g = 0 can be conveniently understood in terms of amplification [34]: starting from a wavepacket in the q quadrature localized on one edge, the wavepacket will be amplified (damped) while propagating to the right (left), and vice-versa for the p's.For PBC, this amplification is unbounded, thus making the system unstable.In contrast, for OBC, the amplification terminates at the edges, so that all average moments of the q's are localized to the right and the p's to the left.Adding a coupling g mixes the quadratures together; as they have opposite directionality at g = 0, this mixing diminishes non-reciprocity.As g increases, this mixing eventually prevents chiral amplification altogether, leading to a sharp transition from a non-reciprocal phase for g < ∆ to a reciprocal phase for g > ∆.
Another manifestation of this transition can be seen from the spectrum of the dynamical matrix iΩh which exhibits the Non-Hermitian Skin Effect (NHSE).This means that in the non-reciprocal phase, g < ∆, for PBC the spectrum of iΩh winds around 0 in the complex plane, (i.e. the system is unstable), whereas the spectrum under OBC collapses onto the real line (i.e., the system is stable).This is accompanied by the localization of all the OBC eigenmodes to the edges [31,45].Conversely, in the reciprocal phase g > ∆, the NHSE is absent -the spectrum is always real, regardless of boundary conditions.We note that the spectral properties of iΩh are in complete analogy with the spectral properties of the non-Hermitian BdG Hamiltonian presented in [27] where an EPT was studied in coupled fermionic Hatano-Nelson chains.However, despite this similarity, we will show that the phenomenology in our Hermitian bosonic model is dramatically different.Given our interest in phenomena induced by the NHSE, we will only consider OBC in the what follows.
The parameters are given by r 0 = 1 2 tanh −1 (g/∆), r = for the non-reciprocal case g < ∆, and r 0 = 1 2 tanh −1 (∆/g), ϕ = arctan (w/ g 2 − ∆ 2 ) for the reciprocal phase g > ∆.After performing this transformation, the system is mapped in both cases to a simple tight-binding Hamiltonian The eigenmodes of Ĥ are standing waves, which we denote by bn such that Ĥ = n ε n b † n bn with ε n := −2 w 2 + g 2 − ∆ 2 cos πn N +1 .A few important points are in order.First, note that as one increases g across the non-reciprocity transition (i.e. from below ∆ to above ∆), there is no signature of a transition in the eigenvalue spectrum: the bandwidth of our system always increases monotonically with g.Hence, the non-reciprocity to reciprocity transition we focus on cannot be simply diagnosed by looking at the spectrum of the OBC system.Second, we stress that for both g > ∆ and g < ∆, our system has propagating quasiparticles.In the non-reciprocal phase g < ∆, it is useful to think of the position-dependent squeezing in Eq. ( 4) in terms of a localization length ξ = a/r, with a the lattice spacing.However, the emergence of this effective localization length does not impede quasiparticle propagation (i.e. the group velocity remains finite).Finally, we note that to compute the EE of a subsystem of size A, one can work either with the âj or the dj operators, as they are related to one another by purely local transformations.

III. ENTANGLEMENT PHASE TRANSITION
We now turn to the study of the entanglement scaling across the different phases.In [27], it was argued that non-reciprocity is detrimental for entanglement generation, as the quasi-particle pairs responsible for entanglement growth [46] propagate in the same direction, thus preventing the generation of long-range correlations, and precluding any volume law scaling of EE in the non-reciprocal phase.Here we show that the BKC, while presenting the main features of a non-Hermitian, non-reciprocal system (e.g. the NHSE and non-reciprocal transport) deviates dramatically from this expectation.While the reciprocal phase indeed presents a volume law as expected, we will show that the non-reciprocal phase on the contrary fulfills a super-volume law as defined in the introduction.
We consider the following setup.We fix the initial state to be the physical vacuum (i.e.∀j, âj |ψ⟩ = 0) and allow it to evolve under the BKC Hamiltonian for a long time.Since the system is dynamically stable, the time averaged EE of a subsystem will eventually converge to some steady state value.For Gaussian states, the covariance matrix σ fully determines the EE of any subsystem σ [47][48][49].The EE of a subsystem A of size l is obtained from the relation where and ν n are the positive eigenvalues of iΩσ| A where | A means that we truncate the support of the matrix to A.
Our quantity of interest will be the long-time averaged quantitiy S A , where we use the overline to denote timeaveraged quantities.
Before examining the EE, we look at the time-averaged position-dependent particle density across the chain.We observe that, as expected, the non-reciprocal phase has particles exponentially localized to the edges (Fig. 2a), while the reciprocal phase does not.Thus, one may expect that, in the reciprocal phase, taking a cut of the system from the left and increasing its size will not lead to a significant increase of the particle number and consequently, no significant increase in the EE.However, this is not the case.To illustrate this fact, we plot on Fig. 2b the so-called Page curve [50], i.e the EE of a subsystem of size l as a function of l while keeping the total system size N fixed.One finds that the reciprocal and non-reciprocal phases yield almost identically-shaped curves.We conclude, perhaps surprisingly, that the BKC does not show any area law phase induced by non-reciprocity.
To more fully understand entanglement properties, we can also consider the scaling of subsystem EE in a slightly different manner.Instead of fixing total system size N and varying subsystem size, we can fix the subsystem size l to be a fraction of total system size, l = z × N , and then vary the total system size while keeping z fixed.For concreteness, we take z = 1/4 in what follows.The results for this scenario are presented in Fig. 3a for w = 1, ∆ = 0.25, and a range of g close to ∆. Usually, one does not have to worry about the total system size as long as it is large enough; however, here, we find drastically different phenomena.When scaling the total system size, we immediately observe the emergence of two phases -a "super-volume" law phase where the EE scales as N 2 corresponding to the non-reciprocal phase g < ∆ and a volume law phase where the EE scales as N corresponding to the reciprocal phase g > ∆.The two are separated by a logarithmic scaling N log N when g = ∆.
From these considerations and following the fitting procedure in e.g.[2], we attempt the following finite size scaling for the EE, ) and find that, fixing the critical exponent ν = 0.5, it yields a good quality collapse of the EE for different system sizes onto the same curve, see Fig. 3b.
We thus have established a key result of our work: despite the lack of measurements, postselection or disorder, our BKC model exhibits a clear entanglement phase transition as a function of g, one that coincides with the transition from a reciprocal to non-reciprocal phase.

IV. ANALYTIC PROOF OF ENTANGLEMENT PHASE TRANSITION FOR A MINIMAL BIPARTITION
Computing the post-quench EE analytically, even for free systems, is in general a formidable task [51][52][53].In this section, we provide analytical insight by studying the EE for the simpler minimal bipartition case, i.e. the case where the subsystem A is composed of a single site.We will see that a transition from a phase where the EE is O(N ) to a phase where the EE is O(1) already occurs for this case when increasing g above ∆.16) and (22).We observe good agreement between the two for all but the smallest chain size.The upper dashed line indicates the N r scaling in the non-reciprocal phase, and the lower dashed line indicates the single-site EE at the critical point g = ∆ which separates the two phases.The two distinct scalings, with N in the non-reciprocal phase, and saturating in the reciprocal phase, are readily apparent.b.Scaling collapse for single site EE with ν = 0.5 in Eq. ( 24).The black dashed line is the expansion near the critical point obtained from the analytical form of the entropy (see Eqs. (18,23), For a single site j, the instantaneous symplectic eigenvalue ν t at time t is given by where we recall that d refers to the tight-binding basis defined in Eq. ( 4).The associated EE is simply S 1 = s(ν t ).
For now, we have not specified a particular site.Given the strong spatial non-uniformity in the density in the non-reciprocal phase (see Fig. 2b), one would naturally expect that S 1 would depend strongly on the choice of site, with large values at the boundary, and small values in the middle of the chain.Surprisingly, this is not the case.Eq. ( 4) already highlights how this can be.It shows succinctly that the entanglement is not simply determined by average photon number alone: one also needs to understand how many of these photons are associated with purely local squeezing correlations, and separate out this contribution.As we will see, despite the average density being highly inhomogeneous, S 1 is largely independent of the position j of the chosen site.
Non-reciprocal phase In the non-reciprocal phase, we expect ν t to be exponentially large with N , so that S 1 grows as N .When ν t is large, the EE takes a simple form: to leading order in ν t , S 1 ≈ 1 2 ln ν 2 t .We will compute the time-average of the entropy, S 1 .In general ln ν 2 t ̸ = ln ν 2 t but it would be desirable to use the latter expression as it is much easier to compute.To quantify the error resulting from this approximation, we Taylor expand ln ν 2 t with respect to ν 2 t around ν 2 t and take the average to get : We show in App.F that the second term is of order O(1) because of cancellation of exponential-in-N contributions to the numerator and denominator.Hence, Finally, to compute ν 2 t , we will make use of the fact that a property that we prove in App.F. Remarkably, the neglect of fluctuations, implicit in this approximation, only holds for this difference, and not for each term individually.More explicitly, one finds Our interpretation of this result is that while individual sites are subject to (exponentially) large temporal fluctuations in density, these fluctuations are almost entirely due to fluctuations in the amount of local pairing correlations.The contributions of these fluctuations (density, local pairing) thus cancel each other to leading order when calculating the symplectic eigenvalue and, consequently, the EE.Hence, to compute ν 2 t , one only needs the average covariance σ.Taken together, all these steps considerably simplify the task of computing the EE and allows us to have quantitative results.
To compute σ, we will consider the continuum limit, which we define as follows : let a be the lattice spacing, we consider the limit N → ∞, a → 0, while keeping fixed the dimensionful quantities ξ := a r (localization length), L := a(N + 1) (system size), x := aj and p := πn a(N +1) .We are free to fix the parameter x 0 = aj 0 .For the particular choice of x 0 = L 2 , the correlations have a compact expression (see App. C): ⟨ dx dx ⟩ = (−1) We see that the time-averaged local density and pairing correlations in the squeezed frame are spatially uniform.
At first glance, this could seem surprising, as in this frame, our initial condition (vacuum in the lab frame) is extremely non-uniform in space, due to the position dependent squeezing transformation in Eq. ( 4).However, the resulting uniformity of the time-averaged state can be understood by the dynamics being equivalent to a simple tight-binding chain.Indeed, for such a model, any spatial product state will lead to an average homogeneous profile in the continuum limit.This in turn means that the EE in the minimal bipartition protocol will yield the same result, independent of our choice of which site to single out.Heuristically, this explains the discrepancy between Fig. 2b and c: while the average density in the lab frame is exponentially localized towards the edges, this excess density can largely be attributed to local squeezing, which does not affect entanglement properties, in line with our interpretation of Eq. ( 13) (see App. E for further discussion about the separation between local squeezing and thermal occupation).Inserting (14,15) in (10) leads to Away from the critical point, the large L limit gives us that To obtain the scaling near the critical point, we consider the limit L/ξ ≪ 1 for which the localization length ξ is far greater that the system size.This (see App. D) leads to: Reciprocal phase In the reciprocal phase, the local correlations are given in the continuum limit by where we defined φ := 1 a (π − 2ϕ).We find that the average density in the squeezed frame is both timeindependent and spatially uniform.This is no surprise: in the reciprocal phase, the transformation to go to the squeezed frame (c.f.Eq. ( 4)) is uniform, hence our initial pre-quench state is also uniform.Such a density profile will not evolve under a tight-binding Hamiltonian.In contrast to this, the time-averaged local squeezing correlators above retain a position dependence in their phase.
Away from the transition, i.e. for φ finite, the local pairing correlations tend to 0 in the large L limit.Thus the EE is simply where s(x) is defined in Eq. (7).Close to the transition, the eigenvalue scales like N and we can apply the same set of approximations used in the reciprocal case.The eigenvalue corresponding to (19,20) is Open circles: S (j) th , the prediction for the single site entropy based using the average photon number alone, c.f. Eq. ( 25).This prediction deviates sharply from the true result.This highlights an important caveat: simply using average particle number as a proxy for entanglement can be extremely misleading.
Close to the critical point, the limit φL ≪ 1 leads to which is consistent with the limit Eq. ( 18) from the nonreciprocal phase.Numerical simulation and scaling collapse Numerical simulations of the EE for a single-site are plotted on Fig. 3a and b alongside analytical estimates.We observe that for large system size, the EE of the non-reciprocal phase always goes to the L/ξ scaling, whereas the EE for the reciprocal phase saturates.The expressions (18,23) suggest the following scaling collapse (24) with ν = 0.5, thus proving Eq. ( 9) for the minimal bipartition.Note that, not only the power laws are in agreement with the numerics but also the non-universal 1/15 prefactor of the second term of (18,23), see Fig. 3.
Before leaving this section, we wish to highlight a crucial fact: entanglement in our system cannot be simply predicted from the behaviour of average photon number.A naive argument would be that the average photon number on each site of the post-quench state determines its effective Hilbert space dimension D j , which would then (assuming thermalization) set its entropy.This line of reasoning would suggest that the entropy of a given site should correspond to the entropy of a single bosonic mode in a thermal state, i.e.

S(j)
1,th = s ⟨â † j âj ⟩ ∼ log ⟨â † j âj ⟩ ∼ log D j (25) where the first approximation holds for large particle numbers.If this reasoning were true, then in the nonreciprocal phase, the entropy versus position curve should be peaked at the edges, reflecting the skin-effect-induced localization of the particle number density.As shown explicitly in Fig. 4, this prediction is manifestly incorrect.The numerics here matche the analytic arguments presented above: the true single-site entropy is almost independent of position, and shows no signature of localization.As such, simply understanding how average particle number depends on system parameters does not immediately let one understand entanglement properties.The discrepancy between particle number and EE is further explored in App.E.

V. GENERALIZED GIBBS ENSEMBLE
In this last section, we show how the previous results can be extended to understand the entanglement entropy of small subsystems for size l satisfying l N ≪ 1.To do this, we make a local thermalization hypothesis that the subsystem will be described by a Generalized Gibbs ensemble (GGE) state [54,55].The GGE ansatz amounts to the assumption that expectation value of local observables can be extracted from where we recall that the bn modes refer to the standing wave basis, Z is a normalization factor, {β n , γ n , γ * n } are thermodynamic variables fixed by the initial condition and n is defined by the relation ε n = −ε n .Note the contribution of the pairing terms b † nb † n , bn bn that do not appear in Ĥ.For free bosonic system, this is simply the time average of ρ : ρGGE = ρ.
In a similar fashion, the GGE ansatz for entanglement is that the EE of a subsystem A can be extracted from ρGGE [56], i.e. that one has : where ν n GGE are the positive eigenvalues associated to σ GGE Ω.Given the breaking of translational invariance in the steady state, as exemplified by the strong inhomogeneity in the local occupation number as seen on Fig. 2 b, and the fact that the GGE is agnostic about the cut chosen for A, one may expect this approach to fail.On the other hand, the fact that the local correlations are spatially uniform in the tight-binding frame suggests that, for the purpose of computing the EE, the GGE might be enough.
In the non-reciprocal phase, fixing our gauge parameter to x 0 = 0, the positive eigenvalue of (D2) in the continuous limit is (see App. D for the derivation) We see that the localization length ξ is the natural scale separating the long and short wavelength physics.The suppression of ν p for pξ ≪ 1 is a direct consequence of low-momentum standing waves having small wavefunction amplitudes near the system boundary.
In the localized regime ξ → 0, the dominant contribution in the above expression is sinh 2 (L/ξ).Since the entanglement is proportional to the log, we have We see that, because of the exponential scaling of the ν p with L, the different modes' contribution to the entanglement becomes independent of p.
Remarkably, a similar statement is true close to the transition where ξ → ∞.In this regime, the momentum dependence in ( 28) cancels out and we are left with which is again independent of p.
Similarly, in the reciprocal phase we have where we defined φ := 1 a (π − 2ϕ).For φ finite, this quantity becomes p independent in the large L limit.Close to the critical point, φ → 0 keeping φL finite while L → ∞ gives which is also p independent.Thus, we see that in all regimes of interest, for the purpose of computing entanglement, the momentum dependence drops out.This in turns implies that the GGE and the minimal bipartition will match in all the limits mentioned above and thus, Interestingly, this means that for computing the EE, the local thermalization assumption yields accurate results, despite the fact the system is both strongly inhomogeneous and subject to exponentially large fluctuations.Our interpretation is that fluctuations in local quantities mainly comes from a variation in the squeezing strength which leaves the EE property unchanged.

VI. CONCLUSION
Our work demonstrates the existence of an entanglement phase transition in a non-disordered bosonic system undergoing purely unitary evolution.When varying the hopping parameter g below a critical value ∆ the system undergoes a transition from a reciprocal to a non-reciprocal phase, accompanied by a transition from a volume law to a super-volume law for the post-quench entanglement entropy of a subsystem.While our system shares many common features with non-Hermitian systems, it does not involve measurement or post-selection in any way.Our study suggests that the breaking of reciprocity can be associated with entanglement transitions even in settings where there is no competition between unitary dynamics and measurement-induced non-unitary evolution.
It is interesting to contrast our results with the related non-Hermitian fermionic model studied in [27], involving two coupled Hatano-Nelson chains.As discussed, that model exhibits identical spectral properties and NHSE as our system.Ref. [27] also found an entanglement transition coinciding with the breaking of reciprocity, but unlike us, found that entanglement generation was greatly suppressed in the non-reciprocal phase (yielding only area law behaviour).In contrast, our nonreciprocal phase exhibits marked directional transport, but no area law entanglement behaviour (and in fact has enhanced entanglement scaling).This suggests that the uni-directional quasi-particle picture proposed in [27] is not applicable to generic entanglement transitions associated with reciprocity breaking.
It is also interesting to note that in contrast to other studies of bosonic systems, we observe the existence of an EPT despite the absence of measurements [57,58] and non-linearities (see e.g.[59]).We also mention that for fermionic systems, entanglement transition were observed for free, unitary, disordered systems; these were directly tied to either Anderson or many-body localization/delocalization transitions, see e.g.[60][61][62].These disorder-driven EPT are also distinct from the phenomenon we describe, as (apart from boundaries) our system is fully translationally invariant.
While our focus in this work was on post-quench entanglement entropy, it is important to note that the reciprocity-breaking transition in our model can also be characterized with other quantities.This comprises the spectrally-heralded reciprocal-to-nonreciprocal transition already pointed out in [34].Another observable that shows clear signatures of the transition is the scaling of the total particle number with N in the post-quench state, a quantity which is linear in ρ.Such signatures of the transition differs markedly from the phenomenology of standard MiPT, where the transition can be a priori only be characterized using quantities non-linear in ρ.We stress that the phase transition in the non-Hermitian model of [27] could also be characterized using a single observable, the total current.Returning to our model, we stress that even though the reciprocal and non-reciprocal phases differ strongly in terms of their average density, this does not by itself let one infer the existence of an EPT.In general, particle number can be made arbitrarily large by means of local squeezing transformations, something that would have no impact on entanglement.The fact that average density and entanglement properties can be extremely different is demonstrated explicitly in Fig. 4 and App.E, where we observe that the entanglement entropy spatial structure is dramatically different from that of the average particle number.
The EPT demonstrated in this work is experimentally appealing for several reasons.First, since the model is a non-disordered closed system, post-selection is a complete non-issue.Second, all the studied dynamics are Gaussian, which for bosonic systems are generally considered much more experimentally tractable.Finally, as we showed in Sec.IV, the entire EPT can be characterized by a single-site covariance matrix.Hence, to detect and characterize the EPT experimentally, one only needs to characterize the correlations of a single site.
In this work, we have demonstrated and characterized an EPT associated with a transition from non-reciprocity to reciprocity in a particular model, namely, the BKC.Future work could investigate the more general relationship between non-reciprocity and entanglement-in particular, how many of the features of this EPT generalize to other models, and what one can say more generally about the entanglement properties of non-reciprocal models?Finally, we note that while entanglement is a quantum property, one could also investigate a classical version of this model and ask whether the non-reciprocal to reciprocal transition there is also heralded by a transition in correlation measures besides entanglement.
We are interested in the long-time average entanglement that results from the quench dynamics described in the main text, and hence the quantity of interest is which is the EE for some subsystem of our 1D BKC lattice for some fixed set of parameters.Without loss of generality, we will fix w = 1, ∆ = 0.25, and vary g, N , where ∆ and g will be written in units of w.We will estimate S t by numerically calculating S t for some discrete set of times and then taking the mean.Since we want the EE in the quasi-steady state, we only need to perform this calculation up to some finite large T for which the estimate of S t converges to some desired accuracy.Since the evolution of S t is deterministic, we can set the accuracy to any level we want.In general, one might not expect to be able to simulate arbitrarily large times accurately for non-reciprocal systems, due to the issue of numerical ill-conditioning.Fortunately, we can avoid this by performing the simulations in the squeezing frame defined by Eq. ((4)), something that is possible whenever g ̸ = ∆.For g = ∆, the squeezing frame is not well-defined.In that case, we simply performed the simulations in the lab frame, and found it to be stable for all chosen parameters in this work.
The following are additional important points about the numerical approach used to calculate the entanglement entropy.
1.The value of the EE at initial small times t are in general not representative of the quasi-steady state of interest.While they get averaged away at long times, including these points slow down the convergence of our calculation.We thus pick an ini- (cosh (r(j − j 0 )) cosh r 0 + i sinh (r(j − j 0 )) sinh r 0 ) âj where j 0 is an arbitrary "gauge factor" and tanh (2r 0 ) = (C7) These expressions simplify in the continuous limit defined as follows.Let a be the lattice spacing.We consider the limit N → ∞, a → 0, while keeping fixed the dimensionful quantities ξ := a r , L := a(N + 1), x := aj and p := πn a(N +1) .To simplify the expressions, we fix the gauge parameter j 0 = 0.This leads to a. Local correlations in the tight-binding basis In this part we derive the average local on-site correlations in a given spatial frame.As discussed in the main text, those are the quantities necessary to characterize entanglement in the minimal bipartition protocol.The spatial frame where the correlations appear in their simplest form is the tight-binding frame with operators { d † j , dj }.The correlations are related to the one in the diagonal basis by a simple OBC Fourier transform : We will make use of the identity Inserting this identity in the previous relations leads to where we recall that j 0 is an arbitrary "gauge factor".Once again, this expression simplifies in the continuous limit: Note that the j dependent term is no longer here in the continuous limit description.Finally, we can fix the gauge parameter x 0 = L/2 to simplify these expressions: For the local pair annihilation correlation, one obtains (sinh (2r(j − j 0 )) + sinh (2r(N + 1 − j − j 0 ))) .

(C17)
Once again, taking the continuous limit and choosing x 0 = L/2, one gets We thus see that, in the continuous limit defined above, the norms of both correlations are independent of x in the tight-binding frame.For the minimal bipartition, this means that the value of the entanglement entropy will be the same, up to finite size corrections, for all the sites.

Reciprocal phase
In the reciprocal phase, the eigenoperators are bn with tanh (2r 0 ) = g ∆ and ϕ = arctan (w/ g 2 − ∆ 2 ).The conserved correlations in this case are given by Performing the sum for ⟨ bn bn ⟩ in the continuous limit leads to where we defined p := π a − p and φ := 1 a (π − 2ϕ).a. Local correlations in the tight-binding basis As for the non-reciprocal case, the local on-site correlations take their simplest form in the tight-binding frame with operators { d † j , dj }.Performing the inverse Fourier transform leads in this case to: Defining φ := θ a , the last expression simplifies once again in the continuous limit: Appendix D: Computation of entanglement entropy In this appendix, we compute the EE of a subsystem of size l in the limit l/L ≪ 1 using the GGE.We begin by showing that, for our model, this is equivalent to the minimal bipartition approach, both in the reciprocal and non-reciprocal phase.

Equivalence between GGE ansatz and minimal bipartition
Recall that in the GGE approach, the stationary entanglement of a subsystem A of size l is simply assumed to be directly proportional to the total EE of the total system, with the proportionality coefficient fixed by l, In this appendix, we show that simply considering particle number is insufficient to fully characterize the spatial profile of entanglement.
Given only access to the density profile, what kind of EE profile might one expect such a system to have?One reasonable approach would be to think about particle number as a proxy for the size of the local Hilbert space, and in general we expect a larger local Hilbert space to indicate that the site is more entangled with the rest of the system.To make this concrete, suppose a site j has occupation ⟨â † j âj ⟩.One can try to associate an entropy with this density by assuming that when the system has thermalized, the density matrix of the site j will be approximately that of a thermal state, for which the entanglement entropy is given by S (j) th = (⟨â † j âj ⟩ + 1) ln ⟨â † j âj ⟩ + 1 − ⟨â † j âj ⟩ ln⟨â † j âj ⟩, (E1) and use this as an estimate for the actual EE S (j) 1 .
We compare these two quantities by plotting their timeaveraged values against each other in both the reciprocal and non-reciprocal phases (Fig. 6).In both cases, true entanglement profile is flat up to finite-sized effects, whereas the expected thermal entropy reflects localization in the non-reciprocal case and periodic spatial oscillations in the reciprocal case.In both cases, the thermal entropy significantly overestimates the true EE.Another way to appreciate the relationship between particle number and entanglement is to explicitly extract the local squeezing and temperature of the time-averaged one-site density matrix.To do so, note that any diagonal single-site covariance matrix σ can be decomposed into a rotations followed by a squeezing operation on a thermal state: where R is some orthogonal matrix and e 2β±2z are the eigenvalues of σ.The symplectic eigenvalue is entirely determined by β and the local squeezing parameter z does not affect entanglement properties at all.The quantities β, z can be easily obtained by diagonalizing σ.In Fig. 7, we plot the values β j , z j extracted from the time-averaged covariance matrix σ j at site j.These quantities display qualitatively similar features to the entropies of the previous plot, with temperature being analogous to EE and local squeezing displaying the same spatial distribution as the particle number.Notably, this allows us to make the more explicit statement that the spatial non-uniformity arising from non-reciprocity can be entirely characterized by local squeezing operations, which do not affect entanglement properties.
Appendix F: Two enabling claims in the nonreciprocal phase 1.The symplectic eigenvalue squared can be approximated using two point correlators In the main text, we assumed that in the nonreciprocal phase, (F1) To show this rigorously requires a long and tedious calculation, which we will outline here.We will show this in two limits: (1) deep in the nonreciprocal phase, where we fix r and take N to be large, and (2) near the critical point, where we fix N and take g → ∆ − .First, let us set up the problem and classify various conserved quantities.For simplicity, we will pick d to be the operator on the first site of the chain, in the frame where j 0 = (N + 1)/2 (picking any other site does not materially change the calculation, as we will see shortly).

FIG. 1 .
FIG.1.a. Schematic of a 1D BKC lattice, with Hermitian hopping (amplitude g + iw), and pairing / two-mode squeezing (amplitude i∆) on each bond (see Eq. (1)).b.Schematic phase diagram depicting super-volume law and volume law phases that arise as a function of g/∆ (with g, ∆ < w).For g < ∆, the BKC exhibits the NHSE, and has nonreciprocal dynamics leading to super-volume law entanglement scaling (see main text), whereas for g > ∆, the BKC is reciprocal, and exhibits the usual volume law of free systems.c.Log-log plot of long-time-averaged EE for a subsystem formed by the leftmost N/4 sites, divided by system size.Circles are values from numerical simulation, dashed lines are simply a guide to the eye.Since the y axis is divided by N , a slope of 0 (greater than 0) indicates volume law (super-volume law).A clear transition is seen as g is increased above ∆.We fix w = 1, ∆ = 0.25, and the values of g for the lines from blue to purple are 0, 0.2, 0.24, 0.245, 0.25, 0.255, 0.26.d.Scaling collapse of the long-time averaged EE for subsystem of size N/4, with ν = 0.5 in Eq. (9).

FIG. 2 .
FIG. 2. Comparison of the average density and Page curves for the reciprocal (with g = 0.3) and non-reciprocal (with g = 0.24) phases, for w = 1, ∆ = 0.25, N = 32.a. Time-averaged position-dependent particle number ⟨â † j âj⟩, in the reciprocal phase.b.Same as (a) but in the non-reciprocal phase.Note the difference in shape and increase in scale in going from panel a to b. (c) EE for a cut of size l, with the subsystem being the left l sites of the chain in both the non-reciprocal and reciprocal phases.The curves are normalized by their maximum value.Despite the dramatic differences between (a) and (b), the Page curves are hardly distinguishable.

FIG. 3 .
FIG.3.a. Log-log plot of long-time-averaged EE for the EE of the left-most site 1.The fixed parameters are w = 1, ∆ = 0.25.The values of g for the lines from blue to purple are 0, 0.2, 0.24, 0.245, 0.25, 0.255, 0.26.Scatter plots indicate values obtained from numerical simulation, whereas the solid lines are the analytic prediction from Eqs. (16) and(22).We observe good agreement between the two for all but the smallest chain size.The upper dashed line indicates the N r scaling in the non-reciprocal phase, and the lower dashed line indicates the single-site EE at the critical point g = ∆ which separates the two phases.The two distinct scalings, with N in the non-reciprocal phase, and saturating in the reciprocal phase, are readily apparent.b.Scaling collapse for single site EE with ν = 0.5 in Eq. (24).The black dashed line is the expansion near the critical point obtained from the analytical form of the entropy (see Eqs.(18,23),

FIG. 4 .
FIG. 4. Filled circles: one-site entanglement entropy S (j) 1 as a function of position j in the non-reciprocal phase, for parameters N = 32, w = 1, g = 0, ∆ = 0.25.Despite being in the skin-effect phase, the entropy is almost completely uniform.

FIG. 5 .
FIG. 5. (a) Full time evolution of the EE of the (N/4 : 3N/4) biparition for a representative time period for w = 1, ∆ = 0.25, g = 0. We observe that visually, the long-time EE appears to oscillate around some mean value, with oscillations small compared to the mean value.(b) Plot of the variance defined in Eq. (B1), calculated over the period of time required for the mean value of the EE to converge according to the criterion outlined in App. A. (c) Same as (a) but for the minimal bipartition (1 : N − 1).(d) Same as (b) but for the minimal bipartition (1 : N − 1).