Analogue quantum simulation of superradiance and subradiance in trapped-ions
R. T. Sutherland





sutherland11@llnl.gov
Physics Division, Physical and Life Sciences, Lawrence Livermore National Laboratory, Livermore, California 94550, USA



 
We discuss a protocol for the analogue quantum simulation of superradiance and subradiance using a linear chain of N trapped qubit ions with a single sympathetic cooling ion. We develop a simple analytic model that shows the dynamics of the qubit subspace converge to those of a cloud undergoing Dicke superradiance and subradiance. We provide numerical simulations of the full ion chain and show that they converge to the dynamics predicted by our analytic model with no fitting parameters. We also map out the parameter regime needed to reach this convergence. 


Analogue quantum simulation of superradiance and subradiance in trapped ions
    R. T. Sutherland
    December 30, 2023
============================================================================





§ INTRODUCTION

Trapped ions are one of the most promising platforms for quantum simulations <cit.>. Using either lasers <cit.> or microwaves <cit.> to generate spin-motion coupling, they offer the ability to create a large range of effective Hamiltonians that have controllable spin-spin interactions <cit.> and environmental coupling <cit.>. This has enabled analogue quantum simulations of important many-body phenomena such as superradiant phase-transitions <cit.> quantum transport <cit.>, and the Dirac equation <cit.>. There has yet to be an analogue quantum simulation, however, of superradiance and subradiance <cit.>, which are essential parts of many physical processes, including atomic, biological, and condensed matter systems <cit.>. While the study of coherent emission from quantum systems has yielded important results in its tenure, the Hilbert space describing a typical system of interest will grow exponentially with its size. This has prompted physicists to create simplified theoretical frameworks by making physical approximations such as assuming a low <cit.> or high <cit.> fraction of the atoms to be excited, or by leveraging symmetries in interparticle interactions <cit.>. In order to explore larger and less constrained systems, a quantum simulation for coherent radiation will likely be needed. In this work, we will show that such a simulation is possible in a trapped ion system, where the photon emission from a single auxiliary sympathetic cooling ion in an N qubit system (see Fig. <ref>a), converges to the exact photon emission of an N atom cloud exhibiting Dicke superradiance <cit.>, with an effective single atom decay rate that depends on (controllable) field strengths.



In its original formulation, Dicke superradiance is represented by a cloud of N two-level atoms confined to a space that is small compared to the atoms' transition wavelengh <cit.>, ignoring the strong dipole-dipole interactions present in this regime <cit.>. If all the two-level atoms are placed in the excited state, the system only decays into symmetric superpositions of every possible state with a given number of excitations. Here, using the parallel between two-level and spin-1/2 systems, we represent these Dicke states as |N/2,m⟩. In this notation, m is the spin projection onto the z-axis of an N spin-1/2 system with total angular momentum N/2. As the system decays, correlations between the atoms build, causing the enhanced decay rate of the system to be Γ^' (N/2 + m)(N/2 - m + 1), as opposed to the Γ^'(N/2+m) rate expected in the absence of correlations, where Γ^' is the single atom decay rate. The increase in emission rate as the system evolves gives the signature ∝ N^2 photon intensity pulse associated with superradiance. This has been shown to be equivalent to an ensemble of two-level atoms symmetrically coupled to a lossy cavity <cit.>. It has also been shown that this lossy cavity effect may be generated by sympathetically cooling a Coulomb mode of a mixed-species Penning trap, where the mechanism may be used to effectively create steady-state superradiance <cit.>. In this work, we will use a similar underlying physical mechanism, to create an analogue quantum simulation of the superradiant cascade effect described in Dicke's original work <cit.> in a linear chain of ions with a single auxiliary ion  (same or different species) used for sympathetic cooling. We will show that, when the motion is cooled significantly faster than the rate of spin-motion coupling, the temporal dependence of the qubit evolution, as well as the scattered radiation of the auxiliary ion, converges to that expected from Dicke superradiance and subradiance. We expect that the work shown here can be expanded to include effective spin-spin interactions, and eventually be used for simulations that are unfeasible on a classical computer.



§ THEORY


We consider a system of trapped ions that is comprised of N qubit ions, a single motional mode, and an auxiliary ion for sympathetic cooling (see Fig. <ref>a). The elements of the initial representation of the system are written as:

    |ψ⟩ = |qubit⟩|aux⟩|n⟩,

