Single photons from vacuum

Heisenberg's uncertainty principle implies that the quantum vacuum is not empty but fluctuates. These fluctuations can be converted into radiation through nonadiabatic changes in the Hamiltonian. Here, we discuss how to \emph{control} this vacuum radiation, engineering a single-photon emitter out of a two-level system (2LS) ultrastrongly coupled to a finite-band waveguide in a vacuum state. More precisely, we show the 2LS nonlinearity shapes the vacuum radiation into a nonGaussian superposition of even and odd cat states. When the 2LS bare frequency lays within the band gaps, this emission can be well approximated by individual photons. This picture is confirmed by a characterization of the ground and bound states, and a study of the dynamics with matrix product states and polaron Hamiltonian methods.

Introduction.-Quantumfluctuations underly many physical phenomena, e.g. the Lamb Shift [1] or a modification of the atomic decay.They also explain [2], in a so far educated guess, the cosmological constant problem [3].Vacuum fluctuations can be converted into radiation by nonadiabatic changes of the electromagnetic environment [4], as in the dynamical Casimir effect [5][6][7][8][9], the Unruh effect [10], and the Hawking radiations [11,12].All these processes are explained with free-field theories-quadratic Hamiltonians of coupled harmonic oscillators-that result in Gaussian states [4].In order to create vacuum radiation with nontrivial statistics, such as Fock state, we need nonlinearities, such as quantum emitters.
In this work we study the conversion of vacuum fluctuations into nonclassical radiation states.We focus on waveguide QED [13], studying a two-level system (2LS) coupled to a finite-bandwidth, one-dimensional environment of bosonic modes.This low-dimensional realization of the spin-boson model [14] leads to enhanced light-matter interactions.We will push these interactions to the ultrastrong coupling regime, where interaction strength approaches the energy of the quantum emitter, hybridizing it with the electromagnetic environment [15][16][17][18][19] and leading to nonlinear optics with a single photon [20].Using this last setup and conditions [20], we show how to convert vacuum fluctuations into a train of individual photons.Our protocol consists in abruptly switching on and off the light-matter coupling constants, as sketched in Fig. 1c).We demonstrate that this nonadiabatic process is mediated by photon bound states, which we characterize numerically and analytically.These bound states, once the emitter approaches the band-gap, allow the emission of individual photons without violating the parity constraints of this model.Finally, we also prove that the two-level system serves not only as a precursor of vacuum radiation but also as a detector of quantum fluctuations.
Our proposal is fundamentally different from the emis-sion that occurs when a 2LS is ultrastrongly coupled to a single cavity mode [21][22][23].Here, apart from the control on the emission, we are dealing with the 2LS coupled to a continuum: this results in a strongly correlated state of light and matter and the appearance of the photon bound states that enable our vacuum-triggered singlephoton emitter.Finally, our theoretical proposal can be realized within current experimental status with superconducting circuits [24].Using ultrastrongly coupled flux or transmon qubits in a continuum [18,19] should allow us to test the generation of nonclassical vacuum radiation, enlarging the family of quantum field theory ideas [7][8][9] that superconducting circuits can emulate.
Model.-We study the spin-boson model, a continuum of bosonic modes coupled to a 2LS [25] The σ ± are ladder operators of the 2LS and appear with the excitation energy of the 2LS ∆.The Pauli matrix σ x couples with strengths g k to the bosonic field operators {a † k , a k }, expressed in momentum space.We consider a dispersion relation ω k with N momenta k and a band edge that allows us to control the vacuum emission.Without loss of generality and to facilitate the numerical simulation, we adopt the band ω k = Ω − 2J cos(2πk/N ) [cf.Fig. 1b)] associated to an array of cavities labeled by positions x [cf.Fig. 1a)], with bosonic operators {a x , a † x }, resonator bare frequency Ω and tunneling constant J between neighbouring cavities.The quantum emitter is coupled to the central cavity x = 0, as in H coupling = gσ x (a 0 + a † 0 ), leading to a uniform distribution of couplings g k = g/ √ N in Eq. ( 1).
If the coupling strength is sufficiently small, g ∆, the rotating-wave approximation (RWA) [26] allows us to replace the interaction term with H coupling ∼ = g(σ + a 0 + H.c.), which conserves the number of excitations.In this limit, the ground state has no excitations |GS RWA = |0; 0 and is a product of the 2LS ground state |0 (|1 is the excited state) and the zero-photon state of the waveguide a k |0 = 0.In other words, within the RWA regime, the emitter is immune to the vacuum fluctuations of the bosonic field.However, if g is strong enough, the RWA breaks down and the ground state of (1) contains excitations: GS|a † x a x |GS = 0, suggesting that the 2LS can extract fluctuations and convert them into radiated light.In this work, we investigate precisely those beyond-RWA vacuum emission of the spin-boson model (1).
Theoretical tools.-Thespin-boson model is not exactly solvable, except for a few points and limits, but matrix-product state (MPS) techniques can be used to obtain numerical results [16,20,27], as explained in Appendix B. In this letter, the numerical simulations are contrasted with analytical approximations based on the polaron transformation [28][29][30].This transformation is a disentangling operation U p that decouples the 2LS from the field simplifying the description of low-energy states.The variational parameters f k are obtained by minimizing the ground-state energy E GS within the polaron ansatz for the g.s.|GS = U p |0; 0 , giving the self-consistent equations The simplified Hamiltonian H p = U † p HU p reads Here, h.o.t.stands for higher-order terms of order O(f 3 ) with two and more excitations.Note how the transformed Hamiltonian conserves the number of excitations and can be treated analytically [30].
The renormalization of the 2LS energy ∆ r [c.f.Appendix C] is a consequence of the entanglement between the emitter and the field.According to the polaron picture, most correlations are captured by the unitary transformation of a product state |GS p = |0; 0 = U † p |GS .In this sense, U p plays here a similar role to the Bogoulibouv transformations [31] used for finding the normal modes which account for the radiation in the Hawking, Unruh or Casimir effects [4].
Spectrum of the spin-boson model.-Thespectrum of the Hamiltonian (1) is essential to understand the dynamics of vacuum-induced photon emission.The photonic band edge causes the appearance of photon bound states: localized excitations in the vicinity of the 2LS [20,[32][33][34].We classify those states according to their parity Π = exp iπ(σ + σ − + k a † k a k ) , which is a symmetry of the Hamiltonian (1), [Π, H] = [Π, H p ] = 0.More precisely, the ground state |GS and the second bound state |E 2 are the first and second eigenstates with even parity Π = +1.The first bound state |Ψ 1 is the lowest eigenstate with odd parity Π = −1.Notice that |E 1 and |E 2 have a well-defined number of particles in the RWA limit (1 and 2 respectively) [34].
We compute these states using both MPS and the polaron Hamiltonian.In the first case, parity can be imposed during the MPS minimization of H; in the second case, we project the polaron Hamiltonian (5) onto a fixed parity and number of excitations and diagonalize it numerically.Fig. 2a) shows the energy of the ground state E GS and of the first two bound states, E 1 and E 2 , as a function of the coupling strength g.Note the excellent agreement between MPS (solid line) and the polaron Hamiltonian calculations (dots).Note also how the first bound state lays just below the one-photon band (gray band) E 1 ≤ k (GS) ≡ E GS + ω k , just as in the RWA model [32][33][34][35][36].The second excited bound E 2 enters the band of propagating single photons.There may be other bound states, but the overlap with propagating photon bands of similar parity turn the RWA bound states into resonances with a finite lifetime, where the excitation is no longer localized [20].Further comparisons between results obtained using MPS and the polaron transforma- We have also analyzed the bound state MPS wavefunctions, |E 1 and |E 2 .These states are localized around the 2LS, as seen in Fig. 2b), which renders the number of photons in real space n x = a † x a x .Interestingly, since the MPS produces wavefunctions in the original frame of reference-i.e. after applying U p onto the polaron states-, we find that these states are actual superpositions of different numbers of photons, as seen in Fig. 2c-d-e).The overall superposition preserves the desired parity of the state but, say, a bound state with two excitations can have a nonzero overlap with a single-photon component.
Vacuum emission .-Inorder to convert vacuum excitations into emitted light, we consider a nonadiabatic protocol where the light-matter coupling strength is rapidly switched on and off, as in Fig. 1c).We begin with an unexcited 2LS with g(t < 0) = 0 and switch on the coupling strength to a value g(t = 0) > 0 beyond the RWA regime.The 2LS immediately begins to emit light to accommodate its new entangled ground state.The emitted photon(s) form a wavepacket that travels with bounded speed max k (∂ k ω k ).After some time the 2LS is no longer emitting and the wavepacket leaves behind a localized cloud around the 2LS.We then suddenly switch off the coupling at t = t off and a second vacuum emission event takes place.
We simulate the dynamics described in the previous paragraph by means of MPS.The initial state is the trivial vacuum |Ψ(t = 0) = |0; 0 , which corresponds to the uncoupled case g = 0, and |Ψ(t) evolves with the full Hamiltonian (1) with g within the ultrastrong coupling regime.In Fig. 3a), we plot the average photon number n x = a † x a x along waveguide, as a function of time t and position x.Note how all perturbations emerge from the 2LS position.We switch off the coupling once the travelling photons are far from the emitter.We choose t off = 350τ , being τ the spontaneous decay rate of the 2LS given by the Fermi's golden rule: τ ≡ J sin(k 0 )/g 2 , with k 0 such that ω k0 = ∆.At this point g(t off ) = 0 and we witness the second photon emission event.
The whole process admits a simple description in the polaron picture.The state just before the quench is This is a superposition of even and odd cat states |α ± ≡ In the limit of weak amplitides, the cats tend to one-and two-photon states [17] and the wavefunction can be written using bound states and propagating photons.Asymptotically we have a state with even parity of the form This wavefunction allows four possible outcomes: the system remains in (i) the ground state or (ii) creates the two-photon bound state with localized photons |E 2 and no emission; (iii) it relaxes to the first odd bound state |E 1 emitting a wavepacket A † 1γ with one photon, or (iv) it relaxes to the ground state emitting two photons A † 2γ .Note that when we write this wavefunction in the laboratory basis |Ψ(t) = U p |Ψ(t) p , the structure of the state is preserved, because the polaron transformation is local in space We have tested numerically that Eq. ( 7) captures the vacuum-triggered emission.First, |Ψ(t) p explains the 2LS oscillations P qb (t) = Ψ(t)|σ + σ − |Ψ(t) , shown in Fig. 3b).These oscillations show the interference of the amplitudes |GS and |E 2 .They settle after the qubit relaxation time, t ≥ τ , and persist until we switch off the coupling at t off when P qb gets frozen.
The MPS simulations also confirm that the system emits photons mainly in two channels: one photon on the first excited odd bound state and two photons on the ground state as predicted by Eq. ( 7).This is shown in Fig. 4a), where we plot the average number of photons n x at time t/τ = 250 and the singlephoton n (1) x |E 1 | 2 and two-photon contributions n (2) As seen, the total average n x is well approximated by the sum of both wavepackets n (1) x .Finally, |Ψ(t) p also explains the second photon emission event.In this case, once we switch off the couplings, the bound states become unstable and decay, releasing their photonic components in the form of propagating photons.The distribution of this emitted light for each bound state matches the statistics in Figs.2c-e).
We can control the vacuum induced emission, selecting the one-photon channel, by playing with the band gap ω k=0 and the bound state energy E 1 − E GS .The energies of the radiating states with one and two flying photons are k (E 1 ) = E 1 (g) + ω k and k1,k2 (GS) = E GS + ω k1 + ω k2 , with respective minima E 1 (g) + ω k=0 and E GS + 2ω k=0 .If we place the emitter in the band gap ω k=0 ∆, the energies k (E 1 ) become closer to the 2LS resonance with respect to k1,k2 (GS), so the two photon component is strongly suppressed (|c 2,0 | 0 in Eq. ( 7)).The selectivity of this process is confirmed in Fig. 4b), where the band gap is five times larger than in Fig. 4a) and all other parameters are equal.The final state has a negligible overlap with |E 2 and the distribution of photons P En contains less than 1% of components with n ph > 2. The state before the second quench is faithfully reconstructed by just its single-photon component A † 1γ |E 1 , and as a result, the second emission is also well approximated by one photon.
Conclusions.-In this work we have studied the dynamics of vacuum fuctuations in ultrastrong waveguide-QED setups.More precisely, we have shown that the nonlinearity of a 2LS, combined with a nonperturbative coupling to a bosonic field, can be used to create a vacuumtriggered single-photon emitter.Our proposal is analogous in spirit to other quantum-field theory inspired proposals, such as the dynamical Casimir effect, which work with nonperturbative and nonadiabatic changes of the theory.In contrast to those experiments, we have shown a minimum setup which extracts nonclassical radiation from vacuum effects, using bound states as mediators of these processes.It is important to remark that this whole study can be repeated using a resonator instead of a 2LS [see Appendix D].In this case, all of the features above disappear, as the emission has a Gaussian statistics that cannot be approximated by one-and two-photon states.
Our proposal and the conditions in this work can be re-  The system emits essentially a single-photon packet.
alized in current circuit-QED devices with superconducting qubits that are ultrastrongly coupled to open transmission lines [18,19].In this exciting platform, stateof-the-art measurement techniques would allow for a detailed reconstruction of the photon wavepackets [37,38].The idealized conditions of instantaneous switching can be significantly relaxed or replaced with an alternative protocol, where the qubit is brought in and out of the photonic band.This simpler protocol gives the same qualitative results, as discussed in Appendix A.
Acknowledgments.-We acknowledge support by the Spanish Ministerio de Ciencia, Innovación y Universidades within project MAT2017-88358-C3-1-R and FIS2015-70856-P.The Aragón Government project Q-MAD and CAM PRICYT Research Network QUITEMAD+ S2013/ICE-2801.EU-QUANTERA project SUMO is also acknowledged.Eduardo Sánchez-Burillo acknowledges ERC Advanced Grant QENO-COBA under the EU Horizon 2020 program (grant agreement 742102).As mentioned in the main text, we propose an alternative protocol to couple and decouple the 2LS.Instead of switching on and off the coupling constant g, we modify the excitation energy of the 2LS ∆: we initially take ∆ = ∆ 0 , with ∆ 0 embedded in the band and, eventually, we strongly detune the excitation energy of the 2LS to ∆ f , with the latter far away from the band (see Fig. 5), so the 2LS is effectively decoupled from the photons.

