Virial Relations for Electrons Coupled to Quantum Field Modes

In this work, we present a set of virial relations for many electron systems coupled to both classical and quantum fields, described by the Pauli–Fierz Hamiltonian in dipole approximation and using length gauge. Currently, there is growing interest in solutions of this Hamiltonian because of its relevance for describing molecular systems strongly coupled to photonic modes in cavities and in the possible modification of chemical properties of such systems compared to the ones in free space. The relevance of such virial relations is demonstrated by showing a connection to mass renormalization and by providing an exact way to obtain total energies from potentials in the framework of quantum electrodynamical density functional theory.


INTRODUCTION
Recent experimental progress allows to manipulate intrinsic properties of molecules by placing them inside photonic structures providing strong and ultrastrong light-matter coupling with a confined electromagnetic vacuum in the optical and infrared regimes. 1−7 Under such conditions, molecular polaritons are formed: hybrid energy states composed of entangled matter-photon degrees of freedom. As a result of the formation of polaritons, the matter system exhibits very different properties than in free space, opening numerous possibilities for manipulating chemical processes. 8−12 These new hybrid polaritonic entities manifest themselves in numerous ways, including modifications of chemical reaction rates, 8 the possibility of long-range energy transfer between matter systems, 6,13 and the suppression of photobleaching, 4

to give just a few examples.
In order to understand and eventually manipulate chemical reactions by strong coupling to confined light modes, the development of theoretical methods that are at the interface between quantum chemistry, solid-state physics, and quantum optics have emerged. 14−42 However, because for large systems these methods always have to rely on some level of approximation, theoretical guidelines that help to assess whether an approximation is valid and thus ultimately predictive are of utmost importance. Our focus will be on such guidelines for methods that treat the matter part of the system such as in electronic-structure theory, 14−35 and not solely as a few level system.
When it comes to electronic-structure problems, the current, most popular approach is density functional theory (DFT). 43 Several exact constraints that functionals of DFT should adhere to are known. 44 Such constraints have served as both, sanity checks for approximate functionals and contributions to the development of new ones. 45 They even led to the improvement of functionals in combination with machinelearning techniques. 46 One of these exact constraints relates to the virial theorem, whose extension to the realm of electronic systems coupled to field modes is the main focus of this work. Recently, a virial-motivated model has been introduced, by which one can obtain accurate singlet−triplet splittings 47 and charge-transfer excitation energies 48 from ground-state DFT calculations. Another recent application of the virial theorem was to provide a way to calculate the "nonadditive" kinetic energy in terms of quantities that can all be obtained through embedding partition-DFT calculations. 49 The virial theorem has been also used to calculate atomic energies in the context of quantum theory of atoms in molecules within the Kohn− Sham DFT formalism. 50 The role of the virial theorem is also important for other classes of quantum-chemistry methods that rely on approximate wave functions, 51 as it is used to assist geometry optimization 52 in electronic-structure codes [e.g., GAMESS (US) 53 ] and as an indicator of the basis-set quality. 54 Importantly, electronic-structure methods are assessed against exact full configuration interaction or other accurate wavefunction method results, which are still numerically feasible for small molecules [e.g., G2/97 and G3/99 theoretical thermochemistry test sets]. These results serve as a consistency check for the reliability of every newly developed method. However, when it comes to many electron systems coupled to photon modes, obtaining such exact reference solutions becomes ever more demanding because of the increased dimensionality of the involved Hilbert spaces. Only very recently, exact calculations became available for HD + and H 2 + molecules as well as for the He-atom. 55 Consequently, theoretical guidelines become even more important.
Our physical starting point is the Pauli−Fierz Hamiltonian, 56 which describes N electrons coupled to quantized field modes. We will restrict ourselves to the dipole approximation, which is valid for a large variety of situations where the spatial extension of the field modes is much larger than the matter system itself, and which is commonly used in cavity QED. Such a model is also useful to capture the effects of the quantized electromagnetic field in free space and thus connects to the old question of how quantum chemistry is affected by the electromagnetic vacuum. Obvious effects are the finite lifetime of excited states or the Lamb shift, a more subtle effect will be highlighted in Section 6 with a brief discussion on mass renormalization. The matter system will be treated within the usual "clamped-nuclei approximation", where the nuclei just contribute to a fixed external scalar potential v(r), as it is the case in most electronic-structure approaches. On the other hand, we do not restrict the number of modes to just one effective mode, as it is commonly done in practice, but allow arbitrarily many modes in order to sample the photon continuum 18,57 and to be able to describe situations where many modes become relevant. 58 Further, our results are not restricted to polaritonic ground states but extend to any eigenstate of the polaritonic system.
Within this setting we present a set of virial relations between the energy expectation values of different parts of the Pauli−Fierz Hamiltonian in an eigenstate. Such relations can be derived by considering the time derivative of expecation values of different time-independent operators, which are zero if an eigenstate of the Hamiltonian is considered. By means of one such operator, consequently called a "virial operator", we obtain a virial relation that is the generalization of the quantum virial theorem for the electronic structure. 59 It will involve two additional terms that stem from the interaction of the electrons with the photons: the dipole coupling energy and the quadratic dipole self-energy. We use the term "virial" for two more identities that we present here, signifying simply that these are the exact relations that eigenstates of the system must satisfy and that involve one or more of the separate energy terms from the Hamiltonian.
Outline of the article: after introducing the physical setup of a system of interacting electrons coupled to photon modes in Section 2, we give a general method of how to derive virial relations from the equations of motion in stationary states in Section 3. Three different virial operators are then used in the following three sections to arrive at different virial relations. First, a basic quantum virial theorem is derived in Section 4, in which the terms from the purely electronic part of the Hamiltonian are combined with the dipole coupling and dipole self-energy. The next virial theorem given in Section 5 is then concerned with photon-field constituents and relates the parts of the photon-field energy to the dipole coupling. Another virial relation that is only given as an estimate is derived in Section 6. It allows to make an interesting connection to mass renormalization because of coupling to the field modes. Finally, in Section 7, the derived virial theorems are employed to map from potential functionals to energy functionals in the Kohn− Sham approach of quantum electrodynamical DFT (QEDFT).
The paper is concluded in Section 8 with a short summary of the main results and a brief outlook.

