Significance

Density-functional theory (DFT) is a well-established method to study many-electron systems. Over the past decades advanced algorithms have been designed that allow even large systems to be solved computationally very efficiently. Only recently, DFT has been generalized to correctly incorporate the quantum nature of light, which becomes important, e.g., in cavity quantum electrodynamics (QED). In general, the accuracy in density-functional theory depends crucially on the capability to construct good approximations that reflect the features of the exact effective potential for the Kohn–Sham system. In this work, we introduce a fixed-point scheme to construct these exact effective potentials for a cavity QED system and demonstrate their features for the ground-state and time-dependent situations.

Abstract

The density-functional approach to quantum electrodynamics extends traditional density-functional theory and opens the possibility to describe electron–photon interactions in terms of effective Kohn–Sham potentials. In this work, we numerically construct the exact electron–photon Kohn–Sham potentials for a prototype system that consists of a trapped electron coupled to a quantized electromagnetic mode in an optical high-Q cavity. Although the effective current that acts on the photons is known explicitly, the exact effective potential that describes the forces exerted by the photons on the electrons is obtained from a fixed-point inversion scheme. This procedure allows us to uncover important beyond-mean-field features of the effective potential that mark the breakdown of classical light–matter interactions. We observe peak and step structures in the effective potentials, which can be attributed solely to the quantum nature of light; i.e., they are real-space signatures of the photons. Our findings show how the ubiquitous dipole interaction with a classical electromagnetic field has to be modified in real space to take the quantum nature of the electromagnetic field fully into account.

In the last decades, the quantum nature of light has inspired many experimental and theoretical developments in physics. In particular, the fields of cavity (1) and circuit (2) quantum electrodynamics (QED) have recently seen exceptional progress. For instance, slow photons in vacuum (3) and two-ion superradiant states (4) have been observed, and only recently the chemical landscape of a molecule has been modified using strong coupling to photons (5, 6), which can be also termed QED chemistry.

However, traditional ab initio approaches developed to investigate large quantum systems (see, e.g., refs. 7⇓⇓⇓⇓–12) are not fully applicable in situations where the quantum nature of light becomes important. These many-body methods ignore typically the quantum-mechanical coupling to photons and usually take only the classical Coulomb interaction into account. Recently, an approach that treats particles (electrons, ions) and the photons on equal footing and closes the gap between traditional many-body and quantum-optical methods has been proposed (13⇓–15). This so-called quantum-electrodynamical density-functional theory (QEDFT) allows the representation of the coupled particle–photon system by two uncoupled, yet nonlinear auxiliary quantum systems. The resulting multicomponent Kohn–Sham systems are subject to effective potentials that take into account the particle–particle (Coulomb) interaction and the particle–photon interaction. If we use approximations to this new type of Kohn–Sham potentials, the resulting equations become numerically feasible and ab initio calculations of large quantum systems (16⇓–18) coupled to photons are possible. Although a wealth of approximations to the particle–particle interaction part of the effective potential are known (see, e.g., refs. 8, 9, 19) at the moment there is only one approximation for the particle–photon part of the effective potential beyond the classical mean-field approximation available (20). Indeed, besides its existence and uniqueness, so far nothing is known about the exact real-space properties of this particle–photon effective potential and how it models the interaction between charged particles and quantized photon fields.

In this article, we present the exact space- and time-resolved Kohn–Sham potential for a coupled multicomponent electron–photon system in an optical cavity. The prototype system that we consider in the present work is a 2D quantum ring containing one electron that is coupled to a single photon mode. We construct the exact effective potential by a fixed-point procedure (21⇓⇓–24) and study ground-state properties as well as the time evolution of the electron–photon system. In the first case we identify electronic states in the weak- and strong-coupling regime, which cannot be generated by any classical light field in dipole approximation. In the time-dependent framework, we analyze the quantum effects of the photon field by putting the field initially into (i) a coherent state and (ii) a linear superposition of two Fock number states and show when the quantum nature of the photon field induces the dominant contribution to the Kohn–Sham electron–photon exchange–correlation (xc) potential. Here we find that the electron–photon interaction is responsible for steps and peaks in the exact Kohn–Sham potential. Similar steps and peaks have been found in purely electronic time-dependent density functional theory (TDDFT) in e.g., charge-transfer processes (25, 26), but these exact features have so far only been observed for the time-dependent case in one-dimensional models (25).

