key: cord-0178158-0o5e1l3w
authors: Li, Tao E.; Nitzan, Abraham; Subotnik, Joseph E.
title: Cavity Molecular Dynamics Simulations of Liquid Water under Vibrational Ultrastrong Coupling
date: 2020-04-10
journal: nan
DOI: nan
sha: d3c6ae54d03099098d9bdf0242a6f745783801b3
doc_id: 178158
cord_uid: 0o5e1l3w

We simulate vibrational strong (VSC) and ultrastrong coupling (V-USC) for liquid water with classical molecular dynamics simulations. When the cavity modes are resonantly coupled to the O-H stretch mode of liquid water, the infrared spectrum shows asymmetric Rabi splitting. The lower polariton (LP) may be suppressed or enhanced relative to the upper polariton (UP) depending on the frequency of the cavity mode. Moreover, although the static properties and the translational diffusion of water are not changed under VSC or V-USC, we do find the modification of the orientational autocorrelation function of H$_2$O molecules especially under V-USC, which could play a role in ground-state chemistry.

Strong light-matter interactions between a vibrational mode of molecules and a cavity mode have attracted great attention of late. 1 The signature of strong interactions is the formation 1 arXiv:2004.04888v1 [physics.chem-ph] 10 Apr 2020 of lower (LP) and upper (UP) polaritons, which are manifested in the Rabi splitting of a vibrational peak in the molecular infrared (IR) spectrum. According to the normalized ratio (η) between the Rabi splitting frequency (Ω N ) and the original vibrational frequency (ω 0 ), or η = Ω N /2ω 0 , one often classifies 0 < η < 0.1 as vibrational strong coupling (VSC) and η > 0.1 as vibrational ultrastrong coupling (V-USC). 2 The investigation of VSC or V-USC in liquid phase was initially suggested by Ebbesen et al , [3] [4] [5] and it was later found experimentally that VSC or V-USC can modify the ground-state chemical reaction rates of molecules even without external pumping. 6 This exotic catalytic effect provides a brand new way to control chemical reactions remotely. As such, there has been a recent push to understand the origins and implications of VSC and V-USC.

While the experimental side has focused on the search for large catalytic effects 7-10 as well as understanding polariton relaxation dynamics through two-dimensional IR (2D-IR) spectroscopy, 11, 12 on the theoretical side, the nature of VSC and V-USC remains obscured.

On the one hand, Rabi splitting can be easily modeled by, e.g., diagonalizing a model Hamiltonian in the singly excited manifold [13] [14] [15] or solving equations of motion classically for a set of one-dimensional (1D) harmonic oscillators. 16, 17 On the other hand, a robust explanation of the catalytic effect of VSC or V-USC remains illusive. [18] [19] [20] [21] For example, as recently shown by us and others, 21-23 the potential of mean force along a reaction pathway is not changed by usual VSC or V-USC setups for standard experiments of interest. Moreover, as demonstrated below, any static equilibrium property of a molecule is not changed under VSC or V-USC.

These findings, unfortunately, show that one cannot explain the observed effect under VSC or V-USC from a static view of point. From such a conclusion, one must hypothesize that the manifestations of VSC or V-USC effect on chemical rates should arise from the modification of non-equilibrium, or dynamical, properties of molecules under VSC or V-USC.

The first step towards proving the above hypothesis is to ascertain whether or not any dynamical property of molecules is actually changed for a realistic experiment, a goal which forms the central objective of this manuscript. In order to investigate whether such modification occurs, below we will model VSC and V-USC using cavity molecular dynamics (MD) simulation, where the nuclei are evolved under a realistic electronic ground-state potential surface. Such an approach is an extension of the usual simplified 1D models where the matter side is evolved as two-level systems [24] [25] [26] or coupled harmonic oscillators. 16, 17, 27, 28 Although such simplified models are adequate enough for studying Rabi splitting qualitatively by fitting experimental parameters, these models usually ignore translation, rotation, collision, as well as the intricate structure of molecular motion, all of which are crucial for determining the dynamic properties of molecules. Therefore, explicit cavity MD simulations become a more appropriate approach for studying all dynamic properties. Moreover, even though one can find a Rabi splitting from 1D models, performing cavity MD simulations is also very helpful for as providing more details about the IR spectrum and this approach can be used to benchmark the validity of 1D models under various conditions.

