Light-Matter Response in Nonrelativistic Quantum Electrodynamics.

We derive the full linear-response theory for nonrelativistic quantum electrodynamics in the long wavelength limit and provide a practical framework to solve the resulting equations by using quantum-electrodynamical density-functional theory. We highlight how the coupling between quantized light and matter changes the usual response functions and introduces cross-correlated light-matter response functions. These cross-correlation responses lead to measurable changes in Maxwell's equations due to the quantum-matter-mediated photon-photon interactions. Key features of treating the combined matter-photon response are that natural lifetimes of excitations become directly accessible from first-principles, changes in the electronic structure due to strong light-matter coupling are treated fully nonperturbatively, and self-consistent solutions of the back-reaction of matter onto the photon vacuum and vice versa are accounted for. By introducing a straightforward extension of the random-phase approximation for the coupled matter-photon problem, we calculate the ab initio spectra for a real molecular system that is coupled to the quantized electromagnetic field. Our approach can be solved numerically very efficiently. The presented framework leads to a shift in paradigm by highlighting how electronically excited states arise as a modification of the photon field and that experimentally observed effects are always due to a complex interplay between light and matter. At the same time the findings provide a route to analyze as well as propose experiments at the interface between quantum chemistry, nanoplasmonics and quantum optics.

Recent years have seen tremendous experimental advances in the nascent field of strongly coupled light–matter systems.[1,2] In particular, new experimental advances have been demonstrated in polaritonic chemistry,[3−5] solid-state physics,[6] biological systems,[7] nanoplasmonics,[8,9] two-dimensional materials,[10,11] or optical waveguides,[12] among others.

In this so-called strong-coupling regime, as a result of mixing matter and photon degrees-of-freedom,[13,14] novel effects emerge such as changes in chemical pathways[15−17] ground-state electroluminescence,[18] cavity-controlled chemistry for molecular ensembles,[19,20] or optomechanical coupling in optical cavities,[21] new topological phases of matter,[22] super-radiance,[23] or superconductivity.[24]

Due to the inherent complexity of such coupled fermion-boson problems described in general by quantum electrodynamics (QED), the theoretical treatment is usually drastically simplified. One common approximation is to restrict the description of the system to simplified effective models that heavily rely on input parameters. Current state of the art in the theoretical description of strong light-matter coupling very often employs a few-level approximation. This approximation leading to the Rabi or Jaynes-Cummings model[25,26] in the single-emitter case, or the Dicke model[27] in the many-emitter case, is however often not sufficient,[28,29] in particular, when observables besides the energy are of interest,[29] such as in experimental setups involving the modification of chemical reactivity.[1]

Alternatively, in linear spectroscopy, the current theoretical description is built on the semiclassical approximation.[30] Herein, the many-particle electronic system is treated quantum mechanically and the electromagnetic field appears as an external perturbation. As an external perturbation, the electromagnetic field probes the quantum system, but is not a dynamical variable of the complete system (see also Supporting Information, S1). Since in the strong-coupling regime light and matter must be on the same level, a semiclassical approximation is not adequate, and the feedback between light and matter has to be considered.

It is, however, long known that the radiative lifetimes are finite. Furthermore, experimentally excited-state properties are usually inferred from (de)excitations of the photon field, which is in stark contrast to the usual semiclassical theoretical description based solely on the electronic subsystem.

In free-space, this mismatch can be circumvented since excited-state properties such as radiative lifetimes of atoms and molecules can be calculated perturbatively using the theory of Wigner-Weisskopf,[31] employing the Markov approximation. However, this perturbative treatment of the coupling of light and matter becomes insufficient in the case that strong light–matter coupling is achieved, for example, due to many emitters or due to reducing the mode volume of a cavity. In such cases, the Markov approximation breaks down and the Wigner-Weisskopf theory is not applicable anymore.[32] Additionally, it is not straightforward how to extend the original formulation of Wigner-Weisskopf to many electronic levels and, hence, to an ab initio treatment of electronic systems.

As a consequence, the current literature shows a large gap for situations, where light and matter is strongly coupled and observables such as excited-state densities, radiative lifetimes, or electron-photon correlated observables of interest. A good example is the control of the radiative lifetimes of single molecules[33,34] by changing the environment. In such cases, the properties of the many-body system are changed, for example, the excitation energies and lifetimes are strongly modified. This happens because certain modes of the photon vacuum field are enhanced which can lead to a strong coupling of light with matter. Alternatively, increasing the number of particles leads to an enhancement of the coupling due to the self-consistent back-reaction of matter onto the photon field and vice versa. It is important to realize that such changes are nonperturbative for the photon field as well as for the matter subsystem and hence need a self-consistent implementation. This fact is most pronounced in the appearance of polaritonic states and their influence on chemical and physical properties of matter.[1,13]

In this paper, we close this gap by presenting a practical and general framework that subsumes electronic-structure theory, nanoplasmonics, and quantum optics. We present a description that challenges our conception of light and matter as distinct entities[35] and that expresses the excited states as modifications of the photon field. We do so by introducing a linear-response formalism for coupled matter-photon systems. This formalism leads naturally to the ability to calculate radiative lifetimes in arbitrary photon environments, including free-space, high-Q optical cavity or nanoplasmonic structures. We make this approach practical by introducing a linear-response framework for quantum-electrodynamical density-functional theory (QEDFT).[13,14,36−38] This development is specifically timely since QEDFT has now been successfully applied to real systems in equilibrium,[39] which demonstrates the feasibility of ab initio strong-coupling calculations, yet an accurate and efficient approach to excited states within QEDFT has been missing. This work therefore furthermore closes a gap within the QEDFT framework. We further want to note that, there have been different studies in literature that are devoted to including the classical feedback of the light field to the matter systems all in dipole approximation, such as for specific systems[40,41] or reduced dimensionality.[42] The presented work not only generalizes these approaches, but also provides a clear path to how to include the quantum effects of the light field for this feedback.

Light–Matter Interaction in the Long Wavelength Limit