Supplemental Material for: Single photons from vacuum
The emitted field and the dynamics of the 2LS population computed with MPS are rendered in Fig. 6.The emitted field behaves qualitatively as in the couplingdecoupling protocol considered in the main text (see Fig. 3a)).The 2LS population here, however, still evolves, after the quench, depicting some oscillations.The amplitude of these oscillations decreases as ∆ f increases (not shown).
We propose this alternative protocol for the sake of generality.We focus on the coupling quench in the main text because the limit with no coupling is well defined (g = 0), whereas we need to take ∆ f → ∞ with the new protocol to totally decouple the 2LS and the photons.

Appendix B: Matrix-product states
As we indicated in the main text, we use the MPS technique to compute the eigenstates and dynamics of the system.Let us justify why we can do it.
Our initial condition is the state with no excitations.After the nonadiabatic driving, the state is a combination of some of the lowest-energy states with a few flying photons (see Eq. ( 7)).As we are in the low-energy  After the quench, the population still evolves (compare to Fig. 3).Apart from the detuning after the final quench ∆ f , the parameters are those of Fig. 3.
sector, we expect our state to fulfill the area law [39], that is, it will be slightly entangled.Therefore, we may use matrix-product states [40][41][42][43][44], since it is valid for 1D systems when the entanglement is small enough.This ansatz has the form This state is constructed from L sets of complex matrices where each set is labeled by the quantum state s i of the corresponding site.The local Hilbert space dimension d i is infinity, since we are dealing with bosonic sites.However, during the dynamics, processes that create multiple photons are still highly off-resonance.Then, we can truncate the bosonic space and consider states with 0 to n max photons per cavity.So, the composite Hilbert space is H = i C di , where the dimension is d i = n max + 1 for the empty resonators and d i0 = 2(n max + 1) for the cavity with the 2LS.We thus expect the state of the photon-2LS system to consist of a superposition with a small number of photons.In our simulations, we checked that n max = 5 is enough for good convergence.
The number of variational parameters is (L − 1)D 2 (n max + 1) + 2D 2 (n max + 1).In general, the ma-trix size D increases exponentially with L for typical states, whereas its dependence is polynomial if the entanglement is small enough, which usually occurs for lowenergy states.Thus, the number of parameters increases polynomially with L for slightly entangled states.In our simulations, D 10 − 20 proved to be enough.
Our work with MPS relies on three different algorithms.(i) The most basic one is to create the trivial initial state |Ψ(t = 0) = |0; 0 , which is actually a product state.This kind of state can be reproduced using matrices of bond dimension D = 1, so each matrix is just a coefficient A si i = δ si1 .(ii) The second algorithm is to compute expectation values from MPS.This amounts to a contraction of tensors that can be performed efficiently [43], and allows us to compute single-site operators a † i a i , σ z , for instance.(iii) Finally, we can also approximate time evolution, both in real and imaginary times, repeatedly contracting the state with an approximation of the unitary operator exp(−iH∆t) for short times, and truncating it to an ansatz with a fixed D. Since our problem does just contain nearest-neighbour interactions, it is sufficient to rely on a third-order Suzuki-Trotter formula [45].Taking imaginary times, we can obtain the ground state and excited states by solving the equation i d dt P |Ψ = P HP |Ψ .Here, P is either the identity (for the ground state) or a projector that either selects a well defined quantum number (e.g. the parity) or projects out already computed states (for instance the ground and first-excited states).In either case, given a suitable initial state, the algorithm converges to the lowest-energy state of the Hamiltonian in the subspace selected by P .

