Quantum Zeno effect in the multimode quantum Rabi model
aa r X i v : . [ qu a n t - ph ] A p r Quantum Zeno effect in the multimode quantum Rabi model
Shu He , , , ∗ Chen Wang , † Xue-Zao Ren , Li-Wei Duan , and Qing-Hu Chen ‡ Department of Physics and Electronic Engineering,Sichuan Normal University, Chengdu 610066, China School of Science, Southwest University of Science and Technology, Mianyang 621010, China Department of Physics, Zhejiang Normal University, Jinhua, Zhejiang 321004, China Department of Physics, Zhejiang University, Hangzhou 310027, P. R. China (Dated: April 9, 2019)We study the quantum Zeno effect (QZE) and quantum anti-Zeno effect (QAZE) of the multimodequantum Rabi model(MQRM). We derive an analytic expression for the decay rate of the survivalprobability where cavity modes are initially prepared as thermal equilibrium states. A crossoverfrom QZE to QAZE is observed due to the energy backflow induced by high frequency cavity modes.In addition, we apply a numerically exact method based on the thermofield dynamics(TFD) theoryand the matrix product states(MPS) to study the effect of squeezing of the cavity modes on theQZE of the MQRM. The influence of the squeezing angle, squeezing strength and temperature onthe decay rate of the survival probability are discussed.
PACS numbers: 03.65.Ge, 02.30.Ik, 42.50.Pq
I. INTRODUCTION
Frequent measurements on a quantum system maycause its dynamical evolution slow down or acceler-ated. This phenomenon, known as the quantum Zenoeffect(QZE) or quantum anti-Zeno effect (QAZE)[1–4],has been observed in experiments of trapped ion[5], cav-ity QED[6] and nuclear spin ensembles[7]. QZE has alsobeen considered as a powerful strategy to implement thequantum control, including quantum communication[8,9], quantum information protection[10], decoherencesuppression[11], purification and cooling[12].Recently, the QZE has been observed in the circuit-QED system where a superconducting flux qubit is cou-pled to a transmission line[13]. Compared to the tradi-tional optical experimental domain, the time scales in-volved in the circuit-QED system are much larger andcan be resolved with ordinary electronics. Moreover,the ultrastrong-coupling strength between the the qubitand the resonator has already been realized in circuit-QED systems[14], making it possible to observe QZE andQAZE in the strong coupling regime. As one of con-sequences in this regime, additional modes of the elec-tromagnetic resonator become increasingly relevant[15].To capture such multiply modes effect, the widely usedquantum Rabi model[16], describing an qubit interactingwith a single electromagnetic mode, has to be generalizedto its multimode version, known as multimode quantumRabi model(MQRM)[17, 18]. Recent studies have shownthat such higher-lying electromagnetic modes of the res-onator has a profound impact on various quantum op-tical phenomenons in the strong coupling regime[17–21]. ∗ Electronic address: [email protected] † Electronic address: [email protected] ‡ Electronic address: [email protected]
Thus, it is an interesting question that what role do mul-tiple modes play in the quantum Zeno dynamics of theMQRM. In addition, previous studies have shown thatthe squeezing of the resonator mode has a significant in-fluence in the quantum Zeno dynamics of a qubit in thequantum Rabi model[22]. It is natural to raise a questionthat how does the squeezing of multiple resonator modesaffect the QZE and QAZE of the MQRM.To address the above problems, we exploit a numer-ical exact method based on the matrix product statesand time-dependent variational principle[23] to study theQZE and QAZE of the MQRM. By restricting the evo-lution in the single excitation subspace under the frame-work of TFD, we derive an analytical expression for thedecay rate of survival probability when multiply modesare initially prepared as a thermal equilibrium states.We observe a crossover between the QZE to the QAZEunder repeated projective measurements. By numeri-cally calculating the energy transport between the TLSand multiply modes, we show that this crossover is at-tributed to the energy back flow from the high frequencymodes to the qubit. Moreover, we generalize the initialstate of multiple modes to a squeezed thermal state andstudy effects of squeezing phase angle and amplitude tothe QZE of the MQRM. We find that the decay of sur-vival probability is accelerated by non-vanishing squeez-ing strengthes. Squeezing angles also significantly affectthe decay rate. Particularly, high frequency modes inthe MQRM cause a positive shift on the critical squeez-ing angles where the decay rate reaches its extremum.This paper is organized as follows. In Sec II, we brieflyintroduce the MQRM and the numerical exact methodwe used to obtain the evolution of the MQRM at finitetemperature. In Sec III, we study QZE the MQRM atfinite temperature. We discuss the crossover of the sur-vival probability decay rate from QZE to QAZE and itsrelation to the energy transport between the qubit andmultiple cavity modes. Then, we extend the initial stateof multiple cavity modes to squeezed thermal states anddiscuss the effect of the squeezing on the QZE in Sec IV.We close this paper with a short summary in Sec V.
II. MODEL AND METHODOLOGYA. The multimode quantum Rabi model
In this paper, we focus on the MQRM which de-scribes a qubit or a two-level system(TLS) interactingwith multiple quantized photonic modes(e.g., cavity orresonator). The Hamiltonian of MQRM can be writtenas[18, 24]( ~ = 1): H = H + H B + H int H = ∆2 σ z , H B = M − X m =0 ω m a † m a m H int = M − X m =0 g m σ x ( a † m + a m ) (1)where σ x , σ z are standard pauli operators and a † m , a m arethe creation and annihilation boson operators for the m thmode with the frequency ω m = ( m + 1) ω and couplingstrength g m = √ m + 1 g to the TLS. M is the total num-ber of modes and depends on the specific physical real-ization. If not mentioned, we set M = 15 throughout thisstudy which is already feasible in the current circuit-QEDexperiment[17]. Since we are interested in the resonantsituation, ∆ = ω is assume in this study.To study the QZE of the MQRM, we focus on thesurvival probability of the TLS under successive idealprojective measurements with operator σ z . The initialstate is assumed as the product state: the TLS is at theexited state | ↑i and cavity modes are squeezed thermalstates: ρ (0) = | ↑ih↑ | ⊗ ρ B , ρ B = S ( ξ ) (cid:0) Z ( β ) − e − βH B (cid:1) S † ( ξ )(2)here Z ( β ) = Tr e − βH B is the partition function. β =1 / ( k B T ) where T is the temperature and k B is the Boltz-mann constant. S ( ξ ) is the squeezing operator: S ( ξ ) = Y m S m ( ξ m ) = Y m exp (cid:18) − ξ ∗ m a m + 12 ξ m ( a † m ) (cid:19) (3)where ξ m = r m e iφ m is the squeezing parameter for the m th mode. For simplicity, we assume a uniform squeez-ing for all M modes, i.e. ξ m = ξ = re iφ . B. MPS-based numerical method
In the following, we briefly introduce the numericallyexact method based on TFD and MPS to simulate the dynamics. We first remove the squeezing operator in theinitial state (2) by applying the unitary transform S ( ξ )to the Hamiltonian (1): H ′ = S † ( ξ ) HS ( ξ ) = H + H ′ B + H ′ int H ′ B = M − X m =0 ω m (cid:0) Aa † m a m + Ba † m + B ∗ a m (cid:1) H ′ int = M − X m =0 g m (cid:0) Ka † m + K ∗ a m (cid:1) σ x (4)where A = cosh r + sinh r, B = e iφ cosh r sinh r, K =cosh r + e iφ sinh r .According to the TFD method[25, 26], the evolution ofthis transformed Hamiltonian (4) from the initial state ofthermalized cavity modes is equivalent to a Schr¨odingerequation of a modified Hamiltonian ˆ H defined in an en-larged Hilbert space(see Appendix):ˆ H = H + ˆ H B + ˆ H int ˆ H B = M − X m =0 ω m (cid:2) A (cid:0) a † m a m − b † m b m (cid:1) + (cid:0) Ba † m − Bb m + H.c (cid:1) (cid:3) ˆ H int = M − X m =0 g m ( Ka † m cosh θ m + Kb sinh θ m ) σ x + H.c (5)The the vacuum initial state is | ˆ ψ (0) i = | ↑i ⊗ M − Y m =0 | i a m | i b m (6)where b † m , b m are boson operators of the fictitious modesand θ m = arctanh( e − βω m / ) . The evolution of the ex-pectation value of an arbitrary operator O that affects inoriginal Hilbert space can be straightforwardly calculated: h O ( t ) i = h ˆ ψ ( t ) | O | ˆ ψ ( t ) i (7)Finally, the TFD Schr¨odinger equation is simulated byMPS-based numerically exact method [23, 27, 28] whichhas been widely used and proved to be highly efficient insolving quantum many body dynamics[29, 30]. C. Survival probability and effective decay rate
Now we turn to the QZE and QAZE of the MQRM.The QZE [1] can be described by the survival probability P sur ( nτ ) which is defined as the probability of finding theinitial state after n successive measurements with equaltime interval τ . The measurement considered in this pa-per is assumed to be an ideal projecting measurementsof operator σ z followed by a post selection regarding tothe positive measurement result.The survival probabilitycan be written as: P sur ( nτ ) = Tr h Pe − i ˆ Hτ ρ (0) e i ˆ Hτ P i n = P n sur ( τ ) (8)where P = | ↑ih↑ | is measurement projecting operator.In the short interval time limit τ →
0, one can furtherwrite P sur ( t = nτ ) in an exponentially decay form[31, 32]: P sur ( nτ ) = exp( − γ ( τ ) t ) (9)where γ ( τ ) is an effective decay rate: γ ( τ ) = − τ ln [ P sur ( τ )] (10)Eq.(10) is valid when the measurement interval τ is rel-atively short where the measurement disturbance to theenvironment(here is cavity modes) can be neglected thusthe state after each projecting measurement collapse tothe identical state as Eq. (6)[33]. Therefore, we restrict gτ < γ ( τ ) is a crucial quantity tocharacterize the QZE and the QAZE[34, 35]: ∂γ ( τ ) /∂τ > ∂γ ( τ ) /∂τ < γ ( τ ) /γ < γ = γ ( τ → ∞ ) is the natural decay rate[34],new definitions through ∂γ/∂τ retains the core physicalpicture of QZE and QAZE without calculating γ thatmay not exist in some models[33, 36–38]. Throughoutthis paper, we use this new criterion to classify QZE andQAZE, the potential crossover point is denoted as τ c , i.e. ∂γ ( τ ) /∂τ (cid:12)(cid:12) τ c = 0 . III. RESULTS AND DISCUSSIONA. Thermal equilibrium initial state
We first consider the QZE when cavity modes are pre-pared as the thermal equilibrium state. The survivalprobability of Hamiltonian (5) with the product initialstate (6) can be written as[22]: P sur ( τ ) = Tr (cid:2) ρ ( τ ) | e ih e | (cid:3) = Tr (cid:2) e − i ˆ Hτ ρ (0) e i ˆ Hτ | e ih e | (cid:3) = Tr " ∞ X n =0 ( − iτ ) n ! [ ˆ H, ρ (0)] n | e ih e | (11)where ρ = | ˆ ψ (0) ih ˆ ψ (0) | and [ ˆ H, ρ (0)] n =[ ˆ H, [ ˆ H, ρ (0)] n − ]. A commonly used approximationto obtain the decay rate of the survival probability is tokeep terms in (11) up to τ , this leads to the decay of P sur ( τ ) in a quadratical form: P e ( τ ) ≈ − (cid:18) ττ z (cid:19) (12) ω τ -1 β ω = ∞ β ω =0.5 β ω =0.2 β ω =0.1 ω τ -6 -4 -2 γ ( τ ) g/ ω = 0.01 M =1g/ ω = 0.01 M =15g/ ω = 0.2 M =1g/ ω = 0.2 M =15 (b)(a) FIG. 1: (Color online) Effective decay rate of the survivalprobability γ ( τ ) for thermalized initial states. Results withlines are obtained from Eq.(15) while circles are results fromnumerical exact MPS-TDVP method. (a) Comparison of de-cay rates between single mode and multiple modes at zerotemperature. (b)effective decay rate at different temperaturefor g/ω = 0 . where τ z is known as ”Zeno time”[39]. Thus for N re-peated projective measurement with equal interval τ , thesurvival probability can be approximated as an exponen-tial decay: P ( t = N τ ) = (cid:20) − (cid:18) ττ z (cid:19) (cid:21) N ≈ e − γ ( τ ) t (13)with the effective decay rate γ ( τ ): γ ( τ ) = τ X m g m (cid:0) cosh ( θ m ) + sinh ( θ m ) (cid:1) (14)Since Eq.(14) is independent of ∆ and ω m , its validityis limited to the case where the measurement interval τ is much shorter than the typical time scale of each mode , that is τ < /ω m . When high frequency modes areincluded, it severely precludes the availability of γ ( τ )compared to the Rabi model where only single resonantmode is considered.Instead, we employ an alternative method based onthe thermofield dynamics(TFD)[25, 40]. This methodtransforms the evolution of a mixed thermalized initialstate into another evolution with a pure initial state inan enlarged Hilbert space whose total number of degreesof freedom is double compared to the original one. Re-cently, this method has received much attention and hasbeen widely used to study finite temperature dynamicsof quantum electron-vibrational systems[41–43] and openquantum systems[44, 45]. Details of TFD can be foundin Ref[26, 41].By restricting the evolution of the TFD Schr¨odingerequation to the single expiation subspace, we can obtainan analytic expression for the decay rate of the survivalprobability(see Appendix): γ th ( τ ) = τ M X m =0 cosh θ m g m sinc (cid:20) τ ( ω k − ∆)2 (cid:21) + sinh θ m g m sinc (cid:20) τ ( ω k + ∆)2 (cid:21) (15)where sinc(x) ≡ sin( x ) /x . Different from the γ inEq(14), γ th shows apparent dependence of frequencies ω k and ∆ and has a leading scaling of g k / (∆ − ω k ) whichis similar to the decay rate of the spontaneous emis-sion in the multimode Percell effect due to non-resonantmodes[46].In Figure(1), the analytical result γ th ( τ ) agree wellwith numerical exact ones in a large range of param-eters including strong coupling and high temperatureregimes. Specifically, for zero temperature case( Fig-ure(1.a)) the decay rate γ ( τ ) for a single mode Rabimodel shows a monotonic increasing with the increase ofmeasurement interval τ for both weak( g/ω = 0 .
01) andstrong( g/ω = 0 .
2) coupling strengthes, clearly demon-strates a pure Zeno effect. However, by including thehigh-lying photonic modes( M = 15), the decay rate ofsurvival probability presents a nonmonotonic behavior,which increases with τ before ω τ < . γ ( τ ) can be re-lated to the population of the TLS in the excited state | ↑i . Since the energy of the TLS is E TLS ( t ) = ∆ h σ z ( t ) i / P sur ( τ ) = ( h σ z ( τ ) i + 1) / E TLS ( τ ) / ∆ + 1 / dγ ( τ ) dτ <
0, this re-quires: dγ ( τ ) dτ < ⇔ − τ (cid:18) d ln( P sur ) dτ − ln( P sur ) τ (cid:19) < gτ < E TLS ( τ ) can be approximatelyexpressed as: E TLS ( τ ) / ∆ = 1 / − a ( gτ ) − a ( gτ ) + O ( gτ ) (18)where a > E TLS is maximized at τ = 0. Inserting Eq.(18)into Eq.(17), we obtain: dγ ( τ ) dτ ≈ dγdτ (cid:12)(cid:12)(cid:12)(cid:12) τ =0 + 43 τ a g (cid:18) dγdτ (cid:12)(cid:12)(cid:12)(cid:12) τ =0 + 12 a g (cid:19) + O ( τ ) < ω t M -0.200.2 -0.2-0.15-0.1-0.050 γ ( τ ) dE TLS dt energy back flow g/ω = 0 . β ω = 0 . (a)(b) FIG. 2: (Color online) Energy transport rate between theTLS and multiple modes. (a) Decay rate γ ( τ ) and total en-ergy transport rate of the TLS at gω = 0 . , βω = 0 .
5. (b)Evolution of mode’s energy operator h E m ( τ ) i with same pa-rameters as (a). For arbitrary system, dγ/dτ | τ =0 > τ →
0, the occurrence of theQAZE requires a <
0, i.e. d E TLS /d τ >
0. As depictedin Figure(2).a, the crossover point τ c is close to the pole ofthe energy transport rate as d E TLS /d τ = 0. As shownin (2), by further calculating the evolution of the energyoperator for each mode E m ( τ ) = ω m h ψ ( τ ) | a † m a m | ψ ( τ ) i ,we observe an energy backflow from the higher frequencyparts of cavity modes to the TLS during the crossoverpoint τ c . It is this energy backflow induced by high fre-quency modes that leads to d E TLS /d τ > and makesit possible to observe QAZE. On the contrary, the lowfrequency part of modes always absorbs energy from theTLS during gτ <
1, leading to the pure QAZE in thesingle mode Rabi model as shown in Figure(1).a.
B. Squeezed thermal initial states
In this section, we generalize the initial state of thecavity modes to the squeezed thermal states .From the view of the Hamiltonian (5), effects of thesqueezing come from two aspects: 1).The renormaliza-tion of mode frequencies and effective coupling strengthesby factors A and K respectively. 2).The nonharmonicterms in ˆ H B with the factor B .The effect of 1) can be included in the analytic expres-sion (15) by replacing g m → | K | g m and ω m → Aω m .By increasing r , the effective coupling strength is in-creased, leading to a larger decay rate. Moreover, due φ / π γ ( τ ) × -3 φ / π γ ( τ ) / g ω τ φ / π × -3 M = 15M = 9M = 5M=1 r = 0r = 0.1r = 0.3r = 0.5 (c)(b)(a)
FIG. 3: (Color online) Relation between squeezing angles andthe decay rates of survival probability γ ( τ ). (a). Decay ratesof different squeezing angles for various numbers of cavitymodes at zero temperature, the maximal and minimal crit-ical angles are marked by solid symbols, other parametersare g/ω = 0 . , ω τ = 0 . , r = 0 .
3. (b). Decay rates as afunction of squeezing angles at zero temperature for differentsqueezing strength r , data denoted by lines are obtained from g/ω = 0 .
01 while hollows are obtained from g/ω = 0 . ω τ = 0 .
1. (c) Decay rates for differentmeasurement interval τ and squeezing angle φ , other param-eters are g/ω = 0 . , r = 0 . to the factor | K | , γ ( τ ) depends on the squeezing angle φ . Since | K | = cosh(2 r ) + cos φ sinh(2 r ), the critical an-gle where γ ( τ ) takes the extreme value may be φ max = 0and φ min = π respectively. Same critical angles havebeen observed on the single-mode rabi model with ini-tial squeezed field[22] and the TLS in a squeezed bathunder RWA[47]. However, this is not the case in thecurrent MQRM where effects of both multimodes andnon-rotating terms are taken into consideration. Specifi-cally, by increasing the number of modes from M = 1 to M = 15, critical angles θ min and θ max gradually increasewith a fixed angle difference θ min − θ max = π as shown inFig(3.a). Interestingly, as depicted Fig(3.b), such criticalangles are independent of squeezing strength r and cou-pling strength g . This independence agrees with resultsof the TLS in squeezed bath[47]. Surprisingly, the scalingof γ ( τ ) ∼ g for thermal equilibrium initial states shownin Eq(15) may still holds for squeezed thermal bath bynoting that γ ( τ ) /g as a function of φ for different cou-pling strengthes g in Fig(3.b) collapses to the same curve.The measurement interval τ also affects such criticalangles, as shown in Fig(3.c). Note when τ →
0, the shiftof critical angles approaches zero. The crossover fromQZE to QAZE can still be observed for the squeezedthermal initial state. In addition, a specific squeezingangle can either highlights or downplay this transition: ω τ -0.500.51 l og ( β ω ) × -3 φ / π γ ( τ ) × -3 β ω = 1.0 β ω = 0.5 β ω = 0.3 (b)(a) FIG. 4: (Color online) Finte temperature effect on the decayrate γ ( τ ). (a) γ ( τ ) as a function of squeezing angle φ at differ-ent temperatures with ω τ = 0 .
1. (b) The decay rate γ ( τ ) forvarious measurement interval τ and temperature, φ = π/ g/ω = 0 .
01 . For both subplots, g/ω = 0 . r = 0 . as shown in Fig(3.c), during transition point τ c , γ ( τ ) aresignificantly increased for φ = φ max and suppressed for φ = φ min , while with the increase of τ , the differencebetween γ ( τ ) at two critical angles is narrowed.Finally, the effect of temperature on the decay rate forthe squeezed thermal initial state is described by θ m inˆ H int where higher temperature leads to larger effectivecoupling strength and consequently larger decay rate, asshown Fig(4.b). Moreover, increasing temperature alsocauses reduction of the shift of critical angles. This isnot a surprising result. As explained above, the shiftof critical angles are related to the Ba † , m + B ∗ a m in H B .However, with the increase of temperature, ˆ H int becomesdominant compared to the nonharmonic term, thus canbe neglected in the high-temperature limit and the re-lation between γ ( τ ) and φ are only determined through | K | . IV. CONCLUSION
In summary, we study the QZE and QAZE of theMQRM where cavity modes are prepared in thermal equi-librium state. We derive an analytic expression of the de-cay rate of the survival probability which shows a transi-tion from QZE to QAZE. By calculating the evolution ofenergy flow between the TLS and cavity modes, we showthat such QZE-QAZE transition is related to the energybackflow from high frequency modes to the TLS. Fur-thermore, we generalize the initial state of cavity modesto squeezed thermal states and study effects of squeez-ing angle and strength by applying a numerically exactmethod. We find that the decay of survival probabil-ity is accelerated by non-vanishing squeezing strengthes.Moreover, squeezing angles also significantly affect thedecay rate: compared to θ max = 0 and θ min = π inwhich decay rates take maximal and minimal values inthe single-mode Rabi model, high frequency modes in theMQRM cause a positive shift on these critical squeezingangles. Finally, with the increase of the temperature, thedecay of the survival probability is accelerated and theshift of critical squeezing angles are reduced. V. ACKNOWLEDGEMENT
Shu He is supported by the National Natural Sci-ence Foundation of China under Grant No. 11804240.Chen Wang is supported by the National Natural Sci-ence Foundation of China under Grant Nos. 11704093and 11547124. Liwei Duan and Qing-Hu Chen are sup-ported by the National Natural Science Foundation ofChina under Grant Nos.11674285 and 11834005.
Appendix A: A Brief introduction to TFD
We briefly introduce the theory of thermofield dynam-ics in this Appendix.The time evolution of an arbitrary Hamiltonian H atfinite temperature can be described by the Liouville-vonNeumann equation( ~ = 1): ∂ t ρ ( t ) = − i [ H, ρ ( t )] (A1)where the initial state is assumed to be at thermal equi-librium: ρ (0) = 1 Z ( β ) exp ( − βH ) (A2)here Z ( β ) = Tr e − βH is the partition function and β =( k B T ) − . k B is the Boltzmann constant and T is thetemperature.According to Ref[48], the time evolution of such den-sity matrix (A1) can be equivalently described by amodified Hamiltonian ˆ H = H − ˜ H and correspond-ing Schr¨odinger equation defined in an enlarged Hilbertspace( ~ = 1): i ∂∂t | ψ ( t ) i = ˆ H | ψ ( t ) i (A3)˜ H represents the fictitious Hamiltonian which can be de-rived from the original Hamiltonian H [25].By defining the vector | I i : | I i = X k | k i| ˜ k i (A4)where | k i ( | ˜ k i ) are arbitrary basis vectors of the origi-nal(fictitious) space, the density matrix ρ ( t ) and wavefunction | ψ ( t ) i have the following relation[49]: | ψ ( t ) i = ρ ( t ) / | I i (A5) Particularly, the initial state at thermal equilibrium is: | ψ (0) i = ρ (0) / | I i = Z ( β ) − / exp (cid:18) − βH (cid:19) | I i (A6)The expectation value for an arbitrary operator A de-fined in the original Hilbert space {| k i} can be calculatedby: h A ( t ) i = h ψ ( t ) | A | ψ ( t ) i = Tr { ρ ( t ) A } (A7)Thus the evaluation of h A ( t ) i by using the TFD wavefunction | ψ ( t ) i is equivalent to the that by using the den-sity matrix ρ ( t ).Now let us turn to the Hamiltonian of the MQRM (1).To describe its evolution from the initial state (2) underthe framework of the TFD, we first remove the squeez-ing operator in the initial state by applying a unitarytransform S ( ξ ) defined in (3) to the Hamiltonian (1): H ′ = S † ( ξ ) HS ( ξ ) = H + H ′ B + H ′ int H ′ B = M − X m =0 ω m (cid:0) Aa † m a m + Ba † m + B ∗ a m (cid:1) H ′ int = M − X m =0 g m (cid:0) Ka † m + K ∗ a m (cid:1) σ x (A8)where A = cosh r + sinh r, B = e iφ cosh r sinh r, K =cosh r + e iφ sinh r , this is the Hamiltonian (4) in themain text.We can choose the fictitious Hamiltonian ˜ H as :˜ H = − X m ω m (cid:0) Ab † m b m + B ∗ b † m + Bb m (cid:1) (A9)where b † m ( b m ) is the creation (annihilation) boson oper-ator defined in the fictitious Hilbert space correspondingto the the mode a m . This leads to the total Hamiltonianfor the TFD schr¨odinger equation:ˆ H = ∆2 σ z + X m ω m (cid:0) Aa † m a m + Ba † m + B ∗ a m (cid:1) − X m ω m (cid:0) Ab † m b m + B ∗ b † m + Bb m (cid:1) + M − X m =0 g m (cid:0) Ka † m + K ∗ a m (cid:1) σ x (A10)The corresponding initial TFD wave function can bewritten as: | ψ (0) i = | ↑i ⊗ | ( β ) i (A11)where | ( β ) i is often referred to as the thermal vacuumstate: | ( β ) i = Z ( β ) − / e − βH B / | I i (A12)This thermal vacuum state can be further transformedto the ground state(vacuum state) of the modes { a m } and { b m } by applying the so-called Bogoliubov thermaltransformation[50]: | ( β ) i = e − iG | i (A13)where | i = Q M − m =0 | i a m | i b m and the transformationoperator G is defined as: G = − i X m θ m ( a m b m − a † m b † m ) (A14)where θ m = arctanh( e − βω m / ) (A15)Instead of solving the TFD schr¨odinger equation (A10)with the initial state (A12) , one can apply the inverseBogoliubov thermal transformation to the Hamiltonianˆ H . By considering the following new relations: e iG a m e − iG = a m cosh( θ m ) + b † m sinh( θ m ) e iG b m e − iG = b m cosh( θ m ) + a † m sinh( θ m ) (A16)We have: | ˆ ψ ( t ) i = e iG | ψ ( t ) i , i ∂∂t | ˆ ψ ( t ) i = ˆ H θ | ˆ ψ ( t ) i (A17)ˆ H θ = e iG ˆ He − iG is the final Hamiltonian (5) in the maintext:ˆ H θ = H + ˆ H B + ˆ H int H = ∆2 σ z ˆ H B = M − X m =0 ω m (cid:2) A (cid:0) a † m a m − b † m b m (cid:1) + (cid:0) Ba † m − Bb m + H.c (cid:1) (cid:3) ˆ H int = M − X m =0 g m ( Ka † m cosh θ m + Kb sinh θ m ) σ x + H.c(A18)and corresponding initial state (6): | ˆ ψ (0) i = | i ⊗ | i = | ↑i ⊗ M − Y m =0 | i a m | i b m (A19) Appendix B: Derivation of Eq.(15)
In this Appendix, we show the detailed derivation tothe analytic expression of the decay rate γ th ( τ ).The evolution of the MQRM from a product initialstate where resonator modes are at thermal equilibriumis govern by the TFD Schr¨odinger equation of the Hamil-tonian ( A
18) with r = 0(without squeezing):ˆ H th = ∆2 σ z + X m ω m ( a † m a m − b † m b m )+ g m (cid:18) cosh θ m ( a † m + a m ) + sinh θ m ( b † m + b m ) (cid:19) σ x (B1) To obtain an analytic expression, we further apply anapproximation by restricting the evolution in the singleexcitation(SE) Hilbert space. This can be achieved bywriting the wave function at arbitrary time t as: | ψ ( t ) i = χ ( t ) | ↑i| i + X m p m ( t ) | ↓i| i a m | i b m + X m q m ( t ) | ↓i| i a m | i b m (B2)Thus the Hamiltonian (B1) can be simplified as:ˆ H th,SE = ∆2 σ z + X m ω m ( a † m a m − b † m b m )+ g m (cid:18) cosh θ m ( σ − a † m + σ + a m ) + sinh θ m ( σ − b † m + σ + b m ) (cid:19) (B3)where σ ± = ( σ x ± iσ y ).Inserting the wave function (B2) into the TFDSchr¨odinger equation, we have: i dχ ( t ) dt = ∆2 χ ( t ) + X m g m cosh θ m p m ( t ) + g m sinh θ m q m ( t ) i dp m ( t ) dt = X m ( ω m − ∆2 ) p m ( t ) + g m cosh θ m χ ( t ) i dq m ( t ) dt = X m ( − ω m − ∆2 ) q m ( t ) + g m sinh θ m χ ( t ) (B4)applying the following transformation: χ ( t ) = ˜ χ ( t ) exp( − i ∆2 t ) p m ( t ) = ˜ p m ( t ) exp (cid:20) − i (cid:18) ω m − ∆2 (cid:19) t (cid:21) q m ( t ) = ˜ q m ( t ) exp (cid:20) − i (cid:18) − ω m − ∆2 (cid:19) t (cid:21) (B5)we have: i dχ ( t ) dt = X m g cosh θ m ˜ p m ( t ) exp[ − i ( ω m − ∆) t ]+ g m sinh θ m ˜ q m ( t ) exp[ − i ( − ω m − ∆) t ] i d ˜ p ( t ) dt = X m g m cosh θ m ˜ χ ( t ) exp[ i ( ω m − ∆) t ] i d ˜ q ( t ) dt = X m g m sinh θ m ˜ χ ( t ) exp[ i ( − ω m − ∆) t ] (B6)Integrating last two equations in (B6) and substitutinginto the first one, we have: d ˜ χ ( t ) dt = − X m cosh θ m g m Z t ˜ χ ( t ′ ) exp[ − i ( ω m − ∆)( t − t ′ )] dt ′ − sinh θ m g m Z t ˜ χ ( t ′ ) exp[ − i ( − ω m − ∆)( t − t ′ )] dt ′ (B7)The solution of the integro-differential equation abovecan be obtained iteratively. Since we are limited to thefrequent measurement limit where gτ <
1, only the shorttime behavior of ˜ χ ( t ) is concerned. This can be obtainedby taking first iteration and inserting the initial condition χ (0) = 1:˜ χ ( t ) ≈ − X m Z t ( t − t ′ ) cos θ m g m exp[ − i ( ω − ∆) t ′ ] dt ′ − X m Z t ( t − t ′ ) sinh θ m g m exp[ − i ( − ω − ∆) t ′ ] dt ′ ≈ exp (cid:26) − t X m g m (cid:2) cosh θ m F (∆ , ω ) + sinh θ m F (∆ , − ω ) (cid:3)(cid:27) (B8)where F (∆ , ω ) = 2 sin ( ω − ∆2 ) tt ( ω − ∆) − i ( ω − ∆) t − sin( ω − ∆) tt ( ω − ∆) (B9)Substituting (B8) into the definition of the decay rateof the survival probability (10) , we finally obtain: γ th ( τ ) = τ M X m =0 cosh θ m g m sinc (cid:20) τ ( ω k − ∆)2 (cid:21) + sinh θ m g m sinc (cid:20) τ ( ω k + ∆)2 (cid:21) (B10)where sinc(x) ≡ sin( x ) /x . Appendix C: MPS-based numerically exact method
As derived in the previous Appendix, the evolution ofthe MQRM from a squeezed thermal states can be de-scribed by the TFD schr¨odinger equation of the Hamil-tonian (5) with the initial state (6). To implementMPS-based numerical method, the Hamiltonian is rewrit-ten in a compact form of the following matrix productoperator(MPO)[28]: ω t h σ + σ - i D max = 5 , N max = 10 D max = 7 , N max = 20 D max = 10 , N max = 30 D max = 14 , N max = 50 r = 0 . β ω = 0 . r = 0 . β ω = ∞ FIG. 5: (Color online) The convergence of the MPS-TDVPmethod for the population evolution of the MQRM. Otherparameters are g/ω = 0 . , φ = π/ , M = 15. ˆ H = M − Y i =0 W [ i ] b ! W s M − Y j =0 W [ j ] a W [ m ] b = I I H b m , G b m , I W s = I σ x I σ z / , σ x , I W [ m ] a = I G a m I H a m , , I (C1)where I is the identity operator of the TLS and the bosonspace H a m = Aω m a † m a m + Ba † m + B ∗ a m H b m = − Aω m b † m b m − Ba m + B ∗ a † m G a m = g m Ka † m cosh θ m + H.c G b m = g m Kb sinh θ m + H.c (C2)Note for the open boundary condition, only first rowof W [1] b and first column of W [ M ] a are used. With thisMPOs, we employ the recently proposed time-dependentvariational principle(TDVP) method[23] to simulate theevolution. For numerical results throughout this study,the convergence can be achieved by setting the maximalbound dimensions of the MPS D max ≤
15 and the cutoffboson number for each cavity mode N max ≤
80. Thetime step of simulations throughout this paper is set to gδt = 10 − . [1] Baidyanath Misra and EC GEORGE Sudarshan. Thezenos paradox in quantum theory. Journal of Mathemat-ical Physics , 18(4):756–763, 1977.[2] Paolo Facchi and Saverino Pascazio. Quantum zeno dy-namics: mathematical and physical aspects.
Journal ofPhysics A: Mathematical and Theoretical , 41(49):493001,2008.[3] F Francica, F Plastina, and S Maniscalco. Quantum zenoand anti-zeno effects on quantum and classical correla-tions.
Physical Review A , 82(5):052118, 2010.[4] Florian Sch¨afer, Ivan Herrera, Shahid Cherukattil,Cosimo Lovecchio, Francesco Saverio Cataliotti, FilippoCaruso, and Augusto Smerzi. Experimental realizationof quantum zeno dynamics.
Nature communications , 5,2014.[5] Wayne M Itano, Daniel J Heinzen, JJ Bollinger, andDJ Wineland. Quantum zeno effect.
Physical ReviewA , 41(5):2295, 1990.[6] Julien Bernu, Samuel Del´eglise, Cl´ement Sayrin, StefanKuhr, Igor Dotsenko, Michel Brune, Jean-Michel Rai-mond, and Serge Haroche. Freezing coherent field growthin a cavity by the quantum zeno effect.
Physical reviewletters , 101(18):180402, 2008.[7] Gonzalo A ´Alvarez, DD Bhaktavatsala Rao, Lucio Fry-dman, and Gershon Kurizki. Zeno and anti-zeno polar-ization control of spin ensembles by induced dephasing.
Physical review letters , 105(16):160401, 2010.[8] Sima Pouyandeh, Farhad Shahbazi, and Abolfazl Bayat.Measurement-induced dynamics for spin-chain quantumcommunication and its application for optical lattices.
Physical Review A , 90(1):012337, 2014.[9] Abolfazl Bayat and Yasser Omar. Measurement-assistedquantum communication in spin channels with dephas-ing.
New Journal of Physics , 17(10):103041, 2015.[10] Adriano Barenco, Andre Berthiaume, David Deutsch,Artur Ekert, Richard Jozsa, and Chiara Macchiavello.Stabilization of quantum computations by symmetriza-tion.
SIAM Journal on Computing , 26(5):1541–1557,1997.[11] Almut Beige, Daniel Braun, Ben Tregenna, and Peter LKnight. Quantum computing using dissipation to remainin a decoherence-free subspace.
Physical Review Letters ,85(8):1762, 2000.[12] Noam Erez, Goren Gordon, Mathias Nest, and GershonKurizki. Thermodynamical control by frequent quantummeasurements. arXiv preprint arXiv:0804.2178 , 2008.[13] K Kakuyanagi, T Baba, Y Matsuzaki, H Nakano, S Saito,and K Semba. Observation of quantum zeno effect ina superconducting flux qubit.
New Journal of Physics ,17(6):063035, 2015.[14] P Forn-D´ıaz, Juan Jos´e Garc´ıa-Ripoll, Borja Peropadre,J-L Orgiazzi, MA Yurtalan, R Belyansky, Christopher MWilson, and A Lupascu. Ultrastrong coupling of a singleartificial atom to an electromagnetic continuum in thenonperturbative regime.
Nature Physics , 13(1):39, 2017.[15] Sal J Bosman, Mario F Gely, Vibhor Singh, AlessandroBruno, Daniel Bothner, and Gary A Steele. Multi-modeultra-strong coupling in circuit quantum electrodynam-ics. npj Quantum Information , 3(1):46, 2017.[16] II Rabi. On the process of space quantization.
PhysicalReview , 49(4):324, 1936. [17] Neereja M Sundaresan, Yanbing Liu, DariusSadri, L´aszl´o J Sz˝ocs, Devin L Underwood, MoeinMalekakhlagh, Hakan E T¨ureci, and Andrew A Houck.Beyond strong coupling in a multimode cavity.
PhysicalReview X , 5(2):021035, 2015.[18] Mario F Gely, Adrian Parra-Rodriguez, Daniel Both-ner, Ya M Blanter, Sal J Bosman, Enrique Solano, andGary A Steele. Convergence of the multimode quantumrabi model of circuit quantum electrodynamics.
PhysicalReview B , 95(24):245115, 2017.[19] Carlos S´anchez Mu˜noz, Franco Nori, and Simone De Lib-erato. Resolution of superluminal signalling in non-perturbative cavity quantum electrodynamics.
NatureCommunications , 9, 2018.[20] Simone De Liberato. Light-matter decoupling in the deepstrong coupling regime: The breakdown of the purcelleffect.
Physical review letters , 112(1):016401, 2014.[21] Juan Jose Garcia-Ripoll, Borja Peropadre, and SimoneDe Liberato. Light-matter decoupling and a term de-tection in superconducting circuits. Scientific reports ,5:16055, 2015.[22] Ion Lizuain, Jorge Casanova, Juan Jos´e Garc´ıa-Ripoll,JG Muga, and Enrique Solano. Zeno physics inultrastrong-coupling circuit qed.
Physical Review A ,81(6):062131, 2010.[23] Jutho Haegeman, Christian Lubich, Ivan Oseledets, BartVandereycken, and Frank Verstraete. Unifying timeevolution and optimization with matrix product states.
Physical Review B , 94(16):165116, 2016.[24] Moein Malekakhlagh, Alexandru Petrescu, and Hakan ET¨ureci. Cutoff-free circuit quantum electrodynamics.
Physical review letters , 119(7):073601, 2017.[25] Masuo Suzuki. Density matrix formalism, double-spaceand thermo field dynamics in non-equilibrium dissipa-tive systems.
International Journal of Modern PhysicsB , 5(11):1821–1842, 1991.[26] Maxim F Gelin and Raffaele Borrelli. Thermalschr¨odinger equation: Efficient tool for simulation ofmany-body quantum dynamics at finite temperature.
Annalen der Physik , 529(12):1700200, 2017.[27] Ulrich Schollw¨ock. The density-matrix renormalizationgroup in the age of matrix product states.
Annals ofPhysics , 326(1):96–192, 2011.[28] Michael L Wall, Arghavan Safavi-Naini, and Ana MariaRey. Simulating generic spin-boson models with matrixproduct states.
Physical Review A , 94(5):053637, 2016.[29] Florian AYN Schr¨oder and Alex W Chin. Simulatingopen quantum dynamics with time-dependent variationalmatrix product states: Towards microscopic correlationof environment dynamics and reduced system evolution.
Physical Review B , 93(7):075105, 2016.[30] Benedikt Kloss, Yevgeny Bar Lev, and David Reichman.Time-dependent variational principle in matrix-productstate manifolds: Pitfalls and potential.
Physical ReviewB , 97(2):024307, 2018.[31] P Facchi, H Nakazato, and S Pascazio. From the quan-tum zeno to the inverse quantum zeno effect.
PhysicalReview Letters , 86(13):2699, 2001.[32] Sabrina Maniscalco, Jyrki Piilo, and Kalle-Antti Suomi-nen. Zeno and anti-zeno effects for quantum brownianmotion.
Physical review letters , 97(13):130402, 2006. [33] Shu He, Qing-Hu Chen, and Hang Zheng. Zeno and anti-zeno effect in an open quantum system in the ultrastrong-coupling regime. Physical Review A , 95(6):062109, 2017.[34] H Zheng, SY Zhu, and MS Zubairy. Quantum zeno andanti-zeno effects: without the rotating-wave approxima-tion.
Physical review letters , 101(20):200404, 2008.[35] Xiufeng Cao, JQ You, H Zheng, AG Kofman, and FrancoNori. Dynamics and quantum zeno effect for a qubit ineither a low-or high-frequency bath beyond the rotating-wave approximation.
Physical Review A , 82(2):022119,2010.[36] Adam Zaman Chaudhry and Jiangbin Gong. Zenoand anti-zeno effects on dephasing.
Physical Review A ,90(1):012101, 2014.[37] Yu-Ran Zhang and Heng Fan. Zeno dynamics in quantumopen systems.
Scientific reports , 5:11509, 2015.[38] Wei Wu and Hai-Qing Lin. Quantum zeno and anti-zenoeffects in quantum dissipative systems.
Physical ReviewA , 95(4):042132, 2017.[39] LS Schulman, A Ranfagni, and D Mugnai. Characteristicscales for dominated time evolution.
Physica Scripta ,49(5):536, 1994.[40] Hiroomi Umezawa, Hiroshi Matsumoto, and MasashiTachiki. Thermo field dynamics and condensed states.1982.[41] Raffaele Borrelli and Maxim F Gelin. Quantum electron-vibrational dynamics at finite temperature: Thermo fielddynamics approach.
The Journal of chemical physics ,145(22):224101, 2016.[42] Raffaele Borrelli and Maxim F Gelin. Simulation of quan-tum dynamics of excitonic systems at finite temperature: an efficient method based on thermo field dynamics.
Sci-entific reports , 7(1):9127, 2017.[43] Lipeng Chen and Yang Zhao. Finite temperature dynam-ics of a holstein polaron: The thermo-field dynamics ap-proach.
The Journal of chemical physics , 147(21):214102,2017.[44] Lu Wang, Yuta Fujihashi, Lipeng Chen, and Yang Zhao.Finite-temperature time-dependent variation with mul-tiple davydov states.
The Journal of chemical physics ,146(12):124127, 2017.[45] Shu He, Chen Wang, Li-Wei Duan, and Qing-Hu Chen.Zeno effect of an open quantum system in the presenceof 1/f noise.
Physical Review A , 97(2):022108, 2018.[46] AA Houck, JA Schreier, BR Johnson, JM Chow, JensKoch, JM Gambetta, DI Schuster, L Frunzio, MH De-voret, SM Girvin, et al. Controlling the spontaneousemission of a superconducting transmon qubit.
Physicalreview letters , 101(8):080502, 2008.[47] DF Mundarain and J Stephany. Zeno and anti-zeno effectfor a two-level system in a squeezed bath.
Physical ReviewA , 73(4):042113, 2006.[48] Masuo Suzuki. Thermo field dynamics of quantum spinsystems.
Journal of statistical physics , 42(5-6):1047–1070, 1986.[49] Masuo Suzuki. Thermo field dynamics in equilibrium andnon-equilibrium interacting quantum systems.
Journal ofthe Physical Society of Japan , 54(12):4483–4485, 1985.[50] Yasushi Takahashi and Hiroomi Umezawa. Thermo fielddynamics.