Dynamical symmetries and crossovers in a three-spin system with collective dissipation

Understanding and controlling the dynamical behavior of quantum systems has seen flourishing interest [1–3] propelled by theoretical and experimental progress that has made it possible to observe and manipulate such systems with unprecedented accuracy. Much attention has also been devoted recently to the notion of dynamical phase transitions in such systems, relating them to the nonanalyticity of, e.g., the Loschmidt echo [4] or the logarithm of a biased partition function in large-deviation (LD) theory [5], which has a natural interpretation in terms of the statistics of rare trajectories observed in experiments. The study of the dynamics of quantum systems through LD methods [6–8] emerged recently both as an extension of the theory as applied to classical systems [9–12] and as a dynamical complement to the standard analysis of equilibrium phase transitions in many-body systems [13]. Here, nonanalyticities in the LD free-energy function of a system, extracted from the equations governing its dynamical behavior, are identified in the literature with dynamical phase transition points [6].

Following [6], in this paper we are interested in studying the statistical properties of rare quantum-jump trajectories [14] of a system that interacts with a heat bath driving the system out of equilibrium. We consider the dynamical LD properties of a simple three-spin quantum open model, which departs from those recently studied in two respects: first, dissipation is due to nonclassical bilinear jump operators; and second, we consider a current-like dynamical order parameter. The two central results in this paper are (a) the observation of intermittency between dynamical phases of distinct activity, itself a consequence of the reducibility of the dynamics in an appropriate limit due to the collective jump operators; and (b) the existence of a Gallavotti-Cohen symmetry in LD functions associated with the time-asymmetric order parameter, analogous to that found in driven classical systems [15], which gives rise to a fluctuation theorem [16] relating to the quantum jump rate. Our system therefore provides a minimal but extendible model that uncovers the effects of thermal baths and the nontrivial interplay between local and global decay channels [17] on the non-equilibrium dynamics of a quantum system.

We start with three spins-1/2, which we label $j=1,2,3$, placed at the vertices of an equilateral triangle (see figure 1), interacting via an Ising-type interaction, and in a uniform magnetic field. Suppose that spin 1 is in a harmonic trap, allowing it to move on an axis parallel to the line joining the two other spins, whereas spins 2 and 3 are held tightly pinned. Since the spin–spin interaction depends on the distance between each pair of spins, the motion of spin 1 will modulate the interactions between them. By damping the harmonic motion, we effectively couple the collective spin degree of freedom to this environment. Consequently, the thermal environment will act as a driving force on the spin system. More precisely, the thermal bath will drive the global degree of freedom of the spin chain. This mechanism can be generalised to other systems with similar symmetry properties with respect to the motion of a particular component. We explore the trajectories of this collective degree of freedom, finding evidence of coexisting dynamical phases and nontrivial dynamical symmetries. It is worth pointing out that individual spins can be coupled to the motion of a harmonic oscillator through the use of trapped ions [18], by embedding solid-state qubits into mechanically compliant structures [19], in nanomechanical resonator arrays [20], on graphene layers [21], or on diamond surfaces [22]. Furthermore, this system lends itself well to being extended by adding more spins, thereby changing the symmetries of the model.

Zoom In Zoom Out Reset image size

We will first derive the master equation for the spin system by adiabatically eliminating the motion of spin 1. The resulting reduced dynamics will be investigated first by means of the quantum version of LD theory, following [6]. This procedure gives access to the statistics of the trajectories of the system, painting a clear picture of its dynamical behavior. Following this, we complement these analytical insights through the use of quantum-jump Monte Carlo simulations [23] that give access to a transparent physical interpretation of the processes occurring. The Hamiltonian describing the three-spin system interacting with a harmonic oscillator as described above can be written $\hat{H}={{\hat{H}}_{{\rm s}}}+{{\hat{H}}_{{\rm m}}}+{{\hat{H}}_{{\rm s}-{\rm m}}}$ with

