Dynamics of large anisotropic spin in a sub-ohmic dissipative environment close to a quantum-phase transition

We investigate the dynamics of a large anisotropic spin whose easy-axis component is coupled to a bosonic bath with a spectral function $J(\w)\propto \omega^s$. Such a spin complex might be realized in a single-molecular magnet. Using the non-perturbative renormalization group, we calculate the line of quantum-phase transitions in the sub-ohmic regime ($s<1$). These quantum-phase transitions only occur for integer spin $J$. For half-integer $J$, the low temperature fixed-point is identical to the fixed-point of the spin-boson model without quantum-tunneling between the two levels. Short-time coherent oscillations in the spin decay prevail even into the localized phase in the sub-ohmic regime. The influence of the reorganization energy and the recurrence time on the decoherence in the absence of quantum-tunneling is discussed.


Introduction
Understanding the influence of the environment on the non-equilibrium dynamics of quantum systems remains one of the challenging questions of theoretical physics. A finite number of quantum mechanical degrees of freedom, a large spin or a qubit interacting with a infinitely large bath of non-interacting bosons with a continuous energy spectrum represents a typical class of model examples for such systems. Its simplest version, the spin-boson model [1], has contributed tremendously to our understanding of dissipation in quantum systems [2].
Single-molecular magnets (SMM) based on Mn 12 -actetate or Fe 8 -complexes behave as large single spins at low temperatures. The magnetic anisotropy of the molecules prevent a simple relaxation of the spin. The resulting hysteretic behavior in the magnetization has been subject to intensive experimental [3,4,5,6] and theoretical investigations [7,8,9,10] (see Ref. [5] for an overview). Using Wilson's numerical renormalization group (NRG) method [11,12], it was shown [13] that the combination of quantum tunneling and antiferromagnetic coupling to a metallic substrate induces a Kondo effect in such SMM. The recently developed time-dependent numerical renormalization group method (TD-NRG) [14,15] has been employed to analyze the real-time dynamics of such large spin [16].
In this paper we will investigate the equilibrium and non-equilibrium spin dynamics of a anisotropic spin coupled to a dissipative bosonic bath. Such a spin complex might be realized in a single-molecular magnet in which the physical properties depend on the type of coupling to the environment. Previously, detailed expansion of the spin-lattice relaxation [17] in powers of the local spin operator have been considered to estimate the spin-relaxation rates [17,8] by Fermi's golden-rule arguments. Vorrath and Brandes have reduced the spin-bath coupling to a linear term proportional to the S z operator [10]. A Weisskopf-Wigner coupling to the environment would describe the interaction with a fluctuating quantized magnetic field. This might be relevant for recent dephasing experiments [6] in SMM. Here, we will restrict ourselves to diagonally coupled spin-boson model interaction [1,10].
The real-time dynamics in the spin-boson model is usually investigated using the noninteracting-bilp approximation [1,18] or by Bloch-Redfield type of approaches [19]. Recently, it has been demonstrated that the non-perturbative numerical renormalization group [11,20,21,12] and its extension to non-equilibrium dynamics, the TD-NRG [14,15], is particularly suitable to access the quantum-critical region of the spin-boson model in a sub-ohmic environment [20,21] and is able to describe the real-time dynamics on short-time scales as well as exponentially long-time scales [22].
In order to gain a better understanding of the non-equilibrium dynamics, we also present the equilibrium properties of the model. We will map out the quantumcritical line separating a delocalized phase in weak coupling and a local phase at strong coupling to the dissipative sub-ohmic bath. The continuous quantum-phase transitions found for 0 < s < 1 and integer spin J resemble the phase diagram previously reported for the spin-boson model by Bulla et al. [20]. No phase transition was found for half integer spin J. Information on the equilibrium dynamics is provided by the equilibrium spin correlation function C(ω). The non-equilibrium spin dynamics of a anisotropic spin complex coupled to a dissipative bosonic bath is investigated in response to a sudden change of magnetic field close to quantum phase transitions.

