Enhanced entanglement in multi-bath spin-boson models

The spin-boson model usually considers a spin coupled to a single bosonic bath. However, some physical situations require coupling of the spin to multiple environments. For example, spins interacting with phonons in three-dimensional magnetic materials. Here, we consider a spin coupled isotropically to three independent baths. We show that coupling to multiple baths can significantly increase entanglement between the spin and its environment at zero temperature. The effect of this is to reduce the spin's expectation values in the mean force equilibrium state. In contrast, the classical three-bath spin equilibrium state turns out to be entirely independent of the environmental coupling. These results reveal purely quantum effects that can arise from multi-bath couplings, with potential applications in a wide range of settings, such as magnetic materials.


Introduction
Entanglement is a uniquely quantum resource that has proved fundamental to the development of quantum information, such as quantum computing and quantum key distribution [1].Understanding how entanglement arises in microscopic systems is, therefore, vital for the development of quantum technologies [2].There is an ever-increasing body of work that investigates the role entanglement plays in thermodynamics [3] and whether it can be exploited as a resource [4][5][6].However, identifying what unique effects Charlie R. Hogg: she/her/hers c.r.hogg@exeter.ac.uk and potential advantages entanglement and other uniquely quantum effects can have on thermodynamic processes in general settings still remains one of the key questions of the field of quantum thermodynamics [3].
In the study of open quantum systems, the spin-boson (SB) model has long served as the standard model for the dynamics of a two-level system [7], analogous to the Caldeira-Leggett model for quantum Brownian motion [8].It describes a two-level quantum system -the spininteracting with an environment modelled as a bosonic bath of harmonic oscillators [9].It has proved highly successful in a variety of settings, from modelling spontaneous emission in two-level systems [10,11] to probing coherences in double quantum dots [12,13].In particular, the presence of entanglement between the spin and its environment is well-documented in single-bath SB systems [14][15][16][17][18][19].
For many practical applications, however, we wish to consider a spin subjected to multiple sources of noise.Indeed, from recovering classical equations of magnetisation dynamics [20] to modelling multiple vibrational interactions in singlemolecule junctions [21], multi-bath couplings naturally arise in a variety of settings.Therefore, it is highly relevant to understand the unique effects that arise due to the presence of multiple baths.However, scaling up to multiple baths is not trivial since the contributions from each bath are generally non-additive [22][23][24][25][26]. Indeed, assuming additivity has been shown to lead to unphysical results.The authors of [25] explain this practically in the context of electron counting statistics.There, they show how, under such an assumption, currents through a double quantum dot coupled to multiple environments can arXiv:2306.11036v3[quant-ph] 20 May 2024 become non-zero, even when there is no lead bias.Furthermore, the effect of multiple baths on the spin-boson model quantum phase transition has been explored [27][28][29][30][31][32][33][34].In particular, the purelyquantum effect of coupling via non-commuting operators has been found to lead to a 'frustration of decoherence' [27,28], which prevents the localisation of the spin into a common bath eigenstate.Moreover, the heat currents through a two-level system arising from coupling to multiple environments held at different temperatures [35] can be altered when the operators that couple the spin to the environment do not commute [36][37][38].Nevertheless, given the increased complexity of multi-bath couplings, many questions still remain about how quantum effects manifest in these models [17].
To explore this, one choice is to consider a spin that couples isotropically to three baths.This is natural because, as an angular momentum operator, spin has three components and each of them can couple to a bath.Three-bath spin-boson models are of direct relevance for the modelling of magnetic materials, where only three-bath models recover the well-known Landau-Lifschitz-Gilbert (LLG) equation [20].Additionally, the presence of three independent bosonic baths is required to described phonons within the underpinning lattices of bulk materials [39].
In this paper, we thus consider an isotropic three-bath SB model, where the spin is coupled arbitrarily strongly to the environment.We focus here on the spin equilibrium states: understanding the static properties of such systems is vital, for example, in magnetic modelling [40][41][42].Specifically, we demonstrate how the entanglement between a spin and its environment at zero temperature can be enhanced in the SB model by coupling to multiple baths.This, in turn, reduces the spin expectation values.We also highlight the inherently quantum nature of this effect by comparison to the classical three-bath model and its spin equilibrium states.