where |qubit⟩≡|↑↑↓...↑⟩ represents the internal qubit states, |aux⟩ represents the internal state of the auxiliary ion, and |n⟩ represents a motional mode containing n phonons. Here, |↑⟩(|↓⟩) is the excited(ground) state of an individual qubit. The auxiliary ion is either in the strongly radiating state, |e⟩, or the non-radiating ground state, |g⟩; this is achieved with a laser-driven cycling transition on an ion that may or may not be the same species as the qubits. We assume a system initially cooled to the motional ground state, n=0, where the auxiliary ion is in the |g⟩ state. The qubits are then symmetrically driven on the red sideband with a continuous field, while a separate red sideband drive is applied to the auxiliary ion. The Hamiltonian for this system is:

    Ĥ = ħΩ_s{Ŝ_+â + Ŝ_-â^†} + ħΩ_a{σ̂_+â + σ̂_-â^†}.

Note that we are working in the rotating frame with respect to the qubit and motional frequency, and have made the rotating wave approximation. Here, â(â^†) are the phonon annihilation(creation) operators, Ŝ_z,+,-≡∑_j^Nσ̂_z,+,-^j represents a collective spin operator on the qubit subspace, where σ̂_z,+,-^j is a Pauli spin operator for the j^th qubit, and σ̂_+(σ̂_-) are Pauli spin operators for the auxiliary ion. The amplitudes in this equation are the Rabi frequencies of the red sidebands acting on the qubits, Ω_s, and the red sideband acting on the auxiliary ion, Ω_a. 

In order to perform the many-ion calculations in Sec. <ref>, we simplify the equations by noting that when the operator, Ŝ_+,-, acts on a superradiant state, given by:

    |N/2,m⟩≡{(N-m)!m!/N!}^1/2∑_qubits|↑↓↓↑...↑⟩,

where the sum in the above equation is over all states with N/2 + m qubits in the excited state, it only couples to |N/2,m±1⟩ with non-zero matrix elements. Note that |N/2,m⟩ is an eigenstate of Ŝ^2 and Ŝ_z. This allows us to greatly reduce the qubit subspace from the 2^N original elements to the N+1 possible |N/2,m⟩ states. We can then rewrite Eq. (<ref>) as:

    Ĥ   =    ħΩ_s{D̂_+â + D̂_-â^†}  + ħΩ_a{σ̂_+â + σ̂_-â^†},

where we have rewritten the qubit operators as:

    D̂_-≡ ∑_m{(N/2 + m)(N/2 - m + 1)}^1/2|N/2,m-1⟩⟨N/2,m|,

and D̂_+ = D̂_-^†, thereby greatly reducing the complexity of the above equations. We assume that the radiation from the qubits themselves is negligible. The auxiliary ion, however, is taken to radiate quickly. We account for this using the Lindblad formalism for a single two-level system's decay, where the photon bath has been traced over. The full master equation for the density operator of the system, ρ̂, is:

    ρ̇̂̇   =    -i/ħ[Ĥ,ρ̂] + Γ{σ̂_-ρ̂σ̂_+ - 1/2σ̂_+σ̂_-ρ̂ - 1/2ρ̂σ̂_+σ̂_-},

where Γ is the decay rate of the auxiliary ion. Note that the full Eq. (<ref>) is used for all the numerical calculations of trapped ions in this work, with no further approximations.



 §.§ Analytic Model


The superradiant cascade is where the Dicke state, |N/2,m⟩, emits a photon into a Markovian bath and is projected onto |N/2,m-1⟩ at a rate of Γ^' (N/2 + m)(N/2 - m + 1) (Fig. <ref>b). We simulate the dynamics of this system using trapped ions with the protocol illustrated in Fig. <ref>c. Initially, the qubit subspace is in the |N/2,m⟩ state, the auxiliary ion is in the |g⟩ state, and the motional mode is in the |0⟩ state, making the overal state |N/2,m⟩|g⟩|0⟩. The red sideband of each of the qubits is then driven with a rate Ω_s so that the system evolves into the |N/2,m-1⟩|g⟩|1⟩ state. Simultaneously, the red sideband of the |g⟩↔|e⟩ transition is driven at Ω_a, a much faster rate than Ω_s, taking the system from |N/2,m-1⟩|g⟩|1⟩ to |N/2,m-1⟩|e⟩|0⟩. Finally, at a rate Γ much faster than Ω_s, the auxiliary ion decays, taking the state from |N/2,m-1⟩|e⟩|0⟩ to |N/2,m-1⟩|g⟩|0⟩, repeating this cycle until |N/2,-N/2⟩|g⟩|0⟩ is reached. Note that we are assuming the transition |N/2,m-1⟩|g⟩|1⟩↔|N/2,m-2⟩|g⟩|2⟩ can be ignored because we are exploring a regime where sympathetic cooling is much faster than this process.

