Photon waiting time distributions: a keyhole into dissipative quantum chaos

TL;DR

This paper investigates how photon waiting time distributions can reveal chaotic versus regular behavior in a driven Kerr-nonlinear cavity, providing a non-invasive method to identify quantum chaos through photon emission statistics.

Contribution

It introduces a novel approach to distinguish chaotic and regular regimes in open quantum systems using photon waiting time distributions and quantum Lyapunov exponents.

Findings

Chaotic regimes show power-law asymptotics in photon waiting times.

Photon emission monitoring can differentiate chaotic and regular states.

Method enables non-invasive chaos detection in quantum optical systems.

Abstract

Open quantum systems can exhibit complex states, which classification and quantification is still not well resolved. The Kerr-nonlinear cavity, periodically modulated in time by coherent pumping of the intra-cavity photonic mode, is one of the examples. Unraveling the corresponding Markovian master equation into an ensemble of quantum trajectories and employing the recently proposed calculation of quantum Lyapunov exponents [I.I. Yusipov {\it et al.}, Chaos {\bf 29}, 063130 (2019)], we identify `chaotic' and `regular' regimes there. In particular, we show that chaotic regimes manifest an intermediate power-law asymptotics in the distribution of photon waiting times. This distribution can be retrieved by monitoring photon emission with a single-photon detector, so that chaotic and regular states can be discriminated without disturbing the intra-cavity dynamics.

Equations

$$H(t)=\frac{1}{2}\chi {\hat{a}}^{†2}{\hat{a}}^2+iF(t)({\hat{a}}^{†}-\hat{a}).$$

$$\dot{\varrho }=\mathcal{L}(\varrho )=-i[H,\varrho ]+\mathcal{D}(\varrho ),$$

$$\mathcal{D}(\varrho )$$

$$i\dot{\psi }=H\psi -\frac{i}{2}{V}^{†}V\psi .$$

$$\frac{d}{dt}||\psi ||=-{\psi }^{\ast }{V}^{†}V\psi ,$$

$$\xi (t)$$

$$\dot{\xi }=-\frac{1}{2}\gamma \xi +F(t)-i\chi |\xi {|}^2\xi .$$

$$\frac{d}{dt}||\psi ||=-\gamma n(t)||\psi ||,$$

Full text

Quantification of regimes emerging in open quantum systems driven out of equilibrium is a problem, interesting in several respects. In particular, it could help bridging non-equilibrium quantum phenomena with manifestations of classical dissipative chaos, such as local instability, bifurcations, strange attractors, etc Ott . Considerable effort was devoted to quantum generalization of Lyapunov exponents (LEs), the most popular and powerful means to quantify classical chaos Toda1987 ; Karol1992 ; Manko2000 . Recently, a new wave of research in this direction has been initiated by the idea of out-of-time correlation functions as possible quantum LE analogs Galitski2017 ; Galitski2018 ; Hirsch2019 . In our recent work Yusipov2019 , we proposed yet another approach to ‘quantization’ of LEs, based on unraveling of the Markovian master equation that describes the evolution of the system’s density matrix into an ensemble of trajectories. However, all the approaches suffer from a common problem: While it is possible to define the corresponding observables and correlators and then treat them analytically or calculate numerically, it is much harder (or not possible at all) to access them in a real-life experiment. Here we demonstrate that ‘chaos-regularity’ transitions can be captured by looking at the statistics of dissipative jump events that are accessible in quantum optics and called there “distribution of photon waiting times” Carmichael1991 ; Flindt2019 .

**Exemplifying in a simple model of an open periodically modulated quantum system, we show that transitions from regular to chaotic regimes (beforehand defined in terms of LEs) are marked by transitions from exponential waiting-time distributions to those with intermediate power-law scaling. Since photon emission events can be resolved in an experiment, e.g. by using single-photon detectors, discrimination between chaotic and regular regimes can therefore be achieved without disturbance of the system dynamics. **

I Introduction

