Squeezed displaced entangled states in the quantum Rabi model

TL;DR

This paper explores analytical solutions of the quantum Rabi model, revealing superpositions of displaced and squeezed states, and introduces a model interpolating between different spectral regimes to analyze spectral statistics.

Contribution

It presents a novel single parameter model for the quantum Rabi model that bridges discrete and continuous spectra, analyzing spectral statistics near spectral collapse.

Findings

Solutions are superpositions of displaced number and squeezed coherent states.

Spectral statistics show interweaving fluctuations in first neighbor differences.

Second neighbor separation remains constant near spectral collapse.

Abstract

The quantum Rabi model accepts analytical solutions in the so-called degenerate qubit and relativistic regimes with discrete and continuous spectrum, in that order.We show that solutions in the laboratory frame are the superposition of even and odd displaced number states, in the former, and infinitely squeezed coherent states, in the latter, of the boson field correlated to the internal states of the qubit. We propose a single parameter model that interpolates between these discrete and continuous spectrum regimes to study the spectral statistics for first and second neighbor differences before the so-called spectral collapse. We find two central first neighbor differences that interweave and fluctuate keeping a constant second neighbor separation.

Equations

$$\hat{H}=\omega {\hat{a}}^{†}\hat{a}+\frac{{\omega }_0}{2}{\hat{\sigma }}_z+g({\hat{a}}^{†}+\hat{a}){\sigma }_x.$$

$$H(\delta )=2\omega (1-\delta ){\hat{a}}^{†}a+\delta {\omega }_0{\hat{\sigma }}_z+g({\hat{a}}^{†}+\hat{a}){\hat{\sigma }}_x,$$

$${\hat{H}}_{FG}(\delta )={\hat{H}}_{+}(\delta )|+⟩⟨+|+{\hat{H}}_{-}(\delta )|-⟩⟨-|,$$

$${\hat{H}}_{\pm }(\delta )=2\omega (1-\delta ){\hat{a}}^{†}a+g({\hat{a}}^{†}+\hat{a})\mp \delta {\omega }_0\hat{\Pi },$$

$$(z+\tilde{g})\frac{d}{dz}{\varphi }_1(z)-x{\varphi }_1(z)=\widetilde{{\omega }_0}{\varphi }_2(z)$$

$$(z-\tilde{g})\frac{d}{dz}{\varphi }_2(z)+(2{\tilde{g}}^2-x-2\tilde{g}z){\varphi }_2(z)=\widetilde{{\omega }_0}{\varphi }_1(z),$$

$$G(x;0)$$

$$j{K}_j(x)=K_{j-1}(x)f_{j-1}(x)-K_{j-2}(x)$$

$$f_j(x)$$

$${\hat{H}}_{I\pm }\equiv \underset{\delta \to 0}{\lim }{\hat{H}}_{\pm }(\delta )=2\omega {\hat{a}}^{†}a+g({\hat{a}}^{†}+\hat{a}),$$

$$\hat{D}(\alpha ){\hat{H}}_{FG}(0){\hat{D}}^{†}(\alpha )=2\omega {\hat{a}}^{†}a-\frac{g^2}{2\omega }.$$

$$|{\psi }_n⟩=\frac{1}{2}\left(|{\varphi }_{+}⟩|+⟩+|{\varphi }_{-}⟩|-⟩\right),$$

$$|{\varphi }_{\pm }⟩=\left[{\hat{D}}^{†}(\alpha )\mp (-1{)}^n\hat{D}(\alpha )\right]|n⟩,$$

$$\text{Q}(\beta )=\frac{e^{-(|\beta {|}^2+|\alpha {|}^2)}}{4\pi }\sum _{r=1,-1}{\left|\sum _{k=0}^{\infty }\frac{{\beta }^{\ast k}}{\sqrt{k!}}P(n,k,\alpha )\left[(-1{)}^{|n-k|}-r(-1{)}^n\right]\right|}^2$$

$$P(n,k,\alpha )=(-1{)}^{|n-k|\frac{(1+\text{sgn}(n-k))}{2}}\frac{{\alpha }^{|n-k|}}{|n-k|!}\sqrt{\frac{\max (n,k)!}{\min (n,k)!}}\text{U}(-\min (n,k);|n-k|+1;{\alpha }^2)$$

