# Nonclassical Radiation from Thermal Cavities in the Ultrastrong Coupling Regime

###### Abstract

Thermal or chaotic light sources emit radiation characterized by a slightly enhanced probability of emitting photons in bunches, described by a zero-delay second-order correlation function . Here we explore photon-coincidence counting statistics of thermal cavities in the ultrastrong coupling regime, where the atom-cavity coupling rate becomes comparable to the cavity resonance frequency. We find that, depending on the system temperature and coupling rate, thermal photons escaping the cavity can display very different statistical behaviors, characterised by second-order correlation functions approaching zero or greatly exceeding two.

###### pacs:

42.50.Pq, 42.50.Ar, 85.25.-j, 03.67.LxThermal radiation has a special place in modern physics. In the search for a solution to the discrepancies between the observed energy spectrum of thermal radiation and the predictions of classical theory, Planck was led to introduce the revolutionary concept of quanta Planck (). Thermal or chaotic light sources emit radiation that is characterized by an enhanced probability of emitting photons in bunches HBT (). In the course of the successful attempt to explain this effect, Glauber established the basis of quantum optics Glauber (); NobelGlauber (). Even more recently the study of thermal emission has continued to provide amazing results. A thermal light-emitting source is often presented as a typical example of an incoherent source. However, it has been shown that the field emitted by a thermal source made of a polar material is enhanced by more than four orders of magnitude and displays first-order coherence in the near-field zone GreffetPRL (). Moreover, by introducing a periodic microstructure into such a polar material, a thermal infrared source can be fabricated that displays first-order coherence over large distances GreffetNature (). While first-order coherence and spectral properties of thermal sources can be manipulated and tailored, their second-order coherence is known to be completely absent NobelGlauber () resulting in the small bunching described by .

Here we investigate photon-coincidence counting statistics of thermal cavities in the ultrastrong coupling regime, where the strength of the interaction between an emitter and the cavity photons becomes comparable to the transition frequency of the emitter or the frequency of the cavity mode . Ultrastrong light-matter interactions have recently been achieved both in semiconductor and superconducting systems where exceptionally high field amplitudes for the photons and emitters with very large dipole moments can be realized guenter09 (); Niemczyk (); Todorov (); Schwartz (); Hoffman (); Scalari (). These structures are attracting increasing interest due to the possibility of manipulating the physical properties of the cavity quantum electrodynamic ground state.

The photon statistics of chaotic sources (like thermal cavities) and lasers can usually be explained classically. In contrast, strongly nonlinear photonic systems can emit individual photons well separated in time from each other when excited coherently or operating very far from thermal equilibrium. For the systems considered so far, such a scenario is however known to not persist when the coupled system is driven by thermal noise induced by reservoirs at finite temperature. Indeed, the standard quantum optics master equation (ME), generally used to study the dynamics of cavity QED systems Carmichael (), predicts for such systems independently of temperature and coupling strength. The interaction between atoms and cavity photons is most often neglected when considering the coupling of this system to an environment. Recently it has however been shown that this simplification, which leads to the standard quantum optics ME, can generate unphysical effects in the ultrastrong coupling regime Blais (). Another key issue is the failure of standard quantum optical normal order correlation functions to describe photodetection experiments for such systems DeLiberato2 (); RidolfoPrl2012 (). The theoretical treatment of this regime thus requires a description that goes beyond the standard techniques of quantum optics.

Exploiting generalized correlation functions as introduced in RidolfoPrl2012 () and a ME that fully takes into account the qubit-resonator coupling Davies1974 (); Blais (), we investigate the photon-coincidence counting statistics of thermal sources for arbitrary light-matter coupling. For this purpose we consider a single mode resonator coupled to a two-level quantum emitter where each subsystem is coupled to independent thermal baths of harmonic oscillators at a common temperature . We concentrate on the zero-detuning and low-temperature cases, where more striking deviations from the standard results appear. We moreover have corroborated the generality of our findings by confirming that similar results are found for one cavity mode coupled to multiple emitters or one emitter coupled to multiple cavity modes, see supplementary material supplement ().

A particularly well suited technology for such an experiment are superconducting circuits wallraff04 (); Norireview () which have recently emerged as an excellent platform for microwave on-chip quantum-optics experiments and where second-order correlation function measurements for quantum menzel10 (); Bozyigit (); Lang () and low temperature thermal fields Mariantoni10 () have been performed using quadrature amplitude detectors.