Reference <cit.> derived the rate that population in |N/2,m⟩ decays into |N/2,m-1⟩ (Fig. <ref>b). Similarly, we wish to determine the effective decay rate, Γ_eff, representing the transfer of population from |N/2,m⟩|g⟩|0⟩ to |N/2,m-1⟩|g⟩|0⟩ (Fig. <ref>c). To do this, we assume an initial population in |N/2,m⟩|g⟩|0⟩, and calculate the rate of decay out of the entire {|N/2,m⟩|g⟩|0⟩,|N/2,m-1⟩|g⟩|1⟩,|N/2,m-1⟩|e⟩|0⟩} subspace. We represent the density matrix for this subspace as ρ̂_sub. 

To begin, we note that Ĥ only couples terms within ρ̂_sub, resulting in no entanglement with the rest of the system. We can therefore focus on the decay of an initial population out of ρ̂_sub in order to determine Γ_eff. The Γσ̂_-ρ̂σ̂_+ term in Eq. (<ref>) represents the transfer of population into ρ̂_sub, which is not a component of the process that we currently wish to analyze; we can thus ignore this term for our analytic model. These considerations allow us to write Eq. (<ref>) for the subspace as:

    ρ̇̂̇_sub = -i/ħ(Ĥ_cρ̂_sub - ρ̂_subĤ_c^†),

where Ĥ_c≡Ĥ -iħΓ/2σ̂_+σ̂_- is a non-Hermitian Hamiltonian acting on ρ̂_sub <cit.>. The above equation is equivalent to a density matrix for a wave function being propagated by Ĥ_c. We use this fact to simplify the analytic calculation of Γ_eff, considering only a wave function for the subspace acting under Ĥ_c. 

Letting c_0, c_1, and c_2 be the probability amplitudes of the |N/2,m⟩|g⟩|0⟩, |N/2,m-1⟩|g⟩|1⟩ and |N/2,m-1⟩|e⟩|0⟩ states, respectively, we can find the equation of motion for c_1:

    ċ_1   =     - iΩ_s{(N/2 + m)(N/2 - m + 1)}^1/2c_0 - iΩ_ac_2.

We now assume that (upon reaching quasi-equilibrium) ċ_1≃ 0 relative to Ω_a and Γ. This allows us to set the left hand side of Eq. (<ref>) to 0 and solve for c_2 in terms of c_0:

    c_2 = -Ω_s{(N/2 + m)(N/2 - m + 1)}^1/2c_0/Ω_a,

where we can now easily solve for Γ_eff, defined by Γ_eff(N/2 + m)(N/2 - m + 1)|c_0|^2 = Γ|c_2|^2. This results in a final equation:

    Γ_eff = (Ω_s/Ω_a)^2Γ.

This is analogous to an N atom system undergoing superradiance with an effective single atom decay rate of Γ^' = Γ_eff. 



§ RESULTS



 §.§ Superradiance



For an ensemble exhibiting the form of superradiance described in Ref. <cit.>, a Dicke state, |N/2,m⟩, decays at a rate given by Γ^' (N/2 + m)(N/2 - m + 1). This indicates that if an ensemble is initialized to the |N/2,N/2⟩ state then the photon emission rate of the ensemble, γ = Γ⟨σ̂_+σ̂_-⟩ = Tr{Γσ̂_+σ̂_-ρ̂}, will initially increase with time, followed by a rapid decrease, when the system runs out of energy, emitting a pulse of photons with intensity ∝ N^2. In this subsection, we demonstrate that γ from a single auxiliary ion in an N qubit ion chain converges to the exact time dynamics expected in a superradiant system with a single atom decay rate of Γ^' = Γ_eff, when Γ≫Ω_s and Ω_a^2/Γ≫Ω_s.

