Dissipative cooling of spin chains by a bath of dipolar particles

We consider a spin chain of fermionic atoms in an optical lattice, interacting with each other by super-exchange interactions. We theoretically investigate the dissipative evolution of the spin chain when it is coupled by magnetic dipole-dipole interaction to a bath consisting of atoms with a strong magnetic moment. Dipolar interactions with the bath allow for a dynamical evolution of the collective spin of the spin chain. Starting from an uncorrelated thermal sample, we demonstrate that the dissipative cooling produces highly entangled low energy spin states of the chain in a timescale of a few seconds. In practice, the lowest energy singlet state driven by super-exchange interactions is efficiently produced. This dissipative approach is a promising alternative to cool spin-full atoms in spin-independent lattices. It provides direct thermalization of the spin degrees of freedom, while traditional approaches are plagued by the inherently long timescale associated to the necessary spatial redistribution of spins under the effect of super-exchange interactions.


Introduction
The coupling of quantum systems to environments can lead to various consequences: typically it produces decoherence towards classical behaviors [1], but specific situations are in the spotlight that would on the contrary produce quantum correlations [2,3,4]. This is an exciting new paradigm, and also opens the fascinating perspective of "environment engineering" to protect or produce entangled states, for quantum information processing and for quantum simulation. Subjecting the various many-body systems explored with cold atoms to wellcontrolled dissipative processes is a promising prospect to explore the diversity of these situations [5,6]. Among the many-body models implemented on ultracold atoms and molecules, lattice spin models have attracted a large attention [7,8,9,10,11,12,13,14,15,16,17,18]. The question of how dissipation affects them and whether it can be tuned to specific purposes emerges spontaneouly in the context of excited Rydberg atoms [19,20,21,22]. It also appears to be an important prospect for ground state atoms in optical lattices [3,23,24,25,26,27,28,29,30], and starts to be experimentally explored [6].
Recent proposals [31,23,24,27] involve the use of light as a bath, and spontaneous emission as the dissipative process. In the present work, we explore binary atomic mixtures [3,28,32], one species acting as many-body quantum system, the other as bath. Namely, the quantum system is a spin chain of fermionic atoms, and it is coupled to a Bose Einstein condensate of a different species. The low-energy many-body states of the spin chain are driven by nearest-neighbour super-exchange interactions. Magnetic dipole interaction between fermions and bath lead to spin flips in the chain, associated with spontaneous phonon emission in the BEC. Thus, spin-thermalization can arise, due to the spin-orbit coupling conveyed by dipole-dipole interactions, an effect which is connected to the Einstein-de Haas effect [33]. As we show here, spin-orbit coupling offers a possibility to directly cool the collective spin degrees of freedom in a spin-chain. Such a collective coupling of the spin chain to phonons in the BEC can be seen as an analog of superradience in quantum optics, with the spontaneous collective emission of phonons instead of photons. We point out that the possibility to couple to the total spin of the chain is particularly relevant to the quantum simulation of the spin-full lattice Hubbard model with ultra-cold atoms. Indeed, in most experiments the collective spin is a conserved quantity which is typically not under control. This can become a limitation to reach the many-body ground state which has a singlet character at half-filling [34]. Our work suggests for the first time to use the intrinsic spin-orbit coupling provided by dipolar interactions, in order to dynamically couple to the collective spin length of a fermionic spin chain and allow reaching its many-body ground state. This approach is complementary to previous studies which used dipolar interactions to cool the mechanical degrees of freedom of a thermal or quantum gas taking benefit of their spin degrees of freedom [35,36].
Formally, as presented in sections 2 and 3, we consider a fermionic gas with two spin states in an optical lattice. Due to repulsive contact interactions a Mott insulating regime is reached, in which density fluctuations are strongly reduced [37,38]. At half filling, the system then consists of a spin chain, described by the Heisenberg spin-spin Hamiltonian. a) The spin chain is composed of fermions, each of which has two quasi-degenerate spin states. ∆ is the difference of single-particle energy between the two fermionic spin states. They interact via dipolar interactions with a bath of bosons with a strong magnetic dipole moment, polarized in a magnetic field. Provided ∆ k B T g B µ B B, inter-species dipolar interactions free the magnetization of the fermions, but energy conservation maintains the magnetization of the bosons. b) The dipolar interaction enables changes in the collective spin chain state, such as changing neighbour correlations from triplet to singlet. For a chain of fermions with spin super-exchange interaction energies ∼ J > k B T , the spin chain can lower its energy by exciting a density mode of the bath.
To complement this spin-chain Hamiltonian, we consider dipole-dipole interactions between atoms in the spin chain and a spatially overlapping BEC. We treat the interaction with this bath using the Fermi golden rule approach. For quantitative predictions, we consider in section 4 a 1D spin chain made of alkali atoms (for which inter-atomic dipole-dipole interactions can safely be neglected compared to super-exchange energies), which interact via dipole-dipole interactions with a bath of a highly dipolar species such as Dy, Er, or Cr [39,40,41,42,43,44]. We find that for realistic experimental parameters, the spin-chain thermalizes close to its highly entangled (singlet) ground state in a few seconds. The final spin-chain temperature is close to the initial temperature of the BEC. Our simulations indicate that the thermalization rate is roughly independent of the length of the spin chain.