Modelling of a large spin in a dissipative environment
The local Hamiltonian modeling single-molecular magnets [5] consists of a single large spin of size J subject to an easy-axis anisotropy energy A and quantum-tunneling terms [5,13] B 2n , which induce transitions between the eigenstates of S z : These quantum-tunneling terms stem from a reduction of the U (1) easy-axis spin symmetry to the discrete symmetry group of the SMM molecule. Assuming A ≫ |B 2n | and a positive A , the energies of states |J z are located on an inverted energy parabola. As depicted in Fig. 1, there is a fundamental difference between integer and halfinteger values of the spin [5,13]: the two ground states | ± J are connected via the quantum-tunneling terms B 2n only for integer values of J, while for half-integer values the ground states remain disconnected. Without further transition mechanism, a halfinteger anisotropic spin cannot relax into thermodynamic ground state after switching off the spin-polarizing external magnetic field. We use a generalization [10] of the spin-boson model (SBM) to describe the coupling to a dissipative bath. Here, b † q creates a boson in the mode q with energy ω q . The coupling between spin and bosonic bath is completely specified by the bath spectral function [1]   The asymptotic low-temperature behavior is determined by the low-energy part of the spectrum. Discarding high-energy details, the standard parametrization is where the dimensionless parameter α characterizes the dissipation strength. ω c is a cutoff energy and is used as the energy scale throughout the paper. The value s = 1 corresponds to the case of ohmic dissipation. In the usual spin-boson model (SBM), the local spin Hamiltonian is parametrized by an energy splitting ǫ and a quantum-tunneling term ∆ where σ i are Pauli matrices [1].

Equilibrium NRG for a bosonic environment
In all our calculations we have used Wilson's numerical renormalization group (NRG) method. Wilson's NRG method is a very powerful tool for accurately calculating equilibrium properties of quantum impurity models [12]. Originally the NRG was invented by Wilson for a fermionic bath to solve the Kondo problem [11,23]. The method was recently extended to treat quantum impurities coupled to a bosonic bath [20,21], combination of fermionic and bosonic baths [24] and real-time dynamics out of equilibrium [25,14,15]. The non-perturbative NRG approach has been successfully applied to arbitrary electron-bath couplings strength [20,21,24,22]. At the heart of this approach is a logarithmic discretization of the continuous bath, controlled by the discretization parameter Λ > 1; the continuum limit is recovered for Λ → 1. Using an appropriate unitary transformation [11], the Hamiltonian is mapped onto a semi-infinite tight-binding chain, defined by a sequence of finite-size Hamiltonians H m with the impurity coupled to the open end. The tight-binding parameters t m linking consecutive sites of the chain m and m+1 falls off exponentially as t m ∼ Λ −m . Each bosonic chain link is viewed as representative of an energy shell since its energy w m also decreases as w m ∼ Λ −m establishing an energy hierarchy. Both ensure that mode coupling can only occur between neighboring energy shells, which is essential for the application of the renormalization group procedure. To this end, the renormalization group transformation where H m is the Hamiltonian of a finite chain up to the site m. The annihilation (creation) operators of site m are denoted by a m (a † m ) and w m represents the energy of the bosonic mode of site m. Note that the rescaling of the Hamiltonian H m by Λ ensures the invariance of the energy spectrum of fixed point Hamiltonians under the RG transformation R[H m ]. For a detailed review on this method see Ref. [12].
The RG transformation (8) is used to set up and diagonalize iteratively the sequence of Hamiltonians H m . In the first step, only the large spin coupling to the single bosonic site m = 0 is considered. It turns out to be sufficient [20,21,12] to include only the N b lowest lying bosonic states, where N b takes typical values of 8−12. The reason for that is a quite subtle one: the coupling between different chain links decays exponentially and is restricted to nearest-neighbor coupling by construction, both essential for the RG procedure. In each successive step (i) a finite number of N b bosonic states of the next site m + 1 are added, (ii) the Hamiltonian matrices are diagonalized and (iii) only the lowest N s states are retained in each iteration. The discarding of high energy states is justified by the Boltzmannian form of the equilibrium density operator when the temperature is lowered simultaneously in each iteration step to the order T m ∝ Λ −m w c .
Denoting the set of low-lying eigenstates by |r N and the corresponding eigenvalues E r (N ) ∝ O(1) at iteration N , the equilibrium density matrix ρ 0 is given [12] byρ where Z N = r e −βE N r andβ are of the order O(1), such that T N = w c Λ −N /β. The thermodynamic expectation value of each local observableÔ is accessible at each temperature T N by the trace The procedure described above turns out to be very accurate because the coupling t m between the bosonic sites along the chain are falling off exponentially so that the rest of the semi-infinite chain contributes only perturbatively [11,12] at each iteration m, while contributions from the discarded high-energy states are exponentially suppressed by the Boltzmann factor.