The convergence of our trapped ion system to Dicke superradiance is shown in Fig. <ref> for systems of N = 1,5,10,20 and 30. This is done using Eq. (<ref>) without approximation. We here show the temporal evolution of γ normalized by Γ^' and N, γ^'≡γ/NΓ^'. We compare γ^' for the entire superradiant cloud with γ^' from the auxiliary ion in an N qubit trapped ion system, and set Γ^' to Γ_eff for comparison. For the trapped ion system, we choose values of Ω_a≃ 81.2Ω_s and Γ = 200Ω_s, chosen to be well-within the regime of convergence (discussed below). At very short times, γ^' from the auxiliary ion in the trapped ion system is close to zero. However, at a timescale ∝ 1/Γ, the system reaches quasi-equilibrium, and, as a result, follows the dynamics discussed in Sec. <ref>. This initial “ramp-up" time, is shown in the inset of Fig. <ref>. Here, after a time, t, that is short compared to Ω_s, each system of N qubits converges to the system they are supposed to simulate. 



In order for γ^' versus t from the auxiliary ion to converge to that of Dicke superradiance, the parameters of the system must be set up such that the system reaches quasi-equilibrium on a time-scale that is fast compared to the system's overall evolution (∝ 1/Γ_eff), and that, at quasi-equilibrium, population is removed from the motional mode at a rate that is fast compared to Ω_s such that there is no back-action into the initial state. As we show in Fig. <ref>, these two conditions are met when Γ≫Ω_s and Ω_a^2/Γ≫Ω_s, respectively. 

The first condition, Γ≫Ω_s, is due to the fact that the damping provided by the photon emission from the auxiliary ion will cause the system to reach quasi-equilibrium at a rate ∝Γ. Since the auxiliary ion will mimic superradiant emission when it has reached quasi-equilibrium (see Sec. <ref>), it must reach quasi-equilibrium on a timescale that is much shorter than the overall evolution of the system. Also, if the system is not damped at a much faster rate than population is put in, back-action will lead to Rabi oscillations; this will significantly complicate the dynamics.  We show this in Fig. <ref>a, a system of N=10 qubits coupled to an auxiliary ion that decays at a rate Γ. Note that we choose a set value of Ω_a^2/Γ = 30Ω_s in this graph, so that—as will be discussed below—population is removed from the motional mode at a rate that is large enough to reach convergence (see Fig. <ref>b). Here it can be seen that for values of Γ≳ 75 Ω_s the evolution converges to that expected from Dicke superradiance, but for smaller values of Γ the system deviates from the desired evolution.

The second condition, Ω_a^2/Γ≫Ω_s, is that once the system has reached quasi-equilibrium, population must be removed from the motional mode (|N/2,m-1⟩|g⟩|1⟩ state) at a rate that is significantly faster than it is put in. We can quantify this rate by using a similar argument as was used to obtain Eq. (<ref>), but for the reduced Hilbert space of just |N/2,m-1⟩|g⟩|1⟩ and |N/2,m-1⟩|e⟩|0⟩. As we originally assumed that the population of |N/2,m⟩|g⟩|0⟩ was fixed in Eq. (<ref>), we now do the same for |N/2,m-1⟩|g⟩|1⟩. This gives:

    ċ_2 = -iΩ_ac_1 - Γ/2c_2.

Setting the left hand side to 0, we get an effective decay rate out of c_1 of 4Ω_a^2/Γ; this means that, at quasi-equilibrium, the value of Ω_a^2/Γ must be large compared to Ω_s, so that the transfer of population from the motional mode to the quickly radiating state is fast enough for the former to act as a lossy cavity. This is shown in Fig. <ref>b, where we can see that our calculations converge to Dicke superradiance only for values of Ω_a^2/Γ≳ 20Ω_s; for larger values of Ω_a^2/Γ, except for Rabi oscillations at small t, we find that the system converges. Figure <ref> also shows that when Γ≫Ω_s, the dynamics of the system is entirely dictated by the value of Ω_a^2/Γ. This is here seen in the fact that calculations for systems such that Γ = 200Ω_s and Γ = 400Ω_s fall on the same lines, reaching convergence when Ω_a^2/Γ≳ 20Ω_s.

In terms of measuring the properties of our superradiant system, one could directly measure the photon emission from the auxiliary ion; for illustrative purposes, we have shown this for the figures in this section. Depending on the particular experiment, however, it may be easier to directly measure the internal states of the qubits using standard techniques <cit.>. This will be particularly useful in when measuring systems that do not radiate strongly, as discussed in the next subsection.



 §.§ Subradiance




