Collective transition quenching in the presence of multiple competing decay channels

Abstract

We present a theoretical framework for ‘collective transition quenching’, a quantum many-body dissipative phenomenon that occurs in systems with multiple collective decay channels. Despite the competition, interactions suppress all but the dominant decay transition, leading to a ‘winner takes all’ dynamic where the system primarily settles into the dominant ground state. We prove that, in the presence of permutation symmetry, this problem is exactly solvable for any number of competing channels. Additionally, we develop an approximate model for the dynamics by mapping the evolution into a continuity equation for a fluid, and show analytically that the dominant transition ratio converges to unity with increasing system size as a power-law, for any branching ratio. This near-deterministic preparation of the dominant ground state has broad applicability. As an example we discuss a protocol for molecular photoassociation where collective dynamics effectively acts as a catalyst, amplifying the yield in a particular final state. Our results open new avenues for many-body strategies in the preparation and control of quantum systems.

The suppression of decay pathways into undesired quantum states is crucial for the control and manipulation of open quantum systems. Simplified theoretical models often rely on ‘closed’ transitions where the excited state decays predominantly to a specific ground state. Yet, in practice, real-world emitters seldom adhere to the idealized paradigm of two-level systems. For instance, highly excited atomic states (such as Rydberg states) have many lower-energy states accessible by spontaneous emission [1]. In photochemistry, decay to a single ground state is often inefficient due to numerous competing pathways arising from electronic, vibrational, and rotational degrees of freedom [2, 3]. Solid-state emitters (such as color centers or dye molecules) also suffer from parasitic decay from phonon sidebands [4, 5]. Achieving closed transitions in experiments is challenging, and often involves isolating two-level systems from more complex internal structures or employing repumping techniques to redirect population back into the desired states. However, these approaches are inherently limited by the natural, single-particle branching ratios of the transitions.

The natural branching ratios and dynamics of decay can, in fact, be modified by engineering the quantum system and its environment. One common approach is to tailor the dielectric environment by placing emitters within optical cavities [6], waveguides [7, 8], or other photonic structures [9, 10] such that the desired decay channel is Purcell-enhanced. Numerical studies with multilevel atoms [11, 12, 13, 14] and molecules [15] have suggested collective emission as an alternative to circumvent limitations from single-particle branching ratios. These proposals rely on many-body, transient superradiance [16, 17], a phenomenon characterized by avalanche-like behavior [18, 19, 20] where decay into a given ground state enhances the probability of subsequent emission into that same state. This process effectively steers the emitters towards a specific ground state. A comprehensive analytical treatment of this physics remains lacking, which is critical for fully understanding and exploiting its potential.

We term the phenomenon where correlated decay suppresses all but the most dominant emission path as collective transition quenching. As shown in Fig. 1(a), we consider an ensemble of emitters coupled to a reservoir. The emitters have multiple decay channels, each leading to different final states. Each decay channel can be collectively enhanced by many-body correlations that emerge dynamically, leading to competition between them. We find approximate steady-state solutions for the populations of the different ground states by modeling the quantum dynamics through a continuity equation for a fluid. We prove that, despite the competition, the population density of the dominant ground state exhibits a power-law convergence to unity for any branching ratio, with the power-law exponent characterized by the ratio between dominant and subdominant decay rates. This is supported by a rigorous analysis of the exact steady-state solution. We apply our framework to the problem of photochemistry, where sub-optimal branching ratios limit the effectiveness of molecule creation, direct laser cooling [21], and optical imaging. Specifically, in molecular photoassociation of strontium dimers, we demonstrate that collective transition quenching greatly enhances sample purity.

(a)

Figure 1: Increased population of the dominant transition (with suppressed fluctuations) due to collective decay in multilevel systems. (a) $N$ multilevel emitters with $d$ ground states and a single excited state decay collectively to an environment (at rates $\Gamma_{1},\ldots,\Gamma_{d}$), such as a “bad” cavity, a waveguide in the “mirror configuration” (where the relative distance between the emitters is a half-integer multiple of the resonance wavelength $\lambda_{0}$) or free space (for the latter, a dense ensemble of subwavelength volume is required to preserve permutational symmetry). (b) Dissipative dynamics can be modeled as a random walk between permutationally symmetric ground states ($d=2$ depicted), labelled by $\ket{n_{1},\ldots,n_{d}}$, where $n_{\mu}$ denotes the population of the ground state $\ket{g_{\mu}}$. (c) Marginal probability distribution $P(n_{1})=\sum_{n_{2}}P(n_{1},n_{2})$ for the population density of the ground state $\ket{g_{1}}$, for $N=150$ and $d=2$ in the steady state, after complete depletion of the fully-inverted initial state. The distributions for collective decay with $r_{2}=1$ (gray), and collective decay with $r_{2}=0.5$ (blue) are obtained via numerical simulations of Eq. (2). The binomial distribution is plotted for independent decay with $r_{2}=0.5$ (green).