Three-bath spin-boson model
In the following, we will consider a single spin coupled to three bosonic heat baths [20], as illustrated in FIG. 1.Using an open system model, we write the total Hamiltonian where H S is the Hamiltonian for the system (spin), H B the Hamiltonian for the baths, and H int the Hamiltonian characterising the interaction between the system and baths [9].For a spin in an external magnetic field B ext , we have where γ is the gyromagnetic ratio of the spin, and S = (S x , S y , S z ) is either a vector of spin components for the classical case or vector of spin operators for the quantum case, both with spin length S 0 = ℏ/2.For the following, we will align B ext along z.The Hamiltonian for the baths is given by [43] where Π ω = (Π ω,x , Π ω,y , Π ω,z ) and X ω = (X ω , Y ω , Z ω ) are the vectors of three-bath momentum and position operators of the oscillator mode at frequency ω, respectively.Again, these will either be components of a vector classically or of a vector of operators in the quantum case obeying the canonical commutation relations, e.g.
In what follows, we will be considering isotropic bath couplings, so the three-bath interaction Hamiltonian takes the form [20] with c ω a scalar that sets the coupling strength between the spin and bath mode at frequency ω.For a study of anisotropic couplings, where c ω is a tensor, see [44].The spectral density J(ω) is related to the coupling function through J(ω) = c 2 ω /2ω [45].Setting c ω to be of Ohmic form (i.e.J(ω) ∝ ω [9]) not only allows for the recovery of the LLG equation [20], but also has justification from condensed matter physics.Specifically, the resulting spectral density is equivalent to using a Debye density of states, which is commonly used to model three-dimensional materials [39].However, recent experiments have revealed the presence of memory effects in magnetic materials at ultrafast timescales [46], which requires one to go beyond a linear spectral density to regimes where these effects naturally arise [20].An ideal choice of spectral density to model these strong-coupling effects is a Lorentzian [20], which we will employ here.
One of the key questions in magnetism is how to model the equilibrium properties of magnetic materials [40][41][42].The equilibrium state of the system interacting with a bath at temperature T is routinely given as the Gibbs state [47] where β = 1/k B T , k B is Boltzmann's constant, and Z S = tr e −βHs is the partition function of the system.Note here how the state only depends on the temperature of the environment and not on the strength or form of the interactions [47].However, on the quantum scale, these interactions naturally become relevant, and so the equilibrium state of the system often differs from the canonical Gibbs state [47].In such regimes, one can instead employ the 'mean force' (MF) state [45,[47][48][49][50][51][52][53], which has been used throughout biology and chemistry to capture the effect of environmental interactions at the smallest of scales [54,55] where Z = tr e −βH is the total partition function.Whilst computing the classical MF state at strong coupling is tractable, tracing out the environmental degrees of freedom in the quantum case is much tricker [45].Indeed, analytical expressions only exist for particular regimes and for a single bath [51].It remains an open question of how to extend these results to general coupling On the left is the original system: a spin coupled strongly to three bosonic heat baths.On the right is the transformed system with three modes at frequency Ω, the reaction coordinates (RCs), extracted from each of the bosonic heat baths, leaving the residual baths (RBs).By appropriately choosing forms of the spectral density defining the coupling between the spin and each RC, and the RCs and RBs, the systems can be shown to generate the same dynamical propagators [56].
strengths and multiple baths.One way to solve this issue numerically is by using the reaction coordinate (RC) mapping technique.