$$\begin{eqnarray}&&{{\hat{H}}_{{\rm s}}}=\alpha \mathop{\sum }\limits_{i=1}^{3}\hat{\sigma }_{{\rm x}}^{i}+\mathop{\sum }\limits_{\langle i,j\rangle }\hat{\sigma }_{{\rm x}}^{i}\hat{\sigma }_{{\rm x}}^{j}-B\mathop{\sum }\limits_{i=1}^{3}\hat{\sigma }_{{\rm z}}^{{\rm i}},\end{eqnarray} \tag{ 1 }$$

being the spin-chain Hamiltonian under uniform magnetic fields α and $B\ll \alpha$, ${{\hat{H}}_{{\rm m}}}={{\omega }_{{\rm m}}}{{\hat{b}}^{\dagger }}\hat{b}$ the harmonic oscillator Hamiltonian, and

$$\begin{eqnarray}&&{{\hat{H}}_{{\rm s}-{\rm m}}}=-g({{\hat{b}}^{\dagger }}+\hat{b})\hat{\sigma }_{{\rm x}}^{1}(\hat{\sigma }_{{\rm x}}^{2}-\hat{\sigma }_{{\rm x}}^{3})\end{eqnarray} \tag{ 2 }$$

the interaction Hamiltonian between the spin chain and the harmonic oscillator; $\langle i,j\rangle$ denotes a sum over nearest neighbours. ${{\hat{H}}_{{\rm s}-{\rm m}}}$ follows from observing that the interaction between any pair of spins depends on the distance between them. In the geometry illustrated in figure 1, when spin 1 moves in a direction parallel to the line joining spins 2 and 3, the distance between particles 1 and 2 decreases (increases) by the same amount that the distance between spins 1 and 3 increases (decreases). Upon identifying $\hat{x}={{\hat{b}}^{\dagger }}+\hat{b}$ as the dimensionless position operator for spin 1, we arrive at the given form for the interaction Hamiltonian. Different geometries or numbers of spins can also give rise to similar Hamiltonians and effects as the ones we discuss below.

We now move into the interaction picture with respect to ${{\hat{H}}_{0}}={{\hat{H}}_{{\rm s}}}+{{\hat{H}}_{{\rm m}}}$, setting $\hat{b}\to \tilde{b}\;{{e}^{-i{{\omega }_{{\rm m}}}t}}$ and $\hat{\sigma }_{{\rm x}}^{j}\to \tilde{\sigma }_{+}^{j}{{e}^{-iBt}}+\tilde{\sigma }_{-}^{j}{{e}^{iBt}}$, where $\tilde{\sigma }_{\pm }^{j}=\tilde{\sigma }_{{\rm x}}^{j}\pm i\tilde{\sigma }_{{\rm y}}^{j}$ is the spin–flip operator for the $j{\rm th}$ spin, and the tilde distinguishes interaction-picture operators from Schrödinger-picture ones. Assuming that ${{\omega }_{{\rm m}}}=2B$ and performing a rotating-wave approximation allows us to consider only the time-independent terms in the interaction-picture interaction Hamiltonian ${{\tilde{H}}_{{\rm s}-{\rm m}}}$, resulting in ${{\tilde{H}}_{{\rm s}-{\rm m}}}\approx -g(\tilde{b}\;{{\tilde{\sigma }}_{+}}+{{\tilde{b}}^{\dagger }}\;{{\tilde{\sigma }}_{-}})$, where the collective operators ${{\tilde{\sigma }}_{\pm }}=\tilde{\sigma }_{\pm }^{1}(\tilde{\sigma }_{\pm }^{2}-\tilde{\sigma }_{\pm }^{3})$. Assuming that $\mid ag\mid b$ sets the longest time-scale of the dynamics, we can adiabatically eliminate the harmonic-oscillator degree of freedom [23]. To do so, we follow the projection operator technique described in detail, e.g., in [24] and the supplemental information for [25]. We write a master equation, valid up to second order in g, that governs the evolution of the reduced interaction-picture density matrix $\tilde{\rho }(t)$ for the spin-only system as