We consider an (undriven) ensemble of $N$ identical emitters with a level structure consisting of a single excited state $\ket{e}_{i}$ and $d$ ground states labelled $\ket{g_{\mu}}_{i}$, with $\mu\in\{1,\ldots,d\}$ and $i\in\{1,\ldots,N\}$. The emitters are symmetrically coupled to a Markovian environment, with separate couplings for each transition $\ket{e}\to\ket{g_{\mu}}$. Physically, this can be realized by coupling the emitters to near-resonant cavity modes (in the “bad cavity” or “weak coupling” limit), to a single-mode waveguide in the “mirror configuration” (where the emitters are separated by a half-integer multiple of the wavelength [22, 23]), or via free-space interactions in a dense ensemble of subwavelength volume, as depicted in Fig. 1(a). We assume that photons from different transitions can be resolved either by polarization or frequency (for the waveguide configuration, we assume that the transition frequencies are approximately equal to fulfil the mirror condition). The system dynamics is modelled by the master equation

$$\dot{\rho}=-i\left[\sum_{\mu=1}^{d}\chi_{\mu}\hat{A}_{\mu}^{{\dagger}}\hat{A}_{\mu},\rho\right]+\sum_{\mu=1}^{d}\Gamma_{\mu}\mathcal{D}[\hat{A}_{\mu}]\rho,$$ (1)

where $\chi_{\mu}$ are the coherent interaction rates, $\hat{A}_{\mu}=\sum_{i=1}^{N}(\ket{g_{\mu}}\bra{e})_{i}$ is the collective lowering operator on the transition $\ket{e}\rightarrow\ket{g_{\mu}}$, and $\mathcal{D}[\hat{A}]\rho\equiv\hat{A}\rho\hat{A}^{\dagger}-\{\hat{A}^{\dagger}\hat{A},\rho\}/2$. Collective dissipation occurs with rates $\Gamma_{\mu}$ (identical to the single-emitter ones), and we assume the ordering $\Gamma_{1}\geq\Gamma_{2}\geq\ldots\geq\Gamma_{d}$, such that $\Gamma_{1}$ is the dominant decay rate. For convenience, we denote the ratio between dominant and subdominant decay rates as $r_{\mu}\equiv\Gamma_{\mu}/\Gamma_{1}$.

(a)

Figure 2: (a) Velocity vector field of the continuum model for $r_{2}=0.5$ and $N=10$. The dashed line represents the flow trajectory in the single-particle approximation starting from the fully excited state $(x_{1},x_{2})=(0,0)$, where $x_{1(2)}$ denotes the population of $\ket{g_{1(2)}}$. Darker color indicates lower velocity. The red solid line shows the neutrally stable ground state attractor. The unphysical region (where the ground state population is larger than $N$) is shown in gray. (b) Population density outside the dominant ground state [$1-\text{Dominant transition ratio (DTR)}$] against number of emitters $N$, for various decay ratios $r_{2}=\Gamma_{2}/\Gamma_{1}$. The points are obtained from numerically solving the rate equation (2), while the dotted lines denote the approximate formula (9). Solid lines denote the exact asymptotic solution (10). (c) Cumulative distribution (i.e., total probability of the dominant ground state population being at most $n_{1}$) $P_{c}(n_{1})=\sum_{n_{1}^{\prime}=0}^{n_{1}}P(n_{1}^{\prime})$. Data obtained from numerical simulation of the rate equation (2), with $r_{2}=0.4$. Darker colors indicate larger number of emitters, with $61\leq N\leq 10^{4}$. All plots are made for $d=2$ ground states.