PHYSICAL SETUP
.., N is coupled to the field modes with fundamental light-matter coupling λ α , which includes the effective coupling strength and polarization direction (for details see (7) in Schafer et al. 60 ), and the mode frequency ω α and mode coordinate q α , α = 1, 2, ..., where an infinite, countable number of modes is permitted. The nature of the field modes is not important here, they could either be due to a cavity environment or a free-field continuum. For the electrons, we will always assume a position basis and thus write simply r i for the respective position operator. The gradient with respect to any coordinate is written ∇ i and the divergence is the inner product with ∇ i , so the momentum operator is −i∇ i and the Laplacian ∇ i 2 . On the side of the bosonic field modes, the coordinates are attached to the operators qα and their conjugates pα, connected through the canonical commutation relation Note that in the literature, q̂α,p̂α sometimes take exactly the opposite roles, but this just amounts to a relabeling. Physically, the mode coordinate corresponds to the displacement field, while its conjugate represents the magnetic field. Atomic units will be employed throughout the article, meaning that the elementary charge, electron mass, ℏ, and 4πε 0 are all set to unity by definition. The system's Hamiltonian is given in length gauge 14,60 (2) and the different parts will be given in the following. We use atomic units throughout the manuscript. The kinetic part in length gauge contains only the canonical momentum −i∇ i and includes no reference to the mode coordinates The potential terms V̂, Ŵare the influence of an external scalar potential v(r) and the usual Coulombic repulsion between the electrons As the last term that depends exclusively on the electron coordinates, we have the quadratic dipole self-energy term with the coupling vectors λ a ∈  d . What follows is the dipolecoupling energy and the purely bosonic field energy in length gauge Journal of Chemical Theory and Computation pubs.acs.org/JCTC Article The last term is the coupling of the mode to an external force where, we have f α as the time-derivative of an external current. 14 Because some of the operators above are squares of self-adjoint operators, or sums of such squares, their expectation value is always positive. This positivity holds for ⟨T⟩, ⟨Ŵ⟩, ⟨Ĥd⟩, ⟨Ĥb⟩ ≥ 0.
Note that only T̂contains electron momentum operators, while V̂, Ŵ, and Ĥd contain only electron positions. On the other hand, Ĥb and Ĥe xt contain only mode operators, and Ĥc is the only part of the Hamiltonian that includes electron positions and pα. This will be important for commutators: if a commutator is between operators where the first includes only electronic coordinates and the second only mode coordinates, then it will automatically be zero. Additionally to the canonical commutation relation (1), we will use the following basic commutators The Hamiltonian Ĥallows a more compact form, where we gather all terms that include mode coordinates, Ĥb, Ĥc, Ĥd, into one field-energy part with a mode momentum pα Although not used further, this form is useful to see that p̂α corresponds to the displacement field of the mode, shifted by the electronic dipole, and that Ĥb + Ĥc + Ĥd, made out only of squares of self-adjoint operators, is a positive operator and thus ⟨Ĥb⟩ + ⟨Ĥc⟩ + ⟨Ĥd⟩ ≥ 0. For more details concerning the importance of including also the terms Ĥc and Ĥd, see refs 60 and 61. For a broad class of external potentials that include the usual molecular potentials, described in Reed and Simon, 62 Section X.2, the whole Ĥis then bounded below too and thus it is able to support a ground state. We will always assume in the following that v(r) is from that class and is such that Ĥhas an eigenstate, ĤΨ = EΨ.
Before we move to the main subject of this work, we want to point out that the following results do also hold when one considers classical instead of quantum fields. The classical field here means that the ground state in the photon sector is fully determined by the ⟨pα⟩. In this case, we substitute the coupling term like where the last term is just a constant shift to not double-count the energy contribution because of the mean−field interactions. This allows for a factorization of the wave function in the matter and the photon sector and we have two coupled, nonlinear equations. In the matter part, we have to keep the dipole self-energy term Ĥd to still allow for eigenstates. 61 The solution in the photon sector is a tensor product of harmonicoscillator eigenstates for ⟨pα⟩ = 0 and coherent states for ⟨pα⟩ ≠ 0. Because only ⟨pα⟩ enters, we can even replace the Schrodinger equation for the photon sector including the substitution (12) by the mode-resolved Maxwell's equation (30). This highlights that by the substitution (12) we actually consider the coupled Maxwell−Schrodinger problem. 12 The virial relations derived in the subsequent work are therefore also valid if we consider classical fields.