$$\begin{eqnarray}&&{{\partial }_{t}}\tilde{\rho }=\int _{0}^{\infty }{\rm d}\tau \;{\rm T}{{{\rm r}}_{{\rm m}}}\left\{ {{\mathcal{P}}_{{\rm m}}}{{\mathcal{L}}_{{\rm s}-{\rm m}}}{{e}^{{{\mathcal{D}}_{{\rm m}}}\tau }}{{\mathcal{L}}_{{\rm s}-{\rm m}}}{{\mathcal{P}}_{{\rm m}}}\tilde{\rho } \right\},\end{eqnarray} \tag{ 3 }$$

where ${{\mathcal{P}}_{{\rm m}}}\tilde{\rho }=\tilde{\rho }\rho _{{\rm m}}^{{\rm ss}}$, with $\rho _{{\rm m}}^{{\rm ss}}$ the stationary solution of the harmonic oscillator, ${\rm T}{{{\rm r}}_{{\rm m}}}\left\{ \bullet \right\}$ denotes the trace over the motional degree of freedom, ${{\mathcal{L}}_{{\rm s}-{\rm m}}}\bullet =-i[{{\tilde{H}}_{{\rm s}-{\rm m}}},\bullet ]$, and

$$\begin{eqnarray}&&{{\mathcal{D}}_{{\rm m}}}\bullet =\frac{\kappa }{2}(\bar{n}+1)(2{{\hat{b}}^{\dagger }}\;\bullet \;\hat{b}-\{\hat{b}{{\hat{b}}^{\dagger }},\bullet \})+\frac{\kappa }{2}\bar{n}(2\hat{b}\;\bullet \;{{\hat{b}}^{\dagger }}-\{{{\hat{b}}^{\dagger }}\hat{b},\bullet \}),\end{eqnarray} \tag{ 4 }$$

is the Lindblad-form dissipator associated with a damped harmonic oscillator connected with a damping rate κ to a thermal bath whose average number of excitations is $\bar{n}$ [23, 26]. The resulting master equation for the spin-system density matrix reads

$$\begin{eqnarray}&&{{\partial }_{t}}\tilde{\rho }=(\bar{n}+1)\Gamma \;{{\tilde{\mathcal{D}}}_{\downarrow }}[\tilde{\rho }]+\bar{n}\Gamma \;{{\tilde{\mathcal{D}}}_{\uparrow }}[\tilde{\rho }],\end{eqnarray} \tag{ 5 }$$

with $\Gamma ={{g}^{2}}/\kappa$. We also have

$$\begin{eqnarray}&&{{\tilde{\mathcal{D}}}_{\downarrow }}[\bullet ]=2{{\tilde{\sigma }}_{-}}\;\bullet \;{{\tilde{\sigma }}_{+}}-\{{{\tilde{\sigma }}_{+}}{{\tilde{\sigma }}_{-}},\bullet \},\end{eqnarray} \tag{ 6 }$$

with ${{\tilde{\mathcal{D}}}_{\uparrow }}[\bullet ]$ that is obtained from ${{\tilde{\mathcal{D}}}_{\downarrow }}[\bullet ]$ by exchanging ${{\tilde{\sigma }}_{-}}$ and ${{\tilde{\sigma }}_{+}}$. Moving out of the interaction picture, we obtain a master equation for the reduced spin density matrix in the Schrödinger picture simply written as ${{\partial }_{t}}\rho =\mathcal{W}[\rho ]$ with the corresponding super-operator $\mathcal{W}[\bullet ]$ defined as

$$\begin{eqnarray}&&\mathcal{W}[\bullet ]=-i[{{\hat{H}}_{{\rm s}}},\bullet ]+\Gamma (\bar{n}+1){{\mathcal{D}}_{\downarrow }}[\bullet ]+\Gamma \bar{n}\;{{\mathcal{D}}_{\uparrow }}[\bullet ],\end{eqnarray} \tag{ 7 }$$

