Entanglement generation between distant spins via quasilocal reservoir engineering

The generation and preservation of entanglement is a central goal in quantum technology. Traditionally, dissipation in quantum systems is thought to be detrimental to entanglement, however dissipation can also be utilised as a means of generating entanglement between quantum spins that are not directly interacting. In particular entanglement can be generated between two qubits, or multi qubit systems via a collective coupling to a reservoir. In this work, we explore multiple spin domains pairwise coupled to different reservoirs and show that entanglement can be generated between spins which are not coupled to each other, or even coupled to the same reservoir.


I. INTRODUCTION
Quantum entanglement enables quantum computation and communication that surpass their classical counterparts in various tasks [1][2][3][4][5][6].Therefore, the generation and preservation of entanglement represents a significant goal within the field of quantum technology.Typically, dissipative effects are known to be detrimental to quantum entanglement [7], and significant efforts are made into isolating the quantum system from its environment.However, the environment can also be used as a means of generating entanglement between different quantum systems [8][9][10][11][12].
With recent advancements in hybrid quantum systems, we can engineer both the system and the environment to perform the task we are interested in [13][14][15][16].Such hybrid systems can involve many different elements, from solid state to atomic, molecular and optical [17][18][19].Hybrid quantum systems naturally provide a rich environment to study complex phenomena in open quantum systems [20].Our focus in this work is the behavoir of many spin ensembles collectively coupled to bosonic reservoirs and the entanglement that can be generated through this system.The collective coupling of multiple spins to the same bath has been shown as a mechanism to transfer excitations between spins [21][22][23] and also as a means to charge quantum batteries [24].Furthermore mutliple spin domains pairwise coupled to mulitple reservoirs can transfer excitations to and from domains which are not directly coupled [25].
Open quantum systems form a complex and dynamic backdrop to study the generation and dynamics of entanglement (see Ref. [26] for a review).It is known that entanglement can be generated between two qubits via the interaction of the qubits with the same reservoir [9][10][11][27][28][29][30][31][32] or the interaction of two oscillators with the * josephine.dias@oist.jpsame reservoir [33].A collective coupling to a common bath can also be used to generate entanglement between multiple spins [12] and spin domains [22], even if the reservoir is in a thermal state [10,11].For instance, a system of two spin domains collectively coupled to the same bosonic reservoir was explored and its steady state was shown to have entanglement between spin-domains [22].
In order to illustrate the mechanism by which initially separable systems can become entangled via a collective coupling to a reservoir, consider an elementary example of two spins A and B collectively coupled to the same zero-temperature reservoir and initialised in the separable state |ψ⟩ AB = |↑⟩ A ⊗|↓⟩ B .This state can also be written as |ψ⟩ AB = 1 √ 2 (|Ψ + ⟩ AB + |Ψ − ⟩ AB ) which is an equal superposition of the symmetric triplet state |Ψ + ⟩ AB and antisymmetric singlet state |Ψ − ⟩ AB .Under the action of the collective dissipator J − A + J − B , the triplet state decays to the ground state |↓⟩ A ⊗ |↓⟩ B , however the singlet state cannot decay.Therefore, the steady state for this system is [9] and this is an entangled state between spins A and B. If we were to quantify the amount of entanglement in this state, we could use the entanglement of formation [34][35][36] (see Sec. III) which for this two-qubit state is E F ≈ 0.354, that is 0.354 ebits.We emphasise here that the spins were not initially entangled, nor did they evolve according to a Hamiltonian coupling individual spins.The entanglement present in the steady state of (1) is solely due to the relaxation of the spins via their collective reservoir coupling.
To go beyond this simple example of two spins, consider the system depicted in Fig. 1 (a) (which is the system considered in Refs.[21,22]).Here we have a system of N A spins in domain A, and a single spin in domain B. With this system initialised in the state of after collective spin relaxation (for a zero-temperature bath), the steady state will be in the limit of N A → ∞.That is, through the collective coupling of the spins to the bosonic reservoir, the steady state of domain B is excited even though it initially started in its ground state [21,22].As can be seen from ( 3), the steady state in the limit of N A → ∞ is separable.However, for finite N A , collective spin relaxation produces entanglement in the steady state between domains A and B. In fact, the logarithmic negativity [37,38] for the steady state peaks for a spin population of N A = 5 and has a value of E N ≈ 0.55 [22,39].
To illustrate how we can extend the model of Fig. 1(a) in Ref. [21,22], consider the system depicted in Fig. 1(b): In addition to the two spin domains of Fig. 1(a) we add another spin domain C collectively coupled with domain B to a second reservoir.That is, domains A and B are collectively coupled to the same zero-temperature reservoir and domains B and C are collectively coupled to a different zero-temperature reservoir.This means we have the simultaneous action of two collective dissipators, J − A + J − B and J − B + J − C .For an initial state of we can expect that for N A sufficiently large, the coupling to the first reservoir will result in the spins in A relaxing to the ground state and the single spin in B will move to its excited state.For large N A , domains A and B will relax at a time that scales with 1/N A due to the superradiant decay [40,41].We can then consider the coupling to the second reservoir.In this case, the single spin in B will be excited (for large N A ) and the single spin in C has already been initialised in its ground state.
With the collective coupling to the second reservoir, we can expect these two spins to end in an entangled state.
In fact, for N A → ∞, we can expect the relaxation of B and C to proceed to the steady state described in equation (1) which is an entangled state with entanglement of formation E F = 0.354, that is 0.354 ebits.
In this work we explore the generation of entanglement through multiple spin domains with multiple reservoir couplings.Unlike previous works investigating entangling qubits through interaction with a common reservoir [9][10][11][27][28][29][30][31][32], we will show here generation of entanglement between two spins in separate domains which have never interacted and are not even coupled to the same reservoir.In particular, we will also see the amount of entanglement (between qubits not coupled to the same reservoir) to be of similar level to that of the previously mentioned example of two qubits collectively coupled to the same zero temperature reservoir.This paper is arranged in the following way: In Sec.II we present details of the model considered in this work.In Sec.III, we use that model to analyse the capabilities of our system to distribute entanglement.We then look at how entanglement distribution is affected by realistic conditions such as individual decay, dephasing and thermal reservoirs in Sec.IV.Finally, we will also show how multiple systemreservoir couplings can be used to generate tripartite entangled states in Sec.V, before presenting our summary and conclusions in Sec.VI.

