A Quantum Chemistry Approach to Linear Vibro-Polaritonic Infrared Spectra with Perturbative Electron–Photon Correlation

In the vibrational strong coupling (VSC) regime, molecular vibrations and resonant low-frequency cavity modes form light–matter hybrid states, vibrational polaritons, with characteristic infrared (IR) spectroscopic signatures. Here, we introduce a molecular quantum chemistry-based computational scheme for linear IR spectra of vibrational polaritons in polyatomic molecules, which perturbatively accounts for nonresonant electron–photon interactions under VSC. Specifically, we formulate a cavity Born–Oppenheimer perturbation theory (CBO-PT) linear response approach, which provides an approximate but systematic description of such electron–photon correlation effects in VSC scenarios while relying on molecular ab initio quantum chemistry methods. We identify relevant electron–photon correlation effects at the second order of CBO-PT, which manifest as static polarizability-dependent Hessian corrections and an emerging polarizability-dependent cavity intensity component providing access to transmission spectra commonly measured in vibro-polaritonic chemistry. Illustratively, we address electron–photon correlation effects perturbatively in IR spectra of CO2 and Fe(CO)5 vibro-polaritonic models in sound agreement with nonperturbative CBO linear response theory.

We explicitly derive approximate mode effective charges and vibro-polaritonic Hessian matrix elements for the CBO-PT (1) and CBO-PT(2) linear response approach.In order to simplify notation, we define polarization-projected permanent dipole derivatives and vertical transition dipole moments as Further, we rewrite the light-matter coupling constant of the kth-cavity mode as which allows us to express all terms via the universal coupling constant, g 0 , which simply depends on the electric permittivity, 0 , and the cavity volume, V cav .Further, g 0 has units of √ E h /(ea 0 ), since [g k ] = E h /(ea 0 ) and Alternatively, one directly obtains (S1.4) * ericwfischer.sci@posteo.de† peter.saalfrank@uni-potsdam.de The second-order energy in Rayleigh-Schrödinger perturbation theory is given by We derive approximate mode effective charges in the normal mode representation.First, we linearize CBO-PT(n) permanent dipole moment components as with a constant term, D (0,n) 00,κ , which will not contribute to vibro-polaritonic IR intensities.
The mth-vibropolaritonic transition in CBO-PT(n) is characterized by an intensity corresponding to a transition between vibro-polaritonic states, |0 and |m (n) , which read in harmonic approximation where |m (n) explicitly depends on the order n of CBO-PT(n).Here, |0 is the global ground state of all molecular normal and cavity modes, which is assumed to be not influenced under VSC.Further, the vibrational polariton, |m (n) , resembles a linear combination of singly excited normal mode, |1 , and cavity mode states, |1 λ k , with expansion coefficients, q (n) mj and c (n) mλ k , where our notation implicitly assume all remaining modes to reside in their respective ground states.
We now introduce the approximate mode effective charge for the linearized permanent dipole moment which turns with the expansion of |m (n) into .
Further, matrix elements between harmonic oscillator states evaluate to which eventually leads to the CBO-PT(n) mode effective charge in normal mode representation In harmonic approximation, the CBO-PT(0) cPES is given by From here, the zeroth-order molecular and cavity block contributions follow immediately as = ω 2 k δ kk δ λλ , (S1.27) QC i,λk = 0 , (S1.28) with molecular normal mode frequencies, ω i , and cavity mode frequencies, ω k .

C. CBO-PT(1) Hessian Matrix Elements
The first-order CBO-PT Hessian is determined by the first-order perturbative energy correction (S1.30) In the second term, we introduced the resolution-of-theidentity (RI), α | = 1 (e) , in the adiabatic electronic subspace between projected dipole moments of the DSE term.We now apply the Franck-Condon (FC) approximation i.e., only vertical transition dipole moments are considered at the molecular reference configuration, R 0 .Then, we have where the last term depends only on the nuclear reference configuration due to the FC approximation and will therefore not contribute to the Hessian.We linearize the polarization-projected permanent dipole moment around the reference configuration (double-harmonic (DH) approximation) with mass-weighted normal mode coordinate, Q i , and In DH-and FC-approximations, Hessian-relevant contributions in W 00 reduce to The molecular normal mode block is only determined by the DSE term (∝ g 2 0 ) due to its quadratic dependence on normal mode coordinates, such that where we obtained the third line via the identity where an additional factor two results from the fact, that both Kronecker delta products provide the same term.
The forth line follows from the assumption of identical polarization vectors for all cavity modes 2Nc λ,k which allows us to perform the sum over mode index k explicitly.Next, the interaction matrix elements follow as where we treat, (λk), as a block-matrix index, such that every k specifies a corresponding 2 × 2-subblock and λ addresses related subblock elements.Finally, we have CC λk,λ k = 0 .(S1.41)

D. CBO-PT(2) Hessian Matrix Elements
The CBO-PT(2) Hessian is determined by the corresponding energy correction where the matrix element, W µ0 , follows in DH-/FCapproximations as Here, all constant or linear terms, which do not contribute to the Hessian, were neglected.Further, we consider only vertical excitation energies in FC approximation Starting with the molecular contribution as determined by the squared second term of Eq.(S1.43),we find , where we took into account the FC approximation in the first line, i.e., derivatives do not act on vertical excitation energies, ∆ 0µ , and used identity Eq.(S1.37) in the third line.
Further, twice the cross term in Eq.(S1.43)provide corrections to the normal mode-cavity coupling block as , (S1.46) and the cavity block finally results from the squared first term of Eq.(S1.43) as , (S1.47) where we used a slightly more general variant of Eq.(S1.37) in the third line The static polarizability tensor, α 0 , is defined by elements [2] , (S1.49) which we can exploit to rewrite CBO-PT(2) Hessian matrix elements above.Since elements, α 0 ab , can be obtained from linear response theory, Eq.(S1.49)allows for circumventing the unfavourable convergence properties of the sum-over-states expression.[2] Further, since polarization-projected transition dipole moments in CBO-PT(2) Hessian matrix elements do not differ for different cavity mode indices, k, we can introduce polarization-projected static polarizability tensor elements as such that CBO-PT(2) Hessian matrix elements can be compactly expressed as , (S1.51)