$$⟨\hat{x}⟩$$

$$⟨\hat{y}⟩$$

$${\hat{H}}_{R\pm }=\underset{\delta \to 1}{\lim }\hat{H}(\delta )=g({\hat{a}}^{†}+\hat{a})\mp {\omega }_0\hat{\Pi },$$

$${\hat{H}}_{R\pm }=2g\hat{x}\mp {\omega }_0\hat{\Pi }.$$

$$|{\psi }_0⟩$$

$$|{\psi }_x⟩$$

$$|{\phi }_{\pm }⟩=\underset{\mu \to 1}{\lim }({\gamma }_{+}\mp {\gamma }_{-})e^{\frac{1}{4}\xi }{e}^{(1-\mu )\frac{{\zeta }^2}{2}}\hat{S}(-r)\left[D(\zeta )\mp D^{†}(\zeta )\right]|0⟩,$$

$$|{\Psi }_0^{\pm }{⟩}_{FG}$$

$$|{\Psi }_x^{\pm }{⟩}_{FG}$$

$$\left\{{\gamma }_{+},{\gamma }_{-}\right\}$$

$$\text{Q}(\beta )=\frac{1}{{\pi }^{5/4}}\underset{\mu \to 1}{\lim }{e}^{\frac{1}{2}\xi -{\zeta }^2\mu }{e}^{-|\beta {|}^2-|\zeta {|}^2}[{\left({\gamma }_{+}-{\gamma }_{-}\right)}^2{\left|\sum _{n,m}G\left(-\mu ,2m+1,2n+1,\frac{3}{2}\right)\right|}^2$$

$$+{\left({\gamma }_{+}+{\gamma }_{-}\right)}^2{\left|\sum _{n,m}G\left(-\mu ,2n,2m,\frac{1}{2}\right)\right|}^2],$$

$$G(\mu ,m,n,\frac{1}{2}+s)$$

$$\times i^{m+n+s(m+1)}{\left(\frac{2}{\mu }\right)}^{s-\frac{m+n}{2}}\frac{\sqrt{{\left({\mu }^2-1\right)}^{s+\frac{1}{2}}}}{\Gamma \left(\frac{n+3}{2}-s\right)\Gamma \left(\frac{m+3}{2}-s\right)}$$

$$\text{leavevmode}\text{leavevmode}{}_2{F}_1\left(-\frac{m+s}{2},-\frac{n+s}{2};\frac{1}{2}+s;\frac{1-{\mu }^2}{{\mu }^2}\right),$$

Full text

I Introduction

The quantum Rabi model (QRM) is the lowest dimensionality Hamiltonian describing light-matter interaction,

[TABLE]

It describes a boson field, with frequency $\omega$ and represented by the annihilation (creation) operator $\hat{a}$ ($\hat{a}^{\dagger}$), interacting with a qubit, with frequency $\omega_{0}$ and represented by Pauli matrices $\sigma_{j}$ with $j=x,y,z$. Trapped ions Lv et al. (2018) and superconducting circuits Mezzacapo et al. (2014) provide highly controllable experimental platforms for the quantum simulation of the model in the different interaction regimes defined by the coupling strength to field frequency ratio Casanova et al. (2010); Yoshihara et al. (2016); A. F. Kockum, A. Miranowicz, S. De Liberato, S. Savasta and F. Nori (2019).

Through its history, the QRM has motivated the development of computational tools for both spectral and dynamic calculations Wolf et al. (2012); Maciejewski et al. (2014). The model conserves parity and is solvable Braak (2011); Moroz (2013). It is possible to diagonalize it in the qubit basis Moroz (2014, 2016). In the so-called adiabatic approximation, it is possible to estimate the spectrum and its eigenstates E. K. Irish (2007); Shen et al. (2016); Xie et al. (2019). Generalizations that account for asymmetry between so-called rotating and counter-rotating terms as well as driving showed the existence of degeneracies in the spectrum Li and Batchelor (2015). Extension for more than one qubit Chilingaryan and Rodríguez-Lara (2013); Wang et al. (2014); B. M. Rodríguez-Lara, S. A. Chilingaryan and H. M. Moya-Cessa (2014); Zhang and Chen (2015); S. A. Chilingaryan and B. M. Rodríguez-Lara (2015) or field Travěnec (2012); Huerta Alderete and Rodríguez-Lara (2016); Cong, Lei and Sun, Xi-Mei and Liu, Maoxin and Ying, Zu-Jian and Luo, Hong-Gang (2017, 2019) have been constructed. The latter was used for the simulation of para-particles in trapped-ion setups Huerta Alderete and Rodríguez-Lara (2017); Huerta Alderete et al. (2017); Huerta Alderete and Rodríguez-Lara (2018).