EQUATIONS OF MOTION FOR VIRIAL RELATIONS
The basic idea to arrive at virial relations is to use the simple result that the time derivative of the expectation value of any time-independent operator Âis zero in an eigenstate of Ĥ. To see this one could choose Ψ such that the Schrodinger equation is 13) and evaluate the equation of motion Here and in the following, we use the short-hand notation ⟨Â⟩ = ⟨Ψ, ÂΨ⟩ for the expectation value of an operator with respect to the state Ψ. It must not be noted separately which state Ψ is meant because we will always only consider one and the same exemplary eigenstate. Because similarly we arrive from (13)  it follows ⟨[Ĥ, Â]⟩ = 0 for eigenstates. This result is known as the "hypervirial theorem". 63 Now [Ĥ, Â] could be all sorts of different operators for which we then know that the expectation value is zero, so if we arrive again at terms from the Hamiltonian (2), we have found a virial relation. The difficulty lies in finding the appropriate operators Â.

ELECTRONIC VIRIAL THEOREM WITH
MODE-COUPLING The well-known quantum virial theorem for a system of electrons interacting via Coulomb repulsion and under the influence of an external scalar potential is 59 Here, r = r 1 , ∇ = ∇ 1 and all the particles yield the same term N times because of the assumed antisymmetry of the wave function. Only for special types of external potentials we also arrive back at the expectation value of V̂. However, take the monomial potential v(r) = r n and the virial theorem is 2⟨T⟩ + ⟨Ŵ⟩ = n⟨V̂⟩. The question that puts itself from the previous section is, which operator can be used such that ⟨[Ĥ, Â]⟩ yields exactly the relation (16). The answer is the so-called "virial operator" 64 Â= ∑ i r i ·∇ i and we will demonstrate this with the full mode-coupled Hamiltonian (2) right away, where all parts except of the mode-only Ĥb contribute. We will treat them one after another. Note that in a double sum over all particles, only the i = j part will remain because r i ·∇ i only affects the i-th particle coordinate. Thus

Journal of Chemical Theory and Computation
Taking all these results together the new electronic virial theorem with mode-coupling is found We see that all terms from the original Hamiltonian are involved, except the field-energy of the cavity modes Ĥb and the energy from the external force on the modes Ĥe xt . Thus, this virial relation connects all the constituents of the system which could be expressed of purely electronic degrees of freedom with the dipole-coupling energy which is the only term of the mixed electron-boson nature.
Some comments are in order here. First, by setting the fundamental light-matter coupling to λ α = 0, we recover the virial for electronic structure (16) as requested. An alternative route to derive the basic virial theorem is from the forcebalance equation for stationary states. 14 If the force totals to f(r) = 0 at every point r, as will be the case for eigenstates, then taking the space integral ∫ f(r)·r dr results in (22). Let us also point out that the implications of the force balance equation for functional construction in the context of time-dependent QEDFT 14,16,19 has been discussed in Tokatly. 65