Time-dependent NRG
While the equilibrium properties are fully determined by the energy spectrum of the Hamiltonian, the non-equilibrium dynamics requires two conditions: the initial condition encoded in the many-body density operatorρ 0 and the Hamiltonian H f which governs its time-evolution. For a time-independent Hamiltonian, the density operator evolves according toρ(t > 0) = e −iH f tρ 0 e iH f t . All time-dependent expectation values Ô (t) are given by where set = 1.
We obtain the density operatorρ 0 from an independent NRG run using a suitable initial Hamiltonian H i . For instance, by choosing a large local magnetic field in H i , we can prepare the system such that S z is fully polarized. To investigate decoherence we choose H i loc such that is an eigenstate of an appropriate H i loc with the lowest eigenenergy. In general, the initial density operatorρ 0 contains states which are most likely superpositions of excited states of H f . For the calculation of the real-time dynamics of the large spin it is therefore not sufficient to take into account only the retained states of the Hamiltonian H f obtained from an NRG procedure. The recently developed time-dependent NRG (TD-NRG) [14,15] circumvents this problem by including contributions from all states. It turns out that the set of all discarded states eliminated during the NRG procedure form a complete basis set [14,15] of the Wilson chain which is also an approximate eigenbasis of the Hamiltonian. Using this complete basis, it was shown [14,15] that Eq. (11) transforms into the central equation of the TD-NRG for the temperature T N where O m r,s = r; m|Ô|s; m are the matrix elements of any operatorÔ of the electronic subsystem at iteration m, and E m r , E m s are the eigenenergies of the eigenstates |r; m and |s; m of H f m . At each iteration m, the chain is formally partitioned into a "system" part on which the Hamiltonian H m acts exclusively and an environment part formed by the bosonic sites m + 1 to N . Tracing out these environmental degrees of freedom e yields the reduced density matrix [14,15] ρ red s,r (m) = e s, e; m|ρ 0 |r, e; m (14) at iteration m, whereρ 0 is given by (9) using H i . The restricted sum trun r,s in Eq. (13) implies that at least one of the states r and s is discarded at iteration m. Excitations involving only kept states contribute at a later iteration and must be excluded from the sum.
As a consequence, all energy shells m contribute to the time evolution: the short time dynamics is governed by the high energy states while the long time behavior is determined by the low lying excitations. Dephasing and dissipation is encoded in the phase factors e i(E m r −E m s )t as well as the reduced density matrix ρ red s,r (m). Discretization of the bath continuum will lead to finite-size oscillations of the realtime dynamics around the continuum solution and deviations of expectation values from the true equilibrium at long time scales. In order to separate the unphysical finite-size oscillations from the true continuum behavior, we average over different bath discretization schemes using Oliveira's z-averaging (for details see Refs. [26,15]). We average over N z = 16 different bath discretizations in our calculation.
Previously, the TD-NRG has been successfully applied to the simple spin-boson model [15,22], the SMM coupled to a fermionic bath [16] and electron-transfer in a dissipative environment [27].

Integer Spin J and quantum-phase transitions
In Fig. 2, we plot the zero-temperature phase diagram in the sub-ohmic regime for a moderate value of an integer spin (J = 3.) Depicted is the rescaled critical coupling strength α c (2J) 2 as a function of the bath power-law exponent in the sub-ohmic regime. The line marks the quantum-phase transition (QPT) between a localized and a delocalized phase. We have obtained the critical coupling strength by investigating the NRG level flow which shows that the QPT is associated with a distincted fixed point. For s > 1, no QPT exist [20,21]. Analyzing the NRG fixed-point spectra [11], the same qualitative picture emerges for all integer values of J. We noted that the critical coupling strength depends (i) on the value of the spin J and (ii) the values of the quantum tunneling rates B 2n . Only for s → 1, the critical coupling becomes independent of B 2n and always approaches which can be understood from the analytical form of H I . The fixed-point spectra always agrees with those of the spin-boson model. This establishes a mapping of our model at low temperatures onto an effective spin-boson model. This mapping becomes intuitively clear by inspecting Fig. 1(b). It shows a coupling between the two lowest states | ± J on the parabola. For T → 0, it is expected that the model can be replaced by an effective spin-boson model with an tunneling rate ∆ which is a complicated function of parameters A and B 2n . Therefore, it is not surprising that our calculations for the anisotropic large spin model reveal the similarities to the results by Bulla et al. [21], who established such a quantum phase transition in the SBM. Since the coupling of the two states | ± J to the bosonic bath is proportional to J z , we find that α c (s)(2J) 2 resembles the phase boundary of spin-boson model. The critical coupling α c is proportional to 1/J 2 .
Such a behaviour is also observed in the equilibrium spin-spin correlation function C(ω) = ℑm ≪ S z |S z ≫ (ω − iδ). The correlation function C(ω) is calculated using the sum-rule conserving algorithm [28,29]. We plot C(ω) on a log-log scale for a series of coupling constants and the sub-ohmic regime s = 0.5 in Fig. 3(a). The high energy feature stemming from the easy-axis energy splitting proportional to A is   clearly visible. At lowest coupling strength α = 10 −4 , the second peak at ω ≈ 3 × 10 −3 is well pronounced, originating from the weak effective quantum tunneling rate which is a function of parameters of H loc : A, B 2n as well as J. For α < α c , C(ω) vanishes according to C(ω) ∝ ω s , while C(ω) diverges as ω −s for ω → 0 and α > α c .