The static and dynamical behavior of the coupled electron–photon systems that we consider in the present work is given by the following Hamiltonian (14, 15),H^(t)=∑i(−ℏ22m∇→i2+vext(rit))+e24πϵ0∑ij,i>j1|ri−rj|+∑α12[p^α2+ωα2(q^α−λαωα⋅R)2]+jextα(t)ωαq^α.[1]The first part of the Hamiltonian describes the electronic subsystem and contains the nonrelativistic kinetic energy, the external-potential energy due to a classical external potential vext(rt) and the Coulomb-interaction energy. The second part of the Hamiltonian accounts for the presence of the photons, where the electron–photon interaction is described in dipole approximation, i.e., the dipole-moment operator R=∑iri couples linearly to the photon displacement coordinate q^α=−ℏ2ωα(a^α†+a^α), which is proportional to the quantized displacement field component of the αth mode, i.e., D^α=−ϵ0ωαλαq^α/e. Although we restrict ourselves here to dipole coupling, all findings shown below also apply for beyond-dipole situations.* The electron–photon coupling is given by λα=λαeα, where eα is the polarization vector of the photon field. The value of the coupling can be changed from the weak- to the ultra-strong-coupling limit in circuit QED experiments (27, 28). (For later reference, we measure λα in units of [meV/nm]). In addition, the Hamiltonian contains a quadratic electron self-interaction ves(R)=∑α(λα⋅R)2/2 and the photons interact furthermore with a classical external current jextα(t). We emphasize that both external potentials, i.e., vext(rt) and jextα(t), can be used to control the quantum system.

In QEDFT the electron–photon system is exactly described by two reduced quantities that couple to the external control fields (13⇓–15). In the case of the Hamiltonian [1] these reduced quantities are the usual electronic density n(rt)=〈n^(r)〉, where n^(r)=∑iδ(r−ri), and the expectation value of the photon coordinates qα(t)=〈q^α〉. In principle, we only need to calculate these expectation values and can then determine (for a fixed initial state) all properties of the electron–photon system. To calculate n(rt) and qα(t) one only needs to solve the corresponding coupled equations of motions for these two basic variables in the system, i.e., the Ehrenfest equations (14, 15).∂t2qα(t)=ωαλα⋅R(t)−ωα2qα(t)−jextα(t)ωα,[2]∂t2n(rt)=−∇→⋅Q(rt)+1m∇→⋅(n(rt)∇→vext(rt))+1m∑αλα⋅∇→〈Ψ|n^(r)(λα⋅R−ωαq^α)|Ψ〉,[3]with the local-force density of the electrons given by Qk(rt)=∑l∂lTkl(rt)+Wk(rt), where the first term describes the momentum-stress forces and the second term is responsible for the forces due to the particle–particle interactions. To solve these implicit equations, we would need to find explicit expressions in terms of n(rt) and qα(t) for the different force densities. We note that, in the equations of motion for the photon coordinates, all terms are explicitly known, and hence the unknown expressions that take care of the proper description of the electron–photon interactions are contained solely in the electronic equation. To make approximations for these unknown quantities easier, one can adopt a Kohn–Sham scheme, such that approximations in terms of the force densities of the uncoupled and noninteracting system become possible. This approach has been applied highly successfully to electronic-structure calculations (see, e.g., refs. 8, 9, 19). In a Kohn–Sham approach to the electron–photon system the missing forces are accounted for by the effective potential vs(rt) that naturally splits into two parts vs(rt)=vext(rt)+vMxc(rt), where vext is the external potential of the original problem, and vMxc (mean-field exchange correlation) denotes the effective potential due to the interaction with the photons. In general, vs(rt) contains both, the contributions due to the Coulombic electron–electron repulsion in the original many-body problem, and in addition the contributions from the electron–photon interaction, here in dipole approximation. However, because our aim is to investigate the effective potentials due to the coupling between photons and electrons, we restrict ourselves here to a single electron in an semiconductor GaAs quantum ring, which is placed in a cavity and is assumed to couple to a specific cavity photon mode, as depicted schematically in Fig. 1.

The figure schematically illustrates a 2D optical cavity containing one atom, with a single electron. The coupling of the electron to the cavity mode at resonance frequency ωα and with electron–photon coupling strength λα modifies the dynamics of the electron density n(rt), which moves in the external potential vext(rt).