Our fundamental description of how the charged constituents of atoms, molecules, and solid-state systems, that is, electrons and positively charged nuclei, interact is based on QED;[13,43−45] thus, the interaction is mediated via the exchange of photons. Adopting the Coulomb gauge for the photon field allows us to single out the longitudinal interaction among the particles, which gives rise to the well-known Coulomb interaction and leaves the photon field purely transversal. Assuming then that the kinetic energies of the nuclei and electrons are relatively small allows us to take the nonrelativistic limit for the matter subsystem of the coupled photon–matter Hamiltonian, which gives rise to the so-called Pauli-Fierz Hamiltonian[13,37,45] of nonrelativistic QED. In a next step, one then usually assumes that the combined matter–photon system is in its ground state such that the transversal charge currents are small and that the coupling to the (transversal) photon field is very weak. Besides the Coulomb interaction, it is then only the physical mass of the charged constituents (bare plus electromagnetic mass[45]) that is a reminder of the photon field in the usual many-body Schrödinger Hamiltonian. In this work, however, we will not disregard the transversal photon field, which makes the presented framework much more versatile and applicable to situations of quantum mechanics and quantum optics at the same time (see also Appendix section, Photonic Observables and Radiative Lifetimes).

Spectroscopy from Quantum Description of Light–Matter Interaction

In the following, we consider cases in which the semiclassical approximation breaks down, as outlined in the introduction. In principle, QEDFT can be formulated for each level of theory of QED as presented in ref.[37] As a consequence, the formalism outlined in this paper can be straightforwardly extended to more general formulations, including full minimal coupling beyond the dipole approximation[46] (in dipole approximation, only energy can be transferred between charged particles and the light field, but not momentum; thus, the dipole approximation is insufficient to describe processes such as the radiation reaction). In this manuscript, to illustrate the concepts, we restrict the discussion in the following to the dipole approximation and the length-gauge.

To this end, from the Pauli-Fierz Hamiltonian, we make the long-wavelength or dipole approximation in the length-gauge[47] since the wavelength of the photon modes are usually much larger than the extent of the electronic subsystem, as well as the Born–Oppenheimer approximations for the nuclei (the inclusion of the nuclei is straightforward;[48] however, the presented formulation is perfectly suited to provide the photon-dressed modified potential-energy surfaces for the nuclei and, hence, access to modifications of chemical reactions in, e.g., optical cavities[15,49]), which leads (in SI units) to[36,37,50]where Ĥe(t) is the standard many-body electronic Hamiltonian.[51] We further restrict ourselves to arbitrarily many but a finite number M of modes α ≡ (k, s), with s being the two transversal polarization directions that are perpendicular to the direction of propagation k. The frequency ωα and polarization ϵα that enter in λα = ϵαλα, with and mode function S(r)[37] define these electromagnetic modes. S(r) is normalized, has the unit 1/ with the volume V, and we choose a reference point r0, where we have placed the matter subsystem to determine the fundamental coupling strength (all results presented in this paper are independent of r0). These photon modes couple via the displacement coordinate , where q̂α is given in terms of photon annihilation âα and creation âα† operators, to the total dipole moment R = ∑er (throughout this paper, we use the implicit definition e = −|e|). The q̂α appears in the contribution of mode α to the displacement field D̂α = ϵ0ωαλαq̂α.[47] Further, the conjugate momentum of the displacement coordinate is given by . Besides a time-dependent external potential v(r, t), we also have an external perturbation jα(t) that acts directly on the mode α of the photon subsystem. Here, jα(t) is connected to a classical external charge current J(r, t) that acts as a source for the inhomogeneous Maxwell’s equation. Formally, however, due to the length-gauge transformations, the jα(t) corresponds to the time-derivative of this (mode-resolved) classical external charge current[36,37] (see also Appendix section, Self-Consistency of the Maxwell’s Equation). Physically the static part jα,0 merely polarizes the vacuum of the photon field and leads to a static electric field.[38,52] The time-dependent part δjα(t) then generates real photons in the mode α. This term is also known as a source term in quantum field theory,[43] where it generates the particles (here the photons) that are studied. From this perspective, it becomes obvious that instead of using δjα(t) one could equivalently slightly change the initial state of the fully coupled system by adding incoming photons that then scatter off the coupled light–matter ground state.[45]

Linear Response in the Length Gauge

With the Hamiltonian of eq in length gauge we can then in principle solve the corresponding time-dependent Schrödinger equation (TDSE) for a given initial state of the coupled matter–photon system Ψ0(r1σ1, ..., rσ, q1, ..., q)where σ corresponds to the spin degrees of freedom. However, instead of trying to solve for the infeasible time-dependent many-body wave function, we restrict ourselves to weak perturbations δv(r,t) and δjα(t) and assume that our system is in the ground state of the coupled matter–photon system initial time (in principle, also other initial states, e.g., an uncorrelated matter–photon state could be chosen). In this case, a first-order time-dependent perturbation theory can be used to approximate the dynamics of the coupled matter–photon system (for details, see Supporting Information, S2). This framework gives us access to linear spectroscopy, for example, the absorption spectrum of a molecule. Traditionally, if we made a decoupling of light and matter, that is, we assumed Ψ0(r1σ1, ..., rσ, q1, ..., q) ≃ ψ0(r1σ1, ..., rσ) ⊗ φ0(q1, ..., q), we would only consider the matter subsystem ψ (the photonic part φ would be completely disregarded). Physically, we would investigate the classical dipole field that the electrons induced due to a classical external perturbation δv(r, t). To determine this induced dipole field we would only consider the linear response of the density operator n̂(r) = ∑δ(r – r), which would be given by the usual density–density response function in terms of the electronic wave function ψ0 only. In the following, we suppress the spin component of the wave function and focus exclusively on the spatial and mode dependence, i.e., Ψ(r1, ..., r, q1, ..., q; t).