where ${{D}_{\uparrow ,\downarrow }}:={{e}^{i{{{\hat{H}}}_{{\rm s}}}t}}{{\tilde{\mathcal{D}}}_{\uparrow ,\downarrow }}{{e}^{-i{{{\hat{H}}}_{{\rm s}}}t}}$. We similarly define ${{\hat{\sigma }}_{\pm }}$ by transforming ${{\tilde{\sigma }}_{\pm }}$ and shall drop the label t when it can be understood from the context. By tracing out the motion of the damped harmonic oscillator, we have obtained an effective damping acting on a collective degree of freedom of the spin system through the jump operators ${{\hat{\sigma }}_{\pm }}$. This collective damping is the source of the interesting behavior we shall explore in the following.

As a first step in exploring the LD behavior of our system, we associate a counting process related to the flow of excitations into or out of the system with the collective jump operators ${{\hat{\sigma }}_{\pm }}$. There are two counting processes ${{K}_{\pm }}$ associated with ${{\hat{\sigma }}_{\pm }}$. ${{K}_{+}}$ counts excitations emitted into the bath and ${{K}_{-}}$ excitations absorbed from it. We then define an overall counting process K, which counts the net number of excitations emitted into the bath due to the collective spin flips $K:={{K}_{+}}-{{K}_{-}}$ (in contrast to the total activity given by ${{K}_{+}}+{{K}_{-}}$ [6]). Next, we can unravel the master equation of the reduced density matrix by projecting it onto a particular number of jump events, i.e., ${{\partial }_{t}}{{\rho }^{K}}={{P}^{K}}\mathcal{W}[\rho ]$, where P^K is a projector over trajectories with K net jump events, and ${{p}_{K}}(t)={\rm Tr}\{{{\rho }^{K}}(t)\}={\rm Tr}\{{{P}^{K}}\rho (t)\}$ represents the probability of observing such a trajectory [14]. The moment-generating function associated with this probability p_K(t) can be written [6]:

$$\begin{eqnarray}&&Z(t,s)=\mathop{\sum }\limits_{K=0}^{\infty }{{e}^{-sK}}{{p}_{K}}(t)={\rm Tr}\{{{\rho }_{s}}(t)\}.\end{eqnarray} \tag{ 8 }$$

with ${{\rho }_{s}}(t)=\sum _{K=0}^{\infty }{{e}^{-sK}}{{\rho }^{K}}(t)$ the Laplace transform of the density matrix with respect to the net excitation exchanges K between the system and the bath; we call s the 'bias parameter.' The Laplace-transformed density matrix evolves according to the modified master equation ${{\partial }_{t}}{{\rho }_{s}}=(\mathcal{W}+{{\mathcal{V}}_{s}})[{{\rho }_{s}}]$, where

$$\begin{eqnarray}&&{{\mathcal{V}}_{s}}[\bullet ]=\Gamma \left[ (\bar{n}+1)({{e}^{-s}}-1){{{\hat{\sigma }}}_{-}}\;\bullet \;{{{\hat{\sigma }}}_{+}}+\bar{n}({{e}^{s}}-1){{{\hat{\sigma }}}_{+}}\;\bullet \;{{{\hat{\sigma }}}_{-}} \right].\end{eqnarray} \tag{ 9 }$$

In the long-time limit, LD theory applies and we can write $Z(t,s)\to {{e}^{t\theta (s)}}$, where $\theta (s)$ represents the systemʼs dynamical free energy [6, 27]. Consequently, we have $\theta (s)={{{\rm lim} }_{t\to \infty }}{\rm ln} \left( {\rm Tr}\{{{\rho }_{s}}\} \right)/t$ ⁴ . $\theta (s)$ is also given by the eigenvalue of ${{\mathcal{W}}_{s}}:=\mathcal{W}+{{\mathcal{V}}_{s}}$ with the largest real part [6] (which can be shown to be real [28]).