General idea
This work focuses on 1D chains of interacting spinful fermionic atoms, prepared in optical lattices. We consider a system with effective spin 1/2 particles -that can be realized by energetically selecting two Zeeman states of an alkali atom in its ground state. These fermions occupy the ground Bloch band of the lattice. Dipole-dipole interactions between fermions are assumed to be negligible. Fermionic Mott insulators are observed when the on-site interaction U between spin states dominates as compared to the Bloch band width, proportional to the hopping term t, and when the temperature k B T U [45,46,37,38,47]. There, for one fermion per site the evolution is well depicted by an effective next-neighbor magnetic interaction Hamiltonian between localized spins [48]. In 1D, with J = (−U + √ 16t 2 + U 2 )/2 4t 2 /U . This involves effective spin half operators Σ ef f , that do not generically coincide with the true spin operators Σ of the atom. For example, if the fermionic species is 40 K, its true spin length is 9/2. Two selected states could be |m F = −9/2 , |m F = −7/2 . Σ ef f are the Pauli spin-1/2 matrices, for example Σ z ef f | − 9/2 = − 1 2 | − 9/2 while Σ z | − 9/2 = − 9 2 | − 9/2 . In this paper we will describe 1D spin chains with lengths up to 7 sites. The eigenstates are generically strongly entangled, favoring next-neighbor singlet correlations at low energy and triplet correlations at high energy. For example, if we consider only two sites, the singlet state |S is the ground state, with energy −3J/4. The three triplet states |T 0 , |T +1 , |T −1 with pseudo-spin projection 0, ±1, share the same energy J/4.
To open up the dissipative evolution between these spin eigenstates, a bath of strongly dipolar bosons is introduced (such as Dy, Er, Cr), which interacts with the fermion spins. Dipolar interactions between one bosonic atom spin S at position r and one fermionic atom spin Σ at position r + dr, with respective Landé factors g B and g F , are anisotropic: using dimensionless (i.e., divided by ) spin operators. Here µ 0 is the vacuum magnetic permeability and µ B the Bohr magneton. This interaction involves the true spin operators ( Σ, S) of both species. This expression can be developed into the various products of the spin operators components Σ ±,z S ±,z , defining S ± = S x ± iS y and Σ ± = Σ x ± iΣ y . Two families of processes appear. First, relaxation and exchange terms that change the bath magnetization, Σ ±,z S ± . Here we consider a bath initially polarized by a magnetic field in the stretched state that minimizes the Zeeman energy. As illustrated in fig. 1, for a strong enough field (g B µ B B k B T, J), these processes will be either forbidden or negligible due to energy conservation. Second, relaxation and Ising terms that conserve the bath magnetization: Σ ±,z S z . These terms are the ones driving the dissipative evolution of the spin chain. We will describe their effect using the Fermi golden rule, to evaluate the rate Γ i→f of transfer between any two collective spin eigenstates |i spin , |f spin , associated with energy transfer to the bath density modes. Writing p i (p f ) the probability that the spin chain occupies an eigenstate |i spin (|f spin ), all occupation probabilities will evolve following a set of equations dp where generally Γ i→f = Γ f →i . We now will evaluate the rate associated with each possible i → f process, based on the strength of the inter-species coupling H int and the bath density of states. Figure 2. System geometry. The quantization axis is defined by the magnetic field along z, and the spin chain is oriented along z . The emission of one bath excitation with wavevector k, with spherical angle coordinates (θ, φ), acts on each spin in the chain with amplitude V dd ( k) exp(−i k · R j ). The angular dependence of V dd ( k) is illustrated as blue surfaces around each of the two spins (here for the Σ z S z component of V dd ): the amplitude of V dd ( k) is given by the distance between surface and origin, along the direction of k. Due to the excitation phases exp(−i k · R j ), emission can change the chain spin correlations, except when k ⊥ z .