A recent progress in such fields of experimental quantum physics as cavity quantum electrodynamics Walter2006 , cavity optomechanics Aspelmeyer2014 , artificial atoms You2011 and polaritonic devices Feurer2003 , gave the theory of open quantum systems book a new impetus. Different from the typical single-particle models of quantum optics Carmichael1991 , such systems are essentially multi-component and thus require many-body models to describe them. These models, when considered out of equilibrium, display asymptotic states that are sometimes interpreted as quantum versions of chaotic (or regular) attractors Ott . Until recently, such classification was performed mostly visually, after projecting quantum states onto a suitable classical phase space Spiller1994 ; Brun1996 ; Hartmann2017 ; Ivanchenko2017 ; Poletti2017 ; Yusipov2019 . The quantitative approach to dissipative Quantum Chaos Haake2010 is possible by means of quantum versions of Lyapunov exponents (LEs) Pikovsky . There is no unique definition of quantum LEs; instead, there are several alternative generalizations ranging from early attempts to define the exponents in terms of quasi-classical distributions Toda1987 ; Karol1992 ; Manko2000 to stochastic Schrödinger equation Persival1992 ; Jacobs2000 ; Pattanayak2008 ; Jacobs2009 ; Pattanayak2009 ; Pattanayak2018 and further on, to a very recent approach based on the concept of out-of-time-correlators Galitski2017 ; Bohrdt2017 ; Galitski2018 ; Hirsch2019 . However, all these generalizations suffer a common problem: the corresponding quantifiers are hard (or not possible at all) to measure in a real-life experiment.

In our recent work Yusipov2019 , we introduced a version of quantum LEs based on the so-called ‘quantum trajectory’ unraveling plenio ; daley (also called ‘Monte Carlo wave-function (MCwf) method’ zoller ; dali ), which replaces the original quantum Markovian master equation book with an ensemble of quantum trajectories. The evolution of the trajectories is governed by a non-Hermitian Hamiltonian in a continuous-time manner, interrupted by random jump-like dissipative events. While this version, in principle, opens a way to theoretical characterization of quantum dissipative chaos in physically relevant setups (e.g., in cavity quantum electrodynamic systems Walter2006 ), the experimental estimation of such quantum LEs is no less challenging than for the prior versions.

Here we propose an experimentally feasible approach to detecting chaotic regimes in cavity-like open quantum systems. Namely, using a periodically modulated Kerr-nonlinear cavity as a model, we demonstrate that the transition to quantum chaos, quantified with the LEs from Ref. Yusipov2019 , is associated with the appearance of an intermediate power-law asymptotics in the distribution of waiting time (that is the time between two consequent photon emissions from the cavity Vyas ; Carmichael1989 ; baran ). This distribution can be sampled in an experiment by using single-photon detection techniques photon1 ; photon2 .

II Model

We consider a photonic mode in a leaky Kerr-nonlinear cavity, periodically modulated by an external coherent EM field Spiller1994 ; Brun1996 . Its unitary dynamics is governed by the Hamiltonian (we set $\hbar=1$)

[TABLE]

Here, $\chi$ is the photon interaction strength, $\hat{a}^{\dagger}$ and $\hat{a}$ are photon creation and annihilation operators, so that $\hat{n}=\hat{a}^{\dagger}\hat{a}$ is the photon number operator. $F(t)=F(t+T)$ describes periodic modulation. We use the two-valued quench-like driving with period $T$, that is $F(t)=A$ within $0<t\leq T/2$ and $F(t)=0$ for the second half period $T/2<t\leq T$.

Photons can be emitted from the cavity. In principle, they can also be pumped in by thermal environment, but here, similar to the setups in Refs. Spiller1994 ; Brun1996 , we work in the zero-temperature limit, assuming that the pumping rate is zero. Evolution of the cavity is described by the Lindblad master equation book ; alicki ,

[TABLE]

where the first term in the r.h.s. captures the unitary dynamics of the system, determined by Hamiltonian (1), while the second term reflects coupling to the environment. The dissipative term involves a single jump operator,

[TABLE]

which describes emission of photons by the cavity into the zero-temperature environment. The dissipative coupling constant $\gamma$ is assumed to be time-independent.