Appendix C: Further polaron tests
We complement the main text with more calculations within the polaron picture.
First, we show the g-dependence of ∆ r , Eq. ( 4).Notice that it is a self consistent equation, so it is solved numerically.The results are given in figure 7a).We also show the probability for the 2LS to be excited, which can be directly computed from the frequency renormalization, namely: The results are plotted in Fig. 7b).We see the reduction of ∆ r , which resembles the spin-boson.The reason, as seen in Eq. (C1), is that the 2LS and the field are hybridized in the ultrastrong regime, so P qb > 0. On top of that, we check that for g ≤ 0.5 both the polaron transformation and MPS results agree.We now compare the ground state wave-function given by the polaron against the numerical MPS calculations.Using the polaron picture, the ground state can be written as: with (Eq.( 3)) Expanding the exponential we get, Therefore, the ground state one-photon coefficients (second term in (C4)) are given by (up to a normalization) f k .These coefficients are compared in Fig. 8 for g = 0.2.The agreement is rather good.
In the main text, we have assumed that the polaron trasnformation is virtually local.This means that far away from the 2LS the transformation is close to the identity.This makes sense, since the ground state is nontrivial only around the 2LS (Cf.Fig. 1b)).In order to verify this guess, we transform the polaron coefficients f k to real space: In Fig. 9 we show that f x is different from zero only around the 2LS, demonstrating the local character for the polaron transformation in our case.