Laying out the Fermi golden rule
The geometry of our system is as follows. The two atomic species are trapped by the same lattice laser beams, but see a different potential depth due to their different electronic structures. We will only consider situations where the fermions are effectively arranged into disconnected 1D chains, mostly out of computational time limits, while the bosons retain 3D coherence (see section 4.3). The formalism is nevertheless generic, and dissipation imposed on 2D or 3D lattice structures would deserve further scrutiny. The quantization axis z is defined along the magnetic field. The spin chain axis z can differ from z (see fig. 2). Several restrictions are made for the sake of simplicity, that are not fundamental requirements. i) We consider that the external magnetic field is strong enough to maintain the bath polarization: g B µ B B > k B T, J. ii) We restrict ourselves to pseudo-spin 1/2 chains, where only two neighboring Zeeman states of the fermions are ever populated, and assume that the Zeeman energy difference between these two states ∆ can be made negligible, ∆ J, or irrelevant (see in section 4.2 considerations on how to impose this). As J U this implies ∆ U : all single-particle energies are smaller than U , such that with the proper chemical potential the occurrence of holes and doublons remains negligible.
We now describe the interaction between the dipolar bosons and the spin chain. In the fermionic Mott regime for large enough U/t, with one atom per site, we make the simplified assumption that k B T U ensures negligible residual number fluctuations [46,47,49, 50] (perfect Fock state on each site). Under this approximation the density and spin operators for the fermions factorize, and we can label fermions by their site. The density distribution of the j-th fermion n j F ( r) is constant, set by the ground band Wannier function on site j. For the bosons, as the BEC is energetically constrained to the polarized state, spin and density also factorize. This enables us to write the coupling between spin chain and condensate as: Here,n B is the bath density operator in second quantization. V j dd is the interaction V dd of eq. 2, acting on the fermionic spin Σ j of the fermion on the j-th site under the influence of a boson of spin S at position r. It includes the added constraint of fixed BEC magnetization: the S + and S − terms are neglected, and the spin operator S z is replaced by the spin length S of this species. We use Fourier transformation and write: As in dipolar relaxation, total intrinsic spin is not conserved. The magnetization changing collisions are tied to the spin orbit coupling inherent to dipolar interactions [51]. We now describe the evolution of the spin chain under the assumption that the BEC behaves as a bath. Proper justification of this will be given once the formalism is entirely laid out, in section 3.6. The processes that transfer energy from the spin chain to the BEC connect an initial spin state |i spin to a final spin state |f spin with energy difference Starting from an initial bath state |i bath , the corresponding rates are: This equation involves a summation over all possible final bath states |f bath , with bath energy difference E bath if . We express all equations for a finite box of volume V with periodic boundary conditions, so that this summation is discrete. Regarding the bath, we start our derivation in the zero-temperature limit, only to generalize to finite T in the later section 3.5. Quantum depletion is neglected (see section 4.3). Thus, the initial bath state is for now described by a coherent state |BEC with average atom number N 0 .

Matrix elements
We now derive an explicit expression of the coupling element of eq. 7. Using eqs 4 and 5, we have: There, the corresponding physical process is a transfer of momentum k 1 from the fermions to the bosons. The lattice potential constrains the fermions to their unchanged Wannier states, whose finite momentum widths allow for this momentum transfer : is the ground band Wannier function centered on the site coordinate R j . We calculate the bath matrix element using the discrete basis of plane waves in cubic volume V with periodic boundary conditions: , we obtain: This expression explicitly relates the change in collective spin state to the interference of phonon radiation patterns from the different fermions, as an operator acting collectively on all the spins of the chain appears: (see fig. 2). The operator V col defines selection rules regarding which collective spin states are coupled. For these rules, we diagonalize the effective Hamiltonian H mag together with the total pseudo-spin length operator Σ 2 ef f,tot = ( j Σ ef f,j ) 2 and with both the total spin and pseudo-spin projection operators, e.g. Σ z tot = j Σ z j . In this eigenstate basis, if |i spin and |f spin differ by up to one in pseudo-spin length Σ ef f,tot and 0 in projection Σ z tot , and the corresponding two-body operator is Σ z S z . If |i spin and |f spin differ by up to one in pseudo-spin length Σ ef f,tot and ±1 in projection Σ z tot , and the corresponding two-body operator is Σ ± S z . Otherwise, the coupling vanishes.

Bath description
The coupling term (eq. 9) is sensitive to the potentials seen by the bosons. Here, we assume that the bosons in the lattice retain large coherence over the 3D sample, which we check independently (see section 4.3). We develop accordingly the BEC matrix element in equation 9 by transforming annihilation and creation operators from the plane wave basis ψ( q) to the Bloch wave basis φ( q), and then to the Bogoliubov excitation basis b( q) and BEC stateĉ BEC . The first transformation, between plane waves and Bloch waves, is based on the decomposition of the Bloch band structure of our lattice, of axes (x , y , z ). We assume a cubic lattice with same period along all three axes, and reciprocal lattice wavevector length K. Then, where q = q x e x + q y e y + q z e z is in the first Brillouin zone, and we define the notation n = n x e x + n y e y + n z e z . The lattice can have different depths along the three axes, such that the a nx , a ny , a nz implicitly differ. We omit in our notations the band indices: in this work only the ground band is ever populated given the relevant energies. We generalize the notation: The second step is the Bogoliubov transformation [52]: where δ K is a Kronecker symbol. The prefactors (u( q), v * (− q)) depend on kinetic and interaction terms (E q , V q respectively) that can be strongly affected by the lattice potential: E q is the non-interacting energy for wavevector q of the Bloch wave. The interaction term V q depends on wavevector due to the dipolar interaction within the polarized dipolar condensate.
In the absence of a lattice potentials for bosons, . We generalize this in the presence of the lattice potential seen by the bosons, using a tight binding approximation: we introduce boson Wannier functions We obtain: whereg is the bath contact interaction parameter renormalized by the lattice. For a 3D lattice confinement such that the bosons experience on-site interaction U B , this interaction parameter readsg = U B (λ/2) 3 , with λ/2 the lattice spacing. The Bogoliubov quasi-particles have an We finally evaluate the coupling element eq. 9, using the basis transformations of bosonic operators, equations 13 and 14. The coupling Hamiltonian H int is of second order in boson field operators, which we just decomposed between Bogoliubov modes and BEC state mode. Consequently, the initial and final bath states can only be coupled if they differ by zero, one or two Bogoliubov annihilation and creation operators. The dominant terms for dissipation are those differing by one excitation, which also contain one BEC operator with macroscopic expectation value √ N 0 . Thus, for T = 0 and neglecting quantum depletion, the state |f bath only differs from |BEC by the creation of a Bogoliubov excitation, parametrized by its wavevector q: We assume the BEC to occupy the zero-momentum Bloch wave state. Then, we get without explicit tight binding approximation The tight binding approximation may also be applied along some dimensions to simplify eq. 18. Applying it to all three dimensions leads to From equation 19 we can see that the concentration of bosons onto Wannier functions by the lattice has an impact on the inter-species dipole coupling ‡, despite the long-range nature of the dipolar interaction. Due to the anisotropic nature of V col , the lattice can be beneficial or harmful to the coupling strength, as the reciprocal space summation n can involve terms of opposite signs, depending on the weights imposed by the Wannier functions.