In this work, however, since we do not assume the decoupling of light and matter, the full density–density response is taken with respect to the combined ground-state wave function Ψ0 and is consequently different to the traditional density–density response. Further, since we can also perturb the photon field in the cavity by δjα(t), which will subsequently induce density fluctuations, the density response δn gets a further contribution leading toHere the response function χ(rt, r′t′) corresponds to the density–density response but with respect to the coupled light–matter ground state and χ(rt, t′) corresponds to the density response induced by changing the photon field. In the standard linear-response formulation, due to the decoupling ansatz, changes in the transversal photon field would not induce any changes in the electronic subsystem. Since obviously we now have a cross-talk between light and matter, we accordingly have also a genuine linear-response of the quantized light fieldwhere χ(t, r′t′) is the full response of the photons due to perturbing the electronic degrees, and χ(t, t′) is the photon–photon response function. For an alternative definition of δqα(t), we also refer the reader to eq S14 in the Supporting Information. The response function χ(t, rt′) is, in general, not trivially connected to χ(rt, t′), due to the different time-ordering of t and t′.

The entire linear-response in nonrelativistic QED for the density and photon coordinate can also be written in matrix form.[53] In this form, we clearly see that the density response of the coupled matter–photon system depends on whether we use a classical field δv(r, t), photons, which are created by δjα(t), or combinations thereof for the perturbation. Furthermore, we can also decide to not consider the classical response of the coupled matter–photon system due to δn(r, t), but rather directly monitor the quantized modes of the photon field δqα(t). This response yet again depends on whether we choose to use a classical field δv(r, t) that induces photons in mode α or whether we directly generate those photons by an external current δjα(t), and we also see that the different modes are coupled, that is, that photons interact. Similarly, as charged particles interact via coupling to photons, also photons interact via coupling to the charged particles. Keeping the coupling to the photon field explicitly therefore, on the one hand, changes the standard spectroscopic observables and, on the other hand, also allows for many more spectroscopic observables than in the standard matter-only theory.

Maxwell–Kohn–Sham Self-Consistent Linear-Response Theory

The problem with this general framework in practice is that already in the simplified matter-only theory, we usually cannot determine the exact response functions of a many-body system. The reason is that the many-body wave functions, which we use to define the response functions, are difficult, if not impossible, to determine beyond simple model systems. So, in practice, we need a different approach that avoids the many-body wave functions. Several approaches exist that employ reduced quantities instead of wave functions.[54−56] The workhorse of these many-body methods is DFT and its time-dependent formulation TDDFT.[57−59] Both theories have been extended to the general coupled matter–photon systems within the framework of QED.[13,36−38,48]

QEDFT allows us to solve instead of the TDSE equivalently a nonlinear fluid equation for the charge density n(r, t) coupled nonlinearly to the mode-resolved inhomogeneous Maxwell’s equation.[36−38,60] While these equations are, in principle, easy to handle numerically, we do not know the forms of all the different terms explicitly in terms of the basic variables of QEDFT, that is, (n(r, t), qα(t)). To find accurate approximations, one then employs the Kohn–Sham (KS) scheme, where we model the unknown terms by a numerically easy to handle auxiliary system in terms of wave functions. The simplest approach is to use noninteracting fermions and bosons that lead to a similar set of equations, which are however uncoupled. Enforcing that both give the same density and displacement field dynamics gives rise to mean-field exchange-correlation (Mxc) potentials and currents.[52,61,62] Formally, this Mxc potential and current is defined as the difference of the potential/current that generate a prescribed internal pair in the auxiliary noninteracting and uncoupled system (vs([n], r, t), jαs([qα], t)); (the subscript “s” is usually not explained in the density-functional literature, but we can assume that it refers to “single particle”, as the potential often appears in effective single-particle equations[62]) and the potential/current that generates the same pair in the physical system defined by eq , which we denote by (v([n, qα], r, t), jα([n, qα], t)), that is,In the time-dependent case, we only have a mean-field contribution to the Mxc current,[36,38] where the total dipole moment is written as R(t) = ∫drern(r, t). Further, we have ignored the so-called initial-state dependence because we assume (for notational simplicity and without loss of generality) in the following that we always start from a ground state[62,63] of the matter–photon coupled system. In this way, we can recast the coupled Maxwell-quantum-fluid equations in terms of coupled nonlinear Maxwell-KS equations for auxiliary electronic orbitals, which sum to the total density ∑|φ(r, t)|2 = n(r, t), and the displacement fields qα(t), that is,Here we use the self-consistent KS potential vKS([v, n, qα], r, t) = v(r, t) + vMxc([n, qα], r, t) that needs to depend on the fixed physical potential v(r, t),[62] and instead of the full bosonic KS equation for the modes α, we just provide the Heisenberg equation for the displacement field. Although the auxiliary bosonic wave functions might be useful for further approximations, it is only qα(t) that is physically relevant, and thus, we get away with merely coupled classical harmonic oscillators, that is, the mode resolved inhomogeneous Maxwell’s equation. To highlight the extra self-consistency due to coupling between light and matter, we contrast the traditional electron-only KS theory with the Maxwell KS theory in Figure . It is then useful to divide the Mxc potential into the usual Hartree-exchange-correlation (Hxc) potential that we know from electronic TDDFT and a correction term that we call photon-exchange-correlation potential (pxc), that is,Clearly, the correction term vpxc will vanish if we take the coupling |λα| to zero and recover the purely electronic case. Since by construction the Maxwell KS system reproduces the exact dynamics, we also recover the exact linear-response of the interacting coupled system (see also Supporting Information, S3). We can express this with the help of the Mxc kernels defined by the functional derivatives of the Mxc quantitiesand use the corresponding definitions for the Hxc kernel (that only for the variation with respect to n has a nonzero contribution) and the pxc kernels. We note that, using eq , we explicitly findand gM(t, t′) vanishes, since jα, in eq has no functional dependency on qα. Via these kernels, we find with χ(rt, r′t′) and χ(t, t′), where χ(t, t′) ≡ 0 for α ≠ α′, the uncoupled and noninteracting response functions thatand accordingly for the mixed matter–photon response functionsHere we employed the formal connection between response functions and functional derivatives χ(rt, r′t′) = δn(r, t)/δv(r′, t′), as well as χ(t, t′) = δqα(t)/δjα′(t′) and accordingly for the auxiliary system. The Mxc kernels correct the unphysical responses of the auxiliary system to match the linear response of the interacting and coupled problem. So, in practice, instead of the full wave function, what we need are approximations to the unknown Mxc kernels. Later we will provide such approximations, show how accurate they perform for a model system and then apply them to real systems. If we decouple light and matter, that is, Ψ0 ≃ ψ0 ⊗ φ0, and disregard the photon part φ0 (as is usually done in many-body physics), we recover the response function of eq with fMxc ≡ 0, and fMxc → fHxc. The response function, which is calculated with the bare matter initial state ψ0, then obeys the usual Dyson-type equation relating the noninteracting and interacting response in TDDFT[64,65] with vMxc([n, qα], r, t) → vHxc([n], r, t).