Appendix D: Linear system
In this appendix, we solve a linear bosonic system analogous to the model of the main text (see Hamiltonian (1)).We compare both spectra as a function of the coupling strength and argue when the spin-boson emission will be nonGaussian.
By linear, we mean that the Heisenberg equations of motion are linear in the bosonic operators, or, equivalently, that the Hamiltonian is quadratic in those operators, that is where α ij and β ij are complex numbers.Besides, α ij = α * ji to ensure H † = H.We might be tempted to substitute the 2LS terms in the original Hamiltonian Eq. (Hamiltonian) by a new resonator: σ − → b and σ + → b † ; however, this does not guarantee the positiveness of the Hamiltonian.In particular, we choose the Caldeira-Leggett model, with Hamiltonian This Hamiltonian is equivalent to the one in Eq. ( 1) by replacing the 2LS by a resonator with operators b, b † and the excitation energy of the 2LS ∆ by a shifted energy ∆, which reads: This shift ensures the positivity of the Hamiltonian.
As argued before, if we apply the nonadiabatical protocol proposed in the main text to this system, it will emit Gaussian light.To see more details on the difference between this and our nonlinear system (Eq.( 1)), we diagonalize this Hamiltonian in the following.Notice that Eq. (D2) can be written in the form of Eq. (D1).We diagonalize Eq. (D2) by means of a Bogoliubov transformation [31].We rewrite the Hamiltonian as that is, ω0 and 2ω 0 , with the gaps of the bound states with respect to the ground state E 1 − E GS and E 2 − E GS (see Fig. 2b)).As seen, both spectra deviate as g increases.Therefore, the behavior of the emitted field for the nonlinear model is expected to deviate from Gaussianity for large enough values of g (roughly at g = 0.2 for the particular parameters of the plot).