Zero-temperature rates
At this stage we conclude our calculation for T = 0. To reveal explicitly the role of box volume, we replace the discrete summation on |f bath ( q) by a continuous one on q.
where the 3D Dirac function results in an integration on the 2D anisotropic wavevector surface q(θ, φ, E if ) defined by conservation of energy: ( q) = E if . In this paper, these energies will always lie in the ground band for the bosons.

Finite temperature rates
At finite temperature, the initial bath state is not a pure state anymore. The theory above can then be generalised. We use a bath eigenstate basis |{n q } = q (b † q ) nq |BEC consisting of Fock states of the Bogoliubov modes labelled q and of energy q . The statistical probability of any such eigenstate |{n q } init is The rate between two spin chain states |i spin , |f spin with respective energies differing by As explained in section 3.3, the dominant terms in the decomposition of H int are those involving one Bogoliubov operator and the BEC mode.
If the spin energy diminishes (E f < E i ), the most contributing |{n q } final differ from |{n q } init by one added excitation (phonon emission). Bose enhancement factors appear as excitation states are not necessarily initially empty: n init q ≥ 0. After statistical averaging over |{n q } init , as explicit in eq. 22, the rate Γ i→f is modified by finite temperature T through an additional factor: where n q = 1 e E if /k B T −1 depends on the spin energy released in the phonon. If the spin energy increases (E f > E i ), the most contributing |{n q } final differ from |{n q } init by one less excitation. The rate Γ i→f , associated to phonon absorption, is now non-zero. It scales with the inverse (phonon emission) process as: This assymetry between the processes pumping energy out (eq. 23) or in (eq. 24) the spin chain leads to thermalization between the spin chain collective degrees of freedom and the bath motional temperature. This is of particular interest for a sufficiently cold BEC, k B T J, such that n q 1.

Validity of the bath approximation and Fermi Golden rule
We conclude this section by defining properly the initial assumption of treating the BEC as a bath, (i) with which energy exchange processes are irreversible (justifying the Fermi golden rule), and (ii) such that its state is negligibly affected by interaction with the chain.
(i) For a given process i → f , dephasing of the populated phonon modes is required to prevent coherent oscillations of the energy between chain and bath. As for spontaneous emission of light, this is satisfied even for an isolated BEC provided the number N modes of excitation modes with energies close enough from E if is much larger than one. We estimate N modes by integrating the shell of wavevectors corresponding to excitations resonant with E if , within an energy bandwidth given by the coupling strength: This amounts to a product of density of state and coupling element, taking into account the anisotropy of the coupling. For fixed BEC mean density N 0 /V , matrix elements of H int scale as 1/ √ V , such that N modes scales as √ V and incoherent dynamics is always justified at the thermodynamic limit. Once in the incoherent regime, the rate is independent on volume. For small clouds that would not satisfy this density-of-state criterion, irreversibility requires dissipation or damping of the phonons by other means, either inherent to finite temperature BECs [53,54,55,56], or external such as evaporation. The latter is experimentally relevant, as BEC temperatures lower than J are typically obtained in traps with low depths.
(ii) Furthermore, throughout the numerical results, the bath excitations will be assumed to be thermally distributed at a fixed temperature. Precisely, we assume that the emissions and absorptions of phonons affect negligibly the populations n q . For a strongly damped bath, this may be ensured by a sufficient heat capacity. For weakly damped baths as most BECs, this will be ensured either by evaporation, or by N modes being larger than the number of phonons emitted within the same energy band. The latter criterion improves with the ratio of boson atom number to fermion atom number.