Half-integer Spin J
A completely different picture emerges for half-integer spins J. Independent of the bath exponent s, no phase transition is found. As seen already in Fig. 1, the quantum tunneling matrix elements B 2 and B 4 cannot connect the two states | ± J . The model can be mapped onto a spin-boson model with a vanishing tunneling rate ∆ at low temperatures. To illustrate this, we compare the NRG level flow [12] of the large spin model with the NRG level flow of a simple spin-boson model and ∆ = 0. The fixed-point spectra for m → ∞ is identical for both models as depicted in Fig. 4. The coupling constant α determines the crossover from high-temperature to the low temperature fixed point: this crossover scale accidentally coincides in both models for our choice of coupling constants. C(ω) depicted in Fig. 3(b) is similar to the results of Fig. 3(a) at high energies: the peak due to the easy-axis splitting A is also present. However, a pronounced peak at the effective quantum-tunnel energy is absent for any half-integer J.

Non-equilibrium dynamics
For the investigation of the non-equilibrium dynamics, the local spin is prepared in a maximal polarized initial state with the expectation value S z (0) = J by applying a strong local magnetic field in z-direction. The external magnetic field is switched off at t = 0. Fig. 5 displays the real-time dynamics of S z (t) for a series of coupling strength α and an ohmic bath (s = 1), J = 4 and the easy-axis splitting A = 0.1. A weak quantum-tunneling was enabled by B 2 = 10 −3 and B 4 = 2 × 10 −3 . As a reference, the dynamics of a decouple spin (α = 0) is included in the graph. For small α, damped coherent oscillations with a superimposed short-time dynamics of small amplitude are clearly visible. The oscillation frequency decreases with increasing α and vanished at a finite coupling α. Above that value -here α ≈ 0.075 -only a spin decay is observed. This corresponds to the regime 0.5 < α < 1 of the SBM. When approaching the critical coupling α c (J = 4) ≈ 0.016, S z (t) only decays insignificantly on a very short time scale of the order O(1) and remains constant for t → ∞, as seen in panel Fig. 5(a). This behaviour occurs in the localized phase of the model where the effective tunneling rate is renormalized to zero.
A different picture emerges in the sub-ohmic regime. The real-time dynamics of a large anisotropic spin with J = 5 is depicted for different coupling strength in coupling regime, coherent high and low frequency oscillations are clearly visible, in addition to a slow decay of the amplitude. At a critical coupling of α = 0.00061, we found a quantum-phase transition in the equilibrium NRG calculations. As discussed in Sec. 4.1, the quantum-critical point separates a delocalized from a localized phase. As seen in Fig. 6, coherent oscillations prevail in the short-time dynamics even in the localized phase close to the quantum-phase transition. A finite time is required before the strong correlations can build up. In the long-time limit, S z (t) approaches a finite value in the localized phase which should be proportional to the order parameter for the QPT. These results are very similar to the dynamics reported in the subohmic spin-boson model [22]. In our case, larger number of eigenstates of H loc yields additional coherent oscillations absent in the simple model. We expect that coherent oscillations prevail to longer time-scales with decreasing exponent s. Simultaneously, the critical coupling α c reduces strongly with deceasing exponent s, as shown in Fig. 2.