This restriction to a single electron has the advantage that we only have contributions of the electron–photon interaction in the effective Kohn–Sham potential vs(rt) and can exclusively study only their behavior.† In addition in single-electron problems, we can incorporate the electron self-interaction ves(r)=∑α(λα⋅r)2/2 exactly, thus the exact Kohn–Sham potential is given by vs(rt)=vext(r)+ves(r)+vMxc(rt). We treat the electron in this semiconductor medium by using an effective mass m=0.067me (29) and we use vext(r)=12mω02r2+V0e−r2/d2, where r2=x2+y2, to describe the quantum ring. This potential effectively confines the electron in a harmonic trap, which contains a Gaussian peak in the center. For the system at hand, we choose the experimental parameter values (29) ℏω0=10meV, V0=200meV, d=10nm, and the effective dielectric constant κ=12.7ϵ0. The 2D electronic system has a nondegenerate ground state and a twofold degeneracy in the excited states (29). We choose the photon frequency ωα in resonance with the transition between the ground state and the first excited state in the electronic system. Thus, ℏωα=1.41meV. For simplicity, we restrict ourselves to one of the two independent polarization directions of the field mode and use eα=(1,1) without loss of generality. In experiments weak- to ultra-strong coupling (27, 28), which corresponds here to λα≥0.134meV/nm, has been realized. To be comparable to such experiments, in this paper we choose three values for the electron–photon interacting strength λα=1.68⋅10−3meV/nm, 3.36⋅10−3meV/nm (weak-coupling limit), and 0.134meV/nm (ultra-strong coupling limit).

The Kohn–Sham scheme then decouples the two subsystems, which leads us to two evolution equations of the formiℏ∂tϕ(rt)=−ℏ22m∇→2ϕ(rt)+vs(rt)ϕ(rt),[4]iℏ∂t|α,t〉=12[p^α2+ωα2q^α2]|α,t〉+jsα(t)ωαq^α|α,t〉,[5]where the Kohn–Sham photon wavefunction is given by |α,t〉=∑ncn|n,t〉 and |n,t〉 are the Fock number states of cavity mode α. The Kohn–Sham construction furthermore requires that the initial state of the Kohn–Sham system has to have the same density n(r,0) and time-derivative n˙(r,0) as the coupled system. The same is required for the basic variable in the photon system, i.e., qα(0) and q˙α(0).

To determine the in general unknown effective potential vs(rt) in terms of n(rt) and qα(t), we use a fixed-point method originally developed for purely electronic TDDFT (21⇓⇓–24). Although naively one could expect that a fixed-point iteration is also needed to determine the effective current jsα(t), from Eq. 2 this current is known explicitly, i.e., jsα(t)=jextα(t)+ωα2λα⋅R(t). Hence we only need to determine vs(rt), for which we use the fixed-point formula−1m∇→⋅(n(rt)∇→vk+1(rt))=∂t2[n([vk],rt)−n(rt)]−1m∇→⋅(n([vk],rt)∇→vk(rt)).[6]To find the fixed points of this equation, we solve first the Schrödinger equation for vk (using zero-boundary conditions) and from the exact many-body solution we determine the corresponding n[vk]. Next, we use a multigrid solver to invert the Sturm–Liouville problem in Eq. 6, which yields vk+1. This procedure is repeated until convergence to the fixed point has been reached. To speed up convergence, we use a direct inversion in the iterative subspace approach (30). We have tested the validity of this approach in the time-independent situation by comparing to the well-known analytic inversion formula for one-electron and two-electron singlet problems (31);vs(0)(r)=ℏ22m∇→2n0(r)n0(r)+E0.[7]Using exact diagonalization (32), we are able to calculate the exact ground state of the correlated electron–photon system. We use a 127 × 127 2D real-space grid for the electron and 40 photon number states. This amounts to a dimensionality of the full problem of 127 × 127 × 40 = 645,160 basis functions. From this we then determine with the above iteration procedure the exact vs. Because we treat the interaction with the photon field in dipole approximation, the mean-field field contribution to the effective potential is vM(rt)=−ωαqα(t)λα⋅r.‡ If the impact of the quantum nature of the cavity light field would be negligible, the classical field would be the only contribution. This is expected to hold in the case of large photon numbers or coherent photon states or equivalently when the electronic system is driven by a strong classical external laser pulse. In Fig. S1 and Movie S1 we provide a supplementary example of a correlated electron–photon propagation with nonfactorizable initial state driven by such an external laser pulse. In this correlated example the electronic dynamics is mainly driven by the external laser pulse and the classical mean-field approximation performs reasonably well. Thus, if we use vMxc=vM+vxc, the nondipole corrections (to all orders in r) due to vxc are a direct measure of the nonclassical light–matter interaction.