Spin chains in K-Dy mixtures
Using one of the most magnetic dipolar species (Cr, Er, Dy) as the bath, dipolar interactions will drive dissipation even on spin chains from the weakly dipolar alkali species. More exotic situations with actual dipole-dipole interactions within the spin chain could be studied, e.g. using an isotope of the strongly dipolar species for the chain as well.
In the prospect of studies of the t-J model, constituting spin chains from fermionic species in the Mott regime is very attractive. Lithium or Potassium are the most prominent in recent experiments, together with alkaline-earth species that however have much too weak Landé factors in the ground state for our present purpose. We will here discuss 40 K in its ground hyperfine state, with spin F = 9/2, and Landé factor g F = 2/9. Fermionic Lithium would be more problematic due to a combination of complications §. A key point is the optical lattice potential seen by the bath species. It generally differs from that seen by the spin chain species. As mentioned earlier (section 3.3), a lattice of opposite sign may often lead to poor interspecies dipolar coupling, concentrating the bosons half-a-space away from the fermions. For Dy or Er in association with K spins, wavelengths > 767 nm would generate lattices of the same sign for both species. However larger J and large interspecies coupling can be obtained at shorter visible wavelengths, using the reddetuned vicinity of "narrow" lines that would dominate over the effect of the wide lines at 421 and 401 nm respectively for Dy and Er. We suggest for example one narrow line of Dy at 626.1 nm (140 kHz) and one narrow line of Er at 582.8 nm (190 kHz) (see figure  3). The intervals 622.5 -624.3 nm for Dy and 578.2 -580.2 nm for Er are attractive, as in these intervals the boson lattice depths can be tuned from 0 to ∼0.4 times the fermion lattice depth while light scattering remains dominated by the fermions. Here we will use 164 Dy, which has the strongest magnetic dipole moment µ = g B Sµ B = 10µ B , but the shorter lattice wavelengths suggested for Er translates in higher super-exchange and fairly similar spin cooling rates despite the lower moment 7µ B .
In the following, we will consider 40 K in the presence of 164 Dy [58], in an optical lattice with wavelength λ L ≈ 623 nm. We effectively create a 2D array of 1D spin chains by using a lattice, as seen by K, of depth 25 E K r along two axes x , y (where E K r = h 2 /(2m K λ 2 L ) is the potassium recoil energy), and 3.5 E K r along one axis z . The spontaneous light emission rate is ∼ 0.2 /s. The on-site interaction energy U = 7.5 t z , where t z is the tunnelling energy along the weak axis z . Fermions can then be prepared in the 1D Mott regime with unit filling [45,59,60]. The super-exchange energy is, along the weak axis, J = h × 632 Hz, while it is along the transverse axes h× 0.07 Hz, thus effectively creating decoupled 1D spin chains. The lattice depth for the bosons, on the other hand, sensitively depends on the lattice wavelength within a few nanometers, and will be used as an independent optimization variable in section combinations differing by their electronic spin state, relaxation of the large Zeeman energy would need to be prevented [51, 57].

Hierarchy of energy scales
The actual implementation requires to respect a certain hierarchy of energy scales, involving the magnetic field B, the BEC temperature T , and the superexchange energy of fermions J (see fig. 1). First, the bath temperature T < J/k B 30 nK is such that low energy spin chain states are favored by coupling to the bath. Second, (J, k B T ) < g B µ B B ensures that the bath remains polarized. This is largely satisfied with fields in the mG range. Although the scheme might work also for a BEC with free magnetization, this study is beyond the scope of our calculations. Third, to study antiferromagnetic models the fermion Zeeman shift ∆ should not favor trivially polarized spin chains. There are several ways to ensure this, with different consequences.
i) The fermionic Zeeman states of interest may be kept quasi-degenerate by g F µ B B J. Combined with our prior assumption of a polarized BEC, this implies g F g B , which is not well satisfied for most of the atomic species combinations. To circumvent this, we propose to add microwave dressing to the spin states [64]: the linear Zeeman effect can be combined with an engineered quadratic shift αm 2 F , such that the detuning between the two lowest Zeeman states is To emulate pseudo-spin 1/2 dynamics in a larger spin species, one further requires E m F =−F +1 − E m F =−F +2 > J, k B T to prevent population of the other Zeeman states. This implies B > J g F µ B 2F −1 2 8 mG. This protocol, and the effect of a residual ∆, are discussed further in Appendix A. The calculations below assume such a setting. Varying ∆ would enable studies of the Heisenberg anti-ferromagnetic model with free magnetization for varying bias field, complementing the more usual fixed magnetization setting in experiments.
ii) Another possibility is to ensure that the energy released by magnetization-changing collisions, ∆, lies in the energy gaps within the boson band structure. Dipolar relaxation is then prevented by energy conservation, as seen in [51,57]. The thermalization of the spin chain would then occur at fixed magnetization (i.e. fixed spin projection along B). This approach requires a much lower degree of field control, and in particular releases the constraint that ∆ J. Cooling of the chain is slightly less efficient, as the magnetizationchanging processes are suppressed. Nevertheless, processes given by eq. 11 are still at play and may thermalize the collective spin at fixed magnetization. The magnetization is set initially by creating the appropriate spin imbalance using standard manipulation techniques. This is within the bounds of state-of-the-art BECs in optical lattices. Already in bulk, 1 nK is achieved in experiments [61]. The lowest temperatures are achieved using shallow traps and eventually adiabatic decompression, which is consistent with the low BEC density specified in the later section 4.3. Furthermore, loading optical lattices can significantly reduce BEC temperatures [62]. The work on thermometry of BECs in lattices [63] shows temperatures down to 3 nK.