In simulations, we limit the number of photons in the cavity mode to $N$, so that the Hilbert space of the numerical model has dimension $N+1$ and can be spanned with the $N+1$ Fock basis vectors, $\{|n+1\rangle\}$, $n=0,...,N$. $N$ is chosen to be large enough so that the average number of photons in the cavity, $\langle N_{\mathrm{ph}}\rangle$, remains substantially smaller. $\langle N_{\mathrm{ph}}\rangle$ depends on parameters of Hamiltonian; yet the main control parameter, which allows to control the mean number of photons is the coupling $\chi$ Spiller1994 ; Brun1996 . Throughout the paper we set $\chi=0.008$ and $\gamma=0.05$. It gives $\langle N_{\mathrm{ph}}\rangle\sim 50$ so we set $N=300$ and verify that it suffices.

We make use of the quantum Monte-Carlo wave function method to unravel the deterministic equation (2) into an ensemble of quantum trajectories zoller ; dali ; plenio ; daley . It allows us to describe evolution of the model system in terms of an ensemble of pure states, $\psi(t)$, governed by the effective non-Hermitian Hamiltonian Spiller1994 ; Brun1996 ,

[TABLE]

The norm of the wave function decays according to

[TABLE]

until it reaches a threshold $\eta$, repeatedly chosen as i.i.d. random number from $[0,1]$. Then a random jump is performed, and the wave function norm is reset to $||\psi(t)||=1$. After that the continuous non-unitary evolution, Eq. (4), is resumed until the next quantum jump, etc. For the model given by Eqs.(1,3), a quantum jump corresponds to emission of a single photon, which can be recorded by a photodetector Carmichael1991 .

The density matrix sampled from a set of $M_{r}$ realizations, that originate from an initial pure state $\psi^{\mathrm{init}}$ for Eq. (4), as $\varrho(t_{\mathrm{p}};M_{\mathrm{r}})=\frac{1}{M_{r}}\sum_{j=1}^{M_{\mathrm{r}}}\ket{\psi_{j}(t_{\mathrm{p}})}\bra{\psi_{j}(t_{\mathrm{p}})}$, converges towards the solution of Eq. (2) at time $t_{\mathrm{p}}$ for the initial density matrix $\varrho^{\mathrm{init}}=\ket{\psi^{\mathrm{init}}}\bra{\psi^{\mathrm{init}}}$.

Following Refs. Spiller1994 ; Brun1996 , we use the complex valued observable of the non-Hermitian photon annihilation operator as a dynamical ‘variable’

[TABLE]

In the mean-field classical limit, $N_{\mathrm{ph}}\rightarrow\infty$, its evolution is given by

[TABLE]

The state of the mean-field system is described by two phase variables, $\{\mathrm{Re}~{}\xi,\mathrm{Im}~{}\xi\}$. This system is essentially non-linear and periodically modulated time and, expectedly, exhibits a spectrum of different asymptotic regimes, from varios periodic orbits to chaotic attractors. The parameter dependence can be visualized with a bifurcation diagram, Fig. 1(a), constructed from stroboscopic values, $\xi_{k}=\xi(t_{0}+kT)$, where $t_{0}$ is a transient time given to the system to land on the corresponding attractor. In particular, one observes that the largest Lyapunov exponent becomes positive upon transition from a fixed point to a chaotic attractor in the stroboscopic map (see Fig. 1b).

III Quantum Lyapunov exponent

To calculate the largest quantum LE, we employ the recently proposed method Yusipov2019 . It is based on the parallel evolution of fiducial and auxiliary trajectories, $\psi_{f}(t)$ and $\psi_{a}(t)$, under Eq. (4), in the spirit of the classical LE calculation Benettin1976 . The auxiliary trajectory is initialized as a normalized perturbed vector $\psi_{a}^{init}=\psi_{f}^{init}+\varepsilon\psi_{r}$, produced with random i.i.d. entries in $\psi_{r}$ and $\varepsilon\ll 1$. The fiducial and perturbed observables, $\xi_{f}(t)$ and $\xi_{a}(t)$, typically remain close to each other over many quantum jump events; as the difference in the observables gets over an upper threshold, $\Delta(t_{k})=|\xi_{f}(t_{k})-\xi_{a}(t_{k})|>\Delta_{max}$, or below a bottom threshold $\Delta(t_{k})<\Delta_{min}$, the perturbed state is set back closer to the fiducial one along the mismatch direction $\psi_{f}(t_{k})-\psi_{a}(t_{k})$, so that $|\xi_{f}(t_{k})-\xi_{a}(t_{k})|=\Delta_{0}$; the wave vector gets normalized and the occurred growth rate recorded, $d_{k}=\Delta(t_{k})/\Delta_{0}$ periodic . The largest LE is estimated following the divergence of the chosen observable as $\lambda(t)=\frac{1}{t}\sum_{k}\ln d_{k}$ Benettin1976 .