The QRM accepts analytic solutions in the so-called degenerate qubit, $\omega_{0}=0$, and relativistic, $\omega=0$, regimes Pedernales et al. (2015). After diagonalization in the qubit basis, the first is reduced to two decoupled harmonic oscillators with discrete spectrum and the second to a Dirac equation in (1+1)D with continuous spectrum. In the following, we introduce a single-parameter QRM that interpolates between these two regimes, diagonalize it in the qubit basis, and write its eigenvalue problem in the Bargmann representation. Our model is engineered to show spectral collapse Ng et al. (1999); M. Hwang, R. Puebla, M. Plenio (2015); Duan et al. (2016); Penna et al. (2017). It is well known that the eigenstates in the degenerate qubit regime are the superposition of even and odd displaced number states correlated to the ground and excited state of the qubit, in that order. For the sake of completeness, we revisit this result. Then, we calculate the expected but unreported eigenstates in the relativistic regime, which are the unbalanced superposition of position states that can be written as infinitely squeezed coherent states correlated to the internal states of the qubit. Afterwards, we numerically explore the statistics of the first and second neighbors spectral differences before the transition to the continuous spectrum regime. Their histograms show the interweaving of two central separations for nearest neighbors that keep a constant second neighbor spectral separation. These results correlate with Braak conjectures regarding the spectrum of the QRM in the Bargmann representation Braak (2011). Finally, we use Husimi Q-function to visually explore the displacement and squeezing of the ground state of our model as the control parameter takes us close to the continuous spectrum regime.

II Model

We propose a single-parameter QRM,

[TABLE]

that interpolates between the so-called degenerate qubit, $\omega_{0}=0$, and relativistic, $\omega=0$, regimes for the extremal values of the control parameter $\delta\in\left[0,1\right]$. We recover the QRM for $\delta=1/2$. A trapped ion quantum simulation Pedernales et al. (2015), $\hat{H}=(\delta_{b}-\delta_{r})\hat{a}^{\dagger}a/2-(\delta_{b}+\delta_{r})\frac{\hat{\sigma}_{z}}{4}+\Omega(\hat{a}^{\dagger}+\hat{a})\hat{\sigma}_{x},$ provides the desired level of control through blue (red) driving sidebands, $\delta_{b}=\omega_{b,r}-\omega_{0}-\nu$ ($\delta_{r}=\omega_{r}-\omega_{0}+\nu$) where the parameter $\omega_{b}$ ($\omega_{r}$) is the frequency of the driving laser detuned to the blue (red) of the ion transition frequency $\omega_{0}$, and the center of mass motion frequency is given by $\nu$. The effective coupling strength $\Omega$ is related to the Lamb-Dicke parameter of the trap and the amplitude of the driving field.

Strictly speaking, our model is solvable Braak (2011); Guan et al. (2018). We focus on the regimes with analytic closed form solution and favour the Fulton-Gouterman procedure to diagonalize it in the qubit basis $\{|+\rangle,|-\rangle\}$ Moroz (2013),

[TABLE]

where the dynamics in the boson sector,

[TABLE]

include the boson parity operator, $\hat{\Pi}=e^{i\pi\hat{a}^{\dagger}a}$. In the following sections, we discuss the closed form analytic solutions in the extremal regimes and provide a numerical study for the transition from discrete to continuous spectrum, Fig. 1.

In this frame, we use Bargmann representation Braak (2011), $a^{\dagger}\to z$ and $a\to d/dz$ with $z\in\mathbb{C}$, to write the eigenvalue problem $\hat{H}_{\pm}\psi_{\pm}=E_{\pm}\psi_{\pm}$ in differential form,

[TABLE]