Dynamical evolution can be understood as a random walk in the subspace of permutationally symmetric states, as shown in Fig. 1(b). Since Eq. (1) preserves permutation symmetry, basis states $\ket{n_{1},\ldots,n_{d}}$ are fully described by occupation numbers, where $n_{\mu}$ is the population of $\ket{g_{\mu}}$. These are typically entangled, as they consist of a symmetric superposition of excitations over $N$ particles. The population of the excited state is $n_{e}=N-\sum_{\mu}n_{\mu}$, with $n_{e}=0$ in the final, steady state. Employing the ansatz $\rho=\sum_{\{n_{\mu}\}}P_{n_{1},\ldots,n_{d}}\ket{n_{1},\ldots,n_{d}}\bra{n_{1},\ldots,n_{d}}$, Eq. (1) reduces to a rate equation [see Supplementary Material (SM) [24]]

$$\begin{split}\dot{P}_{n_{1},\ldots,n_{d}}&=-\left(N-\sum_{\nu=1}^{d}n_{\nu}\right)\sum_{\mu=1}^{d}\Gamma_{\mu}(n_{\mu}+1)P_{n_{1},\ldots,n_{\mu},\ldots,n_{d}}\\ &+\left(N-\sum_{\nu=1}^{d}n_{\nu}+1\right)\sum_{\mu=1}^{d}\Gamma_{\mu}n_{\mu}P_{n_{1},\ldots,n_{\mu}-1,\ldots,n_{d}}\end{split}$$ (2)

for the probabilities $P_{n_{1},\ldots,n_{d}}$ of occupying the state $\ket{n_{1},\ldots,n_{d}}$. This ansatz is justified even if there are coherences in the initial state (i.e., off-diagonal terms in the density matrix), since they are decoupled and hence do not affect population dynamics. While the rate equation holds for any permutationally symmetric initial state, below we choose the fully inverted state (i.e., $\ket{e}^{\otimes N}$).

Although Eq. (2) can be efficiently simulated, it is non-trivial to obtain the steady state analytically. Moreover, the steady state is highly non-unique, since any combination of emitters in any ground state is a possible steady state. The trivial case of just one collective decay channel $(d=1)$ reduces to the problem of Dicke superradiance [16, 17]. We instead focus on the problem of multiple competing transitions.

We quantify collective transition quenching by the dominant transition ratio (DTR), defined as the mean population density of the dominant ground state in the steady state, i.e.,

$$\text{DTR}=\frac{\overline{n_{1}}}{N},$$ (3)

where $\overline{n_{1}}$ is the mean number of emitters in $\ket{g_{1}}$ as $t\to\infty$. For independent emitters, the marginal probability distribution $P(n_{1})=\sum_{n_{2},\ldots,n_{d}}P_{n_{1},\ldots,n_{d}}$ for the dominant ground state is a binomial distribution whose average (normalized by $N$) is the DTR $=(1+\sum_{\mu>1}r_{\mu})^{-1}$, determined solely by the decay ratios and independent of $N$ [see Fig. 1(c)]. Below, we prove that the DTR always converges to 1 as $N\to\infty$ for any number of collective decay channels, assuming $\Gamma_{1}>\Gamma_{\mu>1}$. This effect can be attributed to the superradiant enhancement, which amplifies the dominant transition relative to all subdominant transitions. If $\Gamma_{1}=\Gamma_{2}=\ldots=\Gamma_{d}$, one obtains a uniform distribution for $P_{n_{1},\ldots,n_{d}}$ [13], as shown in Fig. 1(c).

To gain analytical insights from the rate equation, we transform it into a continuity equation of a fluid that flows in Euclidean space. In the limit of large $N$, we make the continuum approximation [25] by setting $n_{\mu}\to x_{\mu}$, $\vec{x}=(x_{1}\ldots x_{d})^{T}\in\mathbb{R}_{+}^{d}$, and $f(n_{1},\ldots,n_{\mu},\ldots,n_{d})-f(n_{1},\ldots,n_{\mu}-1,\ldots,n_{d})\to\frac{\partial}{\partial x_{\mu}}f(\vec{x})$ for an arbitrary differentiable function $f$. The rate equation is then approximated by the continuity equation

$$\frac{\partial}{\partial t}P(\vec{x},t)=-\mathbf{\nabla}\cdot(\vec{v}(\vec{x})P(\vec{x},t)),$$ (4)