There have been a few flavors of cavity MD schemes for electronic strong coupling. [29] [30] [31] For example, Luk et al applied multiscale quantum mechanics/molecular mechanics (QM/MM) simulation for studying the dynamics of electronic polaritons for Rhodamine molecules. 30 By contrast, MD simulations for vibrational strong coupling (VSC and V-USC), to our best knowledge, have not been extensively studied before. Therefore, below we will first establish a framework for cavity MD simulation including implementation details, and second we will investigate the Rabi splitting and the dynamical properties of liquid water.

The motivation for studying liquid water is two-fold: (i) Among common liquids, water shows strong Rabi splitting and strong catalytic effects under VSC or V-USC. 8, 10, 32 More interesting, when the cavity mode is resonantly coupled to the O H stretch mode, experiments 10 have observed that the intensity of the vibrational LP peak is much smaller than the UP peak in the IR spectrum, an observation that cannot be accounted for by standard strong coupling models. (ii) MD simulations of water outside the cavity have been extensively studied and good agreement with experiments can be achieved. [33] [34] [35] Extending such simulations to include coupling to cavity modes is expected to show the cavity-induced spectral changes and provides numbers that are directly comparable to experimental results. This manuscript is organized as follows. Sec. 2 provides theoretical background of VSC and V-USC. Sec. 3 describes how to perform cavity MD simulations for liquid water. Sec. 4 shows simulation results of liquid water. Sec. 5 concludes our numerical findings and provides future directions for research into VSC and V-USC.

The full-quantum Hamiltonian for light-matter interactions reads: 21, 36 

Here,Ĥ M denotes the conventional (kinetic + potential) Hamiltonian for the molecular

where m i ,p i ,r i denote the mass, momentum operator, and position operator for the i-th particle (nucleus or electron), respectively, andV Coul ({r i }) denotes the Coulombic interaction operator between all nuclei and electrons. Under the long-wave approximation, the fieldrelated HamiltonianĤ F readŝ

where ω k,λ ,q k,λ ,p k,λ denote the frequency, position operator, and momentum operator for a photon with wave vector k and polarization direction ξ λ , and the index λ = 1, 2 denotes the two polarization directions which satisfy k · ξ λ = 0. In free space, the dispersion relation

gives ω k,λ = c|k| = ck. 0 and Ω denote the vacuum permittivity and the cavity volume.μ S denotes the dipole operator for the whole molecular system:

where e denotes the electron charge and Z i e denotes the charge for the i-th particle (nucleus or electron).μ S can also be grouped into a summation of molecular dipole moments (indexed by n):μ

Note that the self-dipole term in Eq. (1c) (i.e., theμ 2 S term in the expanded square) is of vital importance in describing USC and is needed to render the nuclear motion stable; see

Refs. [37] [38] [39] for details. Because we will not neglectμ 2 S below, our simulation is valid for both VSC and U-VSC.

When the cavity mode frequency is within the timescale of the nuclear dynamics, the Born-Oppenheimer approximation implies that electrons stay in the ground state. Therefore, we will project the quantum Hamiltonian (1) onto the electronic ground state,Ĥ G QED = Ψ G |Ĥ QED |Ψ G , where |Ψ G denotes the electronic ground state for the whole molecular system. Furthermore, under the Hartree approximation, |Ψ G can be approximated as a product of the electronic ground states for individual molecules: |Ψ G = N n=1 |ψ ng . After such a projection on the electronic ground state, the Hamiltonian (1) reduces tô