where we obviate the subindices and use the auxiliary functions $\phi_{1,2}=e^{\tilde{g}z-g^{2}}\psi(z)$ in terms of the eigenfunction $\psi(z)$. We define scaled coupling and frequency parameters $\tilde{g}=g/\left[2(1-\delta)\omega\right]$ and $\tilde{\omega_{0}}=\delta\omega_{0}/\left[2(1-\delta)\omega\right]$. For the sake of space, we define an auxiliary scaled displaced energy $x=\tilde{g}^{2}+\tilde{E}$ in terms of the scaled eigenenergy $\tilde{E}=E/\left[2(1-\delta)\omega\right]$. In this representation, the parity operator acts as $\hat{\Pi}\psi(z)=\psi(-z)$. A series solution of the auxiliary functions, $\phi_{2}(x;z)=\sum_{j=0}^{\infty}K_{j}(x)(z+\tilde{g})^{j}$ which implies $\phi_{1}(x;z)=\sum_{j=0}^{\infty}K_{j}(x)\tilde{\omega_{0}}(z+\tilde{g})^{j}/{(j-x)}$, allows one to construct an analytic function $G(x;z)=\phi_{2}(-z)-\phi_{1}(z)$. As the auxiliary functions $\phi_{1,2}(\pm z)$ are series expansion of the same function outside the poles, the analytic function $G(x;z)$ must vanish for true solutions for any value of $z$. In particular, when $z=0$ we have

[TABLE]

where the $j$-th coefficient is given by a three-term recurrence relation,

[TABLE]

with initial terms $K_{0}(x)=1$ and $K_{1}(x)=f_{0}(x)$. The auxiliary function is given in the following,

[TABLE]

The solution to the eigenvalue problem reduces to calculate the roots of $G(x;0)$ as Braak demonstrated Braak (2011). This cannot be addressed in a direct way but Braak conjectured that the distance between poles is given by the scaled displaced energy $x$, and that zero, one or two roots can exist between consecutive poles without successive occurrences of two roots. In the following, we will argue that these conjectures are related to the first and second neighbor energy separation.

III Degenerate qubit regime

It is well known that in the degenerate qubit regime, $\delta\rightarrow 0$, the boson sector Hamiltonian reduces to a driven harmonic oscillator,

[TABLE]

that is diagonalized by a displacement,

[TABLE]

The eigenstates are displaced number states, $|n,\alpha\rangle\equiv\hat{D}(\alpha)|n\rangle$, with equally spaced spectrum $E_{n}^{\alpha}=2\omega n-\frac{g^{2}}{2\omega}$. We used the displacement operator $\hat{D}(\alpha)=e^{\alpha\hat{a}^{\dagger}-\alpha^{*}\hat{a}}$ with parameter $\alpha=-g/2\omega$. The eigenstates in the laboratory frame,

[TABLE]

have the familiar form of a maximally entangled boson-qubit state between the unnormalized even (odd) displaced number states,

[TABLE]

and the excited (ground) qubit state. This is not a Schrödinger cat state Schrödinger, E. (1935); Brune et al. (1992); Monroe et al. (1996). The latter needs semi-classical states in the boson sector. We have highly non-classical boson states that belong to different parity sectors. Quantum entangled states with parity properties are useful to detect weak forces Gilchrist et al. (2004).

For the sake of visualization, we calculate Husimi Q-function Husimi (1940) for the reduced boson sector states,

[TABLE]

where the probability of finding $k$ excitations,

[TABLE]

is given in terms of the Tricomi function $\mathrm{U}(n_{1};n_{2};x)$ Villanueva Vergara and Rodríguez-Lara (2015) and the function $\mathrm{max}(n,k)$ ($\mathrm{min}(n,k)$) yields the larger (smaller) value between $n$ and $k$. Figure 2 shows the Husimi Q-function for the ground state, $n=0$, fifth, $n=5$ and tenth, $n=10$, excited states. The mean expectation value for the boson oscillator quadratures in optical phase space, $\hat{x}=\left(\hat{a}^{\dagger}+\hat{a}\right)/2$ and $\hat{y}=i\left(\hat{a}^{\dagger}-\hat{a}\right)/2$, is zero and their squared standard deviations are not:

[TABLE]