Reaction coordinate mapping
This method was first introduced by Garg et al. [56] in 1985 to probe strongly-coupled biological and chemical systems.Since then, it has been employed to study a variety of physical settings [35,38,[57][58][59][60][61][62].This technique involves extracting a mode (the RC) from the bath and coupling the subsystem exclusively to it, and then coupling the RC weakly to a residual bath (RB).As long as the parameters and spectral densities are chosen appropriately, this transformation can be shown to produce identical system dynamics [56].The RC method has been successfully used to study mean force corrections in the standard SB model [45].Here, we extend the RC method introduced in [56] to three baths, in a similar manner to the two-bath model of [38].This is shown Fig. 2, where the spin now couples to three RCs.Firstly, we define a transformed Hamiltonian as follows Here, H ′ S is the new, augmented system Hamiltonian which includes the spin, RCs, and the interaction between them where p, Ω, and x are the momentum, frequency, and position of the three RCs, respectively, and λ is the interaction strength between the spin and each RC.Note here that λ and Ω are the same for each bath, given that in Eq. ( 4) we have considered the same spectral density for each bath.The Hamiltonian H ′ B is now for the residual baths where Π ′ ω and X ′ ω are the three-bath momentum and position operators for the RB modes at frequency ω, respectively.We define the interaction Hamiltonian H ′ int between the augmented system and the RBs as follows where c RC ω is a new coupling function for the interaction between the augmented spin-RC system and the RBs at frequency ω.Finally, the counter term for the RCs is given by This counter term originates from the interaction of the RC harmonic oscillators with the residual baths [9].It is worth noting that no such counter term appears in Eq. (1).In fact, for spin -1/2, the counter term would be proportional to the identity, and so it can be safely ignored.We now set each RC-RB spectral density J RC (ω) = (c RC ω ) 2 /2ω to be Ohmic where Γ is the dissipation strength and Λ the cutoff frequency, which is taken to be larger than all other relevant frequencies.The propagators for the original and transformed systems are identical when the spectral density for Eq.(4) takes a Lorentzian form [56] with resonant frequency Ω and peak width Γ.
The closed form of the mapping here demonstrates the benefit of choosing such spectral densities.Here, we will express the Lorentzian parameters as multiples of the Larmor frequency ω L = γB ext .As in [20], we will write the coupling amplitude λ 2 = ω 2 L α/S 0 , introducing a new parameter α with units of frequency.This allows us to express the spin-RC coupling strength as We now have an expression for the RC Hamiltonian (7), which has been shown to generate the same propagators as the untransformed Hamiltonian (1) [56].For small enough Γ, the augmented spin and RC system is weakly coupled to the RBs.The quantum mean force (QMF) state can then be obtained by taking the Gibbs state of the augmented system and tracing out all three RCs [60] where is the partition function for the spin-RC system.Note that the assumption of small Γ imposes no restriction on the strength of the original coupling between the spin and baths, given that λ can be arbitrarily large [45].We now have all the tools required to calculate the classical and quantum equilibrium spin expectation values s x,y,z = tr[τ MF S x,y,z ]/S 0 for the three-bath model.

Entanglement in the three-bath model
Here, we demonstrate how the presence of multiple baths can have drastic effects on the equilibrium properties of the spin.Using the RC mapping method, we evaluate the expectation values s x,y,z for the QMF for the three-bath SB model (1) -see Eq. (15).Fig. 3a) shows s z as a function of temperature T (solid blue).It is immediately obvious that the QMF s z at T = 0 K starts at a reduced level, compared to the quantum Gibbs state (QG, dashed blue), where the coupling between spin and baths is negligible.For both the QMF and QG, the two spin components s x,y orthogonal to the external field B ext are zero for all temperatures.For the QMF, this is highlighted by the match of the s z expectation value with |s| = s 2 x + s 2 y + s 2 z (grey line).The three-bath interaction Hamiltonian (4) couples the spin to three independent baths, as shown in the inset of Fig. 3a).We now wish to contrast this with the commonly studied case  16).The solid blue line is s z for the quantum mean force (QMF) state, see Eq. (15).The key difference between a) and b) is that the spin components orthogonal to B ext , s x,y , are only non-zero for the single-bath model, as seen from the difference (shaded) between |s| (grey line) and s z .At T = 0 K, this indicates significantly more entanglement between the spin and its environment in the three-bath model.To verify the numerical accuracy of the QMF state, we also plot the quantum dynamical steady states from the RC master equation (QD, blue crosses) [57] and TEMPO (black dots) [26,63,64].The spin length is S 0 = ℏ/2.For the QMF, the parameters are ω 0 = 2ω L , Γ = 0.6ω L , α = 10ω L , and the spin-RC coupling λ is given in (14).For the QD, we further set Λ = 10 10 ω L .The quantum Gibbs state (QG, dashed line, Eq. ( 5)) is shown as a reference.
in which the spin couples to a single bath.Of course, it is worth noting that the three-bath and single-bath models will always be fundamentally different for any choice of single-bath coupling direction.However, as we argue in the following, there is a natural choice of single-bath coupling direction that makes such a comparison fair.We proceed by coupling each spin component in (4) to the same bath with position operator R ω , i.e. the interaction Hamiltonian becomes Here, S R = (S x +S y +S z )/ √ 3 = sin θ 0 cos ϕ 0 S x + sin θ 0 sin ϕ 0 S y + cos θ 0 S z is the spin operator in the direction given by spherical angles θ 0 = arctan √ 2 and ϕ 0 = π/4.Notice here how √ 3 appears naturally as a normalisation to the systembath coupling function c ω .We illustrate this coupling graphically in the inset of Fig. 3b).The interaction terms (4) and (16) permit a fair comparison between the three-and single-bath models in the sense that they do not change the 'weight' of the contribution of each spin component operator to the interaction energy, nor the overall strength of the interaction.In fact, in the threebath case (4), the coupling to the environment must always have both a longitudinal as well as transversal component with respect to the external field.To have a fair comparison, the singlebath case should also have both components and in proportional amount.This fixes the singlebath coupling direction as the one chosen here, see (16).As before, we use the RC mapping method to compute s z for the single-bath QMF (solid blue), which is plotted as a function of T in Fig. 3b).
We will now proceed to quantify the zerotemperature entanglement between the spin and baths for both models.Whilst there are many measures of quantum entanglement, e.g.concurrence and logarithmic negativity [1], many are hard to compute [65].However, in the limit T → 0 K, the global (spin+bath) Gibbs state becomes pure.The reduced state of the spin, τ qu MF , is then mixed whenever the spin is entangled with the baths.The presence of these correlations can be quantified by the entanglement entropy [1] The entanglement-caused mixedness can be further quantified by the spin's purity P. For a spin-1/2, this is [1] The purity is minimal, P = 1/2, when the spinmagnitude is zero, |s| = 0, i.e. all three spin components must be zero.Here, the entanglement is maximal, E max = ln 2. On the other hand, the maximal purity is P = 1, which occurs whenever |s| = 1.Here, the spin state is pure, and hence there is no spin-bath entanglement.This allows us to link the T = 0 K value of |s| from Fig. 3 to the entanglement entropy: namely, the lower the value of |s|, the higher the entanglement entropy between spin and baths.Evaluating Eq. ( 17) numerically for both models, we find i.e. the three-bath entanglement entropy is almost a whole order of magnitude greater than the single-bath entanglement entropy.Note that, whilst the entanglement entropies (19) and (20) are for specific Lorentzian parameters ω 0 = 2ω L , Γ = 0.6ω L and α = 10ω L , this difference is still observed at other strong coupling values.Despite a similar reduction in the value of s z between the single-and three-bath models, see Fig. 3, the origin is fundamentally different.Indeed, the spin in the single-bath model 'compensates' by increasing s x,y .This is clearly demonstrated in Fig. 3b) by |s| (grey line) being much greater than s z for the single-bath QMF state at T ≲ ℏω L /2k B .We can understand this as follows: in the single-bath case, the primary effect of the interaction with the bath is to align the spin towards the coupling axis [45,51].This is in stark contrast to the three-bath model where the spin is aligned parallel to B ext , and the sole impact of the environmental coupling is to introduce disorder into the system via entanglement.This is further illustrated by considering the values of |s| 2 at T = 0 K: |s| 2 ≈ 1 for the singlebath model, corresponding to a nearly pure state [1].Conversely, |s| 2 ≪ 1 in the three-bath case, meaning that the spin state is much more mixed.
To understand this significant difference between the single and three-bath cases, we can consider the ultrastrong coupling regime.For the single-bath model, we know that the ground state of the spin tends towards an eigenstate of the interaction Hamiltonian in this limit [51].The single-bath interaction Hamiltonian (16) is a product operator of spin and bath.Hence, its eigenstates are product states.In contrast, for the multi-bath case, since the S x , S y , and S z coupling operators do not commute, the eigenstates of the interaction Hamiltonian (4) do not admit such a decomposition into product states.This is reminiscent of the 'frustration of decoherence' that has been studied in the context of phase transitions [28], and is unique to multi-bath models.In the single-bath case, there is no such competition between baths, which appears to drastically reduce the entanglement between the spin and its environment.
Up until now, we have focused on the zerotemperature case.As temperature increases, we expect entanglement to decrease [66,67], making the low-temperature limit studied here the most relevant case to quantify the effect of multiple baths.It is also worth noting that, whilst we have discussed in detail here the isotropic model, the enhanced entanglement still manifests in the presence of small anisotropies to the system-bath coupling.
The results presented here demonstrate how the presence of multiple dissipation channels in the SB model can have a significant impact on the equilibrium properties of the spin by substantially increasing entanglement with the baths.Hence, a careful assessment is needed to determine if mapping to a single-bath model is sufficient to describe the physics, or whether a multi-bath treatment is needed.