Figure 1

Schematics of the Maxwell KS approach contrasted with schematics of the usual semiclassical KS theory. While in the semiclassical approach the KS orbitals are used as fixed input into the mode-resolved inhomogeneous Maxwell’s equation in vacuum through the total dipole R(t) = ∫drer∑|φ(r, t)|2 (see also Appendix section, Self-Consistency of the Maxwell’s Equation), in the Maxwell KS framework the induced field acts back on the orbitals, which leads to an extra self-consistency cycle.

Excited States as Properties of the Photon Field

Following the above discussion, the usual response functions will change and response functions are introduced if we keep the matter–photon coupling explicitly. This leads to many exciting consequences. First, we get the completely self-consistent response of the system including all screening, retardation (we note that retardation processes which require the exchange of more than one photon are independent of a dipolar approximation) and other effects that become important when either the matter subsystem is becoming large[66−69] or when strong-coupling situations are considered. Since light and matter influence each other nonperturbatively the usual simplified approximations that only treat one part of the system accurately become unreliable[28,29] (see also the discussion in section ). Second, due to the matter-mediated photon–photon interactions (see Appendix section, Self-Consistency of the Maxwell’s Equation, and Figure ), the Maxwell’s equations become self-consistent. A very interesting consequence is that, in contrast to a purely classical theory, we can theoretically distinguish whether a system is perturbed by a free current (that, in turn, would generate a classical electromagnetic field) or by a free electromagnetic field, for example, a classical laser pulse. Third, the inclusion of the photon modes introduces the missing photon bath that leads to finite lifetimes (see Appendix section, Photonic Observables and Radiative Lifetimes, and section ). In connection to this it becomes important that we suddenly have access to a wealth of observables that describe the photon field. Most importantly, this implies the possibility to completely change our perspective of excited states of atoms and molecules. Indeed, in line with the experimental situation where changes in the photon field give us information on the excited states, we can view excited-state properties as arising from quantum modifications of the Maxwell’s equations in matterThe response of the density is then found with help of the response functions eqs –13. In the usual case of an external classical field δv(r, t) and δjα(t) = 0, we then find the induced field by (suppressing detailed dependencies with ∫ dr → ∫ and ∫dr ∑α → )Here, the first term on the right-hand side corresponds to the noninteracting matter–response. However, due to the electron–electron interaction, we need to take into account also the self-polarization of interacting matter (second term). Finally, the third term describes the matter-mediated photon–photon response. The excited states of the coupled light–matter system are in this description changes in the photon field. That this perspective is actually quite natural becomes apparent if one considers the nature of the emerging resonances for a real system (see Figure ). These resonances are mainly photonic in nature, as they describe the emission/absorption of photons (see Appendix section, Photonic Observables and Radiative Lifetimes). Let us consider now in more detail what the terms on the right-hand side of the modified Maxwell’s equations mean physically. First of all, in a matter-only theory the self-consistent solution of the Maxwell’s equations together with the response of the bare matter-system would correspond approximately to the first two terms on the right-hand side (see Appendix section, Self-Consistency of the Maxwell’s Equation). The photon–photon interaction would not be captured in such an approximate approach. Second, to highlight the physical content of the different terms, we can make the mean-field contributions due toexplicitThe second term on the right-hand side then corresponds to the random-phase approximation (RPA) to the instantaneous matter–matter polarization. Here, a term that corresponds to the dipole self-energy induced by the coupling to the photons arises. The third term on the right-hand side is the RPA approximation to the dipole–dipole mediated photon interaction. To give these terms further physical meaning, note that in the usual perturbative derivation of the van der Waals interaction[44] the first two terms would cancel and leave the photonic dipole–dipole interaction that gives rise to the R–6 for small distances and the R–7 for larger distances. The rest are exchange-correlation (xc) contributions that arise due to more complicated interactions among the electrons and photons. The last term effectively describe photon–photon interactions mediated by matter. In addition, we want to highlight that xc contributions are directly responsible for multiphoton effects, such as two-photon or three-photon processes (see Figure ). If we only keep the mean-field contributions of the coupled problem, we will denote the resulting approximation in the following as photon RPA (pRPA) to distinguish it from the bare RPA of only the Coulomb interaction. We see how the Maxwell’s equations in matter become self-consistent due to bound charges, that is, fields due to the polarization of matter. A new term, the photon–photon interaction, appears. For free charges, that is, due to an external charge current δjα(t), we see similar changes. Clearly, if we had no coupling to matter, then there would be no induced density change and we just find the vacuum Maxwell’s equations coupled to an external current for the electric field. In other terms, the displacement field trivially corresponds to the electric field (see Appendix section, Self-Consistency of the Maxwell’s Equation).

Figure 2

Schematics that contrasts the usual Maxwell’s equation (left) with the fully self-consistent Maxwell’s equation (right). Top: The induced transversal electric field E⊥ as a consequence of the induced polarization P⊥, which can be equivalently expressed in terms of the auxiliary displacement field D⊥. Left: mode-resolved nonself-consistent Maxwell’s equation with no backreaction. The external charge current jα induces the external electric field in Eαtot = Eα + Eαext, which acts as an external perturbation through the dipole. Since the constituents of χ̃ expressed in TDDFT are purely electronic, the induced field does not couple back to the Maxwell field. Right: self-consistent Maxwell’s equation in which jα induces the internal field qα(t) through the electron-photon correlated dipole which has an explicit dependence as seen in the QEDFT form of χ. The self-consistency of the induced field through the dipole introduces nonlinearities in the coupled system and, thus, changes the Maxwell field at the level of linear-response.