which describes a fluid flow in $\mathbb{R}_{+}^{d}$, governed by the position-dependent velocity field $\vec{v}(\vec{x})$ with $\mu$-th component

$$v_{\mu}(\vec{x})=\Gamma_{\mu}\left(N-\sum_{\nu=1}^{d}x_{\nu}\right)(x_{\mu}+1).$$ (5)

The flow comes to a stop as $\sum_{\nu}x_{\nu}\to N$, which physically corresponds to the system approaching its steady state, where no excited state population remains. The velocity field is illustrated for $d=2$ in Fig. 2(a). The fully excited initial state corresponds to the origin of $\mathbb{R}_{+}^{d}$. Due to the non-uniqueness of the steady state, it is not sufficient to simply solve for $\partial_{t}P=0$.

Solving Eq. (4) analytically for arbitrary times is a formidable task. Instead, we employ a single-particle approximation where the fluid is idealized as a point particle initialized at the origin, with velocity dynamics $\text{d}\vec{x}/dt=\vec{v}(\vec{x})$. This reduces the partial differential equation in Eq. (4) to a coupled system of $d$ ordinary differential equations. Dividing the components of $\vec{v}(\vec{x})$, we readily find

$$\frac{\text{d}x_{\mu}}{\text{d}x_{\nu}}=\frac{\Gamma_{\mu}}{\Gamma_{\nu}}\frac{x_{\mu}+1}{x_{\nu}+1}$$ (6)

for any pair of $\mu,\nu\in\{1,\ldots,d\}$. Integrating the above expression yields

$$(x_{\mu}+1)^{\Gamma_{\nu}}=(x_{\nu}+1)^{\Gamma_{\mu}}.$$ (7)

The $d-1$ independent equations of the form (7) define the particle trajectory in $\mathbb{R}_{+}^{d}$ starting from the origin. Generalizing Eq. (7) to other permutationally symmetric initial states is straightforward (see SM [24]).

The steady state solution has a geometrical interpretation as the intersection between the particle trajectory and the hyperplane $\sum_{\mu=1}^{d}x_{\mu}=N$. To find an asymptotic solution, we assume that $x_{\mu}\gg 1$. The trajectory is then described by a simpler set of equations $x_{\mu}^{\Gamma_{\nu}}=x_{\nu}^{\Gamma_{\mu}}$. The constraint $\sum_{\mu=1}^{d}x_{\mu}=N$ can then be rewritten as $\sum_{\mu=1}^{d}x_{d}^{\Gamma_{\mu}/\Gamma_{d}}=N$. For simplicity, we consider the non-degenerate case where all $\Gamma_{\mu}$ are distinct (see SM [24] for the degenerate case). Using the method of dominant balance [26] we obtain $x_{d}\approx N^{r_{d}}$, valid in the regime $1\ll N^{r_{d}}\ll N$. In the regime $1\ll N^{r_{\mu}}\ll N\;\forall\;\mu>1$, an iterative method yields the self-consistent solution

$$x_{1}\approx N-\sum_{\mu=2}^{d}{N^{r_{\mu}}},\quad x_{\mu}=N^{r_{\mu}},\quad\mu>1.$$ (8)

The dominant transition ratio thus reads

$$\text{DTR}\approx\frac{x_{1}}{N}\approx 1-\sum_{\mu=2}^{d}\frac{1}{N^{1-r_{\mu}}},$$ (9)

which converges to unity as $N\to\infty$, with the slowest convergence characterized by the power law $\sim N^{r_{2}-1}$. Equation (9) provides our first theoretical prediction for collective transition quenching. Comparing with numerical simulations of the rate equation (2) for $d=2$ in Fig. 2(b), the approximate formula (9) agrees qualitatively, with higher accuracy attained for smaller $r_{2}$.

Remarkably, we can solve for the DTR exactly in the asymptotic $N\rightarrow\infty$ limit. In the SM [24], we prove that the dynamics governed by the rate equation (2) is integrable, and derive the complete set of $N+1$ independent conserved quantities. We thus overcome the problem of non-unique steady states. However, because of the complexity of the exact solution, we only use it to compute the DTR, which reads

$$\text{DTR}=\frac{\overline{n_{1}}}{N}=1-\sum_{\mu=2}^{d}\frac{\tilde{\Gamma}(1-r_{\mu})}{N^{1-r_{\mu}}},$$ (10)