Decoherence
Decoherence is caused by the interaction of a subsystem with its environment. The concept of decoherence [30,31] provides inside into how a coherence quantum state |s of the subsystem is destroyed by the entanglement with the environment. The reduced density matrix [32] evolves from a matrix describing such a pure state |s into a matrix describing an ensemble. In an appropriate local basis |i , the off-diagonal matrix elements of ρ red i,j must vanish for t → ∞ [30,31]. In a diagonally coupled spin-boson model as described by the interaction H I , Eq. (4), this would be the eigenstates of S z [30]. In the absence of quantum tunneling, S z also commutes with H loc and no energy is exchanged by the coupling H I to the environment. Decoherence is induced by the dephasing with continuum of the bath modes each contributing a different phase shift. For J = 1/2, this is a well studied problem, and it was shown analytically [33,34] that the offdiagonal component of the reduced density matrix for the spin can be written as ρ ↑,↓ (t) = e −Γ(t) ρ ↑,↓ (0), where Γ(t) is given by the exact analytic expression Here T is the temperature. Note that decoherence of the pure quantum state |s does not follow a simple exponential form. Approximating the dephasing by a function ∝ exp(−t/T deph ), which defines a single "dephasing" time-scale T deph , might be insufficient even in very simple models such as the SBM, since Γ(t) cannot be replaced by t/T deph . We initially prepared the spin system in a pure state |s comprising of a linear combination of all S z eigenstates by an appropriate choice of the initial Hamiltonian H i loc . We can gain some information on the decoherence of a large anisotropic spin by measuring the time-dependent expectation value S x (t) = Ŝ x (t), which only depend on the off-diagonal matrix elements of ρ red i,j (t). The quantum-tunneling parameters were set to B 2n = 0 in the calculations. Fig. 7 shows the dynamics of S x (t) on a rescaled dimensionless time axis tα for a integer (J = 5) and a half-integer (J = 9/2) value of the spin in the s = 0.5 sub-ohmic regime. In both cases, we observe an oscillatory decay of S x different from the result for the spin-boson model, Eq. (17). Here, the finite coupling α introduces an effective A via the reorganization energy E 1g The effective energyẼ m ∝ (A−E 1g )m 2 of eigenstates |m of S z induces quantumoscillations whose oscillation frequency is proportional to α for A = 0. This α dependence of the frequency is revealed by plotting S x (t) versus the dimensionless timescale tα so that the oscillation beats coincide and occur on a time scale of tα ≈ O(1). Since all eigenenergiesẼ m are related to each other by integer values proportional to A − E 1g , the local dynamics is characterized by a recurrence time proportional to 1/(A − E 1g ) at which the original local state would reappear in the absence of a spin-environment coupling. However, the envelope function describing the time-dependence of the recurrence amplitude decreases on a time scale proportional to 1/α. We can prolong this recurrence time by an appropriate counter-term in the Hamiltonian. The curves of Fig. 7(c) are obtained for the same parameters as in Fig. 7(a) and by setting A = 2αω c /s. Due to the discretization of the bath continuum, the cancellation is not perfect, but the oscillation frequency is lowered by almost two decades.

Conclusion
The equilibrium and non-equilibrium dynamics of an anisotropic large spin coupled to a dissipative sub-ohmic bosonic environment has been investigated using the nonperturbative NRG. Such a large spin might be realized in a single-molecular magnet. The coupling of the bath modes to the easy-axis component of the spin enhances the tendency for localization. A quantum-phase transition exists only for integer spin values in the sub-ohmic regime, and the quantum-critical point is associated with a new fixed point. These results generalized the work of Bulla et al. [20] for the spinboson model. The spin-spin correlation function C(ω) diverges approximately as |ω| −s at and above the critical coupling α c . The critical coupling for s → 1 is related to the critical coupling of the spin-boson model by α c = α SBM c /(2J) 2 . We have shown that the fixed point of the half-integer model is identical to the S = 1/2 spin-boson model without quantum-tunneling ∆. No quantum-phase transition is found for half-integer J due to the lack of quantum-tunneling between the two states | ± J .
We investigated the spin decay by switching off an external easy-axis magnetic field at t = 0. In the ohmic regime (s = 1), the real-time dynamics of S z (t) is governed by several local frequencies stemming from the easy-axis splitting A and the quantumtunneling matrix elements B 2n . With increasing coupling α, the quantum-tunneling is reduced, but the splitting A is increased by the reorganization energy E 1g . The low-frequency coherent oscillations vanish at some value α but the high-frequency coherent oscillations prevail with a very small amplitude. Further increasing of α beyond the critical coupling α c drives the system in the localized phase. Here, the spin decays weakly on a very short time-scale of the order of the reciprocal cutoff 1/ω c . A completely different picture emerges in a sub-ohmic environment. Oscillatory response is found at short and intermediate time-scales even in the localized phase at coupling strength α > α c .
In the calculation of the decoherence through dephasing, we identified the occurring long-time oscillation frequency with recurrence frequency which is dependent on the reorganization energy. These oscillations can be reduced and suppressed by adding a counter-term in the original Hamiltonian. We have demonstrated that the non-perturbative renormalization group approach can describe such equilibrium and non-equilibrium dynamics close to a quantum-phase transition.