Correlated electron–photon propagation with nonfactorizable initial state and driven by external laser pulse: Initially the system populated the combined electron–photon ground state. A shows the applied external laser pulse, B shows the exact (black) and classical (red) dipole moment. C shows the projections of the time-dependent wave function on the first three electronic eigenstates. D shows the evolution of the Mandel parameter Q(t) and purity γ(t) for the exact (black) and classical approximation (red). E, G, and I show the vMxc(r,t) potential for different times; F, H, and J show the electron density n(r,t) at the corresponding times. The applied laser field v(r,t) = E(t)eα·r was added to Hamiltonian in Eq. 1 with E(t) = E0exp(-(t-t0)2/σ2)sin(ωt) and polarization eα=(1,1). We choose t0 = 0.29 ps, σ= 0.058 ps, ħω = 1.41 meV, and E0 = 0.23 meV/nm.

First, let us investigate the effective potential in the case of the ground state of the multicomponent system. In Fig. 2A we show a ground-state density in the weak-coupling case (λα=1.68⋅10−3 meV1/2/nm). Compared with the cavity-free case, we see a slight prolongation of the density along the x=y axis. If we increase the coupling to strong coupling (λα=0.134 meV1/2/nm), this feature becomes more dominant and the charge density of the electron becomes even separated (Fig. 2B). If we consider the Mxc potential we observe a peak in the middle of the cavity, which models the forces that the photons exert on the electron to elongate its charge distribution (see Fig. 2C for weak coupling and Fig. 2D for strong coupling). In Fig. 2 E and F, we show diagonal cuts through vMxc(r) for the weak (blue) as well as the strong coupling (red) regime. We see how vMxc pushes the density farther apart the stronger the coupling becomes. Such a splitting cannot be generated using a static classical field in dipole coupling. Hence the nondipole contributions of the vMxc mark the nonclassical interaction with the photons, and are a necessary feature to model the exact forces exerted by the photons on the charged particles. To further substantiate our findings we have determined the photon-number expectation value 〈a^α†a^α〉 as well as the purity γ=Tr(ρph2) of the ground state, where ρph is the reduced photon density matrix. Here, a purity value that deviates from 1 indicates that the state is not factorizable into photon and electron wave functions. Therefore, the purity is a measure for electron–photon entanglement. For weak coupling we have 1.18⋅10−3 photons and γ=0.999764 and for strong coupling we find 3.19 photons in the ground state with γ=0.4919. This clearly indicates that the ground state is a hybrid state of the photons and the electron with stronger entanglement in the strong-coupling regime. A further parameter that we consider in our analysis is the Mandel Q parameter (33):Q=〈a^α†a^α†a^αa^α〉−〈a^α†a^α〉2〈a^α†a^α〉,[8]which measures the deviation of the photon statistics from a Poisson distribution and thus is a measure for the quantum nature of the photonic subsystem. If the field is in a quasi-classical state, i.e., in a coherent state, then Q=0. For weak coupling we find Q=3.88⋅10−4, and for strong coupling we have Q=0.4567. This further supports that this model has a highly nonclassical ground state of the coupled matter–photon system.

A shows the ground-state density for a weak-coupling case with λα=1.68⋅10−3 meV1/2/nm and in B we illustrate a strong-coupling case with λα=0.134 meV1/2/nm. The corresponding ground-state Mxc potential for λα=1.68⋅10−3 meV1/2/nm is displayed in C and for λα=0.134 meV1/2/nm displayed in D. In E and F, cuts (blue λα=1.68⋅10−3 meV1/2/nm and red λα=0.134 meV1/2/nm) through vMxc along the diagonal (C) /antidiagonal (D) are shown. The white arrow in A indicates the polarization direction of the field mode.