where $\tilde{\Gamma}(\cdot)$ is the gamma function. The exact solution yields the same scaling as the approximate formula (9), and agrees excellently with numerical simulations, as shown in Fig. 2(b). This agreement explains why the approximation of Eq. (9) becomes more accurate for smaller $r_{\mu}$, since $\tilde{\Gamma}(1-r_{\mu})=1+O(r_{\mu})$. The error in the prefactor of Eq. (9) likely arises from the single-particle approximation.

Our formalism can also be easily extended to include non-collective decay channels, modeled as a leakage at a rate $\Gamma_{\text{leak}}$ [24]. Going back to the fluid model, this adds an extra component $v_{0}=\Gamma_{\text{leak}}\left(N-\sum_{\mu=1}^{d}x_{\mu}-x_{0}\right)$ to the velocity field and modifies $\sum_{\nu=1}^{d}x_{\nu}\to\sum_{\nu=0}^{d}x_{\nu}$ in Eq. (5). A similar analysis for large $N$ yields $x_{0}\approx(\Gamma_{\text{leak}}/\Gamma_{1})\ln N$, in the same regime as the validity of Eq. (9). This justifies the omission of non-collective decay in our model, since $x_{0}\ll x_{\mu}\sim N^{r_{\mu}}$, and the approximate DTR reads

$$\text{DTR}\approx 1-\sum_{\mu=2}^{d}\frac{1}{N^{1-r_{\mu}}}-\frac{\Gamma_{\text{leak}}}{\Gamma_{1}}\frac{\ln N}{N}.$$ (11)

In the non-competing scenario with only one collective decay channel $(d=1)$, the DTR always converges to unity as $\sim\ln N/N$ (even if $\ket{e}\to\ket{g_{1}}$ is not the most dominant transition). This recovers the well-established result of Dicke superradiance in the presence of local decay [15, 27].

The avalanche-like behavior of the dominant transition not only impacts the population of the dominant ground state, but also lowers the fluctuations of the probability distribution dramatically. For $d=2$, we numerically observe (for small $r_{2}$) that the relative fluctuation $\delta n_{1}/\overline{n_{1}}$ in the dominant ground state population vanishes as $N\to\infty$ faster than $\sim N^{-1/2}$ (expected from independent decay). This causes the steady state distribution to become sharply peaked at density $n_{1}/N=1$, as indicated by the cumulative probability distribution in Fig. 2(c) and in the SM [24]. By the Bhatia-Davis inequality [28], the relative fluctuation can be bounded as $\delta n_{1}/\overline{n_{1}}\leq\sqrt{(1-\text{DTR})/\text{DTR}}$ which vanishes like $\sim N^{(r_{2}-1)/2}$. While this bound is not tight, it justifies the single-particle approximation made in our theoretical analysis.

Finally, the convergence time $T$ to the steady state can be estimated within the fluid model under the single-particle approximation. By noting that $\text{d}t=\text{d}x_{1}/v_{1}$, we find (see SM [24])

$$T=\frac{1}{\Gamma_{1}}\int_{\mathcal{C}}\frac{\text{d}x_{1}}{(N-x_{1}-\sum_{\mu>1}x_{\mu})(x_{1}+1)}\approx\frac{2-r_{2}}{\Gamma_{1}}\frac{\ln N}{N},$$ (12)

where $\mathcal{C}$ is the trajectory defined by Eq. (7). The presence of competing collective decay channels affects the well-known superradiant timescale of $\Gamma_{1}T\sim\ln N/N$ [17] only by a constant factor.

(a)

Figure 3: Targeted photoassociation of a molecular dimer. Weakly bound molecules are created by optical excitation and spontaneously decay to various rovibrational ground-electronic states. The dominant transition ratio (DTR) is plotted against the initial number of weakly bound $\text{Sr}_{2}$ molecules in the $(1)1_{u}(\nu^{\prime}=-1,J^{\prime}=1)$ state. Decay is considered to four states, $(\nu=-1,J=0,2)$ and $(\nu=-2,J=0,2)$ with branching ratios 0.54, 0.27, 0.13, and 0.06, respectively. Points show the numerical solution from the Monte-Carlo simulation of the rate equation (2) for the non-competing scenario where only the dominant transition is collectively enhanced ($\square$) and the competing scenario with all four collective channels ($\circ$). The dashed line denotes the analytical prediction (11) for the non-competing scenario with $\Gamma_{\text{leak}}/\Gamma_{1}=0.852$. The solid line denotes the prediction (10) for the competing scenario. Predictions are valid for $N\gg 1$.