Figure 7

First-principles lifetime calculation of the electronic excitation spectrum of the benzene molecule in an quasi one-dimensional cavity: (a) Full spectrum of the benzene molecule, (b) zoom to the Π–Π* transition, where the black arrow indicates the full width at half-maximum (fwhm) ΔE, (c) zoom to a peak contributing to the σ–σ+ transition. The gray spectrum is obtained by Wigner-Weisskopf theory.[31] The dotted spectral data points correspond to many coupled electron–photon excitation energies that together comprise the natural line shape of the excitation. Blue color refers to a more photonic nature of the excitations vs red color to a more electronic nature.

Figure 4

Linear-response spectra for the extended Rabi model (dotted-red) compared to the pRPA (dashed-blue) and RWA (full-orange) approximations and for different coupling strengths λ. (a) Absorption spectra due to matter–matter response, (b) spectra due to photon–photon response, (c) spectra due to matter–photon or photon–matter response. (d) The case for λ = 0.7 shows all excitations that arise in strong coupling. (a–d) Resonant coupling. In (e), the field is halfway detuned from atomic resonance, that is, ω0 = 2 and ωc = 1 with strength and energies shifted to frequencies favoring 2-photon processes. The insets in (d) and (e) zoom into the frequency axis showing a many-photon process.

Illustrative Examples for the Coupled Matter–Photon Response

In this section, we discuss the perspective enabled by the linear response formalism of QEDFT in more detail for a simple and illustrative model system. We discuss a slight generalization of the Rabi model,[70,71] which is the standard model of quantum optics. The Rabi model describes a single electron on two lattice sites/energy levels interacting with a single photon mode. We schematically depict the system in Figure and present all further details of this system in Appendix section, Examples for the Coupled Matter–Photon Response: Details on the Rabi Model.

Figure 3

Two-level system (with excitation ω0) coupled to one mode of the radiation field (with frequency ωc). The matter subsystem is driven by an external classical field v(t) and the photon mode is driven by an external classical current j(t) and both subsystems are coupled with a coupling strength λ.

First, let us analyze the optical spectra for such a system and scrutinize the different approximations to the Mxc kernels. We will compare the numerical exact results, with the mean-field (pRPA) and the rotating-wave approximation (RWA). In Figure a–c we see how the optical spectra of the resonantly coupled system (i.e., δ = ω0 – ωc = 0) change for an increasing electron–photon coupling strength λ. Already for small coupling, the splitting of the electronic state into an upper and lower polariton becomes apparent. Approximately these states are given in terms of the RWA as |+, 0⟩ and |−, 0⟩. The difference in energy between the lower and upper polariton is called the Rabi splitting ΩR and is used to indicate the strength of the matter–photon coupling. In molecular experiments values of up to ΩR/ωc ≃ 0.25 have been measured.[72,73] Up to λ = 0.1 the different spectra for the exact (dotted-red), the pRPA (dashed-blue) as well as the RWA (full-orange) are in close agreement before they start to differ. Already the mean-field treatment is enough to recover the quantized matter-photon responses, even for the coupled matter-photon spectra in Figure c. Consequently, the pRPA seems a reasonable approximation for linear-response spectra even for relatively strong coupling situations. Only upon increasing the coupling strength further and thus going into the ultrastrong coupling regime, the discrepancies becomes large. For ultrastrong coupling (for λ = 0.3 the Rabi splitting is already of the order of 0.5ωc), the approximations do not recover the exact results. Increasing further leads then to not only a disagreement in transition frequencies, but also the weights of the transitions become increasingly different.

Besides a simple check for the approximations to the Mxc kernels, the extended Rabi model also allows us to get some understanding of the response functions χσ, χσ, and χ, where σ is the expectation-value of the corresponding Pauli matrix and describes the density/occupation changes between the two sites/energy levels. This means we consider mixed spectroscopic observables, where we perturb one subsystem and then consider the response in the other. We analogously employ χσ(ω) and χσ(ω), respectively, to determine a “mixed polarizability” (see Supporting Information, S5). If we plot this mixed spectrum (see Figure c displayed in dotted-red for the numerically exact case), we find that we have positive and negative peaks. Indeed, this highlights that excitations due to external perturbations can be exchanged between subsystems, that is, energy absorbed in the electronic subsystem can excite the photonic subsystem and vice versa. The oscillator strength of the photonic spectrum (based on χ) in Figure b provides us with a measure of how strong the displacement field (and with this also the electric field) reacts to an external classical charge current with frequency ω. Similarly, the mixed spectrum (based on χσ or χσ) in Figure c provides us with information on how strong one subsystem of the coupled system reacts upon perturbing the other one. The oscillator strength here is not necessarily positive. What is absorbed by one subsystem can be transferred to the other.

In Figure d,e, we show specifically the absorption spectra of the Rabi model for ultrastrong coupling, that is, λ = 0.7. In this regime, three new peaks arise for the exact case accounting for high-lying excited states with nonvanishing dipole moments due to the strong electron-photon coupling. The new absorption peaks in Figure d, also shown in the inset, describes the resonant coupling case which the RWA and pRPA fail to capture in strong coupling, since processes beyond one-photon are involved. Similarly, Figure e depicts the case where the field is half-detuned from the electronic resonance indicating a two-photon process. Clearly in ultrastrong coupling the absorption peaks are merely shifted close to the bare frequencies of the individual subsystems, but remain dressed by the photon field as new peaks arise due to the coupling. The pRPA and RWA capture the first of the two peaks around ω = 2, which is also the frequency of the atom, but fail to capture higher lying nonvanishing contributions to the spectra. These higher-lying peaks correspond to multiphoton processes. With more accurate approximation for the xc potential results closer to the exact ones can be obtained. We note at this point that the peaks in Figure are artificially broadened and in reality correspond to sharp transitions due to excited states with infinite lifetimes. How to get lifetimes quantitatively will be discussed in the next section.