### Model -

The cavity QED system we explore consists of a single-mode resonator that interacts with a two level system (TLS). This system can be described by the Rabi Hamiltonian (assuming ),

(1) |

where and are the cavity and TLS bare frequencies, is the coupling strength, while denotes the annihilation (creation) operator for the cavity and is the lowering (raising) operator for the TLS. Since we are interested in probing photon-coincidence statistics for arbitrary light-matter interactions where the contribution of counter rotating terms cannot be neglected, we do not make use of the rotating wave approximation (RWA) asquarenote (). Recently, it was shown that in the limit of very large coupling, the ground state and the first excited state of this system can become quasidegenerate Nataf (); noripra2010 (); Emary ().

In order to study thermal emission as well as the statistics of thermal photons we calculate the normal-order correlation functions of the output field. Standard normal order correlation functions were recently shown to not correctly describe the emission properties and photon statistics of systems in the ultrastrong coupling regime RidolfoPrl2012 () as they would, for example, predict an unphysical stream of output photons even for a zero-temperature system, . Following RidolfoPrl2012 (), we here employ correlation functions for the output fields that are valid for an arbitrary coupling strength by expressing the cavity field ( is the rms zero-point field-amplitude) in the atom-cavity dressed basis. In particular, the set-up we have in mind is equivalent to the emission of a thermalized black-box, whose output is coupled to the vacuum of a one-dimensional waveguide. In this case, the output and input operators, and , obey the relation, RidolfoPrl2012 (). Here () is the positive (negative) frequency component of the cavity field , which can be derived by expanding the respective operators in the dressed state basis, namely the eigenstates of as in Eq. (1) ordered according to increasing eigenenergies . Specifically, the time derivative of can be expressed as where and . According to these input-output relations and for input fields in vacuum, the normalized second order correlation function for the output field reads,

(2) |

### Results -

The thermal-equilibrium zero-delay correlation function can be directly calculated from the thermal equilibrium density-operator note (). For a system in thermal equilibrium, statistical properties are related to the density matrix of the canonical ensemble, , that is the most general way to describe such thermalized interacting systems. In the basis where is diagonal, reads,

(3) |

where is the j-th eigenvalue of , the Kronecker delta, and the partition function.

Figure 1 shows the thermal-equilibrium zero-delay correlation function as a function of the effective coupling and temperature for zero detuning ().

The calculated exhibits striking differences from the standard value of . Of particular interest is the region with large effective coupling and low temperature where . This highly nonclassical behavior of thermal photons opens prospects towards the realization of thermal sources of single photons in circuit QED. This anomalous behavior originates from the tendency of the interacting quantum system towards vacuum degeneracy for large couplings. Specifically for increasing coupling the energy of the first excited state converges towards that of the ground state, while the other energy levels remain well separated from that doublet (see Fig. 2b). Hence, at sufficiently low temperature, only the first excited state is significantly populated by thermal noise ( is non-negligible only for ). At the onset of vacuum degeneracy, the ground and first excited state are quantum superpositions with multi-photon components Nataf () and naively one might expect to observe bunching effects. However such photons are mostly virtual and the use of generalized normal order correlation functions shows that the transition can only emit one physical photon at a time.

For large effective couplings and higher temperatures becomes larger than the standard value, . Also this behaviour can be understood from the spectrum of the Hamiltonian in Eq. (1). As conserves the parity of the number of TLS- and cavity-excitations, its eigenstates either contain an odd or an even number of such excitations but not both. The single excitation decays associated with can thus only occur between eigenstate with different parities. As the energy levels that converge to for , c.f. Fig. 2b, contain excitations, the decay process can not occur for . At however a crossing between the second and third excited level occurs and beyond this point the cascaded decay is possible expl (). For temperatures that are sufficiently high to appreciably populate these levels, such cascaded decays can lead to pronounced bunching with , c.f. Fig. 1. Since the energies of excited states decrease as grows, their population increases for a given temperature.