Quantum vs. classical models
We will now compare the equilibrium states of the quantum and classical versions of the three-bath SB model to highlight the quantum nature of the reduction in s z .Fig. 4 shows the expectation values s z for the QMF (solid blue) and the classical mean force state (CMF, solid green, Eq. ( 6)) as a function of temperature T .Remarkably, we find that the environmental corrections to the threebath QMF, discussed before in Section 4.1, vanish entirely for the CMF, which becomes independent of environmental coupling.This is evidenced by the match between expectation values from the CMF and classical Gibbs (CG, dashed  15)) and classical mean force (CMF, solid green, Eq. ( 21)) states.Whilst the CMF and classical Gibbs state (CG, dashed green, Eq. ( 5)) coincide, this is not the case for the quantum model, where there are stark differences between the QMF and quantum Gibbs (QG, dashed blue, Eq. ( 5)) predictions from entanglement.The spin length is S 0 = ℏ/2, and the parameters for the QMF plots are: ω 0 = 2ω L , Γ = 0.6ω L and α = 10ω L .green, Eq. ( 5)) states.
To see why no mean force corrections arise classically, one can explicitly evaluate the CMF.For the three-bath isotropic model, the CMF state (6) takes the form (see Appendix A for details) where Zcl S is the partition function that normalises τ cl MF , and Q the reorganisation energy Since S 2 0 is a constant, the mean-force correction in (21) is cancelled out entirely by the partition function, and we are left with the standard Gibbs state τ cl G = e −βH S /Z cl S -see Eq. (5).This is a unique feature of the highly-symmetric three-bath model considered here [44].We highlight that, despite sharing the same symmetry, such a simplification does not occur in the quantum case.Whilst in the classical isotropic threebath SB model, the spin equilibrium state is environment-independent, we have shown here that entanglement between the spin and baths plays a key role in the quantum model.We demonstrate that such a correction in the quantum case can have a significant impact on the spin equilibrium properties.Indeed, whilst the quantum and classical spins behave in a qualitatively similar manner in the single-bath model, i.e. by aligning along the coupling axis [45], no such analogue exists for the isotropic three-bath case, where strong-coupling effects are due to entanglement with the environment.