We use a recent high-performance realization of the quantum jumps method Volokitin2017 to generate $M_{\mathrm{r}}=10^{2}$ different trajectories. First, we allow each trajectory to evolve up to $t_{0}=2\cdot 10^{3}T$ in order to arrive to the asymptotic regime, and then we follow the dynamics of fiducial and auxiliary trajectories up to $t=10^{3}T$. We determine the dimension of the computational model, $N$, that would be large enough to ensure independence of the results on the chosen value. It has turned out that $N=150$ is already sufficient for the consistent calculation of the quantum LEs almost in the whole range of modulation amplitude, $A$ (cf. Fig.1(b), bottom panel). However, for $A>3.5$ oscillations of the number of photons in the cavity grow, and that leads to the distortion of the dynamics and positive LEs. Ultimately, we find it out that choosing the maximal number of photons $N>250$ provides consistency in the whole range of studied parameters.

Additionally, we verified convergence of the maximal LE to its asymptotic values in the regular and chaotic regimes and confirm that the results are weakly dependent on the choice of observables (e.g., we used photon population number, $n(t)=\langle\psi^{\dagger}(t)|\hat{n}|\psi(t)\rangle$ as an alternative); see Fig. 2.

Depending on parameter values, interaction with the environment can strongly localize quantum trajectories on the classical ones Spiller1994 ; Brun1996 ; Jacobs2000 ; Pattanayak2018 . Our case is notably different, so that the resulting structure of the probability distribution for $\xi$ has only a qualitative resemblance, see Fig. 1(b), top. Nevertheless, working in the essentially quantum regime, we observe that the largest quantum LE becomes positive, see Fig. 1(b). Note, that the chaotic interval in the quantum case is somewhat greater, and fine structure intervals of classical regular dynamics are not reproduced. Such examples, when quantization generates chaotic dynamics, are known in the literature Pattanayak2018 .

To get a more general picture, we performed an extensive round of calculations and obtained the two-parameter dependence of classical and quantum largest LEs on the amplitude and period of modulations; see Fig. 3(a,b). Both the mean-field and quantum models produce visually similar structures of regular and chaotic regimes.

IV Waiting time distribution

As we pointed it out in Introduction, quantum LEs are good objects for theoretical and numerical analysis but their experimental estimation is highly nontrivial.

Can we infer information about intra-cavity dynamics from what is accessible in experiment? Statistics of fluorescent photon emission from cavity is one of the most popular characteristics in quantum optics and appears promising for our purpose Carmichael1989 ; Delteil2014 ; Saez2018 ; Brange2019 . The other reason for the choice is that there is only one dissipative channel in our model system, and a single quantum jump corresponds to a single photon emission, so that the intra-cavity dynamics and radiation are tightly bound.

Evolution of the wave vector norm between consecutive quantum jumps at times $\{t_{k}\}$ is governed by Eq. (5). For the system under study, Eqs.(1,3), it follows

[TABLE]