Optimizing the bath lattice potential
The lattice potential experienced by the bosons has a crucial effect on phonon emission, strongly tied to the anisotropic character of the dipolar interaction. We first illustrate the importance of this anisotropy without the lattice potential. Figure 4 presents phonon emission diagrams for different orientations of a 2-fermions chain in the triplet states. We make two observations on the emission directionality. (a) It is strongly spin dependent, as can be inferred by equations 11 and 12. (b) It is suppressed orthogonally to the spin chain, as the phase terms e −i q· R j in V col are close to 1. Indeed, V col then produces an identical rotation of all spins in the chain, which does not modify the collective spin and therefore is unable to drive any dissipative dynamics. ¶ The bath lattice potential introduces the following effects. (i) It can strongly influence the anisotropy of the phonon dispersion relations, and thus of the resonant phonon wavevector surface (see figure 5). This affects both the density of states and the phase terms e −i q· R j in V col discussed above. (ii) It can generate either constructive or destructive interference through n [..]V col ( q + nK) terms of eqs. 18 and 19, as the dipolar interaction changes sign depending on the orientation of q + nK.
The lattice potential depth for bosons is thus an important control tool, with two main functions: i) tuning the boson dispersion relation, to explore high-density-of-states, and also to stabilize against the bath dynamical instabilities of strongly magnetic species [67,68] (for Dy, dd = µ 0 µ 2 m/12π 2 a sc 1.4); ii) maximizing the constructive interference and thus enhancing the coupling. Optionally, as mentioned in section 4.2, band gaps may also be used to enforce fixed magnetization dynamics with non-zero fermion Zeeman shifts.
In systematic studies of the cooling rates as a function of system parameters, for θ z−z = 0 the rates typically were more sensitive than for θ z−z = 0 to the other parameters, including to the released energy E if . Looking for conclusions that are robust with chain length, we will focus on θ z−z = 0. We find that in such a case where the spin chain is oriented along the magnetic field axis, the lattice transverse to B mostly contributes to enhancing magnetizationconserving couplings, while the longitudinal lattice mostly contributes in stabilizing the dipolar gas, but quickly becomes detrimental to the couplings. Deeper transverse lattices are preferable until one has to consider the deconfinement transition from 1D tubes or 1D Mott states to 3D BEC [69,70,71,72], and the quantum depletion that arises before it [73].
Fully independent optimization of longitudinal (z') and transverse (x',y') boson lattice depths is possible using two laser sources, given the strong wavelength dependence of the boson light shift + . We set our parameters to k B T = 0.3J, a mean density of ρ 0 = 3 · 10 19 /m 3 , ¶ The interaction Hamiltonian may also drive collective precession of the chain under the influence of the magnetic field radiated by the bath dipoles, B int . For a finite size condensate with spherical symmetry, B int = 0 at the center but is not spatially homogeneous. Orders of magnitude for the Larmor frequency of off-centered chains may lie in the 10 Hz range [65,66]. The chain eigenstates with same total spin length but different z projection are coupled and degeneracies are lifted. Considering chains much smaller than the condensate, our expression for rates Γ i→j in this paper remains valid in the new (local) eigenstate basis. The numerical results shown are for a centered chain, with no degeneracy lifting. + Dynamics very similar to those presented below, section 4.4, can also be obtained using a single wavelength.   i.e. a filling of 0.9 bosonic atoms per 3D site, in a transverse lattice depth 12E Dy r , and a longitudinal lattice depth 3.5E Dy r (stabilization against dynamical instabilities operates from ≈ 2.5E Dy r ). The first gap of the band structure is h × 5.2 kHz wide. Based on 1D physics considerations neglecting the longitudinal lattice, the 1D Lieb-Liniger parameter γ = mg 1d / 2 n 1d ≈ 1.3, while the transverse tunnelling is t perp ≈ 6 · 10 −2 µ 1d , which according to [71] and consistently with the experiment [72] would ensure being deep in the 3D coherent regime. From the 3D lattice perspective, we have U/(2t x + 2t y + 2t z ) 0.9, far from the critical value 5.8 where the Mott transition could be reached. Following [73], we estimate quantum depletion based on Bogoliubov theory in the homogeneous system limit to be N qd /N 0 = 1/N 0 d 3 q V /(2π) 3 v( q) 2 5%. The transverse lattice depth is thus kept conservatively low to ensure the BEC 3D coherence.
In this setting, the surface of resonant phonon wavevectors is reduced to two sheets, orthogonal to z (figure 5). Due to the diverging density of states at the edges of the Brillouin zone, the phonon radiation pattern is heavily affected, as visible in fig. 6. As the lattice is more confining in the (x , y ) directions, the terms of n [...]F[|w B | 2 ](− q − nK)V col ( q + nK) in equation 19 bear most weight for n in the (x , y ) plane, in which V col has constant sign. Thus, the summation is constructive.
Finally, we estimate with these settings that for the representative process |T 0 → |S in a length-2 chain, the number of contributing phonon modes is, for V = (20 µm) 3 , N modes 1.7, and for V = (100 µm) 3 , N modes 19. The incoherent description is thus justified without the need to consider evaporation for large but realistic BECs.

Spin chain evolution
We now present the dissipative evolution of spin chains with the above parameters. To get a hint at how our simple picture scales towards the large system limit, we will repeat our calculation for varying 1D chain lengths L. Although our theory is valid for arbitrary chain length, we are limited by the computation time of the radiation pattern and anisotropic density of states, which have to be evaluated for all transition energies between spin chain eigenstates. The evolutions shown below rely on a three step approach. (i) Find all spin chain eigenstates, in a basis that also diagonalizes total pseudo-spin length and total spin projection along the polarizing field. (ii) Compute the rates between all of these states -this being the computationally limiting part as L is increased. We use for efficiency a tight binding approximation in the H int matrix element along the two strong axes (x , y ) (see Appendix B), and along all axes for the V q and thus the Bogoliubov description of the bath. (iii) Given all rates, simulate the dissipative evolution from a statistical distribution of the eigenstate probabilities.
First, we show in figure 7 the evolution with time of the probabilities of the four eigenstates from a length-2 spin chain, for various angles θ z−z of the chain axis with respect to B. The cooling dynamics shows very strong sensitivity with the alignment, as a consequence of the anisotropy of the dipolar interaction. For parallel alignment (θ z−z = 0), the fixedmagnetization transition |T 0 → |S has the fastest rate, with timescale ≈ 0.9 s. As will be demonstrated later, preparing the initial chain at total magnetization close to 0 (i.e., starting from a balanced spin mixture in the chain) may then provide the fastest dissipative route to low energy antiferromagnetic chains.
We demonstrate in Fig. 8 that the spin chain reaches thermal equilibrium with the bath. The plot shows the occupation probabilities of the various spin chain eigenstates after several evolution times, for a chain of length L = 7. We start the evolution from a state with minimal magnetization, i.e. all states with total magnetization ±1/2 have same probability, while all the other states have zero initial probability. This corresponds to a state of infinite spin temperature with the constraint of ±1/2 magnetization. We observe that the evolution converges towards a distribution of probabilities exponential with respect to energy, i.e., a thermal distribution in the many-body states at free magnetization. We furthermore apply an exponential fit to the equilibrated distribution and extract a spin temperature of k B T S = 0.3J which corresponds to the bath temperature. Two seconds of evolution already lead close to a thermal state at spin temperature k B T S 0.45J.
To extend these observations to longer chains, we define global, "thermodynamic" quantities describing the spin chain: the magnetic energy E, rescaled by its initial value, and the global von Neumann spin entropy S = − i p i log(p i ), where p i are the probabilities of the collective eigenstates. In fig. 9, we plot these quantities as a function of time, for chain lengths L ranging from 2 to 7. We actually compare the evolution from samples prepared at minimal initial magnetization, 0 or ±1/2 depending on the chain length parity (left column), and prepared with free initial magnetization (right column). The dissipative evolution itself frees the magnetization, which is responsible for the initial entropy increase for samples prepared at fixed magnetization. Fig. 9(a-b) and (d-e) show that the evolutions of energy and entropy for increasingly large chains tend to converge with L, for both preparation protocols. The equilibration rates Γ equil = d/dt ln E − E(t → ∞) for both preparation protocols ( fig. 9c,f), evaluated at early times, seem to quickly converge to similar and stable (size- independent) values. These observations suggest that long, mesoscopic chains may evolve with similar equilibration rates.
We observe in these plots, already after 2 s, entropies per particle and spin temperatures that compare very favourably with the most recent experiments, e.g. [60,16]. Ultimately, provided the equilibration rate is sufficient to overcome imperfections such as light scattering, a key point is that the outcome is tied to the attainable BEC temperature; and indeed, as mentioned in section 4.2, extremely cold BECs can be produced experimentally, in bulk and in optical lattices [61,62,63]. It thus decouples from the usual challenge to produce low entropy spinful fermions in lattices, tackled e.g. by [74,75].
To conclude, we show in Fig. 10 spin correlations for the largest spin chain (L = 7) before and after the evolution, and compare them to correlations at zero temperature. We use the following correlator definition, normalised to one for full correlation:  The flat negative correlation of the initial state is a finite size bias from the constrained initial magnetization ±1/2. For example, if the first site has magnetization −1/2, the L − 1 other sites have a positive average magnetization per site ∝ 1/(L − 1) that reflects upon the correlator. At the final temperature of 0.3J we observe alternating spin correlations already similar to the ground state, spanning over the whole sample. This correlator is accessible on experiments, especially using quantum gas microscope techniques discriminating two spin states, as in [76,77].

Conclusion and possible extensions
In this paper, we have discussed the possibility to couple the spin degrees of freedom of a spin chain in the Mott insulating regime to the motional temperature of a bath. We have in particular presented a realistic experimental configuration where the spin chain is made of a mixture of two spin states of an alkali atom trapped into an optical lattice, and the bath is provided by a spatially overlapping Bose-Einstein condensate of a strongly magnetic atom. The coupling between the collective spin states of the spin chain and the motional degrees of freedom of the BEC is provided by dipole-dipole interactions, which inherently provide the necessary spin-orbit-coupling. Simulating spin chains of length up to 7, and provided the mechanical temperature of the BEC is low enough, we have found that low temperature (T < J) thermal states of the spin chain can be reached in a few seconds by the dissipative cooling associated to the dipolar interaction of the chain with the bath reservoir. We point out that the nature of the coupling between the spin chain and the BEC is collective. We find that the basic mechanism of cooling is a collective emission of a phonon in the BEC by the spin chain, which has direct analogies with cooperative effects in the context of spontaneous emission. This mechanism may also have interesting consequences for dipolar gas mixtures in new experiments [78], that could be apprehended with our framework.
The efficiency of cooling relies on the existence of phonon modes within the BEC whose energy is comparable to the energy to be dissipated within the spin chain (of the order of the super-exchange energy J), and whose momentum is comparable to the recoil momentum in the lattice K/2. These conditions imply rather strict experimental conditions, in particular on the lattice that controls the BEC phonon spectrum. To circumvent this difficulty, a natural extension would be to use a Fermi gas as a coolant, for which there exists a large density of available states sharing the same energy (close to the Fermi energy), but separated by the Fermi momentum (whose order of magnitude k F is typically comparable to K/2). This situation appears very promising for the release of low energy (i.e. comparable to J) spinexcitations from the spin chain to a Fermi gas.
In this work, we have considered explicitly the case of a coupling which frees the magnetization of the spin chain. To study antiferromagnetism, this implies that two single particle spin states of the species forming the spin chain should be kept almost degenerate, which can constitute an experimental difficulty. We have suggested that the use of lowenough magnetic field associated to micro-wave dressing is a possibility. However, we also pointed out that dissipative cooling at fixed magnetization is possible as well, that circumvents such experimental complication. In this case, it is not necessary to artificially engineer near-degeneracy between two single-particle spin states within the chain. Suppression of magnetization-changing collisions is indeed ensured in a lattice due to energy conservation requirements, when the Zeeman energy falls in a gap of the Bloch band structure.
Finally, we argue that spin-orbit coupling is an essential mechanism to generically provide cooling of a spin-chain down to its lowest energy states. Indeed, in a given realization of an experiment, both the total magnetization and the collective spin of an assembly of spins are not necessarily those of its ground state. While this work focuses on spin-orbit coupling due to dipole-dipole interactions, we are also studying the possibility to artificially engineer spin-orbit coupling using Raman-dressed schemes. associated energy difference between two Zeeman states is ∆ =3 kHz. We consider a microwave which is near resonant with the hyperfine structure of Potassium, at 1.285 GHz, with a detuning δ. Following for example [64] we assume a π polarized microwave field of Rabi frequency Ω/2π = 100 kHz. Due to selection rules, the microwave dressing does not affect the m F =-9/2 state, but introduces an energy shift Ω 2 /4δ on the m F =-7/2 state. In order to bring into degeneracy both states, we postulate Ω 2 /4δ = −∆, which sets the detuning, and the fraction of population off-resonantly excited to the upper hyperfine manyfold F'=7/2, p = Ω 2 /4δ 2 = 4∆ 2 /Ω 2 . For the experimental conditions which we chose, p 0.004. In case this value is too high (i.e. in the case where inelastic collisions between Dy atoms and K atoms in F'=7/2 limit the lifetime), a higher Rabi frequency is needed.
This tuning actually offers the possibility to study the Heisenberg antiferromagnetic model with tunable external polarizing field. The required degree of tunability is to be significantly finer than the superexchange energy J: as shown in Fig. 11, the finite-range 1D antiferromagnetic correlations at our temperature are robust to a bias of 0.2 × J. This bias corresponds to a magnetic field compensated to 0.4 mG. The nearest neighbor antiferromagnetic correlation is robust to a bias of order 1 × J.