Coupled Matter–Photon Response: Real Systems

In this section, we apply the introduced formalism in pRPA approximation to real systems. We make the linear-response formulation practical by reformulating the problem as an eigenvalue equation in the frequency-domain. For electron-only problems this formulation is known as the Casida equation.[65] We refer the reader to Appendix section, Linear-Response Theory as a Pseudoeigenvalue Problem, for a derivation of our extension of the Casida equation, which includes transverse photon fields. For the following discussion, we consider benzene molecules in an optical cavity but the presented approach is not restricted to any specific system.

In Figure , we schematically depict the experimental setup for a photoabsorption experiment under strong light–matter coupling for a single molecule. First we study the prototypical cavity QED setup where a molecule is strongly coupled to a single cavity mode of a high-Q cavity. In the second setup, we lift the restriction of only one mode and instead couple the benzene molecule to many modes that sample the electromagnetic vacuum field without enhancing the coupling to a specific mode by hand. In the third setup, we study the behavior of two molecules in an optical cavity, as well as a dissipative situation, where only a few modes are strongly coupled, embedded in a quasi-continuum of modes. In the last example, we analyze the strong coupling of a single molecule to a continuum of modes. We find a transition from Lorentzian line shape to a Fano line shape[74] for increasing electron-photon coupling strength. These different setups provide us with an ab initio calculation for the spectrum of a real molecule in a high-Q cavity, an ab initio determination of intrinsic lifetimes and an ab initio calculation of the nonperturbative interplay between electronic structure, lifetime, and strong-coupling. The two last situations need a self-consistent treatment of photons and matter alike and cannot be captured by any available electronic-structure or quantum-optical method. All of those examples highlight the rich possibilities and perspectives that the QEDFT framework provides.

Figure 5

Schematic of absorption spectroscopy in optical cavities: Benzene (C6H6) molecule and λα denotes the polarization direction of the photon field.

Strong Light–Matter Coupling

The first results we discuss are a set of calculations, where a benzene molecule is strongly coupled to a single photon mode in an optical high-Q cavity (our approach could also describe strong light–matter coupling for other systems, e.g., nanoplasmonic systems[8] and generalizations to quantum interactions in laser pulses could be done along the lines of ref (75)). We have implemented the linear-response pseudoeigenvalue equation of eq into the real-space code OCTOPUS[76,77] and details of the numerical parameters are given in Appendix section, Numerical Details. The routines used to perform all calculations in this work will be made publicly available. They can be easily transported to any other first-principles code that has the matter linear-response equations implemented to make them ready to describe the complete QED response, i.e., joint matter–photon response, as described in this work.

In the first calculation, we include a single cavity mode in resonance to the Π–Π* transition of the benzene molecule,[76,78] that is, ωα = 6.88 eV. For the light–matter coupling strength λα = |λα|, we choose five different values, that is, λα = (0, 2.77, 5.55, 8.32, 11.09) eV1/2/nm that correspond to a transition from the weak to the strong-coupling limit and the cavity mode is assumed to be polarized along the x-direction.

Since in this manuscript, we focus on electron–photon coupling, we do not consider the coupling to the nuclei. Generalizations are straightforward, for example, along the line of ref (48). In experiment, in particular for molecular systems, the majority of the line-broadening is due to vibrational coupling, see, for example, refs (79 and 80), for the optical spectra of benzene. Strong light–matter coupling for such systems will lead to the splitting of the peak into the lower and upper polariton and both peaks will inherit the vibrational line broadening of the electronic excitation outside the cavity, as has been shown in various experiments, for example, refs (81 and 82).

In Figure , we show the absorption spectra for these different values of λα. We start by discussing the λα = 0 case that is shown in black. This spectrum corresponds to a calculation of the benzene molecule in free space, and the spectrum is within the numerical capabilities identical to ref (76). The spectrum in ref (76) has been obtained using an explicit time-propagation with finite time. In the limit of zero broadening and including all unoccupied states, we would find identical spectra with very long propagated spectra. We stress that here the broadening of the peaks is only done artificially since the photon bath is not included in the calculation. In the examples of sections and 3.4 we include many modes and, hence, sample the photon bath nonperturbatively. We tune the electron–photon coupling strength λα in Figure . We find for increasing coupling strength a Rabi splitting of the Π–Π* peak into two polaritonic branches. The lower polaritonic branch has higher intensity, compared to the upper polaritonic peak. Numerical values for the excitation energy EI, the transition dipole moment xI, and the oscillator strength fI are given in Table in the Appendix. This demonstrates that ab initio theory is able to describe excited-state properties of strong light–matter coupling situations and captures the hybrid character of the combined matter–photon states. Thus, predictive theoretical first-principle calculations for excited-states properties of real systems strongly coupled to the quantized electromagnetic field are now available. This will allow unprecedented insights into coupled light-matter systems, since we have access to many observables that are not (or not well[29]) captured by quantum-optical models.

Figure 6

Absorption spectra for the benzene molecule in free space (black) and under strong light–matter coupling in an optical cavity to ultrastrong coupling (blue). The value for λα is given in units of [eV1/2/nm].

Table 1

Rabi Splitting of the Π–Π* Transition: Electron–Photon Interaction Strength λα = |λα|, Excitation Energy E, Transition Dipole Moment x, and the Oscillator Strength f

λα (eV1/2/nm)	EI (eV)	⟨xI⟩ (A)	fI (a.u.)
0	6.88	0.952	0.546
2.77	6.69	0.721	0.304
2.77	7.03	0.626	0.241
5.55	6.49	0.791	0.355
5.55	7.18	0.550	0.190
8.32	6.28	0.848	0.395
8.32	7.30	0.482	0.149
11.09	6.06	0.896	0.426
11.09	7.41	0.420	0.114

Lifetimes of Electronic Excitations from First-Principles