The trapped ion system we describe in this work can mimic the dynamics of subradiance as well. Since subradiance is a significantly harder effect to probe experimentally (due to its diminished emission) <cit.>, this simulation offers a unique opportunity to probe subradiant effects; due to the exceptional readout capabilities of trapped ions <cit.>, one could generate a system that behaves in a subradiant manner, and probe it by simply measuring the internal states of the qubits; this is likely the optimal way to probe this simulation due to the fact that γ^' is small and thus hard to measure when the system mimics subradiance.

In order to demonstrate the capacity of our system to simulate subradiance, Fig. <ref> shows the fidelity, ℱ≡⟨ψ_0|ρ̂|ψ_0⟩, versus Γ_efft for a 2 qubit system placed in the initial state: 

    |ψ_0⟩≡1/√(2){|↑↓⟩ + e^iϕ|↓↑⟩}|g⟩|0⟩.

This can be generated using a two-qubit entangling gate to create a Bell state followed by a single qubit rotation <cit.>. Note that, here, we simulate the system in the original qubit basis, see Eq. (<ref>), since |ψ_0⟩ is outside of the described |N/2,m⟩ basis. Here we change the value of ϕ to vary the state from superradiant (ϕ = 0), where the state decays at twice the single particle decay rate (ℱ(t) = exp{-2Γ_efft}) to subradiant (ϕ = π), where the initial state does not decay at all ℱ(t) = 1). Even though, in the (idealized) latter case, the auxiliary ion does not emit any radiation, this effect should be observable through the direct measurement of the qubits.



§ EXPERIMENTAL IMPLEMENTATIONS

The two parameter regimes required to generate Dicke superradiance (Γ≳ 75Ω_s and Ω_a^2/Γ≳ 20Ω_s) can, in principle, be met for any choice of Γ and Ω_a by making Ω_s sufficiently weak. This could, of course, result in a small Γ_eff which would lead to a long experimental run time, so it is important to determine what this time would be in a typical experiment. Assuming that Ω_a^2/Γ≃ 20Ω_s, we find Γ_eff≃Ω_a^2/400Γ, and that this condition automatically enforces Γ≳ 75Ω_s for typical values of Γ. As an example, if the auxiliary ion is ^9Be^+ the value of Γ would be approximately 2π× 20 MHz <cit.>, corresponding to the P_3/2→ S_1/2 transition. For a value of Ω_a = 2π× 500 kHz, this would make the convergence criteria for Ω_s≃ 2π× 625  Hz and the timescale of the experiment 1/Γ_eff≃ 5  ms. The speed of this scheme is limited by the fact that a typical sympathetic cooling transition decays quickly and requires a strong laser to sufficiently drive the red sideband. The timescale of a potential experiment as well as the required Ω_a could be decreased if one had access to a cooling transition with a smaller value of Γ. 



§ CONCLUSION


In this work, we showed that a trapped ion chain can be made to exhibit Dicke superradiance and subradiance, in a manner that converges to the exact results from the original system. This offers the opportunity to study the phenomena in a way that allows experimentalists to control the effective single atom decay rate, Γ_eff, by adjusting the field amplitudes of the system. One could even stop the decay of the system at any time and probe any of its particles. For illustrative purposes, we have focused on the simulation of the system originally envisioned by Dicke in Ref. <cit.>, since its symmetries allow us to compare our analogue to the exact answer. Going forward, these ideas could be combined with existing techniques that have used trapped ions to create analogue quantum simulations of other interesting systems <cit.>.  Finally, with various theoretical and technical advances, it could be used to explore the physics of systems not easily probed experimentally or simulated on a classical computer.





§ ACKNOWLEDGEMENTS

The author would like to thank F. Robicheaux, D.H. Slichter, R. Srinivas. M.T. Eiles, S.B. Libby, J. L. DuBois, and S.C. Burd for helpful discussions. Special thanks to F. D. Malone, S.C. Burd, and Y. Rosen for careful reading of the manuscript. Part of this work was performed under the auspices of the US Department of Energy by Lawrence Livermore National Laboratory under Contract DE-AC52-07NA27344. LLNL-JRNL-791496-DRAFT



§ REFERENCES