Figure 1 .
Figure 1.(a) Sketch of the system.The time-dependent 2LSresonator interaction strength is given by g(t).(b) Dispersion relation ω k = Ω − 2J cos k of the model given by Eq. (2).(c) Profile of g(t) for the nonadiabatic protocol to generate fewphoton pulses.

Figure 2 .
Figure 2. Bound states.(a) Eigenenergies as a function of g for ∆ = 0.3.Continuous lines stand for the MPS simulations and the points for the polaron-ansatz computations.(b) Profiles of the bound-state wavefunctions in position space for g = 0.5 and ∆ = 0.3.(c), (d), and (e) Histograms with the weights in the n ph -photon sector for the ground state, the first excited bound state, and the second one.Same parameters as in panel (b).The parameters defining the photonic waveguide are Ω = 1.0 and J = 0.4.The lattice length is N = 400.

Figure 3 .
Figure 3. Dynamics of the system for the quenching protocol described in the main text: the trivial vacuum is the initial state |Ψ(t = 0) = |0; 0 and we switch the coupling on at t = 0. We finally switch g off at t off /τ = 350.(a) Number of photons as a function of time and position.The system emits a wavepacket at t = 0.At t = t off it radiates again.(b) 2LS-excited-state population.It converges to an oscillatory behaviour with frequency E2 − EGS (red-dashed line).When we decouple the 2LS and the cavities, the population remains constant.In both panels, g = 0.5 after the initial quench and ∆ = 0.3.The other parameters are those of Fig. 2