As a possible application of collective transition quenching, we analyze the creation of ultracold diatomic molecules via photoassociation. In one-color photoassociation, laser-cooled atoms form weakly bound molecules that spontaneously decay into various more tightly-bound states [29]. The vibrational and rotational branching ratios are dictated by Franck-Condon factors and angular-momentum Hönl-London factors, respectively. Specifically, we examine ultracold strontium dimers that lack Feshbach resonances and must be produced optically [30, 31]. Strontium dimers have narrow optical transitions and a structureless ground state, making them well-suited for metrology [32, 33, 34, 35]. We consider photoassociation via the state $(1)1_{u}(\nu^{\prime}=-1,J^{\prime}=1)$ [33], where $\nu^{\prime}$ and $J^{\prime}$ are the excited state vibrational and rotational quantum numbers, respectively, and negative vibrational numbers count down from the dissociation threshold. This state has a dominant decay path to the $\nu=-1,J=0$ ground state, with an overall branching ratio of $\sim 54\%$ (molecular parameters are detailed in the SM [24]). Our protocol uses a short, strong pulse to create $N$ weakly bound molecules from a dense sample of ultracold atoms. Alternatively, effective decay directly from unbound atoms into ground-state molecules can be engineered via a continuous off-resonant drive [36].

The avalanche-like behavior dramatically enhances the fraction of molecules in the dominant ground state, as shown in Fig. 3. If only one transition is collective, the dominant transition ratio rapidly approaches unity, as shown in Ref. [15]. This can be realized by engineering the dielectric environment to be frequency-selective (either via a cavity [37] or a photonic crystal [38]). We note that collective transition quenching occurs in the weak coupling limit of cavity QED, distinguishing it from the so-called “polaritonic chemistry” in strong coupling [39, 40]. Figure 3 demonstrates that a nearly-pure sample of molecules is created even under maximal competition, where all decay channels are collective. Collective enhancement of all transitions can occur in broadband cavities (or for molecules with hyperfine structure, which yields smaller frequency differences between levels) or in single-mode fibers [41] within the “mirror configuration” [23]. This effect can also occur in free space with dense, subwavelength molecular clouds. However, in this scenario, collisions can lead to significant losses. For polar molecules, collisions can be prevented via electric field [42] or microwave shielding [43, 44, 45, 46]. Alternatively, they can be suppressed by trapping molecules in optical tweezers to form ordered arrays [47, 48]. As these systems are extended, their dynamics are not constrained to the permutationally symmetric subspace. Nevertheless, we expect collective transition quenching to occur for molecules placed close to waveguides (at arbitrary distances) or arranged in ordered two- and three-dimensional arrays of subwavelength lattice constant [49, 50], albeit with reduced scaling.

In summary, collective decay holds promise for a wide range of quantum systems, from closing open transitions in atoms [12, 14] and preventing parasitic decay in solid-state emitters [51, 52, 5] to directing emission into dielectric nanostructures instead of unwanted modes [53, 20]. While we have primarily focused on an initially inverted ensemble, our treatment captures collective transition quenching from a large range of multiply excited states. Future research directions include studying the potential for many-body-enhanced metrology [54], due to the sensitivity of the ground state populations to the decay ratios. Other interesting avenues include the observation of symmetry breaking, which can occur if multiple dominant transitions have identical decay rates, akin to mirror symmetry breaking predicted in waveguide QED [20].

Acknowledgements.— We thank Kon Leung, Luis A. Orozco, and Maximilian Schemmer for helpful discussions. We acknowledge support by the National Science Foundation through the CAREER Award (No. 2047380) and the QLCI program (grant No. OMA-2016245), the Air Force Office of Scientific Research through their Young Investigator Prize (grant No. 21RT0751) and through grant No. FA9550-1910328, ARO through the MURI program (Grant No. W911NF-20-1-0136), DARPA (Grant No. W911NF2010090), as well as by the David and Lucile Packard Foundation and the Brown Science Foundation. The Institute for Quantum Information and Matter is an NSF Physics Frontiers Center.