One dispersion is always larger, $\Delta x>\Delta y$, as the displacement parameter $\alpha$ is never zero.

IV Relativistic regime

Here the QRM reduces to a $(1+1)$D Dirac equation Gerritsma et al. (2010); Noh et al. (2013); Rodríguez-Lara and Moya-Cessa (2014); Gutiérrez-Jáuregui and Carmichael (2018) and yields an effective boson sector Hamiltonian,

[TABLE]

that can be mapped into a simpler form in optical phase space,

[TABLE]

It is not bounded from below, has a continuous spectrum, and the solution is Dirac delta normalizable. Its eigenstates in the laboratory frame,

[TABLE]

include boson sector eigenstates,

[TABLE]

that are the superposition of squeezed coherent states Satyanarayana (1985); Nieto and Truax (1993); Bishop and Vourdas (1994) with symmetric displacement with respect to the origin, infinite squeezing, but different probability amplitude, $\xi=\ln{(1-\mu^{2})}$, $\zeta=2xe^{\frac{\xi}{2}{}}$ and $\mu=\tanh{(r)}$. Squeezed displaced states are interesting from a fundamental point of view Satyanarayana (1985); Bishop and Vourdas (1994); Nieto and Truax (1993). They are a resource for high precision measurements Huang et al. (2015); Kienzler et al. (2015); Lo et al. (2015) as they increase the signal-to-noise ratio Korobko et al. (2017).

In order to obtain this result, we start from the effective boson sector Hamiltonian and construct the eigenvalues for its square, $E^{2}(x)=4g^{2}x^{2}+\omega_{0}^{2}$. The square root of these eigenvalues provide a guide to construct the eigenstates using the invariant subspace spanned by the position states $\{|\pm x\rangle\}$,

[TABLE]

We use the map of Fock states to optical phase space, $\langle x|n\rangle=\pi^{-1/4}e^{-x^{2}/2}/(\sqrt{2^{n}n!})\mathrm{H}_{n}(x)$, in terms of Hermite polynomials Lebedev and Silverman (1972), and translate them into operator form, $e^{-\hat{a}^{\dagger 2}+2x\hat{a}^{\dagger}}|0\rangle=\sum_{n=0}^{\infty}\mathrm{H}_{n}(x)/\sqrt{n!}\leavevmode\nobreak\ |n\rangle$,

[TABLE]

Considering the canonical pair, $\hat{q}=\sqrt{2\hbar/(m\nu)}\leavevmode\nobreak\ \hat{x}$ and $\hat{p}=\sqrt{2\hbar m\nu}\leavevmode\nobreak\ \hat{y}$, where the parameters of the physical oscillator are the mass of the ion $m$ and the natural frequency $\nu$ of the ion trap, provides all the scaling factors. The nature of the solution does not allow us to provide reliable numerics for the exact parameter value $\delta=1$.

We can write Husimi Q-function for the reduced density matrix of the boson field,

[TABLE]

where we construct a closed form for the weight factor Satyanarayana (1985),

[TABLE]

in terms of Gauss hypergeometric function $\,{}_{2}F_{1}\left(l,k;r;x\right)$ Lebedev and Silverman (1972). The only eigenstate that can be calculated exactly is the bosonic vacuum, $x=0$ and $\delta=1$.

V Numerical analysis

Calculating the spectrum is an important task. Once the spectrum is known, the intricacies of the system are revealed; for example, it is straightforward to calculate its time evolution, find its phase-space configurations Rodriguez et al. (2018), or even model its interaction with an environment González-Gutiérrez et al., 2017. The QRM is one of a few systems that can be solved Braak (2011). Here, we use numerical diagonalization to observe the changes in the spectrum of our single parameter model as it gets closer to the continuous spectrum regime. In our simulations we use a Fock basis of dimension $30\,000$ for the bosons on resonance with the qubit, $\omega=\omega_{0}$, an ultra-strong coupling parameter, $g=0.99\leavevmode\nobreak\ \omega$, and work only in the bosonic sector related to the excited qubit state in the Foulton-Gouterman reference frame. This is equivalent to working in the positive parity subspace of the model in the laboratory reference frame. Similar results are obtained for the negative parity subspace and different coupling parameter values.