Figure 4 .
Figure 4. (a) Profile of number of photons at time t/τ = 250 after the quantum quench described in 3.As seen, we can approximate the emitted field with one-and two-photon components.All the parameters are those of the previous figures (see Figs. 2 and 3).(b) Same as before, increasing the energy of the resonators Ω such that the band gap is now five times larger: Ω = 1.8 (remind that the band gap is Ω − 2J).The system emits essentially a single-photon packet.
Appendix A: Detuning the excitation energy of the 2LS and the band

Figure 5 .
Figure 5. Sketch of the alternative protocol.The coupling constant is not modified, but the excitation energy of the 2LS shifts at time t off away from the pohoton band.The blue shaded region stands for the single-photon band of the photonic model.

Figure 6 .
Figure 6.Dynamics analogous to Fig. 3 for the detuning protocol (see discussion in App.A).(a) Number of photons as a function of both time and position.At t off /τ = 350, we set ∆ f = 10, well outside the photonic band and the bound states radiate photons.(b) Population of the excited state of the 2LS.After the quench, the population still evolves (compare to Fig.3).Apart from the detuning after the final quench ∆ f , the parameters are those of Fig.3.

Figure 7 .
Figure 7. (a) The renormalized frequency ∆r as given by Eq. (4) in main text.(b) Populaton of the excited state of the 2LS in the ground state.In both plots, the continuous lines render results from the MPS simulations, whereas the dots represent the results obtained with the polaron ansatz.The parameters Ω, J, and ∆ are those of Fig. 3.

2 Figure 8 .Figure 9 .
Figure 8.We compare the single-photon component of the ground state computed with MPS and the polaron ansatz.The parameters are those of Fig. 3.

Figure 10 .
Figure 10.Comparison of the gap between the first-and second-excited states of the nonlinear model, E1 − EGS and E2 − E − GS (solid red and black lines respectively), with ω0 and 2 ω0 (dashed red and black, respectively) of the linear model (Eq.(D8)).The parameters are those of Fig.2a).