II. OUR MODEL
The system we examine in this work consists of M separate spin domains, each consisting of N m identical spin-1/2 particles with frequency ω 0 /2π.Although in our system, we assume each spin has the same energy, entanglement can still be generated through a common environment when this is not the case [42].The spins in two neighbouring spin domains are collectively coupled to a bosonic reservoir.Additionally, our system is symmetric under exchange of two spin-1/2 particles in a given spin domain.The Hamiltonian for our system of M spin domains and M − 1 reservoirs is given by: where the operators a † km (a km ) are the creation (annihilation) operators of the k th m bosonic mode of the reservoirs with frequency ω km /2π.The operators J α i are the collective spin operators for the m th spin domain and are given by J α m = 1/2 Nm nm=1 σ α nm where σ α nm are the n th m Pauli spin operators with α = x, y, z.The collective spin raising and lowering operators for the m th spin domain are J ± m = J x m ± iJ y m .In (5), the first term represents the Hamiltonian of the spin domains, the second term represents the Hamiltonian of the bosonic reservoirs and the third term represents the interaction between the spin domains and the reservoirs.In the third term, each term in the sum over M − 1 reservoirs describes the coupling between the m th and (m + 1) th spin domain which are collectively coupled to the m th bosonic reservoir, with t km (t * km ) representing the emission (absorption) amplitudes that fix the spectral density of the reservoirs, Γ m (ω) = 2π km |t km | 2 δ (ω − ω km ).Within the Born-Markov-secular approximation, the master equation in the rotating frame can be written as [20,[43][44][45]: for bosonic reservoirs at zero temperature with Previous work has experimentally implemented a double spin domain system coupled to a common resonator where each spin domain consisted of nitrogen-vacancy (NV) ensembles in diamonds [46].Additionally, the work of Ref. [47] implemented a double nuclear spin-domain system coupled to the same bosonic mode.The nuclear spin ensemble size is on the order of N > 10 12 which means that the time scale on which collective effects like superradiant decay occur is extremely fast (< 10 −10 s) while the associated individual dephasing and thermalisation times are much slower (on the order of T * 2 ∼ 1 ms and T 1 ∼ 40 s respectively).In this work, we primarily consider the case of three spin domains coupled to two reservoirs and therefore our master equation is given specifically by: where we have used the same reservoir coupling constant γ for all system-reservoir couplings.Further, we can introduce a scaled time γt/2 to simplify our considerations.  .Three spin domains pairwise coupled to two separate reservoirs at zero temperature.Collective relaxation of the spins produces entanglement between the spins.Entanglement between spins in outer domains is maximised for a configuration with a large number of spins in the middle domain.

III. ENTANGLEMENT DYNAMICS
In this work, we extend the results of Ref. [22] to more than two spin domains with more than one reservoir coupling.We will show here that entanglement can be generated not just through collective coupling to the same reservoir but also through a chain setup where entanglement will be generated between spins which are not coupled to the same reservoir.
To begin, we will conisder here the situation where domains A and C consist of a single spin in each domain, while the central domain B has N B spins.This is depicted in Fig. 2. Our goal is to entangle the qubits (in domains A and C) which are not coupled to the same reservoir.After this sytem evolves in time according to equation (7), through collective spin relaxation, entanglement is generated between all the spins, however, here we primarily focus on the entanglement between the two single spins in domains A and C. In the remainder of this paper, we will examine how to maximise the entanglement between spins in A and C as well as how robust this entanglement is to realistic imperfections.
Our initial state configuration is: which is depicted in the top part of Fig. 2.However, different initial state configurations can still produce entanglement (see Appendix A).After our initial state (8), evolves in time we then trace out the spins in domain B leaving us with the reduced density matrix ρ AC = Tr B [ρ] of the state of the two spins, one in domain A and one in domain C. To measure the amount of entanglement present in the state ρ AC , we will use the entanglement of formation, which quantifies the minimum amount of pure state entanglement needed to prepare the state [34].For an arbitrary state of two qubits, this quantity is given by [36]: where E is given by where λ i are the square roots of the eigenvalues of the matrix ρρ with ρ = (σ y ⊗ σ y ) ρ * (σ y ⊗ σ y ) [35].
In Fig. 3(a) we plot the entanglement generation dynamics between the two spins in domains A and C as well as the relaxation dynamics of the initially excited domain B.Here we can see that as the number of spins in domain B increases, the decay time decreases.This is indicative of the collective effect of superradiance [40,41] where the radiation is enhanced and sped up in the presence of mutliple emitters.The increase in spin domain size also contributes to enhanced entanglement dynamics -that is the speed of entanglement generation is increased as is the amount of entanglement present between qubits in A and C. Note that only collective dissipation into a zero temperatature reservoir is govering the dynamics, and for this case, the entanglement present in the steady state remains and does not decay.We also emphasise here that entanglement is generated between spins which are not coupled to the same reservoir.For the system size of N B = 12, the entanglement generated between spins in domains A and C is E F ≈ 0.315 ebits.
Next in Fig. 3(b), we plot the maximum value of the entanglement of formation E F (reached at steady state) as well as the value of γt E /2 where t E is the time for entanglement between A and C to reach half of the maximum entanglement present during the system evolution.It can clearly be seen that the entanglement distribution is enhanced (more entanglement and faster distribution) for a large spin domain size in B. Results here are limited to small total system sizes (N A + N B + N C ) due to the memory required to simulate the combined system.As can be seen from Fig 3(b) as we increase N B , the increase in entanglement is rapid at first, however the growth steadily slows.Thus, for even bigger N B the payoff with increased entanglement between A and C will not be significant.We find that for N B → ∞, the steady state of the reduced system ρ AC goes to where , and thus the entanglement of formation of this state is the same as the first example in the introduction, that is E ≈ 0.354 ebits.In Appendix B, we also illustrate the steady state solution for any spin population N B with initial state While the previous results focus on distributing entanglement between spins in the outer domains of a three spin domain chain, entanglement generation through collective spin relaxation is not limited to three domains.In Fig. 4 we present results for the entanglement between the single spins in the outer domains of spin domain chains consisting of four (pink, dashed line) and five (blue, dotted line) spin domains.These dynamics use the initial state depicted in the inset diagrams in Fig. 4 respectively, that is all spins initialised in the ground state with the exception of domain B which has all spins initialised in the excited state.While we only show entanglement for a single spin domain population configura- tion for four and five domains, as before, entanglement between outer qubits is maximised for larger spin popula-tions of the middle domains.Here, we have constrained the total number of spins in all middle domains to be the same (12 spins), however the top line uses 12 spins split between two spin-domains and the bottom line uses 12 spins split between three spin-domains.These results are limited due to the total size of the system however, they show that entanglement can be distributed between qubits separated by two or three spin domains which are in turn connected by collective coupling to three or four reservoirs.For this result, as is the case with the previous results, the only effect taken into account in this model is collective decay into zero-temperature reservoirs.In this case, we observe constant entanglement which does not decay in the steady state.

IV. INDIVIDUAL DECAY AND DEPHASING
The results contained in previous section are interesting but represent the highly specific and idealized system of collective relaxation with a zero-temperature bath and with no other individual effects present.It is important to consider more realistic models, which includes effects such as dephasing, individiual spin relaxation and non-zero temperature reservoirs, and examine how these aspects may alter the previously observed entanglement generation.
For the three domain system depicted in Fig. 2, we now incorporate finite temperature reservoirs as well as dephasing and individual relaxation.This yields the following master equation, again in the rotating frame of the system [47]: In (11), the first and second terms correspond to the collective thermalisation between the pairs of spin domains (A and B, B and C respectively) and the thermal reservoir.Associated with that reservoir is the Bose-Einstein distribution n = 1/ e ℏω0/k B T − 1 for a given temperature T and with Boltzman constant k B .Next, the third and fourth terms in (11) correspond to individual thermalisation of the spins.Finally, the last term corresponds to dephasing.To summarise, the first line in (11) corresponds to collective effects while the second line corresponds to individual effects associated with thermalization and dephasing.Note that the individual and collective thermalisation occurs with the same system-reservoir coupling constant γ.
In Fig. 5(a), we plot the entanglement dynamics of the original system with no dephasing (pink, dashed line) as well as the entanglement between spins in domains A and C with various dephasing rates of γ dep /γ (broken lines).We observe a reduction in peak entanglement generated through the system-bath evolution in time but also a gradual reduction in the entanglement present in the steady state.For the results in Fig. 5(a), the dephasing is the only additional effect present, the reservoirs remain at zero temeperature (n = 0) and there is no individual decay.
We then examine the entanglement generation capabilities of our system under the effects of individual spin coupling to the reservoir.In this case, the spins all have individual as well as the collective dissipation.Again, we shall keep the bosonic reservoir to zero temperature to isolate the effect of individual decay on entanglement dynamics.We find that individual coupling has a signficant detrimental effect on the entanglement present.Shown in Fig. 5(b) is the entanglement dynamics between single spins in domains A and C. Note the significantly smaller scale of the y-axis in Fig. 5(b) relative to the other results in this work.The maximum entanglement of formation produced in this case is E F = 0.0039.Note that as the system size increases (number of spins), the superradiant collective effects occur on a much faster time scale than the individual one and most likely dephasing.For example in NV centers in diamond, the collective relaxation can be of order of microseconds, while the individual one is 10 orders of magnitude larger [48].
Lastly, the temperature of the thermal reservoir also substantially affects the entanglement generated through this system.Shown in Fig. 5(c) is the entanglement dynamics between qubits in A and C for thermal reservoirs of various temperatures for spin domain population of N B = 11.For this result, we keep only collective spin-reservoir coupling to the bath.As can be seen in Fig. 5(c), the entanglement generated through our system is highly sensitive to the temperature of the thermal reservoirs.As the temperature increases, less entanglement is generated and in the steady state entanglement is also reduced.
We have seen in this section the impact of realistic effects on entanglement distribution via collective dissipation with multiple spin domains and reservoirs.While it is clear that entanglement can still be generated between spins in domains A and C, it significantly decreases if individual decay or finite temperatures are taken into account.Further, the constant, non-decaying entanglement in the steady state observed in Sec.III is no longer present once any of these effects are taken into account.In the results of Figs.5(a)-(c) we observe a decline in the entanglement of formation over time.However, individual process are always detrimental to entanglement between spins, and therefore the decay of the entanglement in the results of Fig. 5 is to be expected.

V. TRIPARTITE NEGATIVITY
The previous results in this work focused solely on entangling two qubits through collective spin relaxation with multiple collective spin-reservoir couplings.In this section, we will show how this same process can be utilised to generate tripartite entanglement between three spins.Consider the setup shown in Fig. 6(a).This set up is reminiscent of a spin-star configuration [49,50] but we emphasise again the spins are not directly coupled, merely collectively coupled to the same reservoir as other spins.Here, a central spin-domain (labelled D) contains many spins and is collectively coupled to three bosonic reservoirs with three other spin-domains each containing a single spin only (labelled A, B and C).The system is initialised in the separable state: Under the process of collective dissipation with the reservoir setup illustrated by the diagram in Fig. 6(a), the steady state of A, B and C (with domain D traced out) is given by: ) that is, the steady state of the reduced system of three qubits A, B and C is a mixture of the ground state and the W-state which is |W ⟩ = 1/ √ 3 (|↑↓↓⟩ + |↓↑↓⟩ + |↓↓↑⟩).Therefore, the steady state of this reduced system has tripartite entanglement between the three remaining qubits.In order to demonstrate the presence of tripartite entanglement in the steady state of this system, we calculate the tripartite negativitity, defined as the geometric mean of the negativities of all bipartitions [51]: (14) where each of the terms N i|jk are the bipartite negativities between subsystems i and jk given by N (ρ) = ||ρ Γi || − 1 /2 with ρ Γi corresponding to the partial transpose of the density matrix with respect to subsystem i and ||X|| is the trace norm of the operator X given by ||X|| = Tr √ XX † [37,52].Note that with this definition of the bipartite negativity, maximally entangled Bell states have the negativity of 0.5.Further, for three qubit entangled states, GHZ states yield a tripartite negativity of 0.5 but W states have a tripartite negativity of 0.47.In Fig. 6(b), we give the tripartite negativity generated via the configuration shown in Fig. 6(a) with reservoirs at zero-temperature and no individual effects.The dynamics show clear presence of tripartite negativity in the steady state.For the given spin-domain population of N D = 11, the steady state tripartite negativity is N ABC ≈ 0.043.Additionally, the inset shows the steady state tripartite entanglement for various values of N D , the population of the central domain.As the spin-domain population of D grows, the steady state contribution of the W-state term [c W in (13)] grows larger, and thus the tripartite negativity is maximised for large N D .However, while this increase is rapid at first, it slows as N D continues to increase indicating less payoff in tripartite correlations for further increasing the size of the central domain D. Nevertheless, the results here show the generation of tripartite correlations between three qubits not coupled to each other or even the same reservoir.