FIELD-MODE VIRIAL THEOREM
Next, we derive an analogous equation of motion on the side of the modes, which means that we combine the mode operators into an mode virial operator ∑ α qαpα. In the equation of motion ∑ α ⟨[Ĥ,qαpα]⟩, only three terms contribute and we use the commutators (1), (8), (9), and (10) The relation (23) yields, as one would expect, just the virial theorem for the harmonic oscillator. Together, with the other relations, that just give back the coupling and external-force energies, we arrive at the field-mode virial theorem We see that the field-mode virial theorem has a similar structure than the electronic virial theorem from (22), with an external influence on the right and a connection between purely bosonic parts of the system with the dipole-coupling energy on the left. The field-mode virial theorem (26) is also a consequence of the equation of motion that involves the second time-derivative of Ĥb.
Another useful relation arises if we choose the equation of motion for the much simpler operator pα. Here, only Ĥb contributes and, using (10) (29) Summing those terms and taking the expectation value that must be zero, we have This relation corresponds to the mode-resolved Maxwell equation for the displacement field, see for example Flick et al. 17 In contrast to the connection between the electronic force-balance equation and the electronic virial theorem in the previous section, it does not give rise to a virial theorem by itself, but it will be employed in the derivation of the next virial relation.

MODE-COUPLED VIRIAL ESTIMATE AND
CONNECTION TO MASS RENORMALIZATION In the search for another virial relation that will relate the electronic and field parts of the Hamiltonian, we consider the equation of motion of an operator that includes both electron coordinates and mode coordinates. A promising candidate, as we will see, is the mixed virial operator ∑ i (λ α ·r i )qα. We derive Journal of Chemical Theory and Computation pubs.acs.org/JCTC Article the non-zero commutators with the individual parts of Ĥusing the basic commutation relations (1), (8), and (9).
This already looks promising, because the constituents of Ĥc and Ĥd reappear. To make this connection explicit, we add the terms above, remember that the expectation value of their sum in an eigenstate is always zero, divide by iω α and sum over α. The resulting expectation value must still be zero The rightmost term can be re-expressed with (30) and the defining equation (7) as The term ⟨Ĥc⟩ + 2⟨Ĥd⟩, relating to dipole coupling and dipole self-energy, also already appears in the basic virial theorem (22), just the remaining first term of (35) appears a bit unwieldy if we want to re-express it by parts of the lengthgauge Hamiltonian. Indeed, it corresponds exactly to the coupling term in velocity gauge, cf. Rokaj et al., 61 (10) or, in other words, to the Fourier transform of Ĥc in both, electronic and mode coordinates. We will handle it further by deriving an estimate. Because −iω α −1 λ α ·∇ j and qα are both self-adjoint, and so is their sum, their sum squared is a positive operator, and because they commute we get This is the point where an estimate enters the virial relation that we will derive in this section, needed to handle the first term in (35) in terms of expressions that come directly from parts of the Hamiltonian Hî To proceed and derive an expression for the first sum above, we have to choose a certain distribution of modes. This could be the one that privileges a specific direction by the choice of λ α , like in cavity QED, but let us assume that the modes are those of free space, evenly distributed in k-space in d = 3 dimensions. To every mode we assign a vector λ α ∈ R 3 , , see Schafer et al., 60 and let always two modes occupy a k-space volume of (ω 0 /c) 3 , where ω 0 = 2πc/L is the lowest frequency in a cube with length L. We take two modes per k-space box because of the two different polarization directions of the photons. If we now choose a radius k = ω α /c, then the number of modes within a thin spherical shell with thickness dk is 2·(ω 0 /c) −3 4πk 2 dk. For any mode in the shell with λ α there are two other modes with λ β , λ γ such that they are all orthogonal to each other and thus give The specific orientations of λ α , λ β , λ γ do not matter at this point because the Laplacian ∇ j 2 is spherical symmetric. The above also tells us that for any three modes within the considered shell volume we get one contribution like above. Thus summing over all modes up to a cutoff k ≤ Λ we have for the complete sum μ = ∑ α ω α −2 λ 2 when passing to an integral This value corresponds to the mass renormalization discussed in Hainzl and Seiringer, 66 (3.41). Some comments on the cutoff Λ introduced here are in order. Although somehow arbitrary, it is a standard ingredient of nonrelativistic QED theory, as the bare mass of the electron is given with respect to such a cutoff. Among other choices, one could choose, for example, to allow for field modes with energy smaller than mc 2 , because for larger energies the nonrelativistic theory is anyway invalid. Such an ultraviolet cutoff is not expected to affect the low-energy processes considered in the setting discussed here. Rewriting the virial relation (35) (41) Usually the mass renormalization is derived by requiring the same dispersion relation for the bare electron and an electron coupled to the free-space modes. Here, we have a field-mode virial relation where the kinetic energy of the electrons enters with a small prefactor that relates to change of electron mass due to the coupling to the free-space modes. A different distribution of modes would just lead to a different mass renormalization prefactor μ, while the relation (41) would stay exactly the same.

QEDFT KOHN−SHAM POTENTIALS TO ENERGIES
VIA THE VIRIAL THEOREM As an application which highlights the usefulness of the virial relations, we show in this section how they can be used to recover energies from potentials in the ground-state Kohn− Sham approach to QEDFT. 15,20 This has special significance in cases where the energy functional is previously unknown