Numerical accuracy
In general, computing the QMF state (6) proves to be a difficult task, and, as such, we evaluate it numerically using the RC mapping framework.To verify the accuracy of these QMF states (see Figs. 3 and 4), we also calculate the steady states from both the RC master equation (QD, blue crosses) [57] and, for the three-bath model, multibath TEMPO (TEMPO, black dots) [26,63,64].At sufficiently large times, convergence of both these states to the QMF is observed; see Fig. 3.The agreement of the TEMPO state is of particular importance as it is entirely independent of the RC mapping and its assumptions.
Aside from verifying our results, this convergence also provides strong numerical evidence that the dynamical steady state of the system is indeed the QMF.Whilst, at least for the singlebath case, this is proven to be true in both the weak [68], and ultrastrong [69] limits, there are no such proofs in general for the case of intermediate coupling, especially for multiple baths.
Analytical expressions for the QMF are also limited, although they do exist in particular regimes such as the weak and ultrastrong coupling limits [51], and recently an approach has been proposed to tackle the intermediate regime [70].To further establish the numerical accuracy of the QMF, we compute the three-bath weak limit by extending the calculation detailed in [51] to multiple baths, of which details are given in Appendix B. Fig. 5 shows s z calculated using this weak expansion (WK, light grey dots) alongside the zero-temperature limit (WK T = 0 K, light grey line).As before, we also plot the s z obtained from the QG (dashed blue line) and QMF (solid blue line) states.We see excellent agreement between the predictions of the QMF and our weak expansion.

