Discontinuities in driven spin-boson systems due to coherent destruction of tunneling: breakdown of the Floquet-Gibbs distribution

TL;DR

This paper reveals that in driven spin-boson systems, the stationary state distribution exhibits discontinuities at points of coherent destruction of tunneling, affecting observable emission spectra and explained via Floquet-Redfield theory.

Contribution

It demonstrates the occurrence of distribution discontinuities in dissipative driven two-level systems at coherent destruction points, linking them to observable spectral jumps and providing analytical and numerical analysis.

Findings

Discontinuities in stationary state distributions at tunneling destruction points

Observable jumps in Mollow triplet emission spectra

Analytical and numerical confirmation using Floquet-Redfield formalism

Abstract

We show that the probability distribution of the stationary state of a dissipative ac-driven two-level system exhibits discontinuities, i.e. jumps, for parameters at which coherent destruction of tunneling takes place. These discontinuities can be observed as jumps in the emission of the Mollow triplet. The jumps are the consequence of discontinuities in the transition rates, which we calculate numerically and analytically based on the secular Floquet-Redfield formalism.

Figures

Tables

Table 1. Table 1: Coupling coefficients a μ , λ ( n ) superscript subscript 𝑎 𝜇 𝜆 𝑛 a_{\mu,\lambda}^{(n)} of all possible transitions between different Floquet states | φ λ ​ ( 0 ) ⟩ ket subscript 𝜑 𝜆 0 \left|\varphi_{\lambda}(0)\right> .

$n=0$	${|-1⟩}_x$	${|+1⟩}_x$
${|-1⟩}_x$	${𝒮}_x^{\left(0\right)}$	${𝒮}_z^{\left(0\right)}$
${|+1⟩}_x$	${𝒮}_z^{\left(0\right)}$	$-{𝒮}_x^{\left(0\right)}$

Equations

$$H(t)=\frac{h_{\mathrm{x}}}{2}{\sigma }_x+\frac{h_{\mathrm{z}}(t)}{2}{\sigma }_{\mathrm{z}}+{\hat{\sigma }}_{\theta }\sum _k{V}_k\left(b_k+b_k^{†}\right)+H_{\mathrm{B}},$$

$$\frac{d}{dt}{p}_0$$

$$\frac{d}{dt}{p}_1$$

$$A_{\lambda \leftarrow \mu }^{(n)}$$

$$a_{\lambda \leftarrow \mu }^{(n)}$$

$$|{\phi }_{\lambda }(t)⟩=e^{-i\hat{\Lambda }(t)}|{\phi }_{\lambda }(0)⟩,$$

$$e^{-i\hat{\Lambda }(t)}$$

$$a_{\mu ,\lambda }^{(n)}={\mathcal{S}}_{\mathrm{x}}^n{⟨{\sigma }_{\mathrm{x}}⟩}_{\mu \lambda }+i\cdot {\mathcal{S}}_{\mathrm{y}}^n{⟨{\sigma }_{\mathrm{y}}⟩}_{\mu \lambda }+{\mathcal{S}}_{\mathrm{z}}^n{⟨{\sigma }_{\mathrm{z}}⟩}_{\mu \lambda },$$

$${\mathcal{S}}_{\mathrm{x}}^n$$

$$+g_{\mathrm{z}}{\delta }_{n\text{mod}2,1}{l}_{|n|},$$

$${\mathcal{S}}_{\mathrm{y}}^n$$

$${\mathcal{S}}_{\mathrm{z}}^n$$

$$+g_{\mathrm{x}}\frac{1}{2}{\delta }_{n\text{mod}2,1}{\alpha }_{\mathrm{x}}^{(n)},$$

$${\alpha }_{\mathrm{x}}^{(n)}$$

$$\frac{I_{\mathrm{red}}}{I_{\mathrm{blue}}}\approx \frac{{\left({\mathcal{J}}_{-1}\left(z_0\right)-\frac{h_{\mathrm{x}}}{\Omega }{\alpha }_{\mathrm{x}}^{(-1)}(z_0)\right)}^2}{{\left({\mathcal{J}}_{-1}\left(z_0\right)+\frac{h_x}{\Omega }{\alpha }_{\mathrm{x}}^{(-1)}(z_0)\right)}^2}\approx 1\pm 2\frac{h_x}{\Omega }{\alpha }_{\mathrm{x}}^{(-1)}(z_0)$$