Next, we turn our attention to the time-dependent situation. As initial states of the combined matter–photon system (as well as for the corresponding Kohn–Sham system) we consider two different cases. In both cases we choose factorizable initial states, which consist of the electronic ground state of the unperturbed quantum ring and the photon field in (i) a coherent state with 〈a^α†a^α〉=4 and in (ii) a superposition of the vacuum state and the one-photon state with 〈a^α†a^α〉=0.5. In both examples, we choose the electron–photon coupling strength λα=3.36⋅10−3 meV1/2/nm. To numerically propagate the system, we use a Lanczos scheme and propagate the initial state in (i) 160,000 [(ii) 360,000] time steps with Δt=0.146 fs. We start with the analysis of example 1. Here, we compare the exact dynamics to the dynamics induced by the (self-consistent) classical mean-field approximation. In Fig. 3A we display the time evolution of the exact experimentally accessible dipole moment (black) and contrast it to the classical approximation (red). Because the field is initially in a coherent state and resembles a classical field, the evolution of the classical approximation is for a short time similar to the exact one. Between t=0 and t=6 ps, the classical approximation is the dominant part for most of the time in vMxc(rt) (Fig. 3C). Nevertheless, also in this time interval we find large beyond dipole corrections to vMxc(rt), which appear at the turning points of the dipole evolution and vanish afterward. Precisely at these times, we observe several peaks and steps in the 2D surface plot of the effective potential (Fig. 3 F and I) that describe the nonclassical forces due to the interaction with the photon mode. After t=6 ps, the beyond-dipole correction becomes the dominant part in vMxc(rt) and we find a dominating peak in the 2D surface plot (Fig. 3L). The peak structure of vMxc becomes clearly visible in Fig. 3 G, J, and M, where we plot the diagonal of vMxc for the different time steps. To further analyze this time-dependent system, we computed the (now time dependent) Mandel Q(t) parameter and the purity γ(t).

Coherent state as initial state for the photon mode: In A we display the dipole of the exact (black) and mean-field (red) time evolution. In B we contrast the exact Mandel parameter Q(t) (see main text for definition) (solid black) and γ(t) (dotted black) with the corresponding mean-field values (red). C, F, I, and L show the corresponding Mxc potentials at different times (t=0,3.67,4.53,and7.29ps); D, G, J, and M show the corresponding diagonal cuts of the Mxc potentials; and in E, H, K, and N we present the corresponding densities. The Inset in B shows Q(t) between t = 3 and 5 ps. The negative Q(t) in the exact solution indicates nonclassical behavior in the photon mode. Movie S2 shows the full time evolution from 0 to 23 ps.

In Fig. 3B we contrast the exact results (black) to those found from the mean-field calculation (red). The purity (dotted black) as well as the Mandel Q(t) parameter (solid black) are in agreement with our previous observations, namely that around t=6 ps, where these parameters start to deviate more strongly from the mean-field values, the classical description breaks down. We point out, that in our (decoupled) Kohn–Sham system these parameters are by construction constant and equivalent to the mean-field values, and the Kohn–Sham photon field only changes the number of photons in the coherent state. Hence, the values of these parameters become nontrivial functionals of the initial state as well as n and qα. In particular the assessment of the purity allows us to conclude that the peaks in vMxc are associated with how close to a factorizable (electronic) state the many-body system is. For small times, the system remains close to a factorizable state (purity value close to 1) and we find peaks and steps only at the turning point of the dipole moment, while later in time memory effects become dominant and cause permanent peaks and steps. Finally we note that although we have termed all beyond-dipole contributions to vMxc(rt) as nonclassical (because they come solely from the quantum nature of light), the nonclassicality of the photon field alone is often associated with a negative Q(t). In Fig. 3B we have thus provided an inset to highlight that (up to t=5 ps; Fig. 3B, Inset) such sub-Poissonian statistics, which cannot be described by any probability distribution in phase space, are also present in our prototype system.