VI. SUMMARY AND CONCLUSIONS
We have shown here how, by utilising reservoir engineering, and with a completely separable initial state, entanglement can be generated between two spins which are not coupled to each other or even coupled to the same reservoir.We have focused mostly here on the case where the two spins collectively decay into different reservoirs with a single common spin domain collectively decaying into both reservoirs.For this case, the entanglement is maximised for a large number of spins in the central domain.The speed of entanglement generation is also maximised for a larger spin population in the central domain.Additionally, the entanglement betweeen the qubits not coupled to the same reservoir, can approach the level of entanglement generated when two qubits are collectively coupled to the same zero-temperature reservoir.Entanglement generation via this mechanism is not limited to the three spin-domain case, but can also occur between qubits in the outer domains of four and five domain systems as well as entangling three qubits using a large central domain of many spins.
The results in Sec.IV show that, as expected, when realistic effects are taken into account the entanglement generation between spins suffers significantly.With dephasing or individual spin-reservoir coupling, not only is the maximum entanglement generated reduced, the entanglement in the steady state also decays over time.However, we emphasise that for large enough ensemble size, the collective effects will occur on a timescale that is much faster than that of the individual effects.Additionally, finite temperature reservoirs also reduce the maximum entanglement and cause the steady state entanglement to decay.While the results here predict rather small entanglement generation, further optimisation may significantly increase the maxium entanglement of qubits and make it more robust against imperfections.
. Entanglement of formation between single spins in A and C generated via collective spin relaxation (7) for different initial states.The various initial states are labelled in the figure .Here we have used NB = 10.
MXS0118069605 and the Moonshot R&D Program Grant JPMJMS226C.We are grateful for the help and support provided by the Scientific Computing and Data Analysis section of the Research Support Division at OIST.