Derivatives of this dynamical free energy with respect to s can be used to obtain the activity (net count rate) $k(s)=-{{\partial }_{s}}\theta (s)$ of the system. This quantity is represented in figure 2 for different values of the bath population $\bar{n}$. As is clearly visible from the upper plot, for the value of the bias parameter s = 0, we have a nonanalytic point in $\theta (s)$ for any value of $\bar{n}$. The lower curves show that this point presents two distinct values of the activity k(s). Unbiased dynamics takes place at s = 0; from this we can conclude the existence of two dynamical phases [6, 27]. These two phases have $k({{0}^{-}})\gt 0$ and $k({{0}^{+}})=0$, i.e., one phase is active in the sense that the net rate of excitations exchanged between the bath and the system is nonzero, whereas the other has an exact balance between excitations emitted and absorbed. The two terms composing the dissipative part of the master equation (7) act at different rates; we thus deduce that this balanced phase is inactive and no emission or absorption events occur. Focussing now on the active phase corresponding to $k({{0}^{-}})\gt 0$, we also notice that the activity seems not to depend on the thermal population of the bath $\bar{n}$, since all the curves converge to the same point as $s\to 0$. This effect, seen for small $\bar{n}$, stems from the weak coupling between the bath and the system ($\Gamma \ll \alpha$).

Zoom In Zoom Out Reset image size

A feature of figure 2 that is not seen when considering a 'symmetric' dynamical observable such as the total activity ${{K}_{+}}+{{K}_{-}}$, as in [6, 27, 29–31], is the second point of nonanalyticity in $\theta (s)$ that occurs for $s\gt 0$ when $\bar{n}\gt 0$. Counting processes of the type we consider, unlike 'symmetric' ones, are odd with respect to time reversal. This is related to a Gallavotti–Cohen symmetry [15] due to the driven nature of the systemʼs dynamics [32]. In contrast to most studied examples of systems presenting such dynamical properties, the dynamics of a global, rather than a local, degree of freedo are considered here. Based on the detailed balance exhibited by equation (5), we have that for ${{e}^{{{s}_{0}}}}=(\bar{n}+1)/\bar{n}$ we find $\theta (s)=\theta ({{s}_{0}}-s)$, which yields a fluctuation theorem of the form

$$\begin{eqnarray}&&p_{K}^{\infty }/p_{-K}^{\infty }={{e}^{{{s}_{0}}K}}\;\mathop{\longrightarrow }\limits^{\bar{n}\gg 1}\;{{e}^{K/\bar{n}}},\end{eqnarray} \tag{ 10 }$$

where we have defined $p_{K}^{\infty }:={{{\rm lim} }_{t\to \infty }}{{p}_{K}}(t)$, that relates the infinite-time probability of observing a trajectory with a net count of K to one with a net count of $-K$. This ratio approaches unity as $\bar{n}\to \infty$ in which case the rates for the collective jump process balance (${{s}_{0}}\to 0$), such that $k({{0}^{\pm }})=0$ ⁵ . Conversely, as $\bar{n}\to 0$, the probability for observing negative K goes to zero and the above ratio diverges.

To make the model more realistic, we now add independent damping channels acting on each spin and explore the consequences of these channels on the dynamical behavior of the system. For simplicity, the single-spin baths are also taken to have occupation $\bar{n}$ and be coupled with the rate γ. In equation (7) we set $\mathcal{W}\to {{\mathcal{W}}^{\prime }}$:

$$\begin{eqnarray}&&{{\mathcal{W}}^{\prime }}[\bullet ]\;:=\mathcal{W}[\bullet ]+\gamma (\bar{n}+1)\mathop{\sum }\limits_{i=1}^{3}\mathcal{D}_{\downarrow }^{i}[\bullet ]+\gamma \bar{n}\mathop{\sum }\limits_{i=1}^{3}\mathcal{D}_{\uparrow }^{i}[\bullet ].\end{eqnarray} \tag{ 11 }$$