An overview of the behaviors of is shown in the contour plot of Fig. 2, where four different regions are displayed: i) a region for small (green) where the standard thermal result is recovered (here ) and which broadens for increasing temperatures; ii) a sub-Poissonian or nonclassical region with (blue) in the ultrastrong coupling regime and for sufficiently low temperatures; iii) a region with an intermediate regime (gray) with ; and iv) a super-Poissonian region (red) beyond the standard value of for very large coupling and higher temperatures. We note here that a calculation of within the RWA and using the standard quantum optics ME Carmichael () would always yield the value , regardless of the coupling strength and the temperature.

While zero-delay correlation function can be directly inferred from the thermal density operator, a description of the time dependent dynamics of the open quantum system is required to calculate the time-delayed second order correlation function . A viable description of system-bath interactions typically requires an expansion in the system bath coupling. A suitable way to perform this perturbative expansion consists in writing the Hamiltonian in the basis of its eigenstates Davies1974 (). In this way we obtain the following ME breuer (); Blais (),

(4) |

where and are Liouvillian superoperators describing the losses and the thermal feeding of the system. They read, for with . The relaxation coefficients depend on the spectral density of the baths, , and the system-bath coupling strength, , at the respective transition frequency as well as on the transition coefficients (). is the thermal population at frequency and temperature . In all examples reported in this work, we consider a cavity that couples to the momentum quadratures of fields in one-dimensional output waveguides, assuming that the spectral density is constant and . Hence the relaxation coefficients reduce to , where are the standard damping rates. These assumptions correspond to typical experimental settings e.g. in circuit-QED. In the ME (4) we neglect contributions of dephasing noise and Lamb shifts as they do not affect significantly our findings. A first result worth mentioning is that the steady state solution of Eq. (4) reproduces the thermal state of Eq. (3) that is independent of the damping rates. This confirms that the present approach is able to correctly describe the system’s thermalization and thus corroborates the accuracy of Eq. (4). We solve Eq. (4) using as initial condition to calculate via the quantum regression theorem GardinerZoller ().

Figure 3 displays , for the temperatures and couplings corresponding to the reported markers in Fig. 2 and . Notably, oscillations appear in the gray and green regions. These arise from interferences of the possible decay channels. In fact, for these values of the coupling strength, the separation between the energy levels is not high enough to suppress thermal occupation of higher excited states. Hence decays from upper energy manifolds into lower excited states are possible. In this way, e.g. the third excited energy level can decay into the second or into the first, and interferences of these possible decays result in the observed oscillating behavior. This explanation is corroborated by the fact that the frequency of the oscillations corresponds exactly to , i.e. the energy difference between the first and the second excited state of . Yet, in the blue region, is almost zero and behaves as for a usual TLS, moreover the oscillations in are suppressed due to the low thermal noise and great separation between the energy eigenvalues. Instead, in the red region, shows a slow decay for large which is due to the fact that the lifetime of an excitation in the state is a lot longer than of an excitation in the state . This difference in lifetimes is obvious from widths of the corresponding lines in the spectrum, see Fig. 4. Hence after the emission of a photon from the decay , the probability of detecting a second photon that originates from the decay is substantial for a long delay range . For this reason, in the red region remains bunched for such a long time. This behavior is a consequence of the spectral density of the bath that, being a linear function of the frequencies, tends to narrow the spectral linewidth of the lower resonances.

Deeper insights into the cascade processes involving transitions with different frequencies can be obtained from the frequency filtered second-order cross-correlation functions Moreau (); delValle2012 (). Here we calculate them for the circle-marker case considering the frequencies and ,

where is the time-ordering operator that positions the operators at earlier times on the outer places in the expectation value and with . For , shows pronounced bunching due to the cascaded decay process . For on the other hand, is anti-bunched as after detecting a photon with frequency , the system needs to be re-excited before a photon at frequency can be emitted and detected step ().

Thermal emission is also characterized by its power spectrum, i.e. the Fourier transform of the two time correlation .

We exploit again the quantum regression theorem to calculate the relevant two-time correlation function. Fig. 4 shows for the same parameters as used in Fig. 3 and the corresponding values of temperature and couplings for the different markers as introduced in the caption of Fig. 2. As expected the heights of the spectra increase for increasing temperature. Moreover the resonances have i) different linewidths, as one can see from the definition of the damping rates and ii) different heights. The latter is mainly a consequence of the thermal feeding. For a fixed temperature, the thermal occupation of a spectral resonance (i.e. its height) is determined by that increases as the frequency of the resonance decreases. The presence of two different decay times in the example with the dot marker becomes apparent via the different linewidths of the resonances contributing to the signal, see Fig. 4. Finally, the effects described here could also be measured by populating the system with a quasithermal field distribution realized by mixing a fixed frequency microwave tone with noise sources of different bandwidths Lang ().