$$\hat{U}(t,0)=e^{-i\hat{\Lambda }(t)}{e}^{-i{H}_{\mathrm{eff}}(t)}$$

$$H_{\mathrm{s}}=\frac{h_{\mathrm{x}}}{2}{\sigma }_x+\frac{h_{\mathrm{z},1}}{2}\cos (\Omega t){\sigma }_z,$$

$${\hat{U}}_{\mathrm{rot}}(t)$$

$$\theta (t)$$

$$U_{\mathrm{rot}}(t)$$

$$U_{\mathrm{rot}}^{†}(t){\sigma }_{\mathrm{z}}{U}_{\mathrm{rot}}(t)$$

$$U_{\mathrm{rot}}^{†}(t){\sigma }_{\mathrm{x}}{U}_{\mathrm{rot}}(t)$$

$$=\cos (2\theta ){\sigma }_{\mathrm{x}}-\sin (2\theta ){\sigma }_{\mathrm{y}}$$

$$U_{\mathrm{rot}}^{†}(t){\sigma }_{\mathrm{y}}{U}_{rot}(t)$$

$$H_{\mathrm{s}}(t)=\sum _{n=-\infty }^{\infty }{H}_{\mathrm{s}}^{(n)}{e}^{in\Omega t},$$

$$H_{\mathrm{s}}^{(n)}=\frac{h_x}{2}{\mathcal{J}}_n\left(\frac{h_{\mathrm{z},1}}{\Omega }\right)\{\begin{array}{ll}{\sigma }_{\mathrm{x}}& n\text{even}\\ i{\sigma }_{\mathrm{y}}& n\text{odd}\end{array},$$

$$J_n(x)=\frac{1}{2\pi }{\int }_{-\pi }^{\pi }{e}^{x\sin \tau +n\tau }d\tau$$

$$H_{\mathrm{eff}}$$

$$\hat{\Lambda }(t)$$

$$H_{\mathrm{eff}}=\frac{h_{\mathrm{x}}}{2}{\mathcal{J}}_0(h_{\mathrm{z},1}/\Omega ){\sigma }_x-2\frac{h_x^2}{\Omega }\sum _{k>0}\frac{{\mathcal{J}}_{2k-1}(h_{\mathrm{z},1}/\Omega )}{2k-1}{\sigma }_{\mathrm{z}}+\mathcal{O}\left(\frac{1}{\Omega }\right)$$

$$\hat{\Lambda }(t)$$

$$=\frac{h_{\mathrm{x}}}{\Omega }\sum _{k=1}^{\infty }{l}_{2k}\sin (2k\Omega t){\sigma }_{\mathrm{x}}+\frac{h_{\mathrm{x}}}{\Omega }\sum _{k=1}^{\infty }{l}_{2k-1}(\cos ((2k-1)\Omega t)-1){\sigma }_{\mathrm{y}}+\mathcal{O}\left(\frac{1}{{\Omega }^2}\right)$$

$$=\frac{h_{\mathrm{x}}}{\Omega }{l}_{\mathrm{x}}(t){\sigma }_{\mathrm{x}}+\frac{h_{\mathrm{x}}}{\Omega }{l}_{\mathrm{y}}(t){\sigma }_{\mathrm{y}}+\mathcal{O}\left(\frac{1}{{\Omega }^2}\right),$$

$$a_{\mu ,\lambda }^{(n)}$$

$$=\frac{1}{\tau }{\int }_0^{\tau }⟨u_{\mu }(0)|e^{i\hat{\Lambda }(t)}{U}_{\mathrm{rot}}^{†}\left[g_x{\sigma }_x+g_z{\sigma }_z\right]U_{\mathrm{rot}}{e}^{-i\hat{\Lambda }(t)}|u_{\lambda }(0)⟩e^{-in\Omega t}dt$$

$$e^{-i\hat{\Lambda }(t)}\approx 1-i\hat{\Lambda }(t),$$

$$a_{\mu ,\lambda }^{(n)}=g_x{a}_{\mathrm{x},\mu ,\lambda }^{(n)}+g_z{a}_{\mathrm{z},\mu ,\lambda }^{(n)},$$

$$a_{x,\mu ,\lambda }^{(n)}$$

Full text

In a recent article, Shirai et al. have discussed under which special conditions effective Floquet-Gibbs states arise Shirai et al. (2015, 2016). The probabilities of these periodic Floquet states, characteristic states of periodically-driven systems, are determined by their corresponding quasienergies $\epsilon_{\lambda}$ in a Gibbs-like fashion, thus $p_{\lambda}\propto e^{-\beta\epsilon_{\lambda}}$. However, given the richness of quantum effects in periodically-driven systems and their experimental control, the resulting stationary states are not necessary Gibbs-Floquet states and can possibly exhibit intriguing features.

In isolated ac-driven systems, CDT is defined by an exact degeneracy of two quasienergies. Consequently, the dynamics in these levels is frozen as the effective Hamiltonian is zero Frasca (2003); Stockburger (1999); Barata and Wreszinski (2000). The emission of an ac-driven two-level system ($\left|0\right>,\left|1\right>$) coupled to an environment (sketched in Fig. 1(a)) can exhibit a Mollow triplet as sketched in Fig. 1(b) Mollow (1969); Pigeau et al. (2015); Yan et al. (2016). Thereby, the driving field with frequency $\Omega$ induces transitions between the states which are accompanied by emitting a photon (or phonon, which we focus here) with frequencies $\Omega-\Delta$ ($\left|0\right>\rightarrow\left|1\right>$), $\Omega$ ($\left|0\right>\rightarrow\left|0\right>$ and $\left|1\right>\rightarrow\left|1\right>$), and $\Omega+\Delta$ ($\left|1\right>\rightarrow\left|0\right>$). Here, $\Delta$ denotes the difference between two quasienergies. The emitted phonons can be blue shifted $I_{\rm{b(lue)}}$, unshifted or red shifted $I_{\rm{r(ed)}}$, respectively. Clearly, the unshifted transitions do no change the system state, but the shifted transitions can be used to control it.

Our main findings are summarized in Fig. 1 (c) and (d), where we depict the probability of the Floquet state with the lowest quasienergy $p_{0}$ and the difference of emitted phonons $I_{\rm{b}}-I_{\rm{r}}$ as a function of driving amplitude $h_{z,1}$ and a system-environment coupling angle $\theta$. Remarkably, we observe jumps in both observables at specific values of the driving amplitude, marked with arrows in Fig. 1 (c) and (d). As we will explain in detail, they appear at parameters, at which CDT takes place.

Model system. We consider a generic two-level system, which describes, e.g., the low-energy physics of a double-well potential or a superconducting qubit. It is externally driven and coupled to a thermal environment as sketched in Fig. 1(a). The Hamiltonian reads

[TABLE]

where $\sigma_{\alpha}$ with $\alpha=\left\{\rm x,y,z\right\}$ denote the Pauli matrices, $h_{\rm x}$ denotes the tunneling amplitude and $h_{\rm z}(t)=h_{\rm z,0}+h_{\rm z,1}\cos(\Omega t)$ is the time-dependent energy splitting, where $h_{\rm z,0}$ is the offset, $h_{\rm z,1}$ is the driving amplitude and $\Omega$ is the driving frequency. The bath is quadratic in bosonic operators $b_{k}$ and is coupled via the system operator $\hat{\sigma}_{\theta}=\sin\theta\sigma_{\rm x}+\cos\theta\sigma_{\rm z}$ with strengths $V_{k}$. Depending on the coupling angle $\theta$, it is known that undriven systems can give rise to diverse physical behavior Duan et al. ; Castro Neto et al. (2003); Guo et al. (2012); Kohler et al. (2013); Bruognolo et al. (2014).

Floquet theory describes the dynamics of periodically-driven systems Shirley (1965); Grifoni and Hänggi (1998). Due to the periodic system dynamics with frequency $\Omega=2\pi/\tau$, there are characteristic states of the system which fulfill $\left|\Phi_{n,\lambda}(t)\right>=e^{-i\epsilon_{n,\lambda}t}\left|\varphi_{n,\lambda}(t)\right>$, with quasienergy $\epsilon_{n,\lambda}$ and periodic Floquet state $\left|\varphi_{n,\lambda}(t)\right>=\left|\varphi_{n,\lambda}(t+\tau)\right>$. These states are the analogue to the eigenstates in time-independent systems. Importantly, the Floquet states $\lambda$ are not uniquely defined due to the Brillouin zone index $n$: a state with index $n$ can be related to the $n=0$ state by $\epsilon_{n,\lambda}=\epsilon_{0,\lambda}+n\Omega$ and $\left|\varphi_{n,\lambda}(t)\right>=e^{-in\Omega t}\left|\varphi_{0,\lambda}(t)\right>$. The stroboscopic Floquet states are the eigenstates of the time-evolution operator after one period $\hat{U}_{\rm s}(\tau)\left|\varphi_{n,\lambda}(0)\right>=e^{-i\epsilon_{n,\lambda}\tau}\left|\varphi_{n,\lambda}(0)\right>$.

In Fig. 2(a), we depict the numerically calculated quasienergies $\epsilon_{\lambda}=\epsilon_{0,\lambda}$ as function of $h_{\rm z,1}/\Omega$. Here, the index $\lambda$ can take the values $\lambda=0,1$ and we consider $n=0$. The stroboscopic dynamics follows the effective Hamiltonian $H_{\rm eff}=\frac{h_{\rm z,0}}{2}\sigma_{\rm z}+\frac{h_{\rm x}}{2}\mathcal{J}_{0}(h_{\rm z,1}/\Omega)\sigma_{\rm x}+\mathcal{O}(\frac{1}{\Omega})$. For $h_{\rm z,0}=0$, we find $\epsilon_{\lambda}=\pm h_{\rm x}\mathcal{J}_{0}(h_{\rm z,1}/\Omega)/2$, so that there are degeneracies at the roots of the Bessel function $\mathcal{J}_{0}(h_{\rm z,1}/\Omega)=0$. This is the celebrated CDT effect, as the dynamics at these parameters is frozen. The stroboscopic Floquet states read $\left|\varphi_{\lambda}(0)\right>\approx\left|\text{sign}\left[\mathcal{J}_{0}(h_{\rm z,1}/\Omega)\right](-1)^{\lambda}\right>_{\rm x}$. Accordingly, there is a non analytic crossover of the Floquet state, e.g. $\left|\varphi_{0}(0)\right>=\left|-1\right>_{\rm x}$ to $\left|\varphi_{0}(0)\right>=\left|+1\right>_{\rm x}$, at the roots of the Bessel function. As we explain, the crossover of the Floquet states is the origin of the probability jumps observed in Fig. 1.

Rate equations. An important point to realize is that though the states $\left|\varphi_{n,\lambda}(t)\right>$ of the spin system are equivalent for different $n$ in a closed system, $n$ becomes physically relevant when the system is coupled to a thermal bath $H_{\rm B}$. In this situation, the bath can trigger transitions between different Brillouin zones $n,n^{\prime}$. In Fig. 1(b) we illustrate transitions associated with $\Delta n=n-n^{\prime}=-1$. Using the secular Floquet-Redfield formalism Grifoni and Hänggi (1998); Carreno et al. (2017), one can derive the rate equations

[TABLE]

where $p_{\lambda}$ denotes the probability to be in Floquet state $\lambda$ and $A_{\lambda\leftarrow\mu}^{(n)}$ is the transition probability between Floquet states

[TABLE]

Here, $n_{\rm B}(\omega)$ denotes the Bose distribution and $\Gamma(\omega)=\sum_{k}V_{k}^{2}\delta(\omega-\omega_{k})=\gamma\omega/(\omega^{2}+\omega_{c}^{2})$ (coupling strength $\gamma$, cut-off frequency $\omega_{c}$) denotes the coupling density, which we define for negative frequencies by $\Gamma(\omega)=-\Gamma(-\omega)$, and $\Delta_{\lambda\mu}^{n}=\epsilon_{\mu}-\epsilon_{\lambda}-n\Omega$. The time-dependent operator reads $\hat{\sigma}_{\theta}(t)=e^{i\hat{\Lambda}(t)}\hat{\sigma}_{\theta}e^{-i\hat{\Lambda}(t)}$, where $\hat{\Lambda}(t)$ is defined by

[TABLE]

at which the operator $e^{-i\hat{\Lambda}(t)}$ propagates the Floquet states.

It is easy to show that the coefficients fulfill $a_{\lambda\leftarrow\mu}^{(n)}=\left(a_{\mu\leftarrow\lambda}^{(-n)}\right)^{*}$. As a consequence, the related rates obey the detailed balance condition $A_{\lambda\leftarrow\mu}^{(n)}=A_{\mu\leftarrow\lambda}^{(-n)}e^{\Delta_{\lambda\mu}^{(n)}/T}$, where $T$ is the temperature of the environment. Yet, in general $\left|a_{\mu\leftarrow\lambda}^{(n)}\right|\neq\left|a_{\lambda\leftarrow\mu}^{(n)}\right|$, which gives rise to a break down of the detailed balance relation in the stationary state, thus $p_{0}/p_{1}\neq e^{-\Delta/T}$ with $\Delta=\epsilon_{1}-\epsilon_{0}$.

Stationary state. Figures 2(b) and (c) depict the stationary state probabilities for the coupling $\hat{\sigma}_{\theta}$ with $\theta=\pi/2,\pi/4$, respectively. For $\theta=0,\pi$, the system approaches a Floquet-Gibbs state according to Ref. Shirai et al. (2016). In (b), we find a probability inversion for small $h_{\rm z,1}/\Omega$ and probability jumps at the roots of the Bessel function $\mathcal{J}_{0}(h_{\rm z,1}/\Omega)=0$. These effects can be explained by analyzing the rates in Eq. (3). Due to a generalized parity symmetry, $A_{\lambda\leftarrow\mu}^{(0)}$ vanishes exactly par . Consequently, the rate equations are dominated by the transitions $\Delta n=-1$, as $n_{\rm B}(|\Delta^{n>0}_{\lambda\mu}|)\ll 1$ for $\Delta n>0$ due to a large energy difference $\Delta_{\lambda,\mu}^{n}\approx-n\Omega$, i.e., the transitions marked by the red and blue arrows in Fig. 1(b) (green transitions do not change the state). The corresponding coefficients $\left|a_{\mu\leftarrow\lambda}^{(-1)}\right|^{2}$ are depicted in Fig. 3(a), where we observe that $\left|a_{1\leftarrow 0}^{(-1)}\right|^{2}>\left|a_{0\leftarrow 1}^{(-1)}\right|^{2}$ for, e.g., $h_{z,1}/\Omega<z_{0}$, with $z_{0}$ denoting the first root of the zeroth-order Bessel function. As $n_{\rm B}(\Delta^{n}_{\lambda\mu})\ll 1$ and $\Gamma(\Delta^{(-1)}_{10})\approx\Gamma(\Delta^{(-1)}_{01})$, we find from Eq. (2) that $p_{1}/p_{0}\approx\left|a_{1\leftarrow 0}^{(-1)}\right|^{2}/\left|a_{0\leftarrow 1}^{(-1)}\right|^{2}$ which explains the probability inversion. The rates explain the jump in the probability distribution. In Fig. 3(a) we observe jumps at the CDT positions. The jump is magnified in Fig. 3(b). The non-continuous behavior becomes more clear when considering Eq. (3). At the CDT, the Floquet states switch, thus $\left|\varphi_{\lambda}(0)\right>\leftrightarrow\left|\varphi_{\mu}(0)\right>$, which gives rise to the non-analytic behavior.

A similar reasoning can be applied to the $\hat{\sigma}_{\pi/4}$ coupling depicted in Fig. 2(c). Away from the CDT, the probability distribution mainly corresponds to the Floquet-Gibbs state. This appears as the coupling $\hat{\sigma}_{\pi/4}$ has a $\sigma_{z}$ contribution so that the rates $A_{\lambda\leftarrow\mu}^{(0)}\gg A_{\lambda\leftarrow\mu}^{(-1)}$ result mainly in a Gibbs state. However, the coefficients $a_{\lambda\leftarrow\mu}^{(-1)}$ are almost equal to those of the $\hat{\sigma}_{\pi/2}$ case and thus give rise to a probability jump, yet, to a very small extend.

Importantly, although there is a jump discontinuity, the density matrix remains continuous as a function of $h_{\rm z,1}$. In the Floquet-Redfield formalism, the system density matrix reads $\rho_{\rm s}(t)=\sum_{\lambda}p_{\lambda}(t)\left|\varphi_{\lambda}(t)\right>\left<\varphi_{\lambda}(t)\right|$. As there is a simultaneous switch of $p_{\lambda}$ and $\left|\varphi_{\lambda}(t)\right>$, the density matrix remains continuous. Yet, the probability jump does not depend on how the states $\lambda=0,1$ are labeled, as the labeling can be uniquely defined. Let us consider the system with $h_{\rm z,0}\neq 0$. A corresponding quasienergy spectrum is depicted in Fig. 2(a) with thin lines. For finite but small $h_{\rm z,0}$ the gap closing is released as depicted in Fig. 2(a). Accordingly, the states depend smoothly on $h_{\rm z,1}$, so that the $p_{\lambda}$ are also uniquely defined as can be observed in Fig. 2(b). The ordering of the states $\lambda=0,1$ can be thus uniquely defined in terms of the limit $h_{\rm z,0}\rightarrow 0$.

The interplay of $\theta$ and $h_{\rm z,1}$ can be fully analyzed in Fig. 1(c), where the stationary state strongly depends on $\theta$, though, the system is only weakly coupled to the environment. In particular, for $\theta\approx 0.5\pi$ and small $h_{\rm z,1}$ we find inversion, $p_{0}<0.5$, thus, there is a strong deviation from the Floquet-Gibbs state, which can be found for $\theta=0$.

Phonon emission. Due to the continuous stationary state, the jump cannot be observed in system observables. However, the nonanalytic behavior can be observed in the emission. Every transition $A_{\lambda\leftarrow\mu}^{(-1)}$ is related to the emission of phonons with either energy $\Omega$,$\Omega\pm\Delta$. The corresponding intensities $I_{\rm b,r}=I(\Omega\pm\Delta)$ are given by $I_{\rm b}=\left(\Omega+\epsilon_{1}-\epsilon_{0}\right)A_{0\leftarrow 1}^{(-1)}p_{1},I_{\rm r}=\left(\Omega-\epsilon_{1}+\epsilon_{0}\right)A_{1\leftarrow 0}^{(-1)}p_{0}.$ We depict the blue and red shifted intensities in Fig. 3 (c) and (d). For $\hat{\sigma}_{\pi/4}$, we observe that the two intensities are almost equal. As the stationary state is governed by the rates $A_{\lambda\leftarrow\nu}^{(-1)}$, we find $p_{0}\propto A_{0\leftarrow 1}^{(-1)}$ and $p_{1}\propto A_{1\leftarrow 0}^{(-1)}$, so that $I_{\rm r}\approx I_{\rm b}$ for high frequency $\Omega\gg\epsilon_{1}-\epsilon_{0}$.

For $\hat{\sigma}_{0}=\sigma_{z}$, at which the system approaches a Floquet-Gibbs state, it is known that the rates $A_{\mu\leftarrow\nu}^{n\neq 0}\approx 0$ Shirai et al. (2015). Consequently, here both $I_{\rm r/b}$ vanish. Considering the difference $I_{\rm b}-I_{\rm r}$ in Fig. 1(d), we consequently find that the difference is smooth in $h_{\rm z,1}$ for both limiting cases $\theta=0,\pi/2$. However, in between there is a significant jump, which is strongest for about $\theta\approx 0.3\pi$.

Let us consider the $\hat{\sigma}_{\pi/4}$ coupling. Due to the $\sigma_{z}$ coupling component, the $A_{\lambda\leftarrow\nu}^{(0)}$ rates are dominant, which leads to an (almost) thermalization of the system with its environment. However, the rates $A_{\lambda\leftarrow\nu}^{(-1)}$ still exhibit a jump at the CDT, so that we find jumps in $I_{\rm r}$ and $I_{\rm b}$, as can be observed in Fig. 3(d).

Magnus expansion. To find an appropriate approximation of the numerical results, we employ a rotating wave approximation defined by a unitary transformation Gong et al. (2009) $U_{\rm rot}(t)=\exp\left[-i\sigma_{z}\theta(t)\right]$ with $\theta(t)=\frac{h_{\rm z,1}}{\Omega}\sin(\Omega t)$. The resulting Hamiltonian $H_{r}(t)=U_{\rm rot}^{\dagger}\hat{H}(t)U_{\rm rot}$ also exhibits $\tau$ periodicity. The Floquet state propagator can be written as

[TABLE]

with $e^{-i\hat{\Lambda}_{r}(t)}$ being the Floquet state propator in the rotating frame. The explicit form of $\hat{\Lambda}_{r}$ can be calculated using a standard high frequency expansion Mananga and Charpentier (2011). However, using Eq. (5) without approximation, it is not possible to evaluate the integral in Eq. (3). As $\hat{\Lambda}_{r}\propto\mathcal{O}(\frac{1}{\Omega})$ is small in the high frequency limit, the expansion in Eq. (5) is justified. For consistency, we also expand the Hamiltonian up to the same order, thus $H_{\rm eff}=H_{\rm eff}^{(0)}+\frac{1}{\Omega}H_{\rm eff}^{(1)}+\mathcal{O}\left(\frac{1}{\Omega^{2}}\right)$. Evaluating Eq. (3), we obtain

[TABLE]

where we have defined $\left<\sigma_{\alpha}\right>_{\mu\lambda}=\left<u_{\mu}(0)\right|\sigma_{\alpha}\left|u_{\lambda}(0)\right>$ and

[TABLE]

with $g_{\rm x}=\sin\theta$, $g_{\rm z}=\cos\theta$, $l_{m}=\frac{h_{x}}{m\Omega}\mathcal{J}_{m}\left(h_{z}/\Omega\right)$ for $m>0$ and $l_{m}=0$ for $m=0$. Due to the factor $h_{\rm x}/\Omega$ the coefficients $l_{n}$ can be considered to be small.

The transition coefficients $a_{\mu,\lambda}^{(n)}$ in Eq. (7) are evaluated in Tab. 1. For $n=0$, we find a hermitian structure for the eigenstates $\left|\varphi_{\lambda}(0)\right>\approx\left|-1\right>_{\rm x},\left|+1\right>_{\rm x}$, which are mainly determined by $H_{\rm eff}^{(0)}\propto\sigma_{\rm x}$. For $n=-1$, the transition coefficients are dominated by $\mathcal{S}_{\rm y}^{(-1)}\propto\mathcal{J}_{-1}(h_{z,1}/\Omega)$, which explains the oscillations in Fig. 3(a). Importantly, they do not exhibit a hermitian structure. This leads to a breakdown of the detailed balance relation, and gives rise to the jump in the probability distribution in Fig. 2(b) and (c). This appears as $\left|\varphi_{0}(0)\right>$ switches from $\left|-1\right>_{\rm x}$ to $\left|+1\right>_{\rm x}$ and, simultaneously, $\left|\varphi_{1}(0)\right>$ switches from $\left|+1\right>_{\rm x}$ to $\left|-1\right>_{\rm x}$, causing a jump of $a_{\mu,\lambda}^{(-1)}$.

We can use Eq. (6) to understand the intensity jump in the high-frequency regime. Assuming that $p_{\lambda}\approx 0.5$ as in Fig. 2(b), we find for $h_{\rm z,1}/\Omega\approx z_{0}$

[TABLE]

with a constant $\alpha_{\rm x}^{(-1)}(z_{0})$, and ’$\pm$’ for $h_{z}/\Omega\lessgtr z_{0}$. Thus, interestingly, this ratio (Fig. 3(e)) close to the CDT and consequently the jump magnitude $4\frac{h_{x}}{\Omega}\alpha_{\rm x}^{(-1)}(z_{0})$ scales as $1/\Omega$.

Discussion. The CDT in a driven dissipative system gives rise to surprising effects. Besides the well-investigated freezing of the internal system dynamics at the CDT, the presence of the thermal environment can give rise to counter intuitive jumps in the probability distribution of the Floquet state. Yet, as there is a simultaneous switch of Floquet states and probability, the reduced density matrix remains continuous while crossing the CDT. Consequently, the probability jumps can not be observed in system observables. Moreover, the jump behavior has a drastic consequence on the Mollow triplet, such that the blue and red shifted intensities both exhibit a discontinuity at the CDT. This can be directly measured by phonon spectroscopy. The ratio of both shifted intensities exhibits a jump of about $10\%$.

The underlying physical reason for these discontinuities is a switch of the Floquet states at the CDT. This causes a jump in the system-bath coupling coefficients. Consequently, the effect does not depend on the details of the thermal bath. Yet, the system-bath coupling operator is important. For a pure $\hat{\sigma}_{0}=\sigma_{z}$ coupling as investigated in Ref. Shirai et al. (2015), the stationary density matrix of the system recovers an effective Floquet-Gibbs states. For $\hat{\sigma}_{\pi/2}=\sigma_{x}$, we find extreme deviations from the Floquet-Gibbs state with occupation inversion even for small driving amplitude as observed in Fig. 2(b). With this coupling, the probability jump at the CDT turns out to be most significant. However, there is no signature in the emitted phonons. For $\hat{\sigma}_{\pi/4}$ coupling, though the probability jumps are very small, there is a clear jump discontinuity in the blue and red shifted intensity.

The noncontinuous behaviour is not a consequence of the high-frequency regime, which we considered here to explain the numerical results with analytical calculations. The CDT appears due to an exact degeneracy of the quasienergyies which is persistant even for very low driving-frequencies Grifoni and Hänggi (1998). Consequently, the jump behavior will remain when lowering the driving frequency. Our findings are not restricted to the dissipative two-level system. Similar probability jumps could be also observed in dissipative driven Lipikin-Meshkov-Glick model, which gives rise to many-body CDT Gong et al. (2009). Furthermore, these findings will be important for electronic transport Sánchez et al. (2008); Brandes et al. (2004); Sánchez et al. (2007). It will be interesting to explore the fate of the jumps for stronger environmental coupling using methods such as Restrepo et al. (2016); Magazzù et al. (2018); Xu and Cao (2016); Duan et al. (2017); Lee et al. (2012)

Finally, the analytical treatment based on the Magnus expansion extended with our approximation Eq. (5) is the key tool to understand the discontinuities. The analytical treatment can generalized to other configurations of the driving which could give rise to even more flexible mechanisms of the stationary state, and other exotic stationary states in driven systems.

Acknowledgments. G. Engelhardt gratefully acknowledges financial support from the China Postdoc Science Foundation (Grant No.: 2018M640054) , J. Cao acknowledges support from the NSF (Grant No.: CHE 1800301 and CHE 1836913). Both acknowledge the Natural Science Foundation of China (under Grant No.:U1530401). G. Platero acknowledges MINECO (Grant No.: MAT2017-86717-P). The authors thank Hui Dong and Tilen Čadež for inspiring discussions

I Transformation into a rotating frame

The system Hamiltonian describing the isolated spin system reads

[TABLE]

where we consider the case $h_{\rm z,0}=0$ for simplicity. The transformation into a rotating frame is defined by the unitary operator

[TABLE]

which transforms the system wave function $\left|\Psi(t)\right>=\hat{U}_{\rm rot}^{\dagger}(t)\left|\tilde{\Psi}(t)\right>$, with $\left|\tilde{\Psi}(t)\right>$ being the wave function in the rotating frame.

Evaluating the expressions of the unitary operators, we find

[TABLE]

which gives rise to the transformation of the Pauli operators

[TABLE]

II Rotating wave approximation

This can be expanded as a Fourier decomposition

[TABLE]

at which each Fourier component reads

[TABLE]

where we have used the integral representation of the Bessel functions

[TABLE]

Following Ref. Mananga and Charpentier (2011) , a high-frequency expansion of the system time-evolution operator in Eq. (8) reads

[TABLE]

where we consider the lowest non-vanishing terms in order of $\frac{1}{\Omega}$. For the system under consideration, we find the well-known result

[TABLE]

and

[TABLE]

where $l_{n}(z)=\frac{1}{n\Omega}\mathcal{J}_{n}\left(z\right)$. We use this to calculate the terms $a_{\mu,\lambda}^{(n)}$

[TABLE]

Using the expansion

[TABLE]

we can explicitly solve the integral and find

[TABLE]

with

[TABLE]

This integrals can be solved analytically, giving

[TABLE]

and

[TABLE]

where we have defined

[TABLE]

Evaluating these terms we obtain

[TABLE]

[TABLE]

[TABLE]

[TABLE]

[TABLE]