Appendix A: Alternate initial states
In Sec.III, we showed the dynamics of entanglement generation in our three spin-domain system with the initial state |ψ in ⟩ = |↓⟩ A ⊗|↑ • • • ↑⟩ B ⊗|↓⟩ C where the middle and most populated spin domain B is initialised with all spins in the excited state.This state can also be written in the notation |- . However, this is not the only configuration that can produce entanglement through collective spin relaxation.In the following, we will explore other initial states and their effect on the steady state entanglement produced via this mechanism.
As we consider the collective coupling of many spins to the same reservoir, a necessary condition of the collective coupling is that the spins are indistinguishable from the perspective of the reservoir.Thus we limit ourselves to initial conditions that are permutationally symmetric within each spin ensemble.These states of permutational symmetry are the so-called Dicke states, and apart from the lowest and highest energy Dicke states (ground or fully excited states respectively), the Dicke states have significant entanglement [53,54].In this work, our goal is to start with an initially separable state, and via reservoir engineering, result in the system in a highly entangled steady state.Thus, in exploring alternate initial states, we will limit ourselves to only the separable permutationally symmetric states, which are the fully excited state or ground state of each spin ensemble.
In considering a three spin-ensemble system, this leaves us with two possible initial states for each of the three spin-ensembles, thus leaving eight possible configura-tions.We will briefly go over the outcomes for each of the eight configurations here.Firstly, with the initial state as the ground state of the entire system, this trivially produces no dynamics as it does not decay and thus remains in the fully separable ground state.In Fig. 7, we have plotted the entanglement dynamics for three different intitial states.Again, this is the entanglement of formation between the single spins in domains A and C after the spins in B have been traced out.The dark green, dashed line shows the entanglement for the case studied in the main text of this paper, where domains A and C, of spin populations N A = N C = 1, are initialised in the ground state and the central domain, of spin population N B , is initialised in the excited state.This is indicated by the line labelled On this plot, we also show the entanglement dynamics for the inital condition which is given by the solid orange line.This solid orange line also represents the entanglement dynamics of the equivalent initial condition of Thus, three different initial conditions all produce strikingly similar entanglement levels in the steady state.
The strong entanglement in the steady state reached from the initial condition blue, dotted line in Fig. 7) highlights the novelty of the system considered in this work.To explain this, consider the following: If the three spin ensembles A, B and C were all coupled to a single reservoir, under the collective dissipator J − A + J − B + J − C , the system would relax to the fully separable ground state of the entire system.As can be seen on Fig. 7, the steady state from this initial condition is highly entangled.This means that, as the system relaxes, a significant proportion of the state remains in an entangled dark state and this is purely due to the collective coupling to two zero-temperature reservoirs, that is the simultaneous action of the dissipators produce very little entanglement in the steady state.All three initial states have less than 0.006 ebits of entanglement of formation.Thus, we can say that, for three spin ensembles, of the eight possible initial state configurations, four of these configurations produce significant steady state entaglement.
In Fig. 8, we take a closer look at the entanglement dynamics for the initial state and use this initial state to reproduce the result of Fig. 3(a).As was already observed in Fig. 7, the maximum entanglement generated from intial state (A1) is remarkably similar albeit slightly lower than that generated from initial state (8).This can be seen in Fig. 8 by comparing the dotted, dashed and dash-dotted lines  (8) (that is the result shown in Fig. 3) and these results are given for the same NB which is via the corresponding color.
with their corresponding solid, transparent line of the same color.The system reaches the maximum entanglement slightly faster for initial condition (A1) than (8).Additionally, we find again here, an enhancement in entanglement dynamics, that is faster and more entanglement in the steady state for a larger number of spins N B .However, it is worth emphasising that the number of excitations in the initial state for all lines in Fig. 8 is fixed at 1.This indicates that it is not the collective, superradiant decay of many spins that enhances the entanglement dynamics, but simply the presence of many spins interacting collectively with the same reservoirs.Finally, we ask what is the effect of mixed initial state preparation on the steady state entanglement produced in this scheme.As all other results in this work use pure initial states, it is wothwhile to examine whether imperfect state preparation can be detrimental to the steady state entanglement.To answer this question, we return to the main initial condition examined in this work C and consider the case where the preparation of the central domain B is imperfect and a mixed state is prepared instead.We can model this by the following state that is, a mixture of the target state |↑↑ • • • ↑⟩ and the maximally mixed state, given by the identity matrix I.
We assume that preparation of the ground state can be done with high fidelity, thus the intial state of the total system is given by: with ρ B defined as in (A2).In Fig. 9, we give the steady state entanglement of formation between the single spins in A and C for a given initial state fidelity F 0 .To be precise, F 0 is the fidelity of the initial state ρ B with the tar- On Fig. 9, we can see that while a decrease in initial state fidelity does have a detrimental impact on the steady state entanglement, the scheme still produces significant entanglement even when the initial state is not completely pure. of this particular dark state is always: where N is needed for normalisation and is given by We stress that all numerical simulations (which have been performed for 1 ≥ N B ≥ 12) have yielded the steady state as (B2) with the given dark state (B3).Thus, while we know the dark state for any N B , to fully characterise the steady state we still need to find the dependence of x d on N B , that is the proportion of the state that remains entangled and in the dark state as the system relaxes to its steady state.
To find this dependence of x d on N B we start by examining a simple case of N A = 1, N B = 2 and N C = 1.This system is initialised in the state and will relax to a steady state of the form where for N B = 2, this dark state is: B7) We now define two new states, and, as the initial state (A1) only contains a single excitation, we only need to be concerned with the single excitation subspace of the total system.For convenience, we choose the two states to be the following: In a similar fasion as the example in the introduction, we can now write down the form of the steady state because we know that the first two terms in (B9) will relax to the ground state while the last term is a dark state and does not decay.Thus, the steady state of the total system is: Therefore, for N B = 2, we have x d = 2/5.We are interested in the entanglement between spins in A and C and from (B10), we can trace out the two middle spins belonging to domain B and be left with the reduced denisty matrix of the spins in A and C.This is given by: In fact, the form of the reduced density matrix for this problem is always: While the parameter x d in (B2) describes the proportion of the state that remains in the dark state, the parameter x in (B12) describes the proportion of the reduced state that is in the entangled Bell state, specifically, the triplet state.As shown above, for N B = 2, x = 8/25.This analysis can be repeated for N B > 2 with the crucial step being defining the dark state (B3) and the two other states |ψ 1 ⟩ and |ψ 2 ⟩ and then using such states to express the initial state.While the formula for the dark state is given in (B3), the states |ψ 1 ⟩ and |ψ 2 ⟩ can be generalised for higher spin populations N B in a straightforward way.
For |ψ 1 ⟩, we have all permutations of a single spin up in B |↓⟩ C (B13) while for |ψ 2 ⟩ we have:  know that the dark state component does not decay while the |ψ 1 ⟩ and |ψ 2 ⟩ components decay fully to the ground state.Thus we can find the proportion of the state that remains in the dark state, and subsequently identify the dependence of x d on N B .We shall summarise the findings here: The steady state of the total system is given by (B2), with the dark state defined as in (B3), and with For the reduced system of two single spins in domains A and C, the steady state is defined as in (B12) with To quantify the entanglement of the two qubit reduced state, we again use the concurrence, which for this state (B12) is simply: C = x and thus the concurrence of the two qubit state for a given N B is and the entanglement of formation can be calculated via: In Fig. 10, we show the behavoir of the entanglement of formation for the range of 1 ≤ N B ≤ 10 3 .As the system size N B increases, so too does the proportion of the state that remains in the Bell state.For N B → ∞, this value goes to x ∞ = 1/2.This means that the steady state of the two qubits in this limit is given by |↓↓⟩ ⟨↓↓| , (B19) and the entanglement of formation for this state is E ∞ ≈ 0.354 . (B20)