as $V^{\dagger}V=\gamma\hat{n}$. The norm decay from $||\psi(t_{k-1})||=1$ to the random jump threshold $||\psi(t_{k})||=\eta_{k}$ can be approximated by an exponent with an average rate proportional to an effective number of photons within $t_{k-1}<t<t_{k}$, that is $s_{k}=\gamma n^{(eff)}_{k}$. In other words, one can write $\tau_{k}=t_{k}-t_{k-1}=-\ln(\eta_{k})/s_{k}$. Since $\eta_{k}$ are random i.i.d. numbers on $[0,1]$, the variable $\zeta=-\ln(\eta_{k})$ has the probability density distribution $W_{\zeta}(\zeta)=\exp(-\zeta)$. If an asymptotic density matrix $\rho$ has a regular structure (e.g. unity matrix or unimodal distribution) then $n\approx const$ and, hence, $s\approx const$. In this case, the intervals between jumps also follow Poisson distribution, $W_{\tau}(\tau)=\exp(-\tau)$. Another example is a bimodal asymptotic solution, a result of the quantum analogue of period doubling bifurcation Ivanchenko2017 ; Poletti2017 , also a regular one, that would produce a superposition of two exponents.

When the system has a chaotic attractor, the dynamics of observables, in particular, $n(t)$, becomes irregular. An example for a quantum trajectory deep in the chaotic regime, $A=4.0,T=20$, is presented in Fig.4(a). It is seen that the norm decay rates between each jump differ substantially. Therefore, the probability distribution $W_{s}(s)$ cannot be deduced from simple argumetns any longer. At the same time, it can be estimated numerically, see Fig.4(b), upper panel. The results indicate that the distribution becomes much broader for $A>0.4$, with transition to chaos (Fig.4(b), lower panel). Correspondingly, the resulting distribution of interjump intervals, $W_{\tau}(\tau)$, will no longer be exponential. In the complete absence of a theoretical background in this case, we turn to numerical simulations.

Our main finding is presented in Fig. 5. The key observation is that the waiting time distribution becomes non-Poissonian and acquires a power-law intermediate scaling when the largest quantum LE becomes positive, Fig. 5, see also Fig.4(b), lower panel.

We performed massive sampling over the parameter plane $\{A,T\}$. At every parameter point, the power-law fitting of the sampled PDF was obtained by the least squares linear regression for the log-log scaled distributions. The quality of fit was characterized by the coefficient of determination Regression , $R^{2}\in[0,1]$, with larger values corresponding to better fit. We vary the interval of durations to search for the best fit, the power-law exponent $\alpha$, and request that it spans over at least one decade for the horizontal axis and has $R^{2}>0.98$ note2 . Otherwise, the power-law hypothesis is rejected (black zones in Fig. 5(b)).

It is noteworthy that the zones, where power-law asymptotics in the waiting time distribution is present, are well correlated with the zones, where the corresponding maximal quantum LE is positive (cf. Fig.3(b)). Moreover, there are concordant variations of the values of power-law exponents and of the largest LE.

V Conclusions

By using an open ac-driven Kerr cavity as a model, we found that the photon waiting time statistics can serve a good diagnostic tool to detect dissipative quantum chaos by the appearance of power-law intermediate asymptotics (which can be alternatively characterized by the positive largest quantum LE). The theoretical foundation of this observation remains a challenge. We conjecture that the first step would be to assume the statistical independence of the consecutive times $t_{j-1},t_{j},t_{j+1}$ (which is not guaranteed though), and then to implement the machinery used by Scherr and Montroll to derive power-law statistics of photocurrent in amorphous materials Scher1975 . Implementation of this program is the subject of future study.

Our results open a new perspective for quantification of regimes emerging in open quantum systems far from equilibrium, especially in such fields as quantum electrodynamics, quantum optics, and polaritonic devices, where photon emission statistics is an established and conventional tool Carmichael1989 ; Delteil2014 ; Saez2018 ; Brange2019 . Recently, a substantial progress in a non-demolition detection of individual emitted microwave photons has been achieved (see, e.g., nondem2019 ) so that the field of microwave quantum optics represents the most appealing opportunity to explore dissipative quantum chaos by waiting time statistics.

Potential links to self-organized criticality Markovic2014 ; Munoz2018 and Lévy flights (recently used to model power-law flip statistics of open spin systems Imamoglu2016 ) are another issues worth of consideration.

The authors acknowledge support of Basis Foundation grant No. 17-12-279-1 and President of Russian Federation grant No. MD-6653.2018.2. Numerical simulations were carried out at the Lobachevsky University and Moscow State University supercomputers.