S2. DETAILS ON FRANCK-CONDON APPROXIMATION
We demonstrate that the Franck-Condon approximation is equivalent to neglecting corrections of CBO-PT(2) Hessian molecular (QQ) and light-matter interaction (QC) blocks, which scale as (∆ 0 , to depend on normal mode coordinates.The corresponding first derivative follows as and the second derivative takes the form with ∂ ∂Qi = ∂ Qi .In the FC approximation, we approximate excitation energies to be independent of normal mode coordinates, i.e., being only vertical, with such that all terms in Eq.(S2.2) besides the first one vanish identically.Thus, contributions proportional to (∆  In DH-and FC-approximations, the given matrix elements read W 0µ , (S3.6) . (S3.7) In the following, we consider only leading order corrected IR intensities obtain from Eqs.(S3.5) and (S3.6) and neglect higher-order contributions of Eq.(S3.7),such that The cross term in Eq.(S3.6) is explicitly given by (S3.9) where we used the identity Thus, the linearized CBO-PT(2) dipole moment can be written as , (S3.12) and contains a molecular term with cartesian components (κ = x, y, z) as well as a cavity component The normal-mode derivative of the molecular contribution follows as and the cavity coordinate derivative is given by S4. COMPUTATIONAL AND VIBRO-POLARITONIC DETAILS A. Computational Details Molecular equilibrium structures of CO 2 and Fe(CO) 5 have been optimized via TPSSh/Def2TZVP with a subsequent normal mode analysis providing normal mode frequencies and permanent dipole derivatives.Static polarizability tensor elements were obtained via CAM-B3LYP/Def2TZVP for both molecules.All calculations were performed with Gaussian16.[4] Optimized equilibrium structures are provided in Tab.I.Dipole derivatives for the asymmetric CO 2 stretching mode (TPSSh/Def2TZVP) and non-vanishing static-polarizability tensor elements for CO 2 (CAM-B3LYP/Def2TZVP) are given in Tab.II.Axis are given with respect to optimized equilibrium structures.Same properties for equatorial e -and axial a 2 -modes of Fe(CO) 5 are given in Tab.III.
Vibro-polaritonic IR spectra with Lorentzian peak broadening are evaluated as for a frequency interval, ω ∈ [0, 3600]cm −1 , discretized by N ω = 4000 equidistant grid points and κ given in the main text.In harmonic oscillator systems, upper (+) and lower (-) vibro-polaritonic states resemble symmetric and antisymmetric linear combinations of singly excited molecular vibrational or cavity mode states as given by  In Tab.IV, CBO-PT(n) expansion coefficients are given for n = 1, 2 at light-matter coupling g 0 = 0.03 √ E h /ea 0 .For CBO-PT(1), lower and upper vibrational polariton states are nearly equally determined by light and matter contributions.This is contrasted by CBO-PT(2) results, where the lower polariton state is dominantly photonic and the upper polariton dominantly molecular in character in line with non-perturbative CBO results [3].In Fig. S1, σ (2) IR ( ω) and σ (C,2) IR ( ω) for antisymmetric CO 2stretching mode under VSC are shown for the full frequency range, which reveals a small additional peak in σ ( ω) as they do not carry a photonic component.

S5. GENERALIZATION TO MOLECULAR ENSEMBLES
We generalize both CBO-PT(2) Hessians and IR intensities to ensembles of M aligned, identical molecules.For M molecules, we have M N vib molecular normal modes, which is the new molecular dimension.Specifically, we can simply understand the ensemble like a large molecule with distinct but degenerate sets of normal modes, such that both normal mode and light-matter interaction block entries, QQ and QC, remain in CBO-PT(1) as derived above.Further, for CBO-PT(2) we consider the ensemble correction E and we obtain the same light-matter interaction and molecular corrections, i.e., QQ and QC, as in the singlemolecule limit, since differentiation with respect to ∂Q ia will remove M a in both cases.However, we have to account for the ensemble in the cavity contribution, i.e., the first term of the squared matrix element, W µa0a , which can be now written as i.e., the cavity intensity depends herein on the ensemble size due to its relation to static polarizabilities.

( 2 )
IR ( ω) related to bare molecular transitions of uncoupled CO 2 bending modes, which are absent in σ (C,2) IR a sum with index, a, running over all molecules in the ensemble.Now, the single-molecule matrix element in Eq.(S1.43)simply generalizes toW µa0a ≈ g λk Q ia , (S5.2)

4 )
The second line follows from the assumption of identical molecules in the ensemble.Thus, for a single-cavity mode model, the respective dressed frequency is then given to leading order in g 0 by ωc = ω 2 c to the main text up to a factor M .In line with the Hessian cavity block, we also have to account for ensemble effects in the cavity intensity.The ensemble cavity component of the generalized dipole moment reads here D the same identical-molecule argument as for the Hessian correction.Then the related derivative is simply given by ∂D

TABLE I .
Cartesian coordinates for CO2 and Fe(CO)5 equilibrium structures optimized on TPSSh/Def2TZVP level of theory.
.2)with single excited normal and cavity mode states, |1 as , 0 c and |0 as , 1 c , here illustratively for the CO 2 model.