Figure 1 .
Figure 1.Collective spin relaxation via a collective coupling to a bosonic reservoir can be used to transfer excitations between spins.(a) The double domain system of Ref. [21, 22] (b) Three spin domains pairwise coupled to two reservoirs. entanglement

Figure 2
Figure2.Three spin domains pairwise coupled to two separate reservoirs at zero temperature.Collective relaxation of the spins produces entanglement between the spins.Entanglement between spins in outer domains is maximised for a configuration with a large number of spins in the middle domain.

Figure 3 .
Figure 3. Entanglement generated between spins in A and C using the configuration shown in Fig. 2 .(a) The entanglement of formation EF versus the scaled time γt/2 for various values of NB (the spin-population in domain B) which clearly shows the entanglement generation is higher and faster for large NB.Further, shown in the inset is the normalised collective spin relaxation J z B /NB illustrating the enhanced decay dynamics due to superradiance [40, 41].(b) Steady state entanglement EF and speed of entanglement generation γtE/2 as a function of the spin-population NB.The left axis shows the steady state entanglement of formation EF generated between spins in A and C (green, dotted line).The right axis shows the value of γtE/2 where tE is the time it takes for entanglement between A and C to reach half of its maximum value (purple, dashed line).

FFigure 4 .
Figure 4.The entanglement of formation EF generated between spins in the outer domain of a spin domain chain versus the scaled time γt/2.The pink (dashed) line shows EF for spins in the outer domains of a four spin domain chain connected with three reservoirs, depicted in the inset shown in the upper right corner.The spin domain population configuration for this result is {NA, NB, NC , ND} = {1, 6, 6, 1}.The blue (dotted) line shows the EF for spins in the outer domains of a five spin domain chain with four reservoirs, depicted in the inset diagram below the blue (dashed) line.The spin domain population configuration for this result is {NA, NB, NC , ND, NE} = {1, 4, 4, 4, 1}.