Next, we consider how to obtain lifetimes from QEDFT linear-response theory. In this example, we explicitly couple the benzene molecule to a wide range of photon modes similar as in the spontaneous emission calculation of ref (83) While in ref (83) the system was simulated with 200 photon modes, we choose here now 80 000 photon modes. The energies of the sampled photon modes cover densely a range from 0.19 meV, for the smallest energy up to 30.51 eV for the largest one with a spacing of Δω = 0.38 meV. However, we do not sample the full three-dimensional mode space together with the two polarization possibilities per mode but rather consider a one-dimensional slice in mode space. This one-dimensional sampling of mode frequencies will change the actual three-dimensional lifetimes, but for demonstrating the possibilities of obtaining lifetimes, this is sufficient (a detailed analysis of real lifetimes would include, besides a proper sampling of the mode space, considerations with respect to the bare mass of the particles). The sampling of the photon modes corresponds to the modes of a quasi-one-dimensional cavity. We choose a cavity of length L(83) in the x-direction, with a finite width in the other two directions that are much more confined. Thus, we employ ωα = αcπ/L and , where x0 = L/2 is the position of the molecule in the x-direction. While we have a sine mode function in the x-direction, we assume a constant mode function in the other directions. For this example, we choose a cavity of length L = 3250 μm in the x-direction, L = 10.58 Å in the y-direction, and L = 2.65 Å in the z-direction.

The results of this calculation are shown in Figure . In Figure a, we show the full spectrum. The electron–photon absorption function that has been obtained by coupling the benzene molecule to the quasi one-dimensional cavity with 80 000 cavity modes is plotted in blue. Since we have sampled the photon part densely, we do not need to artificially broaden the peaks anymore. Formulated differently, we can directly plot the oscillator strength and the excitation energies of our resulting eigenvalue equation and do not need to employ the Lorentzian broadening anymore. In Figure , from blue (more photonic) to red (more electronic) for the electron–photon absorption spectrum we plot the different contributions of each pole in the response function. These results confirm our intuition that resonances are mainly photonic in nature and that a Maxwell’s perspective of excited states is quite natural. In (b) we zoom to the Π–Π* transition. Due to quasi one-dimensional nature of the quantization volume, we find a broadening of the peak that is larger than it is for the case of a three-dimensional cavity due to the sampling of the electromagnetic vacuum. This is similar to changing the vacuum of the electromagnetic field. Accordingly the lifetimes of the electronic states are shorter if the electromagnetic field is confined to one dimension and we will discuss this in the next section.

Connection to the Standard Wigner-Weisskopf Theory

If the coupling between light and matter is very weak and neither subsystem gets appreciably modified due to the other, in contrast to the previous strong light–matter coupling case, the radiative lifetimes of atoms and molecules can be calculated using the perturbative Wigner-Weisskopf theory[31] in single excitation approximation, as well as under the assumption of the Markov approximation. These approximations are justified in the usual free-space case, where the results of Wigner and Weisskopf reproduce the prior results of Einstein based on the ad-hoc A and B coefficients. However, it does not include the treatment of ensembles of molecules that effectively enhance the matter–photon coupling strength, as shown below. Under the assumption of the Wigner-Weisskopf theory, the radiative decay rate is given byFor a one-dimensional cavity in x-dimension the results change to[32]For comparison, we show in Figure the peaks in gray that are predicted by the Wigner-Weisskopf theory. Since our sampling is very dense, we find for both peaks shown in the bottom a good agreement with eq .

In fact, if we take the continuum limit for the photon modes, we recover in our framework the lifetimes predicted by the Wigner-Weisskopf theory, including the diverging energy shifts,[84] that is, the Lamb shift. Due to the Lamb shift, our resulting peaks are slightly shifted, due to the divergencies. These divergencies can be handled by renormalization theory. The lifetimes can now be obtained the following way: We measure the full width at half-maximum (fwhm), indicated by the black arrow in (b). In this case, we find ΔEfwhm = 0.0204 eV and the corresponding lifetime τΠ–Π* follows by τΠ–Π* = ℏ/ΔEfwhm = 32.27 fs. Using the Wigner-Weisskopf formula from eq and the dipole moments and energies from the LDA calculation without a photon field, we find a lifetime of 32.21 fs. As a side remark, the same transition using eq has a free-space lifetime of 0.89 ns, roughly in the range of the 2p-1s lifetime of the hydrogen atom of 1.6 ns.

In Figure c, we finally show the ab initio peak of the σ–σ+ transition. We find a narrow ab initio peak that is not as well sampled as the Π–Π*. We note in passing that we find an ionization energy of 9.30 eV using Δ-SCF in the benzene molecule with the LDA exchange-correlation functional. In our simulation, coupling to peaks higher than the ionization energy are broadened by continuum (box) states.

Beyond the Single Molecule Limit and Dissipation in QEDFT

In contrast to the free-space result, where weak coupling as well as the assumption of a dilute gas of molecules are implied, in the case of single-molecule strong coupling[8] or when nearby molecules or an ensemble of interacting molecules modify the vacuum, the usual perturbative theories break down. Changes in the electronic and the photonic subsystem become self-consistent, and the usual distinction of light and matter becomes less clear. In such situations, the linear-response formulation of QEDFT as well as the Maxwell’s perspective of excited-state properties becomes most powerful. Consider, for instance, two benzene molecules weakly coupled to a one-dimensional continuum of photon modes. If the molecules are far apart, we just find the usual Wigner-Weisskopf result. But if we bring the molecules closer (see Figure a), we see that the combined resonance shifts and the combined line width becomes broader, implying a shortened lifetime. In Figure b, we consider the case of single-molecule strong coupling, where a few out of the 80 000 modes have an enhanced coupling strength. In red, we show the spectrum where the molecule is coupled to the continuum, as is also shown in Figure . We then introduce a single strongly coupled mode at the Π–Π* transition energy, and the resulting spectra is shown in green. We note that, in the figure, the cavity frequencies are plotted in dashed lines. The single mode introduces the expected Rabi splitting into the upper and lower polariton and the peaks of the upper and lower polariton become broadened due to the interaction with the continuum. Interestingly, we find a different line broadening for the lower and the upper polaritonic peak, since only the sum of both has to be conserved. The smaller broadening for these two lower polaritonic states implies that the radiative lifetimes of the lower and upper polaritonic states are longer than the lifetime of the excitation in weakly coupled free-space. In blue, we show the spectra, where we have introduced three strongly coupled modes in addition to the cavity 80 000 modes of the continuum. We tune the two additional cavity modes in resonance to the lower and upper polariton peaks of the green plot. We find additional peak splitting, but also a shifting of peak positions, at 7.8 eV.