Here, the ground-state molecular HamiltonianĤ G M = Ψ G |Ĥ M |Ψ G depends on the nuclear degrees of freedom only, and can be expressed aŝ

where the capital lettersP nj ,R nj , and M nj denote the momentum operator, position operation, and mass for the j-th nuclus in molecule n,V inter denotes the intermolecular interactions between molecule n and l.

The field-related Hamiltonian becomes 21

where we defined ng,λ ≡ ψ ng |μ n |ψ ng · ξ λ (5a) Similarly, on the last line in Eq. (4c), the self-dipole fluctuation term 1 2Ω 0 ψ ng |δd 2 ng,λ |ψ ng , which denotes the cavity modification of the single-molecule potential, should also be very small for standard VSC setups where micron-length cavities are used. Therefore, in what follows, we will assume thatV 

Here, to be compatible with standard MD simulations (which requires the information of mass for particles), an auxiliary mass m k,λ for each photon is also introduced:

Note that the auxiliary mass of photon does not alter any dynamics and serves only as a convenient notation for further MD treatment.

The above quantum Hamiltonian, although depending only on the nuclear and photonic degrees of freedom, is still too expensive to evolve exactly. The simplest approximation we can make is the classical approximation, i.e., all quantum operators are mapped to the corresponding classical observables, which leads to the following classical Hamiltonian:

Eq. (8) serves as the starting point of this work. We note that one can go beyond the treatment here by propagating the quantum Hamiltonian (4) using the path-integral technique 43,44 and evolve the ring polymer Hamiltonian with n copies of coupled classical trajectories (aka n beads). In the present manuscript, we focus on the classical system, deferring the path-integral calculation to a later study.

In our classical MD simulations, the simulated system is represented by particles that obey the Newtonian equations of motion:

where the cavity-free force F

inter /∂R nj , and the coupling between particles representing photons and nuclear degrees of freedom is

In order to perform a realistic simulation for VSC or V-USC, we need a macroscopic number (say, 10 9 ∼ 10 11 ) of molecules, 6-9 which is far beyond our computational power if we simulate Eq. (9) directly. To proceed, we assume that the whole molecular ensemble can be divided into N cell periodic cells, in which the molecules evolve identically, i.e., we can approximate the second term on the right of Eq. (9b) by N n=1 d ng,λ = N cell N sub n=1 d ng,λ , where N sub = N/N cell denotes the number of molecules in a single cell. By further denoting

we can rewrite the equations of motion in Eq. (9) in a symmetric form

The form of Eq. (10) has several advantages. First, we simulate the VSC of a macroscopic number of molecules by evolving molecules in a single cell plus the few photon modes that we are interested in. Second, when considering the dependence of Rabi splitting on molecular numbers, we can fix the number of molecules in a single cell (N sub ) and vary only the coupling constant ε k,λ = √ N cell ε k,λ . Such a change is very easy to implement in practice and has the physical interpretation of increasing the number of cells while leaving the number of molecules per cell and the size of the simulation cell fixed.

The question remains as to exactly how we will calculate the ground-state quantities F (0) nj , d ng,λ , and ∂d ng,λ /∂R nj . In general, these properties can be calculated by classical empirical force field or ab initio electronic structure theory. For this initial publication of liquid water, we will use a classical force field -the q-TIP4P/F water model 34 -which provides the simplest description of both the equilibrium and dynamic properties of liquid water.

In the q-TIP4P/F model, the pairwise intermolecular potential is characterized by the Lennard-Jones potential between oxygen atoms plus the Coulombic interactions between partial changes:

where R OO nl denotes the distance between the oxygen atoms and R ij denotes the distance between the partial charge sites in molecules n and l. Within a single H 2 O molecule, two positively partial charge with magnitude Q M /2 are assigned to the hydrogen atoms, and the negatively charge site with magnitude −Q M is placed at R M :