Next, we analyze in detail our example (ii) in Fig. 4. We again compare the exact dynamics to the time evolution induced by the (self-consistent) classical mean-field approximation. For the photon mode we choose in this example as initial state a superposition of the lowest two Fock number states. Although the field is in this case initially in a state that does not resemble a classical situation as in our first example, the evolution of the classical approximation is for small times similar to the exact one. Between t=0 ps and t=2 ps, the classical approximation is the dominant part for most of the time in vMxc(rt) (Fig. 4C). Nevertheless, also in this time interval we find large beyond-dipole corrections to vMxc(rt). After the first turning points of the dipole evolution, a dominant peak in vMxc appears (Fig. 4F). Later, we observe again several peaks and steps appearing in the 2D surface plot of the effective potential (Fig. 4 F and I) that describe the nonclassical forces due to the interaction with the photon mode. After t=2 ps, the beyond-dipole correction becomes the dominant part in vMxc(rt) and we find a dominating peak in the 2D surface plot. The peak structure of vMxc becomes clearly visible in Fig. 4 G, J, and M, where we plot the diagonal of vMxc for the different time steps. The purity (dotted black) as well as the Mandel Q(t) parameter (solid black) are in agreement with our previous observations in the coherent state example, that around t=2 ps, where these parameters start to deviate more strongly from the mean-field values, the classical description breaks down. For small times, the system remains close to a factorizable state, while later memory effects become dominant and cause permanent peaks and steps. Both examples shown in Figs. 3 and 4 illustrate how a Kohn–Sham approach can exactly describe the different regimes of quantized photon fields that interact with matter.

(Color online) Superposition of Fock number states as initial state for the photon mode: In A we display the dipole of the exact (black) and mean-field (red) time evolution. In B we contrast the exact Q(t) (solid black) and γ(t) (dotted black) with the corresponding mean-field values (red). In C, F, I, and L, we show the corresponding Mxc potentials at different times (t=0,1.56,13.92,and37.83ps); in D, G, J, and M, we show the corresponding diagonal cuts of the Mxc potentials; and in E, H, K, and N, we present the corresponding densities. Movie S3 shows the full time evolution from 0 to 50 ps.

In conclusion, we have presented the real-space signatures of the exact effective potentials for a Kohn–Sham approach to cavity QED. We have identified step and peak structures that are reminiscent of the steps and peaks in the exact Kohn–Sham potential of traditional TDDFT, but arise here solely due to the coupling to quantized photon fields. These effective potentials account for the forces that the photons and the electrons exert on each other if we use an uncoupled Kohn–Sham system to describe the coupled matter–photon system. Provided we have a good approximation to the effective potential (20), which includes the peak and step structures observed here that are caused by the nonclassical light–matter interaction, the Kohn–Sham approach can be used to perform ab initio calculations of large quantum systems interacting with photons in a high-Q cavity. In this case we have a valuable computational tool for QED chemistry (5, 6), which would open up a new field of research for the electronic structure community. Besides developing approximations to the Mxc potential, a further important line of research is the extension of the current work to cavities with loss (15), which is the standard situation in most cavity-QED experiments. Work along these lines is currently in progress in our group.

↵*We note that in the dipole approximation the wavelength of the quantized electromagnetic field is assumed to be much larger than the electronic system and hence only the lowest order in a spatial Taylor expansion of the field is kept. In this approximation the spin-dependent magnetization density of the electronic system, which couples to the curl of the photon field (15), has no contribution to the electron–photon interaction. Thus the matter–photon coupling does not depend on spin but still this approximation includes the dominant photon quantum effects, e.g., spontaneous emission. However, if the electron–photon interaction is taken beyond the dipole approximation then the coupling to the photons depends explicitly on the spin of the electrons.

↵†In the case of many electron systems the effective field changes due to nonzero interaction forces Wk(rt) in Qk(rt) of Eq. 3. Additionally the momentum-stress forces ∑l∂lTkl(rt) and the electron self-interaction λα〈Ψ|n^(r)(λα⋅R)|Ψ〉 are changed. We point out that both, the interaction contribution Wk(rt) as well as the explicit electron–photon contribution in Eq. 3 come from the coupling to the photon field in Coulomb gauge. The electron–electron interaction can also be interpreted in terms of the exchange of (virtual) photons and both contributions merge for different gauges (e.g. Lorentz gauge).

↵‡We note that in many-electron problems additional contributions from the electron self-interaction ves(r) also appear in the electron–photon mean-field potential vM(rt) (14).

Bacteria could help tackle the growing mountains of e-waste that plague the planet. Although researchers are a long way from optimizing the approach, some are already confident enough to pursue commercial ventures.

Holographic acoustic tweezers, in which ultrasonic waves produced by arrays of sound emitters are used to individually manipulate up to 25 millimeter-sized particles in three dimensions, could be used to create 3D displays consisting of levitating physical voxels.