Tight binding approximation
We use twice in this paper a tight binding approximation, introducing the bosons Wannier functions and assuming w B j, * ( r)w B k ( r) ≡ δ K (j, k)|w B j ( r)| 2 . When lowering the bath lattice depth along some direction of space, this becomes improper, such that the bath interaction with the spin chain, eq. 18, may be poorly approximated by eq. 19. In the conditions of section 4, we find that for a length-2 chain, the |T 0 → |S transition rate differs between the two expressions by 12 %. The agreement is improved to 1% if we apply the tight binding approximation only along the transverse axes (x , y ): where we use a 3D Wannier wavefunction for the fermions w F (x , y , z ) but a 2D Wannier wavefunction for the bosons along the two tight axes, w B ⊥ (x , y ).
Furthermore, to reduce computation time and study long chains, we calibrate on the length-2 chain calculation the acceptable truncation to Wannier functions in Fourier space, that reduces the reciprocal space summation nz m while not affecting the rates within a 1% tolerance.

Phonon direction sampling
Each rate Γ i→f requires to numerically sample the spherical coordinates (θ, φ) of the emitted excitation wavevectors. We design an approximately homogeneous, isotropic sampling around any point on the sphere as follows.
i) The polar angle θ is varied by a fixed step : θ n = π N θ (n + 1 2 ), with n ∈ [0, N θ − 1]. ii) All perimeters on the unit sphere, defined by a polar angle θ n , are sampled with same angular spacing. This means that the azimuthal angle φ sampling depends on θ n . We furthermore sample each perimeter from a different, random starting point. This sums up as, for each θ n , an azimuthal sampling φ k = random[0, 2π] + 2π N φ (θn) k, with N φ (θ n ) = Integer N 0 φ sin(θ n ) and k ∈ [0, N φ (θ n ) − 1]. We observe that approximately isotropic sampling around any point on the sphere is obtained for N 0 φ 2N θ . Then, we find that the cooling rates have converged to a 6 · 10 −2 stability using N θ = 50. Figures 4 to 6 are constructed with N θ varying from 50 to 200, figures 7 to 10 with N θ = 50.