For parameters, = 0.1852 kcal mol −1 , σ = 3.1589 Å, Q M = 1.1128 |e| (where e denotes the charge of the electron), and γ = 0.73612.

The intramolecular interaction is characterized by

where

Here, R n1 and R n2 denote the lengths of two O H bonds, θ n and θ eq denote the H O H angle and the equilibrium angle. For parameters, D r = 116.09 kcal mol −1 , α r = 2.287Å −1 , r eq = 0.9419Å, k θ = 87.85 kcal mol −1 rad −2 , and θ eq = 107.4 deg.

Given the q-TIP4P/F force field, one can easily calculate the cavity-free force F

nj as a function of the nuclear configurations by standard molecular dynamics packages. The dipole moment can be calculated by

and the derivative ∂d ng,λ /∂R nj is trivial. Our modification is illustrated as the green region in Fig. 1 (10)).

After calculating the forces, we use the interface of I-PI to update momenta and positions.

Due to the user-friendly structure of I-PI, the current cavity MD code should be easily generalized to the cases of ab initio calculation and path-integral cavity MD simulations, results which will be reported in a separate publication. We consider the following scenario for simulation. As shown in Fig. 2 The simulation step is set as 0.5 fs and we store the snapshots of trajectories every 2 fs.

The signature of VSC is the collective Rabi splitting in the IR spectrum. In our MD simulations, the IR spectrum is calculated by linear response theory. For isotropic liquids, the absorption coefficient α(ω) is expressed as the Fourier transform of the autocorrelation function of the total dipole moment µ S : 48-52

Here, n(ω) denotes the refractive index and V denotes the volume of the system (i.e., the simulation cell). The factor ω 2 arises from the absorbed photon energy by the liquid. For Note that, as ε increases, the LP peak is suppressed and the UP peak is enhanced.

is along the z-axis), 32 we need to modify the above equation to

where e i denotes the unit vector along direction i = x, y. Eq. (17) states that the average is performed only along the polarization directions of the detecting signal (i.e., the x and y directions here). When the incident light is unpolarized these two directions are of course equivalent. and UP peaks. More interestingly, our simulation results also suggest that the UP and LP peaks can be largely asymmetric especially when ε is large, which agrees with experimental findings at least qualitatively. 10 In Fig. 4a we plot the Rabi splitting frequency (the difference between the UP and LP frequencies, or ω + − ω − ) as a function of ε. The simulation data (black triangles) can be fit with a linear ansatz (gray line) very well. As mentioned above, because ε = √ N cell ε ∝ √ N , Fig. 4a demonstrates that the Rabi splitting is proportional to the square root of the total number of molecules, which agrees with theoretical expectation and experimental observation: 32,53

where g 0 denotes the coupling constant between a single molecule and the photon mode. Fig. 3 . In Fig. 4a the simulation data (black triangles) are fit linearly (gray line). In Fig. 4b -c the simulation data (blue stars for LP and red circles for UP) are compared with the analytical expressions from a simplified 1D model (Eqs. (19) and (20)), where the parameters are given as ω 0 = ω c = 3550 cm −1 and Ω N is taken as the values in Fig. 4a .

Of particular interest is the asymmetric nature of the LP and UP: this asymmetry is manifest in two aspects. As shown in Fig. 4b-c , both the polariton frequencies and the integrated peak areas of the LP (blue stars) and UP (red circles) show asymmetric scalings as a function of the normalized Rabi frequency (Ω N /2ω 0 , where Ω N is taken from Fig. 4a) , especially in the V-USC limit (the red-shadowed region). Note that the standard treatment of collective Rabi splitting does not account for this asymmetry and the observation of the suppression (or enhancement) of the LP (or the UP) in Ref. 10 was explained by the higher absorption of water and gold cavity mirrors in the LP region. Some insight into the origin of this asymmetry can be obtained from a simple 1D model where N independent harmonic oscillators interact with a single photon mode. As calculated in the Appendix, by taking the self-dipole term into account (to describe V-USC), we obtain