Figure 5 .
Figure 5.The entanglement of formation EF generated between spins in domains A and C versus the scaled time γt/2 with a spin population of NB = 11 in domain B. In (a), we show various dephasing rates γ dep .The pink (long dashed) line shows the effect of dephasing with γ dep /γ = 0.2, the blue (dotted) line uses γ dep /γ = 0.1, the orange (dash-dot) line uses γ dep /γ = 0.05, the green (dashed) line uses γ dep /γ = 0.02, and the purple (solid) line shows the ideal case of γ dep /γ = 0.In (b), we show the entanglement dynamics with both collective and individual decay.In (c), we use a thermal reservoir of temperature T with the following parameters chosen: ω0/2π = 10 GHz and γ = 0.01 Hz.

Figure 6 .
Figure 6.Tripartite quantum correlations generated from collective decay with multiple reservoirs.(a) Initially, all spins are in a separable state (12) with their configuration depicted here.After collective spin relaxation, the reduced state of the three single qubits in spin-domains A, B and C are entangled in the steady state.(b) Tripartite negativity generated between spins A, B and C during spin relaxation with the configuration pictured in (a) for a central domain populated with ND = 11 spins.Shown in the inset is the steady state tripartite negativity for spin population ND ranging from 1 to 11.

Figure 9 .
Figure 9. Steady state entanglement of formation between single spins in A and C for a mixed initial state.The initial state used here is given in (A3).Shown here is the steady state entanglement for a range of spin-ensemble B population 2 ≤ NB ≤ 10 (indicated via color) and initial state fidelities F0.
B14)With such definitions, one needs to express the initial state in terms of |ψ d ⟩, |ψ 1 ⟩ and |ψ 2 ⟩ and from this we