We focus on the statistic of the normalized energy separation between $k$-th nearest neighbors, $s_{k}=(E_{n+k}-E_{n})/\hbar\omega$. Figure 3 shows the histogram for the first nearest neighbor probability distribution, $P(s_{1})$. In the degenerate qubit regime, the energy levels are equidistant and a peak at the value $s_{1}=2$ appears, Fig. 3(a). As the control parameter increases, the histogram displaces towards the origin, two well-defined peaks appear, and the height decreases. For the standard QRM on resonance, $\delta=1/2$ and $\omega=\omega_{0}$, the histogram is centered at the value $s_{1}=1$, Fig. 3(b). As the control parameter increases, the double peak displaces towards the origin. At some critical value of the control parameter, the peak separation starts decreasing and the height increases, Fig. 3(c). As we get closer to the relativistic regime, where the spectrum is continuous, we expect a single peak at the value $s_{1}=0$, Fig. 3(d).

We can correlate this result to Braak conjecture regarding the distance between poles and the distribution of zero, one, and two roots between consecutive poles. First, let us discuss the separation between poles. It is controlled by an integer auxiliary scaled displaced energy $x_{k}=\left[2(1-\delta)E_{k}+g^{2}\right]/\left[2(1-\delta)\right]=k$. In the degenerate case, the energy separation between nearest spectral neighbors is two as the auxiliary scaled energy is $x_{k}=(2\omega E_{k}+g^{2})/(4\omega^{2})=k$. As the control parameter $\delta$ increases, the separation between consecutive nearest spectral neighbors diminishes until, in the relativistic regime, it becomes zero. This correlates with the transition from discrete to continuous spectrum and the displacement of the histogram towards the origin, Fig. 3. Now, the second conjecture, that there exists zero, one, or two roots between consecutive poles, can be related to the double peaked histogram. Sections between consecutive poles that possess two roots will define short nearest spectral neighbor separation. Most probably, the long separation is related to the fact that consecutive two roots occurrences are not feasible.

An interesting feature arises for second nearest neighbors, $s_{2}$, where the spectral separation seems constant. It starts in the degenerate qubit regime as a single peak at the value $s_{2}=4$, Fig. 4(a), that displaces towards the origin at almost constant height, Fig. 4(b) to Fig. 4(c), that should become a single peak at $s_{2}=0$ as we get closer to the relativistic regime, Fig. 4(d). It is straightforward to relate the displacement of the histogram with the separation between poles but we cannot provide a viable conjecture for the constant separation at the moment. As we said before, first neighbor separation points to the existence of two central separations, one short and one large, that interweave and fluctuate, keeping a constant second neighbor separation.

Finally, we follow the ground state of our model to see how the incoherent superposition of odd and even coherent states from the degenerate qubit regime, Fig. 5(a), starts displacing in the $x$-quadrature without noticeable changes in the $y$-quadrature, Fig. 5(b) and Fig. 5(c).

VI conclusions

We proposed a single parameter QRM that interpolates between the degenerate qubit and the relativistic regimes. We provided the analytic solutions in these regimes. The spectrum of the former is discrete and its eigenstates are odd and even displaced number states of the boson field correlated to the excited and ground state of the qubit, in that order. The latter has continuum spectrum and its eigenstates are the unbalanced superposition of infinitely squeezed coherent states correlated to the qubit states.

We conducted a numerical statistical analysis below the transition from discrete to continuous spectrum regimes. We focused on the energy gap between first and second neighbors of the bosonic sector after Foulton-Guterman diagonalization in the qubit basis. We found a structure that points to a distribution of first neighbors separation centered around two peaks with seemingly constant spectral separation for second neighbors. Our statistical results seem to correlate with Braak conjecture regarding the pole and roots distribution for the full QRM. We used Husimi Q-function to follow the process that displaces and squeezes the ground state from the degenerate qubit into the relativistic regime for reasonable control parameter values.

Acknowledgements.

F.H.M.-V. acknowledges funding from CONACYT Cátedra Grupal No. 551. C.H.A. acknowledges funding from CONACYT doctoral grant No. 455378. B.M.R.-L. acknowledges funding from CONACYT CB-2015-01-255230 grant and the Marcos Moshinsky Foundation Research Chair 2018.