where ω 0 and ω c denote the frequencies of the harmonic oscillators and the photon mode.

Given ω 0 = ω c = 3550 cm −1 and Ω N in Fig. 4a Fig. 4b . We see that this analytical result already shows some asymmetry in the positions of the polariton peaks when plotted versus Ω N . While Eq. (19) agrees with our simulation data very well in the VSC limit (the green-shadowed region), the simulation data seem to be more asymmetric than Eq. (19) in the V-USC limit. Such disagreement may arise from the strong intermolecular interactions between H 2 O molecules, which is completely ignored in the simplified 1D model of the Appendix.

Likewise, the simplified 1D model in the Appendix also suggests that the integrated peak areas of the LP and UP are

. Again, as shown in Fig. 4c , Eq. (20) (black dashed lines) matches the simulation data roughly but not quantitatively, which may come from ignoring all the intermolecular interactions in the 1D model. Nevertheless, from Eq.

(20), we find that the asymmetry in the IR spectrum comes from two factors: (i) the factor ω 2 ± and (ii) the angular part sin 2 θ 2 or cos 2 θ 2 . While the first part originates from the absorbed photon energies associated with the vibration modes and is universal for all IR spectrum (so that it is trivial), the second factor is quite nontrivial: at resonance (ω 0 = ω c ) one would naively assume that sin 2 θ 2 = cos 2 θ 2 and this is true if one ignores the selfdipole term (which means ignoring the Ω 2 N term in tan (θ); see the Appendix for details).

However, when the self-dipole term is considered, one finds sin 2 θ 2 < cos 2 θ 2 , which leads to an additional suppression of the LP and the enhancement of the UP.

For liquid water in the cavity, in Fig. 5 , we further investigate how (a) the polariton frequencies and (b) the integrated peak areas of polaritons depend on the cavity mode frequency for ε = 5 × 10 −4 a.u., which is well in the USC regime. The simulation data (scatter points) agree well with the analytical result (dashed black lines) for the simplified 1D model (Eqs. (19) and (20)). As shown in Fig. 5a , the energy difference between the polaritons is minimal at resonance (∼ 3550 cm −1 ), in which the uncoupled O H stretch mode frequency crosses with the cavity mode frequency; see gray solid lines. Such a cross corresponds to the maximally hybridized light-matter state. By contrast, when the cavity mode frequency is larger (smaller) than the molecular frequency, the LP (UP) becomes increasingly dominated by the O H stretch mode (as evident from the uncoupled case for which this mode is represented by the gray horizontal line).

Our model implies that for the uncoupled molecule-cavity case, only the molecular optical transition is coupled to the far field. This suggests that in contrast to the resonance case when the UP peak is larger than the LP peak, when the cavity mode frequency becomes sufficiently large (i.e., the LP is mostly constituted by the matter side), the LP should have a larger peak size than the UP. This finding is confirmed by Fig. 5b. More interestingly, Fig.   5b also shows the symmetric peak size of polaritons occurs when the cavity mode frequency is ∼ 4250 cm −1 , which is far beyond the O H stretch frequency of liquid water (∼ 3550 cm −1 ). Therefore, in principle, from this fact, one would predict that one can engineer the relative ) and (20)). For simulation parameters, ε = 5 × 10 −4 a.u. and all other parameters are the same as in Fig. 4 . For parameters of the analytical expressions, we take ω 0 = 3550 cm −1 and Ω N = 937 cm −1 , which corresponds to the resonant Rabi frequency when ε = 5×10 −4 a.u. (see Fig. 4a ). The gray solid lines in Fig.  5a represents the uncoupled O H stretch mode frequency and the cavity mode frequency. The insert in Fig. 5b plots the cavity mode frequency corresponding to the case of symmetric polaritons (i.e., the crossing point frequency in Fig. 5b ) as a function of Ω N /2ω 0 . strength of polaritons by tuning the cavity mode frequency. Furthermore, the inset of Fig.   5b plots the cavity mode frequency (for which the polariton intensities become symmetric) as a function of Ω N /2ω 0 . Again, we find that for large Ω N /2ω 0 , detecting polaritons with symmetric intensities requires a very large off-resonant cavity mode frequency. (g(r) ) of oxygen atoms in liquid water. The result outside the cavity (solid black) is compared with that inside the cavity (with effective coupling strength ε = 4 × 10 −4 a.u., cyan dashed). All other parameters are set the same as Fig. 3 . Note that g(r) is not changed by VSC or V-USC.