FFigure 10 .
Figure 10.Steady state entanglement of formation between single spins in A and C for the initial condition (A1) evolving under master equation(7) .The purple diamonds indicate data points acquired by numerically solving the master equation which can be done for 1 ≤ NB ≤ 12, while the blue dashed line is (B18) with the concurrence (B17).
Figure 8. Entanglement of formation generated via collective spin relaxation using the initial state where the single spin in domain A is in its excited state while the spins in domains B and C are in their ground state.This state given in (A1) is depicted in the inset.The generated entanglement of formation EF versus scaled time γt/2 is given for NB = 2, 7, 12.We clearly observe an increase in EF increases as NB increases.The transparent, solid lines show the entanglement dynamics of the initial state (7)These states were chosen as it can be easily verified that a system starting in initial state |ψ 1 ⟩ or |ψ 2 ⟩ and evolving under the master equation(7), will relax to the ground state of the total system |↓↓↓↓⟩.Note that the states |ψ need |ψ 1 ⟩, |ψ 2 ⟩ and |ψ d ⟩ as the initial state (B5) can be written in terms of these states: 1 ⟩, |ψ 2 ⟩ and |ψ d ⟩ are all mutually orthogonal and if we wanted to define a basis to span the entire single excitation subspace, the remaining basis state would be |ψ d2 ⟩ = 1/ √ 2 (|↓↑↓↓⟩ − |↓↓↑↓⟩) which is also a dark state under both dissipators J − A + J − B or J − B + J − C .However, we only