Page-curve-like entanglement dynamics in open quantum systems

The entanglement entropy of a black hole, and that of its Hawking radiation, are expected to follow the so-called Page curve: After an increase in line with Hawking's calculation, it is expected to decrease back to zero once the black hole has fully evaporated, as demanded by unitarity. Recently, a simple system-plus-bath model has been proposed which shows a similar behaviour. Here, we make a general argument as to why such a Page-curve-like entanglement dynamics should be expected to hold generally for system-plus-bath models at small coupling and low temperatures, when the system is initialised in a pure state far from equilibrium. The interaction with the bath will then generate entanglement entropy, but it eventually has to decrease to the value prescribed by the corresponding mean-force Gibbs state. Under those conditions, it is close to the system ground state. We illustrate this on two paradigmatic open-quantum-system models, the exactly solvable harmonic quantum Brownian motion and the spin-boson model, which we study numerically. In the first example we find that the intermediate entropy of an initially localised impurity is higher for more localised initial states. In the second example, for an impurity initialised in the excited state, the Page time--when the entropy reaches its maximum--occurs when the excitation has half decayed.

Introduction.Entanglement has become an object of intense study in diverse fields from quantum information to quantum many-body and black-hole physics.A central result is that the von Neumann entropy of subregions in the ground state of many-body Hamiltonians with local interactions grows like the surface area of the subregion and not its volume [1].Computational matrix-product approaches [2] use such entanglement properties for efficient simulation of quantum systems, from quantum many-body systems [3] to non-Markovian open quantum systems [4].Another important result is that the entanglement entropy of ergodic systems out of equilibrium generically grows linearly until it saturates at a value given by the volume law for excited states [5][6][7].From this perspective, the behaviour of the entropy of black holes is unusual.
By considering quantum fields in curved spacetimes, Hawking found that black holes can evaporate by emitting thermal radiation, thus associating a temperature proportional to the inverse mass to black holes [8].This, together with an associated entropy proportional to the area of the event horizon [9], makes certain laws of black hole mechanics look like the laws of thermodynamics [10].By considering the situation in which matter in a pure state undergoes gravitational collapse to form a black hole and then evaporates according to Hawkings semiclassical calculation, one finds that a pure state evolves into a mixed state in a closed system, which stands in stark contrast to the unitarity of time evolution in quantum mechanics.This is the so-called black hole information paradox [11].If unitarity were to hold for quantum gravity, instead of a linear increase of entropy according to Hawkings calculation, the curve describing the entropy should bend down and decrease back to zero once the black hole has fully evaporated, so that the final state remains pure [12].The resulting curve is called the Page curve.Recently [13,14], tremendous progress has been made calculating this behaviour semiclassically, using so-called quantum extremal surfaces [15,16], first developed in the context of the AdS/CFT correspondence [17].
Recently [18], a relatively simple solvable system-plus-bath model [19,20] has been proposed, which shows Pagecurve-like entanglement dynamics.The entanglement entropy is examined as function of time and not as function of subsystem size, as often studied in quantum many-body systems [12,[21][22][23][24]. Similar observations have been made for the system-bath mutual information [25] and the entanglement negativity [26] of open quantum systems.Furthermore, in [18] a more general argument was made, that such an entanglement dynamics is generally expected when the resulting open dynamics forces the system into a low dimensional subspace of its Hilbert space.In the example given, the fermionic system empties out.
Here we make an argument that the condition suggested in [18] can be generically realised for open quantum systems weakly coupled to an environment at low temperature.When initialised in a pure far-from-equilibrium state, the ensuing non-equilibrium dynamics will be accompanied by a high entanglement entropy production.Under our conditions, however, one can expect, that the impurity will eventually settle down to a state close to its ground state, which carries a low entanglement entropy.The entanglement entropy as a function of time will therefore qualitatively look like the Page curve.We will corroborate and illustrate this on paradigmatic open-system models; namely, the exactly solvable harmonic quantum Brownian motion [27][28][29][30][31] and the spin-boson model [32][33][34][35][36][37], which we study numerically.
The main ingredient of our argument is the so-called mean-force Gibbs state [38][39][40], which amounts to the reduced state of the global system-bath canonical equilibrium state.It is conjectured [40] that open quantum systems in large enouth thermal baths generically approach this steady state.The physical picture behind this is the following.As the environment is much larger than the impurity, the state of the impurity is just a perturbation to the global equilibrium and the dynamics essentially consists of the return to equilibrium of the global system, in line with the eigenstate thermalisation hypothesis [41][42][43][44].Under our conditions, that is, weak coupling and low temperature, the mean-force Gibbs state is close to the local ground state and thus carries low von Neumann entropy.
Entanglement dynamics.We study the entanglement dynamics for an impurity, the system, embedded in an environment, the bath.The total Hamiltonian for such a problem takes the form commonly studied in the open quantum systems literature where the subscripts S, B and I stand for system, bath and interaction.Here, and in what follows, operators will be denoted by boldface symbols.Considering the situation where the probe is initialised in a pure state ρ ρ ρ 0 and the environment is in its ground state (which we assume to be non-degenerate), the global product state is also pure.Despite the global state remaining pure under the dynamics generated by the Hamiltonian in Eq. ( 1), the subsystems (impurity and environment) get entangled as they interact.We quantify this by the entanglement entropy, i.e., the von Neumann entropy S of the subsystems given by The entanglement entropy of the system and that of its complement (the bath) are the same.Generically, for large systems avoiding recurrence, is it expected that the entanglement entropy rises with time and reaches a plateau asymptotically [5][6][7].
On the other hand, when the bath is much larger than the system, it is generically expected that the open system reaches a steady state given by the so-called mean-force Gibbs state [38][39][40] where β is the inverse temperature of the bath.This is the expectation as long as no symmetries lead to conserved quantities constraining the system dynamics.Here, and in what follows, we work in natural units, i.e. h = k B = 1.This means that the system reaches a state which looks like the reduction of the system-bath equilibrium at the temperature of the bath.Our environment being in its ground state, this means that we expect a steady state for the impurity taking the form τ τ τ MF (T = 0) = tr B |Ω⟩ ⟨Ω|, where |Ω⟩ is the global ground state.At very weak coupling between system and bath, as commonly assumed in thermodynamics and, more generally, whenever the coupling is given by a surface term, while the system and bath energies are given by large volume terms, the mean force ground state can be approximated well by the system ground state |0⟩, i.e. ρ ρ ρ S (t → ∞) = tr B |Ω⟩ ⟨Ω| ≃ |0⟩ ⟨0|.The entanglement entropy, therefore, must be close to the one of the local ground state, which we assume to be nondegenerate.This means the entanglement entropy must be close to zero.Importantly, this has to hold for arbitrary initial states, in particular states far from equilibrium whose dissipative dynamics result in large entropy production.
Putting our pieces together we have to following picture: Impurity and environment start at a pure state with zero entanglement entropy.The entanglement entropy can increase substantially as they interact during the non-equilibrium dynamics.At long times the entanglement entropy needs to reach the value of the mean-force state, which is close to zero.In between it has to decay, therefore qualitatively following the Page curve.
The above is very similar to evaporation of AdS black holes.This is usually facilitated by coupling them, via a surface term, to an auxiliary CFT in the vacuum state [13].Those are essentially the conditions we have described here.
In the following we study the entropy dynamics quantitatively on concrete models.Despite working with small coupling, we do this without using weak coupling master equations [45][46][47] prevalent in the open quantum systems literature, as they are not applicable at such low temperatures [19].In particular we do not assume Markovianity.Instead we work with exact solutions.
Harmonic quantum Brownian motion.The first example is the exactly solvable harmonic quantum Brownian motion, i.e. a harmonic oscillator linearly coupled to a continuum of harmonic oscillators.This model can, for example, describe a Bose polaron in a condensate [48].The system Hamilton reads where we have set the mass to 1.The bath Hamiltonian is given by the most common model for the environment in the open quantum systems literature.Such an environment could represent the electromagnetic field or phonons in a crystal.Finally the coupling between system and bath is given by The effect of the bath on the system can be encoded via the spectral density J(ω) = π ∑ µ g 2 µ /(2m µ ω µ ) δ (ω − ω µ ).We choose the common Ohmic spectral density with Lorentz-Drude cutoff given by This spectral density is linear for frequencies much smaller than the cutoff Λ and decays for larger frequencies.Other choices of spectral densities would lead to qualitatively similar results.We choose the frequency ω R in Eq. ( 4) to cancel The entanglement entropy S of the oscillator impurity as a function of time starts at zero, as the impurity is initialised in a pure state, for a localised wave packet with width δ = 1/100 (solid blue) and δ = 1/1000 (dashed orange).The interaction entangles the impurity with the environment, which is initially at zero temperature.Their entanglement entropy peaks at an intermediate time.The more localised the initial state, the higher the intermediate peak of the entanglement entropy.Regardless of initial state, the mean-force Gibbs state is approached at long times, which, for weak coupling, is close to the local ground state.As in the Page curve, the entanglement entropy thus has to decrease; here, it saturates close to zero.In the inset the initial state is the local ground state, the maximum of the entanglement entropy (dot-dashed green) still reaches an intermediate maximum, but it is just slightly higher than the final value, as the dynamics is not very far from equilibrium, seen by the fidelity F with the local ground state staying close to one (dotted red).The parameters are ω 0 = 1, γ = 0.001 and Λ = 10.out the distortion due to the bath with regard to the oscillator frequency ω 0 by setting where the counter-term frequency is given by ∆ω 2 = 2 π ∞ 0 dω J(ω)/ω.Such a counter-term naturally arises in many systems of interest [28].For our spectral density it evaluates to ∆ω 2 = γΛ.
Initialising the bath in a thermal state, which is Gaussian, we exploit the fact that the Hamiltonian is quadratic in positions and momenta, thus generating a Gaussianity-preserving dynamics [49].That is, for Gaussian initial conditions with vanishing first moments, the covariances σ xx = ⟨x x x 2 ⟩, σ xp = 1 2 ⟨x x xp p p + p p px x x⟩ and σ pp = ⟨p p p 2 ⟩ fully characterise the state of the impurity.We collect the covariances in the matrix In the Supplemental Material [50] we give details on the exact dynamics of the model.
For Gaussian states we can express the von Neumann entropy, Eq. ( 2), as [51] The entropy of the oscillator impurity, initialised in a localised wave packet with with δ = 1/100, reaches an intermediate maximum not just at absolute zero, i.e.T = 0 (solid blue), but also when the bath is initialised at a finite, but low temperature, here T = 0.5 (dashed orange) and T = 1 (dot-dashed green).(bottom) Dependence of the entropy dynamics on the strength of the coupling, for γ = 0.001 (solid blue), γ = 0.05 (dashed orange) and γ = 0.1 (dot-dashed green).The overall qualitative features of the curve do not change.The value of the entropy at the intermediate peak is of the same order for the different coupling strengths.In contrast, the weaker the coupling, the closer to zero is the asymptotic value.The other parameters are as in Fig. 1.
where λ is the symplectic eigenvalue of the state.This is calculated as the the absolute value of the eigenvalues {iλ , −iλ } of the matrix product ΣΩ, with the symplectic matrix We consider impurity initial states which are Gaussian wave packets highly localised at the origin.We will quantify the position variance as σ xx = δ and the momentum variance σ pp = 0.5 2 /δ and take σ xp = 0.Those are pure states and saturate the Heisenberg inequality.The environment is initialised in its ground state.The global state is pure and therefore the von Neumann entropy is just the entanglement entropy.In Fig. 1 we show the entanglement entropy, calculated by Eq. ( 10), as a function of time.We observe that the entropy, starting from zero as the impurity is initialised in a pure sate, reaches a maximum at intermediate times due to the interaction and then decays back to the entropy of the mean-force Gibbs state, a value close to zero.The intermediate peak of the entanglement is higher, the more localised the initial state is.Here, even if initialised in the local ground state and staying close to it, the intermediate peak is slightly FIG. 3. The entropy S (solid blue) of a two-level impurity as a function of time.The impurity is initialised in the excited state, the bath is at a temperature close to zero.No off-diagonals of the reduced state in the impurity energy basis are generated in the dynamics and the excited state population P1 (dashed orange) monotonously decays.Therefore the Page time corresponds to the time at which the ground state and the excited state are equally populated (dotted red lines).The entropy of this state is the maximum possible value for a two-level system, i.e. ln(2).In the long time limit the entropy converges to the entropy of the mean-force Gibbs state (dot-dashed green), which is close to zero, thus behaving qualitatively like the Page curve.If initialised in the impurity ground state, the entanglement entropy approaches is final value without an intermediate peak (not shown).The parameters are ε = 1, γ = 0.001, Λ = 10 and T = 0.2.higher than the final value.The closeness is calculated by the Uhlmann fidelity, which can easily be calculated for Gaussian states [52].
In Fig. 2 we study the robustness of the entropy dynamics of the impurity against changes of parameters, by varying the environment temperature and coupling strength.In the case of non-zero temperature, the global state is not pure and the von Neumann entropy is therefor not just the entanglement entropy.We observe that for low temperatures the qualitative behaviour of the entropy as a function of time stays the same, but changes when the temperature becomes larger.At high enough temperature, the entropy of the impurity does not decay to a value close to zero, but continues to grow and saturate at a higher value.Again for an environment initially in its ground state, changing the strength of the interaction leaves the qualitative picture unchanged, as long as it stays small.The lower the coupling, the lower the final value of the entropy.Reducing the strength of the interaction, the final value of the entanglement entropy can get arbitrarily close to zero.
Spin-boson model.We now turn our attention to a case study for which no exact solution is available-the spin-boson model, i.e. a two-level system coupled to the continuum of harmonic oscillators.This model can, for example, describe an exciton interacting with phonons [53].The system Hamil- The entropy of the two-level impurity, again initialised in the excited state, as a function of time at different environment temperatures, i.e.T = 0.2 (solid blue), T = 0.25 (dashed orange) and T = 0.3 (dot-dashed green) behaves qualitatively the same and reaches the maximum entropy as long as the temperature is low enough, only the asymptotic values are different.The lower the temperature the closer to zero the asymptotic value.(bottom) Dependence of the entropy dynamics on the overall strength of the coupling to the environment, here by varying the parameter γ, for γ = 0.001 (solid blue), γ = 0.0015 (dashed orange) and γ = 0.002 (dot-dashed green).For smaller coupling the excited state decays slower, leading to a later peak in entropy and the asymptotic value gets closer to zero, as shown in the inset.The other parameters are as in Fig. 3. tonian in this case reads and the coupling between system and bath is given by where the σ α are the Pauli matrices.Like in the previous example the bath is given by Eq. ( 5) and the spectral density by Eq. (7).As the impurity couples to the environment via σ x , for which σ 2 x = 1 1 1, no counter-term is needed as there is no distortion of the system's potential [38,39,54,55].
We turn to numerical methods to solve the dynamics; namely, we use the numerically exact hierarchical equations of motion (HEOM) approach [56,57].We use the QuTiP-BoFiN implementation [58], which is integrated in the QuTiP platform [59,60].To numerically calculate the dynamics using HEOM we have to truncate two expansions.One is the Matsubara expansion of the bath correlation function in an exponential series, we call the number of retained terms N k .The other is the level of the hierarchy for the auxiliary density operators, which we call N C .To speed up convergence, the Tanimura terminator has been used [58].Our calculations have converged for N k = 30 and N C = 2, which we confirmed by comparing against calculations with higher N k and N C .
The two-level impurity is initialised in the pure excited state with density matrix ρ ρ ρ 0 = |1⟩ ⟨1| and the environment in a thermal state at a temperature close to zero.The latter is done for computational reasons, as for lower temperatures the computations are numerically more demanding.The von Neumann entropy is therefore not just the entanglement entropy, as the global state is not pure.
In Fig. 3 we show the von Neuman entropy, Eq. ( 2), as a function of time.We observe that the entropy, starting from zero as the impurity is initialised in a pure sate, reaches a maximum at intermediate times and then decays back to a value close to zero.No off-diagonals in the impurity energy basis are generated during the time evolution.The population of the excited state, also shown in the Figure, monotonously decays to its final value.The time at which the impurity has an equal population in the excited state and ground state is when the entropy reaches its maximum and starts to decay, i.e., the Page time.The entropy at this time is the maximum possible entropy for a two-level system, i.e. ln(2).
In Fig. 4 we study the parameter dependence of the entropy dynamics of the impurity by varying the environment temperature and the coupling strength.We observe that for lower temperatures the final value of the entropy is closer to zero.As the global state is also closer to being a pure state, the von Neumann entropy is closer to the entanglement entropy.Similarly, the lower the coupling between impurity and environment, the lower the final value of the entropy.At the Page time, the entropy reaches the maximum allowed value for a two-level system for all coupling strengths.

Conclusion.
In this letter we studied the entropy dynamics of quantum impurities weakly coupled to an environment close to the absolute zero.We gave a general argument why the entropy as a function of time should qualitatively look like the Page curve when the impurity is initialised in a pure state far from equilibrium.For two paradigmatic open-system models, we quantitatively illustrated this effect without making any approximations.We found that for a localised oscillator impurity, the entanglement entropy show a Page-curve-like behaviour, with the peak value of the entropy depending on how localised the impurity initially was.For a two-level impurity initialised in the excited state we found that the entropy at the Page time reaches the maximally possible value independently of coupling strength, when the excitation has half decayed.Further, we found the effect robust to changes in coupling strength as long as the overall coupling stays small and even for small non-zero temperatures of the environment.
The systems considered here are simple enough to study in detail.In particular, a discretised approximation of the environment could be followed for the whole evolution.This gives a possibility to study of how the environment becomes more pure after the Page time and retains a memory of the system initial state.The study of such simple toy models might there-fore open up new paths to study analogues to the purification of the Hawking radiation after the Page time.
Furthermore, the conditions under which our findings hold are very general and should be suitable for experimental realisation.In particular, the studied models are commonly realised in cold atoms and solid state physics.

SUPPLEMENTAL MATERIAL
The Heisenberg equations of motion for the oscillator impurity can be written as the quantum Langevin equation [20] ẍ x x(t) + ω 2 R x x x(t) − t 0 dτ χ(t − τ)x x x(τ) = F F F(t), (14) where χ(t) is the 'dissipation kernel' and F F F(t) is the so-called quantum stochastic force, which takes the form Introducing the notation Z Z Z = (x x x, p p p) T allows to write the exact solution of Eq. ( 14) as [55] Z Z Z(t) = G(t)Z Z Z(0) with entries of G(t) [G(t)] 12 = L −1 [G(t)] 21 = L −1 and F F F(t) = (0,F F F(t)) T .Here, f (s) := ∞ 0 dt e −st f (t) stands for the Laplace transform and L −1 t [ g ] for the inverse transform.For our spectral density, we have χ(s) = γΛ 2 s+Λ .With this, we can calculate the time evolution of the covariance matrix as Σ(t) = Re ⟨Z Z Z(t)Z Z Z(t) T ⟩ − ⟨Z Z Z(t)⟩⟨Z Z Z(t) T ⟩.
In the long-time limit we get FIG.1.The entanglement entropy S of the oscillator impurity as a function of time starts at zero, as the impurity is initialised in a pure state, for a localised wave packet with width δ = 1/100 (solid blue) and δ = 1/1000 (dashed orange).The interaction entangles the impurity with the environment, which is initially at zero temperature.Their entanglement entropy peaks at an intermediate time.The more localised the initial state, the higher the intermediate peak of the entanglement entropy.Regardless of initial state, the mean-force Gibbs state is approached at long times, which, for weak coupling, is close to the local ground state.As in the Page curve, the entanglement entropy thus has to decrease; here, it saturates close to zero.In the inset the initial state is the local ground state, the maximum of the entanglement entropy (dot-dashed green) still reaches an intermediate maximum, but it is just slightly higher than the final value, as the dynamics is not very far from equilibrium, seen by the fidelity F with the local ground state staying close to one (dotted red).The parameters are ω 0 = 1, γ = 0.001 and Λ = 10.

FIG. 4 .
FIG. 4. (top)The entropy of the two-level impurity, again initialised in the excited state, as a function of time at different environment temperatures, i.e.T = 0.2 (solid blue), T = 0.25 (dashed orange) and T = 0.3 (dot-dashed green) behaves qualitatively the same and reaches the maximum entropy as long as the temperature is low enough, only the asymptotic values are different.The lower the temperature the closer to zero the asymptotic value.(bottom) Dependence of the entropy dynamics on the overall strength of the coupling to the environment, here by varying the parameter γ, for γ = 0.001 (solid blue), γ = 0.0015 (dashed orange) and γ = 0.002 (dot-dashed green).For smaller coupling the excited state decays slower, leading to a later peak in entropy and the asymptotic value gets closer to zero, as shown in the inset.The other parameters are as in Fig.3.