Rabi splitting represents the collective optical response of liquid water. As shown above, although MD simulations can obtain the IR spectrum of the polaritons in a straightforward way, one can argue that since most important features of the IR spectrum can be qualitatively described by the 1D harmonic model (see the Appendix), there is little advantage to perform expensive MD simulations. As has been argued above, the real advantage of the MD simulations is that one can simultaneously obtain many other physical properties of molecules alongside with the IR spectrum. Below we will investigate whether any property of individual H 2 O molecules can be changed under VSC or U-VSC. According to linear response theory, the translational diffusion of H 2 O can be described by the VACF (C vv (t)) of the center of mass of each molecule:

One can calculate the diffusion constant D from C vv (t) by D = 1 3 +∞ 0 C vv (t)dt. C vv (ω), which is shown in Fig. 8b Fig. 9b so we do not report it here.

As for the rotational behavior, according to linear response theory, one must compute the orientational autocorrelation function (OACF, denoted by C l (t)), [54] [55] [56] which is defined

where u n (t) denotes the three principal inertial axes of molecule n at time t, and P l denotes the Legendre polynomial of index l. For simplicity, we will study only the first order of OACF, which means P 1 [u n (0) · u n (t)] = u n (0) · u n (t).

For H 2 O, the z axis of the principal axes coincides with the dipole moment direction. In respectively]. Fig. 9b plots the corresponding spectrum I z 1 (ω), which is defined as