${{\rho }_{s}}$ evolves according to ${{\partial }_{t}}{{\rho }_{s}}=({{\mathcal{W}}^{\prime }}+{{\mathcal{V}}_{s}})[{{\rho }_{s}}]$, where ${{\mathcal{V}}_{s}}[{{\rho }_{s}}]$ stays unchanged from its definition equation (9) when $\gamma /3\ll \Gamma$ and neglecting correlation effects between the various damping channels. $\mathcal{D}_{\uparrow ,\downarrow }^{i}$ are the spin–flip Liouvillians for spin i.

As figure 3 illustrates, introducing extra damping at the level of the individual spins has a strong effect on the dynamical properties of the system. Concentrating on the dynamical free energy, we see that individual spin damping will smooth out the nonanalyticity at s = 0. The same holds for the second nonanalyticity at $s\gt 0$. This smoothing effect is also visible in the activity, which becomes well-defined everywhere. Meanwhile, as is clearly visible in figure 3, for small values of the single-spin damping and low thermal population, the activity remains approximately constant. Conversely, for high thermal population and strong single-spin damping, the activity can switch sign and become negative. Physically, this corresponds to the case where the single-spin damping channel upsets the balance, making it more likely that excitations enter (${{K}_{-}}$) than leave the system (${{K}_{+}}$) through the collective channel, leading to a thermally driven system.

Zoom In Zoom Out Reset image size

The LD approach to dynamical phase transitions yields a transparent physical interpretation based on the statistics of ensembles of trajectories of the system. To explore these statistics, we now conduct an analysis based on a Monte Carlo wave-function (MCWF) simulation (also known as quantum jump Monte Carlo). For classical systems, Monte Carlo simulations or other numerical methods are currently used to obtain the LD function [10, 12]; such methods have also occasionally been applied to a quantum system [27]. The MCWF technique is a well-established method to simulate open system dynamics, such as the one we are interested in, following the ideas set forth in [23, 33]. Using this technique and starting off from a randomly chosen initial condition, we simulate trajectories in equation (11) of jump events related to the operators ${{\hat{\sigma }}_{\pm }}(t)$. For each trajectory generated, we estimate the activity $k(0)$ by calculating the net rate of jump events. We consider a set of 2000 trajectories, each having approximatively 10⁴ jump events, from which we obtain probability distributions for $k(0)$, as represented in figure 4. This figure shows the probability distributions for thermal populations of $\bar{n}=0$ and 5, with (upper panel) our without (lower panel) single-spin damping. Samples of the typical trajectories obtained corresponding to different parts of the distributions are shown in the inset; each vertical line corresponds to a jump event, with upper ones representing an emission from the system to the bath (${{K}_{+}}$) and lower ones the opposite (${{K}_{-}}$). It is clear that when $\gamma =0$ (upper panel), the probability distribution is bimodal, with one peak centred at an activity equal to $k({{0}^{+}})=0$ and the other at $k({{0}^{-}})\gt 0$. In the former case (see inset) the corresponding trajectories have no jump events.

Zoom In Zoom Out Reset image size

The second peak is centred about the same value for both values of $\bar{n}$, in agreement with what is expected from the lower panel of figure 2, where we saw that the activity is independent of the thermal population in the low-temperature limit. Corresponding trajectories are shown in the inset, and as expected demonstrate jumps associated only with ${{\hat{\sigma }}_{-}}$ (i.e., the system losing excitations to the bath) for $\bar{n}=0$; for $\bar{n}=5$ jumps are observed in both directions. It can be see in figure 4 that this second peak of the distribution broadens when the temperature increases. All this yields the interpretation that the Liouville space accessible to the system consists of two disconnected subspaces, one active and one inactive, to which the two peaks in the activity distribution are related. Based on the fraction of K = 0 trajectories in our simulation, we can determine that 25% of the Liouville space is inactive.