### Acknowledgements -

SS acknowledges useful discussions with Omar Di Stefano and Roberto Stassi. This work is part of the Emmy Noether project HA 5593/1-1 and the CRC 631, both funded by the German Research Foundation, DFG.

## References

- (1) M. Planck, in Nobel Lectures, Physics 1901-1921, Elsevier Publishing Company, (1967).
- (2) R. Hanbury Brown and R. Q. Twiss, Nature 177, 27 (1956); Proc. Roy. Soc. (London) A242, 300 (1957); A243, 291 (1957).
- (3) R. J. Glauber, Phys. Rev. 130, 2529 (1963).
- (4) R. J. Glauber, Rev. Mod. Phys 78, 1267 (2006).
- (5) A. V. Shchegrov, K. Joulain, R. Carminati, and J.-J. Greffet Phys. Rev. Lett. 85, 1548 (2000).
- (6) J.-J. Greffet, R. Carminati, K. Joulain, J-P Mulet, S. Mainguy, and Y. Chen, Nature 416, 61 (2002).
- (7) G. Günter, A.A. Anappara, J. Hees, A. Sell, G. Biasiol, L. Sorba, S. De Liberato, C. Ciuti, A. Tredicucci, A. Leitenstorfer, and R. Huber, Nature 458, 178 (2009).
- (8) T. Niemczyk,F. Deppe, H. Huebl, E.P. Menzel, F. Hocke, M.J. Schwarz, J.J. Garcia-Ripoll, D. Zueco, T. Hümmer, E. Solano, A. Marx, and R. Gross, Nat. Phys. 6, 772 (2010).
- (9) Y. Todorov, A.M. Andrews, R. Colombelli, S. De Liberato, C. Ciuti, P. Klang, G. Strasser, and C. Sirtori, Phys. Rev. Lett. 105, 196402 (2010).
- (10) T. Schwartz, J. A. Hutchison, C. Genet, and T. W. Ebbesen, Phys. Rev. Lett. 106, 196405 (2011).
- (11) A.J. Hoffman, S.J. Srinivasan, S. Schmidt, L. Spietz, J. Aumentado, H.E. Türeci, A.A. Houck, Phys. Rev. Lett. 107, 053602 (2011).
- (12) G. Scalari, C. Maissen, D. Turcinkova, D. Hagenmuller, S. De Liberato, C. Ciuti, C. Reichl, D. Schuh, W. Wegscheider, M. Beck, and J. Faist, Science 16, 1323 (2012).
- (13) H. J. Carmichael, Statistical Methods in Quantum Optics 1, 2nd ed., Springer, (2000); Statistical Methods in Quantum Optics 2, Springer, (2008).
- (14) F. Beaudoin, J. M. Gambetta, and A. Blais, Phys. Rev. A 84, 043832 (2011).
- (15) S. De Liberato, D. Gerace, I. Carusotto, and C. Ciuti, Phys. Rev. A 80 053810 (2009).
- (16) A.Ridolfo, M. Leib, S. Savasta and M. J. Hartmann, Phys. Rev. Lett. 109, 193602 (2012).
- (17) E. B. Davies, Commun. Math. Phys. 39, 91 (1974)
- (18) see supplementary material
- (19) A. Wallraff, D.I. Schuster, A. Blais, L. Frunzio, R.-S. Huang, J. Majer, S. Kumar, S.M. Girvin and R.J. Schoelkopf, Nature 431, 162 (2004).
- (20) J. Q. You and F. Nori, Nature 474, 589 (2011).
- (21) E.P. Menzel, F. Deppe, M. Mariantoni, M. Á. Araque Caballero, A. Baust, T. Niemczyk, E. Hoffmann, A. Marx, E. Solano, and R. Gross, Phys. Rev. Lett. 105, 100401 (2010).
- (22) D. Bozyigit, C. Lang, L. Steffen, J. M. Fink, C. Eichler, M. Baur, R. Bianchetti, P. J. Leek, S. Filipp, M. P. da Silva, A. Blais, and A. Wallraff, Nature Phys. 7, 154 (2011).
- (23) C. Lang, D. Bozyigit, C. Eichler, L. Steffen, J.M. Fink, A.A. Abdumalikov, Jr., M. Baur, S. Filipp, M.P. da Silva, A. Blais, and A. Wallraff, Phys. Rev. Lett. 106, 243601 (2011).
- (24) M. Mariantoni, E.P. Menzel, F. Deppe, M.Á. Araque Caballero, A. Baust, T. Niemczyk, E. Hoffmann, E. Solano, A. Marx, and R. Gross, Phys. Rev. Lett. 105, 133601 (2010).
- (25) Notice that Eq. (1) does not contain the so called term, that is the term associated to the squared electromagnetic vector potential term. This term is often neglected, however its absence in Eq. (1) is not an approximation when considering a system that contains an artificial atom (fluxonium Manu ()) made of a Josephson junction coupled inductively to a resonator.
- (26) S. Ashhab and F. Nori, Phys. Rev. A 81, 042311 (2010).
- (27) P. Nataf, and C. Ciuti, Phys. Rev. Lett. 104, 023601 (2010).
- (28) C. Emary and T. Brandes, Phys. Rev. A 69, 053804 (2004).
- (29) The expectation value of a generic one-time operator can be calculated as .
- (30) Note that, motivated by our photodetection theory, we always number the dressed states according to increasing energies.
- (31) H.-P. Breuer and F.Petruccione, The Theory of Open Quantum Systems, Oxford University Press (2006).
- (32) C. W. Gardiner and P. Zoller, Quantum Noise, Springer-Verlag, (2000).
- (33) E. Moreau, I. Robert, L. Manin, V. Thierry-Mieg, J. M. Gérard, and I. Abram, Phys. Rev. Lett. 87, 183601 (2001).
- (34) E. del Valle, A. Gonzalez-Tudela, F. P. Laussy, C. Tejedor, and M. J. Hartmann, Phys. Rev. Lett. 109, 183601 (2012).
- (35) Our calculations correspond to the regime where the delay is large compared to the detector response time. Longer detector response times would smoothen the jump of at , c.f. GonzalezTudela ().
- (36) V. E. Manucharyan, J. Koch, L. I. Glazman, and M. H. Devoret, Science 326, 113 (2009).
- (37) A. Gonzalez-Tudela, F. P. Laussy, C. Tejedor, M. J. Hartmann, and E. del Valle, arXiv:1211.5592, New J. Phys., to appear (2013).