I z 1 (ω) can be regarded as the single-molecule IR spectroscopy along the dipole-motion direction, which describes how a single molecule rotates in the environment. As clearly shown in the zoom-in inset, for large enough ε (in the V-USC limit, or ε ≥ 4×10 −4 a.u.), an additional small peak emerges with intensities 2% ∼ 8% of the peak from a bare molecule. Compared with the IR spectrum of the liquid water in Fig. 3 , these additional small peaks have the same frequencies as the UP peaks, demonstrating the modification of single-molecule rotation under V-USC. Note that for smaller ε (i.e., in the VSC limit), the additional peak will be covered by the large bare-molecule peak and is hardly identifiable. The change of the rotational behavior of individual molecules may possibly change the ground-state chemistry for many scenarios, which should be extensively studied in the future. Lastly, we emphasize that apart from these additional peaks, the width of the bare-molecule peaks is mostly unchanged. (ω) = ω 2 C z 1 (ω)). A zoom-in inset is also plotted in each subplot. The results outside the cavity (black dashed) is compared with those inside the cavity (with effective coupling strength ε as 4 × 10 −4 (cyan solid), 6 × 10 −4 (red dashed), and 8 × 10 −4 a.u. (blue dash-dotted), respectively. All other parameters are set the same as Fig. 3 . Note that OACF is changed by V-USC.

In conclusion, we have performed classical cavity MD simulations under VSC or V-USC.

With liquid water as an example, when the cavity modes are resonantly coupled to the O H stretch mode, we have found asymmetric Rabi splitting of the O H stretch peak in the IR spectrum where the LP is suppressed and the UP is enhanced. Such asymmetry can be inverted (i.e., the LP is enhanced and the UP is suppressed) by increasing the cavity mode frequency. Moreover, while we have found no modification of the static equilibrium properties as well as the translational diffusion of liquid water, we have observed that the 

Starting from the classical Hamiltonian for V-USC (Eq. (8)), let us assume that (i) the molecules are non-interacting 1D harmonic oscillators, (ii) the dipole moment for a single molecule is linear (i.e., d ng,λ = d 0 x n ), and (iii) only a single cavity mode is considered. With these simplifications, the Hamiltonian can be written as:

Here, g 0 ≡ d 0 /2 √ Ω 0 , and we have assumed all masses to be 1. Note that the self-dipole term (the N n=1 x n 2 term in the expanded square above) is necessary for studying V-USC.

With the Hamiltonian in Eq. (A1), the equations of motion now read

Let us define the bright mode as

x n , so that the equations of motion for the bright mode and the cavity mode becomë

where Ω N = 2 √ N g 0 is the usual Rabi frequency. In the matrix form, the above equations can be written as¨

Note that the Ω 2 N term above comes from the self-dipole term.

The polariton frequencies (ω ± ) can be determined by solving the eigenvalues of the matrix K:

At resonance (ω c = ω 0 ), the polariton frequencies are reduced to

In the VSC limit (Ω N ω 0 ), Eq. (A7) can be further simplified as

which is the usual strong-coupling result.

The IR spectrum of molecules is calculated by Eq. (16) . With our 1D model, the IR spectrum is expressed as

where we have neglected all prefactors (including the temperature as we take room temperature throughout this manuscript). According to Eq. (A3), the solution of x B (t) is

where tan (θ) = 2ω c Ω N / ω 2 0 + Ω 2 N − ω 2 c (A10b) By substituting Eq. (A10) into Eq. (A9) and using x B (0)x c (0) = 0, we obtain n(ω)α(ω) ∝ ω 2 cos 2 θ 2 δ(ω − ω + ) + sin 2 θ 2 δ(ω − ω − ) (A11)

The integrated peak areas for LP and UP are I LP ∝ ω 2 − sin 2 θ 2 (A12a)

From the above, we find that the asymmetric peaks come from two origins: (i) the prefactor ω 2 ± and (ii) the self-dipole term in the dipole-gauge Hamiltonian. While the first origin is trivial, the second origin can be understood as follows. If we had neglected the self-dipole term, we would naively take ω 2 0 +Ω 2 N → ω 2 0 in Eq. (A5) and obtain a different expression for θ, tan (θ) = 2ω c Ω N / (ω 2 0 − ω 2 c ) (where the Ω 2 N term now vanishes compared to Eq. (A10b) ). At resonance, we would obtain that tan (θ) = ∞, i.e., θ = π/2 and cos 2 (θ/2) = sin 2 (θ/2) = 1/2.

However, because of the Ω 2 N term, θ < π/2, which implies sin 2 (θ/2) < cos 2 (θ/2). In other words, when ω c = ω 0 , the self-dipole term forces the LP to be further suppressed and the UP be further enhanced.

(42) In principle, when one considers the exact quantum Hamiltonian for systems with lightmatter interactions, all (i) instantaneous interactions between molecules are canceled exactly by the presence of terms that involve (ii) the non-local (and also instantaneous) self-interaction of delocalized photon modes. This exact cancellation allows for causality to be enforced, such that all meaningful intermolecular interactions are carried exclusively the transverse photon field at the speed of light. In the present paper, we do not worry about causality and so we have ignored the details of the cancellation alluded to above, i.e. we do not address how this cancellation is affected by the cavity and the presence of a finite number of cavity modes. In principle, the presence of a cavity leads to a dressedV 

Molecular polaritons for controlling chemistry with quantum optics

Ultrastrong coupling between light and matter

Coherent coupling of molecular resonators with a microcavity mode

Liquid-Phase Vibrational Strong Coupling

Multiple Rabi Splittings under Ultrastrong Vibrational Coupling

Ground-State Chemical Reactivity under Vibrational Coupling to the Vacuum Electromagnetic Field. Angew. Chemie Int

Cavity Catalysis by Cooperative Vibrational Strong Coupling of Reactant and Solvent Molecules

Cavity Catalysis -Accelerating Reactions under Vibrational Strong Coupling

Tilting a ground-state reactivity landscape by vibrational strong coupling

Modification of Enzyme Activity by Vibrational Strong Coupling of Water. Angew. Chemie Int

Two-dimensional infrared spectroscopy of vibrational polaritons

State-Selective Polariton to Dark State Relaxation Dynamics

Theory of the Contribution of Excitons to the Complex Dielectric Constant of Crystals

Multi-level quantum Rabi model for anharmonic vibrational polaritons

Theory for Polariton-Assisted Remote Energy Transfer

Oscillator model for vacuum Rabi splitting in microcavities

Theory for Nonlinear Spectroscopy of Vibrational Polaritons

Cavity Casimir-Polder Forces and Their Effects in Ground-State Chemical Reactivity

Resonant catalysis of thermally activated chemical reactions with vibrational polaritons

A Reaction Kinetic Model for Vacuum-Field Catalysis Based on Vibrational Light-Matter Coupling

On the Origin of Ground-State Vacuum-Field Catalysis: Equilibrium Consideration

Polaritonic normal modes in Transition State Theory

Vacuum field in a cavity, light-mediated vibrational coupling, and chemical reactivity

Quantum trajectory simulation of controlled phase-flip gates using the vacuum Rabi splitting

Quasiclassical modeling of cavity quantum electrodynamics

Benchmarking semiclassical and perturbative methods for real-time simulations of cavity-bound emission and interference

Vacuum Rabi splitting in a plasmonic cavity at the single quantum emitter limit

Effects of exciton-plasmon strong coupling on third harmonic generation by two-dimensional WS 2 at periodic plasmonic interfaces

Atoms and Molecules in Cavities, from Weak to Strong Coupling in Quantum-Electrodynamics (QED) Chemistry

Multiscale Molecular Dynamics Simulations of Polaritonic Chemistry

Tracking Polariton Relaxation with Multiscale Molecular Dynamics Simulations

Vibrational Ultra Strong Coupling of Water and Ice

A general purpose model for the condensed phases of water: TIP4P

Zero point energy leakage in condensed phase dynamics: An assessment of quantum simulation methods for liquid water

Combined electronic structure/molecular dynamics approach for ultrafast infrared spectroscopy of dilute HOD in liquid H2O and D2O

Photons and Atoms: Introduction to Quantum Electrodynamics

LightâĂŞmatter interaction in the long-wavelength limit: no ground-state without dipole self-energy

Relevance of the Quadratic Diamagnetic and Self-Polarization Terms in Cavity Quantum Electrodynamics

Effect of Many Modes on Self-Polarization and Photochemical Suppression in Cavities

Applying electric field to charged and polar particles between metallic plates: Extension of the Ewald method

Cavity quantum electrodynamics in the nonperturbative regime

Fast Parallel Algorithms for Short-Range Molecular Dynamics

CP2K: atomistic simulations of condensed matter systems

Ab Initio Molecular Dynamics Computation of the Infrared Spectrum of Aqueous Uracil

Comparison of path integral molecular dynamics methods for the infrared absorption spectrum of liquid water

Chemical Dynamics in Condensed Phases: Relaxation, Transfer and Reactions in Condensed Molecular Systems

This expression reflects one of several suggestions that were made for a correction factor which relates the quantum time-correlation function to its classical counterparts

Elements of Quantum Optics

Reorientational correlation functions for computersimulated liquids of tetrahedral molecules

Spectroscopic and transport properties of water

Quantum diffusion in liquid water from ring polymer molecular dynamics

This material is based upon work supported by the U.S. Department of Energy, Office of