Conclusions and outlook
In this paper, we explored the equilibrium states of a quantum spin coupled isotropically to three bosonic baths.We found that the environmental coupling induces significant entanglement between the spin and baths at zero temperature.We compared these three-bath results to those of a comparable single-bath model, finding substantially more entanglement in the three-bath case.
In fact, the strong coupling effects of a singlebath manifest largely as coherences, which have a classical analogue [45,71], whereas there exists no such analogue for entanglement.The enhanced effect of entanglement in the three-bath model leads to marked corrections to the canonical Gibbs state in the form of reduced spin expectation values at zero temperature.This contrasts with the equilibrium state of a classical spin vector coupled isotropically to three baths, for which we proved that no such correction occurs, even for large environmental couplings.Finally, we provided an analytical expression for the threebath weak coupling state.However, it remains an open question to find a general analytical ex-pression for the spin equilibrium states in other regimes.
The results shown here give insight into the impact of multiple dissipation channels upon spin equilibrium properties of a spin.The direct application is to magnetic modelling, where such multi-bath Hamiltonians arise naturally [20].Additionally, we expect this work to find use in thermodynamic processes that rely upon the presence of multi-baths, such as thermocurrents [36][37][38] and quantum heat engines [59,72].Finally, we anticipate that our choice of environmental coupling will add to the interest on non-commuting operators in quantum thermodynamics, for example, in recent work on generalised Gibbs states with non-commuting conserved quantities [73,74].

B Weak-coupling limit
We will write the interaction Hamiltonian (4) as where The vector operator B can be written in terms of its vector components as B = i B i e i , i = x, y, z.
Thus the spin is coupled to three baths, which will be assumed to be independent in the sense that B i , B j = 0.Each vector operator component can be expressed as where For three independent bosonic baths, the annihilation and creation operators a w,i and a † ω ′ ,j will satisfy the bosonic commutation rules a ω,i , a † ω ′ ,j = δ ij δ(ω − ω ′ ).The bath Hamiltonian will then be with [H B,i , H B,j ] = 0.The Gibbs state of the bath with all components at the common inverse temperature β will be τ We will now suppose that the combined system comes to equilibrium at the same inverse temperature β as the initial temperature of the bath.The partition function for the combined system will then be where ρS = tr B e −βH is the unnormalised reduced state of the system.Expanding e −βH to second order in the interaction H int , i.e. the Kubo expansion, gives this reduced state as where, for any Schrödinger operator O, the 'thermal interaction' picture operator O(β) is given by O(β) = e βH 0 Oe −βH 0 , with H 0 = H S + H B .We can then write, and noting that [B i (β), S j (β)] = 0, We then find that  The system correlation function i tr S [τ S (β)S i (β ′ )S i (β ′ − β ′′ )] can be readily evaluated using H S = −ω L S z and τ S (β)S ± = S ± τ S (β)e ∓βℏω L .Further simplification follows by making use of the spin operator commutation rules.After carrying out the double integral and taking the trace over the system, we find that the reduced partition function for the system will then be, with Z with We can use the partition function to calculate the expectation value of S z as follows )