Figure 8

(a) Two molecules of benzene strongly coupled to 80 000 cavity modes of an one-dimensional cavity. The further apart the molecules are, the closer the peak gets to the single molecule peak. Also, we notice the doubled peak broadening (shorter lifetime). The gray spectrum is obtained by the Wigner-Weisskopf theory.[31] (b) We show the Rabi splitting in a situation of a single strongly coupled mode with 80 000 cavity modes (green) and three strongly coupled modes with 80 000 cavity modes (blue). The red lines correspond to the same setup as in (a). The dashed lines refer to the frequency of the cavity modes. The peaks become broadened due to the interaction with the continuum.

In the last numerical example, we study the strong coupling to the continuum for the case of a single molecule. The results are shown in Figure . Here, we effectively enhance the light–matter coupling strength by reducing the volume of the cavity along the y- and z-direction. For comparison, we show in red the setup that is also shown in Figure , where the excitations have a Lorentzian line shape consistent with the Wigner-Weisskopf theory, as discussed in the previous section. By gradually reducing the dimensions along the y- and z-direction, we find drastic changes in the line shape of the excitations. These changes lead to the transition of the line shape from a Lorentzian to a Fano line shape, as becomes clearly visible for LL = 0.28 Å.

Figure 9

Ab intio lifetime calculation of the electronic excitation spectrum of the benzene molecule in a one-dimensional cavity along the x-direction with different lengths in the L and L directions. The red spectra refer to the same setup as in Figure . Effectively the electron–photon strength increases with smaller L and L lengths leading to a transition from a Lorentzian line shape to a Fano line shape.

As a summary, we have presented in this section, that lineshapes, as well as lifetimes can be inferred directly from first principle calculations. In the case of a Lorentzian line shape, we find that the width of the calculated peaks (no need to introduce any artificial broadening as commonly done) correspond to the lifetimes. These calculations demonstrate that the ab initio theory is able to capture the true nature of excitations, that is, resonances with finite intrinsic lifetimes, without the need of an artificial bath or postprocessing. Furthermore, we find that the excitations measured in absorption/emission experiments are mainly photonic in nature, and it is only the peak position that is dominated by the matter constituents. This is, of course, very physical, since what we see is the absorption/emission of a photon, not of the matter constituents. Further, since we describe the photon vacuum on the same theoretical footing as the matter subsystem, we have full control over the photon field, making it straightforward to simulate very intricate changes, for example, changing the character of a specific mode out of basically arbitrarily many, and investigating its influence on excited-state properties such as the radiative lifetime. This allows predictive first-principle calculations for intricate experimental situations similar to the ones encountered in refs (33 and 34).

Summary and Outlook

In this work we have introduced a linear-response theory for nonrelativistic quantum-electrodynamics in the long wavelength limit, which can be straightforwardly extended to the full minimal coupling case. Compared to the conventional matter-only response approaches, we have highlighted how in the coupled matter-photon case the usual response functions change, how photon–photon and matter–photon response functions are introduced, how these response functions provide a photonic perspective on excited state properties, how the results lead to self-consistent Maxwell’s equation in matter, and how we can efficiently calculate all these response functions in the framework of QEDFT. By investigating a simple model system, we have shown how the spectrum of the matter subsystem is changed upon coupling to the photon field. Further, we have demonstrated the range of validity of a simple yet reliable approximation to the, in general, unknown mean-field exchange-correlation kernels. Using this approximation, we have presented the first ab initio calculations of the spectrum of real systems (benzene molecules) coupled to the modes of the quantized electromagnetic field. In one example we have calculated the change upon strong coupling to a single mode of a high-Q cavity, which leads to a large Rabi splitting. In the second example we have calculated from first-principles the natural line widths of benzene coupled to a specific sampling of the vacuum field. In the last examples, we demonstrated the abilities to calculate many-molecule systems, as well as dissipative strong-coupling situations, as well as strong coupling to the continuum, where we find a transition from Lorentzian line shape to Fano line shape, where the usual (perturbative) approaches to light-matter coupling fail. These results demonstrate the versatility and possibilities of QEDFT, where light and matter are treated on equal quantized footing. In the context of strong light–matter coupling, for example, in polaritonic chemistry, the presented linear-response formulation allows now to determine polaritonically modified spectra from first principles. Together with ab initio ground-state calculations,[39] QEDFT now provides a workable first-principle description to analyze and predict photon-dressed chemistry and material sciences. In particular, our approach provides a unique practical computational scheme to compute photon-dressed excited-state potential-energy surfaces and nonadiabatic coupling elements[49] that are required for ab initio calculations in the emerging field of polaritonic chemistry. Further, in the context of standard ab initio theory, the linear-response formulation of QEDFT now allows the calculation of intrinsic lifetimes and provides access to quantum-optical observables. Specifically, due to the nonperturbative nature of the approach, quantum-optical problems where the self-consistent feedback between light and matter has to be taken into account, for example, that many molecules change the photon vacuum and hence the Markov approximation breaks down, becoming feasible. For optical physics, the presented linear-response framework presents an interesting opportunity to study the self-consistency of the Maxwell’s equations in matter from first principles. Finally, we want to highlight that although the QEDFT linear-response framework includes the coupling of light and matter, its similarity to the usual matter-only linear-response formulation in terms of a pseudoeigenvalue problem makes it very easy to include in already existing first-principle codes. This, together with the above-discussed possibilities in different fields of physics, shows that there are many interesting cases that can be studied with the presented method.