Consider now the case with damping on the individual spins, corresponding to the one shown in figure 4 in the lower panel. In contrast to the $\gamma =0$ case, we immediately see that the distribution becomes unimodal, with the mean activity decreasing by almost 25% compared to the active trajectories of the $\gamma =0$ situation. This fraction corresponds to the fraction of inactive trajectories observed when $\gamma =0$ (upper panel). The unimodal behavior and the mentioned reduction lead to the interpretation that the single-spin damping channel connects the two previously disconnected parts of the Liouville space. This interpretation is supported by the sample trajectories shown in the inset of the lower panel of figure 4, where one can observe the trajectory 'blinking,' i.e., spontaneously switching between active and inactive behavior. In analogy with what occurs with the total activity (see [27]), the intermittency in trajectories is a form of 'mesoscopic' (i.e., finite time) dynamical phase coexistence, consistent with the fact that the dynamical free-energy is analytic in this case, and the transition therefore becomes a crossover. A simple interpretation of this behaviour can be deduced by looking at the eigenspace of ${{\hat{H}}_{{\rm s}}}$ and the collective spin–flip operators. We start by writing the Schrödinger-picture Hamiltonian ${{\hat{H}}_{{\rm s}}}$ and collective operators ${{\hat{\sigma }}_{\pm }}=\hat{\sigma }_{\pm }^{1}(\hat{\sigma }_{\pm }^{2}-\hat{\sigma }_{\pm }^{3})$ in the computational basis. The resulting 8 × 8 matrices are not trivial to diagonalize analytically. However, it is straightforward to find an eigenvector $\mid \psi \rangle$ such that ${{\hat{\sigma }}_{\pm }}\mid \psi \rangle =0$ and ${{\hat{H}}_{{\rm s}}}\mid \psi \rangle =\epsilon \mid \psi \rangle$ for some nonzero real number . The first pair of conditions render $\mid \psi \rangle$ a dark state of the collective spin–flip operators, which are the operators that enter the dissipative part of the reduced master equation, and the last condition assures that $\mid \psi \rangle$ is an eigenstate of the dynamics, i.e., that a system in $\mid \psi \rangle$ will remain in this inactive subspace. It can similarly be shown that no such $\mid \psi \rangle$ also obeys $\hat{\sigma }_{-}^{1}\mid \psi \rangle =0$, and that no state exists such that $\hat{\sigma }_{\pm }^{1}\mid \psi \rangle =0$; in other words, any inactive state is coupled through the dynamics to the active subspace when single-spin damping is introduced. To sum up, one can find a subspace of the Liouville space that is both inactive and isolated, in the sense that ${{\hat{H}}_{{\rm s}}}$ does not couple it to the rest of the space, and a subspace that is active. When $\gamma \ne 0$, these two partitions are no longer isolated, and the system can switch dynamically between the active and inactive subspaces.

To understand how closely our MCWF results agree with the LD analysis, we present in figure 5 a quantitative comparison between the two, where we plot the activity of the system at s = 0 as a function of the thermal populations, for both low and high temperatures. This plot shows that the MCWF results (data points) agree closely with the results from the LD theory (solid curves). As visible in figure 4, increasing $\bar{n}$ tends to broaden the distribution of the activity. Consequently, the error bars shown in figure 5 grow quickly with $\bar{n}$. As discussed previously, we clearly see that both curves tend to zero as $\bar{n}\to \infty$, where the rates of the two counting processes ${{K}_{\pm }}$ balance such that the net rate of jump events is zero. Note also that the decrease in activity matches the predicted 25%, independently of $\bar{n}$.

Zoom In Zoom Out Reset image size

We sum up by recalling our main results. We have explored a simple, yet intriguing, system consisting of three equidistant spins interacting pairwise, one of which moves in a harmonic trap. This motion gives rise to collective spin dynamics, which are dissipated through the mechanical decay channel. Adopting an LD approach to analyze this system, we observe that its dynamics consist of two distinct dynamical phases, one active and one inactive, possessing different emission statistics. We observed that these two phases can be mixed by introducing damping on the individual spins. All our observations were confirmed through Monte Carlo simulations, which lend themselves to a natural interpretation in terms of ensembles of trajectories observed in repeated experimental runs. The system we explore is not overly complex, but yields a surprisingly rich behavior.

Acknowledgments