Figure 1 :
Figure 1: Three-bath spin-boson model.A single spin coupled to three bosonic heat baths via the S x , S y and S z components of the spin vector (classical) or vector of operators (quantum) S [20].Here, the external field B ext is aligned in z.

Figure 2 :
Figure 2: Three-bath reaction coordinate mapping.On the left is the original system: a spin coupled strongly to three bosonic heat baths.On the right is the transformed system with three modes at frequency Ω, the reaction coordinates (RCs), extracted from each of the bosonic heat baths, leaving the residual baths (RBs).By appropriately choosing forms of the spectral density defining the coupling between the spin and each RC, and the RCs and RBs, the systems can be shown to generate the same dynamical propagators[56].

Figure 3 :
Figure 3: Impact of multi-bath coupling on spin-boson equilibrium states.Panel a) shows the normalised equilibrium expectation s z = ⟨S z ⟩/S 0 for a spin coupled isotropically to three baths X, Y and Z via the S x , S y and S z operators, respectively -see Eq. (4), whilst b) shows s z for a spin coupled to a single bath R through the S R = (S x + S y + S z )/ √ 3 operator -see Eq. (16).The solid blue line is s z for the quantum mean force (QMF) state, see Eq.(15).The key difference between a) and b) is that the spin components orthogonal to B ext , s x,y , are only non-zero for the single-bath model, as seen from the difference (shaded) between |s| (grey line) and s z .At T = 0 K, this indicates significantly more entanglement between the spin and its environment in the three-bath model.To verify the numerical accuracy of the QMF state, we also plot the quantum dynamical steady states from the RC master equation (QD, blue crosses)[57] and TEMPO (black dots)[26,63,64].The spin length is S 0 = ℏ/2.For the QMF, the parameters are ω 0 = 2ω L , Γ = 0.6ω L , α = 10ω L , and the spin-RC coupling λ is given in(14).For the QD, we further set Λ = 10 10 ω L .The quantum Gibbs state (QG, dashed line, Eq. (5)) is shown as a reference.

Figure 4 :
Figure 4: Quantum vs. classical isotropic three-bath spin-boson model.Plotted are the normalised equilibrium spin expectation values s z = ⟨S z ⟩/S 0 over temperature for the quantum mean force (QMF, solid blue, Eq. (15)) and classical mean force (CMF, solid green, Eq. (21)) states.Whilst the CMF and classical Gibbs state (CG, dashed green, Eq. (5)) coincide, this is not the case for the quantum model, where there are stark differences between the QMF and quantum Gibbs (QG, dashed blue, Eq. (5)) predictions from entanglement.The spin length is S 0 = ℏ/2, and the parameters for the QMF plots are: ω 0 = 2ω L , Γ = 0.6ω L and α = 10ω L .

Figure 5 :
Figure 5: Weak coupling limit for a quantum isotropic three-bath spin-boson model.Plotted are the normalised equilibrium spin expectation values s z =⟨S z ⟩/S 0 as a function of temperature for a coupling of α = 0.1ω L .We observe excellent agreement between the predictions of the mean force Gibbs state (QMF, solid blue line), see Eq.(15), and the analytical weak expansion (WK, grey dots) and its zero-temperature limit (WK T = 0 K, light grey line), see Appendix B. The spin length S 0 = ℏ/2, and the parameters for the QMF plots are: ω 0 = 2ω L and Γ = 0.6ω L .The Gibbs state (QG, dashed blue line, Eq. (5)) is shown here for reference.