## Appendix A Supplementary Material

Here we report calculations of for many TLSs coupled to a single cavity mode, and for a single TLS coupled to two cavity modes. The latter case is of particular interest because of the unavoidable coupling of the TLS to higher modes that is often present in a real experimental setup. In each graph we plot as function of the temperature and of coupling strength, imposing and . We show moreover in all the graphs the plane to highlight the differences to the standard ME calculation which applies a RWA. The generalization of our model to many TLSs is straightforward, and described by the Hamiltonian .

Figs. 5, 6 and 7 show the feasibility to reach a robust quantum regime, characterized by the strong antibunching, that persists even in presence of many TLSs. It is worth to note the onset of a significant superbunching for higher coupling strenghts. This is due to the degeneracy of the involved energy levels, that increases with the number of TLSs.

Fig. 8 in turn shows the thermal calculated for a single TLS coupled with two cavity modes as described by the Hamiltonian . In particular, for this latter calculation we used a second bare cavity mode with frequency and a coupling with strength . Such ratios of coupling strength appear for example in circuit QED experiments Niemczyk (). There is an enhancement of the antibunching even though we consider a second bosonic mode. The explanation of such a scenario is that the presence of the second cavity mode is compensated by the stronger coupling that effectively increase the nonlinarity of the whole system.

### Acknowledgements -

This work is part of the Emmy Noether project HA 5593/1-1 and the CRC 631, both funded by the German Research Foundation, DFG.

## References

- (38) T. Niemczyk,F. Deppe, H. Huebl, E.P. Menzel, F. Hocke, M.J. Schwarz, J.J. Garcia-Ripoll, D. Zueco, T. Hümmer, E. Solano, A. Marx, and R. Gross, Nat. Phys. 6, 772 (2010).