# Lorentz symmetry violating Lifshitz-type field theories.

Emiliano Rizza<EMAIL_ADDRESS>Dipartimento di Fisica, Università
di Catania,
via S. Sofia 64, I-95123, Catania, Italy  Dario Zappalà
<EMAIL_ADDRESS>INFN, Sezione di Catania, via S. Sofia 64, I-95123,
Catania, Italy

###### Abstract

ABSTRACT

We discuss the ultraviolet sector of 3+1 dimensional Lifshitz-type anisotropic
higher derivative scalar, fermion and gauge field theories, with anisotropy
exponent z=3 and with explicit breaking of Lorentz symmetry. By discarding
from the action all momentum dependent vertex operators, which is essential to
avoid phenomenologically unacceptable deformations of the light cone, we find
that renormalizable scalar self-interaction and Yukawa-like couplings are, in
general, asymptotically free. However, the requirement of cancelling momentum
dependent vertex operators is incompatible with gauge symmetry and, therefore,
for this kind of theories, gauge symmetry as well as Lorentz symmetry are
recovered only as emergent properties below some energy scale $M$, that must
be constrained from experiments. The quantum corrections to the scalar mass
and their impact on the hierarchy problem are also analyzed.

##  I Introduction

The smoothening of ultraviolet (UV) divergences with the consequent
improvement of the renormalizability of a quantum field theory, induced by an
enhanced number of space and time derivatives of the field, is a long studied
mechanism, especially in the context of field theoretical formulation of
gravity. thirring ; Pais:1950za ; Stelle:1976gc ; deUrries:1998obu ;
Hawking:2001yt ; Collins1

As might be expected, the improved renormalization can be traced back to the
presence of a novel Renormalization Group fixed point and, in fact, fixed
points of higher derivative field theories are known in condensed matter as
Lifshitz points. Horn ; hornrev ; selke1988 ; Diehl ; Diehl3  They can be
classified as isotropic, if the number of derivatives is the same for every
coordinate of the Euclidean space under investigation, or as anisotropic
points, if there is an anisotropy between two sets of coordinates, consisting
in a different number of derivatives with respect to the coordinates of one
set and the number of derivatives of the other set.

When (a Wick-rotated version of) a Minkowskian quantum field is under
investigation, a manifestly Lorentz invariant action is typically expected,
mainly because astrophysical observations pose the energy scale of Lorentz
violating effects above $10^{10}$ GeVellis:2003 ; chenhuang  (or even above
$10^{17}$ GeV,Ellis:2018lca  if the specific theory investigated predicts a
correction to the dispersion relation that grows linearly with the energy).

Actually, Lorentz symmetric, isotropic Lifshitz points, possess a rich phase
structure and are potentially useful for the investigation of some high energy
physics issues. Bonanno:2014yia ; Zappala:2017vjf ; Zappala:2018khg ; Zapp ;
Defenu ; Buccio  However, isotropic higher derivative theories that contain
more than two time derivatives are affected by the Ostrogradski instability,
associated with violation of unitarity. deUrries:1998obu ; Woodard:2006nt
Therefore, such isotropic higher derivative theories can only be regarded as
low energy effective theories of some more fundamental high energy theory that
preserves unitarity, while the only Lifshitz points entitled to describe the
UV completion of a fundamental field theory are those with anisotropic
structure, generated by actions with two derivatives with respect to the time
coordinate and $2z$ derivatives with respect to the space coordinates (below,
we refer to $z$ as the anisotropy exponent). Clearly, this kind of actions do
not respect Lorentz symmetry and one must be careful that their predictions
are not in contradiction with the observational constraints and the symmetry
is recovered at low energy as an emergent property.

The original formulation of a renormalizable anisotropic gravitational theory,
known as Hořava-Lifshitz gravity horava , is constructed on the basis of these
motivations (for a recent version of this formulation see Barv1 ; Barv2 ;
Barv3 ; Barv4  ), and also anisotropic field theories on flat spacetime were
formulated on the same grounds. Anselmi:2007ri ; Iengo ; Dhar:2009dx ;
horava:ym ; eune ; alexandre ; Kikuchi ; Chao:2009tqf ; Solomon:2017nlh

In this Letter, we focus on $3+1$-dimensional field theories, to investigate
on the possibility of constructing a renormalizable, consistent and UV safe
completion of scalar, as well as fermion and gauge theories on flat spacetime.
In particular, we focus on the case with anisotropy exponent $z=3$, rather
than $z=2$, in order to achieve an enhanced suppression of the UV loop
divergences from the modified propagators. As discussed below, for a
$3+1$-dimensional system with $z=3$, the scaling dimension of the scalar field
vanishes and becomes negative for $z>3$. Therefore, we select $z=3$ to avoid
theories with negative scaling dimension of the field.

More in detail, we generalize to non-zero spin fields the approach developed
in zappalast  for scalar fields. In fact, a particular restriction on the
scalar action is enforced in zappalast  to overcome a serious inconvenience,
outlined in Iengo ; alexandre , generated by a class of operators that produce
dissimilar corrections to the effective light cones seen by different
particles, when these interact. This unacceptable effect is actually avoided
by the restriction discussed in zappalast  and, below, we review its
application to the case of scalar fields in Section II, while the case of
fermion and gauge fields is discussed respectively in Sections III and IV.
Conclusions are reported in Section V.

##  II Scalar Fields

The scaling dimension of a scalar field in a 3+1 dimensional Lifshitz-type
theory and generic $z$, is easily read from the derivative sector of the
action

$S_{D}=\int\;d^{3}x\,dt\left\\{\frac{1}{2}\left(\partial_{t}\phi\right)^{2}-\sum_{k=1}^{z}\,\sum_{i=1}^{3}\frac{a_{k}}{2\,M^{2k-2}}\left(\,[\partial_{i}]^{k}\phi\right)^{2}\right\\}\,,$
(1)

by requiring that, under scale transformations with scale factor $b>1$, the
space-time coordinates transform as $x^{i}\to bx^{i}$ and $t\to b^{z}t$, which
defines the following scaling dimensions $[x^{i}]_{s}=-1$ and $[t]_{s}=-z$
and, accordingly, the scaling dimension of the field,
$[\phi]_{s}=(3-z+\eta)/2$, and that of the coefficient of the highest
derivative operator, $[a_{z}]_{s}=-\eta$ (in these relations we also took into
account a non-zero anomalous dimension of the field, $\eta$). Suitable powers
of a mass scale, $M$ are inserted in Eq. (1) to keep track of the canonical
dimensions of the various operators and couplings, so that all $a_{i}$ are
dimensionless pure numbers.

Clearly, if $z=3$ and $\eta$ is neglected, one finds $[\phi]_{s}=0$ and this
determines the structure of the interaction sector, which, added to $S_{D}$,
represents the complete renormalizable action:

$S_{I}=-\,\int\;d^{3}x\,dt\left\\{\sum_{n=2}^{\infty}\frac{g_{n}\,\phi^{n}}{n!\,M^{(n-4)}}+\sum_{k=1}^{3}\,\sum_{i=1}^{3}\left[\sum_{m=1}^{\infty}\frac{w_{m,k}\,\phi^{m}}{M^{2k-2+m}}\right]\;\left(\,[\partial_{i}]^{k}\phi\right)^{2}\right\\}\,.$
(2)

Since $[\phi]_{s}=0$ , the potential term, which corresponds to the first sum
in curly brackets in (2), must include all powers of $\phi$. For the same
reason, the complete renormalizable action must contain the second block in
curly brackets, which is constructed by multiplying the derivative part sector
in $S_{D}$ by any power of $\phi$, while any other term containing higher
powers of the derivatives is, by construction, non-renormalizable. The powers
of the scale $M$ are again selected so that all couplings $g_{n}$ and
$w_{m,k}$ have zero canonical dimension.

As pointed out in Iengo  and, later, analyzed again in alexandre , loop
integrals containing momentum dependent vertices proportional to $w_{m,k}$,
produce corrections to the coupling $a_{1}$, which is related to the ratio of
space and time coordinates and therefore defines the value of the speed of
light. However, in the case of two interacting fields, these corrections to
$a_{1}$ turn out to be different for each field, and one concludes that the
excitations of the two fields have different light cones.

To avoid this flaw, instead of resorting to a quite unnatural ”ad hoc” tuning
of the various $w_{m,k}$, it is suggested to restrict the form of the action
and discard all momentum dependent vertices, by taking $w_{m,k}=0$ for any $m$
and $k$, even though renormalizable.zappalast  In fact, not only the
suppression of all momentum dependent vertices avoids the dangerous
corrections to the speed of light, but one also finds that such vertices are
not generated by quantum corrections associated to the other remaining
vertices in the action. To illustrate this point we analyze the divergences of
the quantum corrections produced by the action $S_{D}+S_{I}$ with $z=3$ and
all $w_{m,k}=0$.

Figure 1: List of divergent diagrams: (a) scalar tadpole; (b) fermionic
tadpole ; (c) two-vertex fermionic loop, dressed by an arbitrary number of
internal scalar propagators.

The degree of divergence of a loop diagram $D_{\Lambda}$ is easily computed:
zappalast

$D_{\Lambda}=6\,\left(1-\sum_{n}V_{n}\right)\,,$ (3)

where $V_{n}$ is the number of vertices with $n$ legs (proportional to
$g_{n}$). We notice that $D_{\Lambda}$ does not depend on the number of
external legs and it is always negative unless the diagram has one vertex
only, $V_{\widehat{n}}=1$ and $V_{n}=0$ for $n\neq\widehat{n}$, regardless the
specific value of $\widehat{n}$. Then, one concludes that the only primitively
divergent scalar diagrams are those with one vertex, an arbitrary number of
external legs, and an arbitrary number of tadpoles, all insisting on the same
vertex. As an example, the diagram with one tadpole and five external legs is
displayed in Fig. 1 and labelled with (a).

The first consequence of this analysis is that loop diagrams with two or more
vertices are not divergent and, since only diagrams with more than one vertex
can carry dependence on some external momentum, it follows that quantum
corrections depending on external momenta must be finite. In particular,
corrections to the parameters $a_{k}$ in Eq. (1) come from external momentum
dependent two-point diagrams and, therefore, they are finite. An explicit
calculation in zappalast , shows that these corrections, at least at energy
scales above $M$, are certainly negligible.

The second consequence comes from the observation that also all operators
proportional to $w_{m,k}$ get quantum corrections from diagrams that depend on
some external momentum and therefore, when the constraint $w_{m,k}=0$ is
enforced in (2), no significant divergent correction to this constraint is
generated. Then, the dangerous corrections to the speed of light are avoided
in a consistent way.

The tadpole (a) in Fig. 1, gives the elementary UV divergence. It is a one
loop diagram that, following zappalast , can be regularized by adopting the
non-Lorentz invariant cutoff on the modulus of the tri-momentum, which is
consistent with the non-Lorentz invariant structure of the action in Eq. (1).
Since we are interested in the UV properties of the theory, we integrate the
tadpole between two momentum scales $\Lambda_{1}$ and $\Lambda_{2}$, with
$\Lambda_{1}>>\Lambda_{2}>>M$ and we assume $M$ to be much larger than any
scale appearing in Eqs. (1) and (2) that, consequently, will be neglected.

Then, for the tadpole (a), we get:

$I_{1}(\Lambda_{1},\Lambda_{2},M)=\frac{M^{2}}{\sqrt{a_{3}}\,(2\pi)^{2}}\ln\left(\frac{\Lambda_{1}}{\Lambda_{2}}\right)+O\left(\frac{M^{4}}{\Lambda_{2}^{2}},\,\frac{M^{2}\Lambda_{2}^{2}}{\Lambda_{1}^{2}}\right)\,.$
(4)

where, without any loss of generality, we can put $a_{3}=1$. zappalast

If we now compute the quantum corrections of a generic coupling $g_{n}$, the
combinatorial factors of each contributing diagram always yield the following
series (which corresponds to the exponential series, if the various couplings
are neglected):

$g_{n}(\Lambda_{2})=\sum_{m=0}^{\infty}\frac{g_{{}_{n+2m}}(\Lambda_{1})}{m!}\;\left[\frac{I_{1}(\Lambda_{1},\Lambda_{2},M)}{2M^{2}}\right]^{m}$
(5)

The series in (5) is simplified if one assumes that the theory has only one
coupling, i.e. $g_{n}=g$ for all $n$. In fact, in this case by using Eq. (4),
one finds

$g(\Lambda_{2})=g(\Lambda_{1})\;\exp\left[\frac{I_{1}(\Lambda_{1},\Lambda_{2},M)}{2M^{2}}\right]\simeq
g(\Lambda_{1})\;\left(\frac{\Lambda_{1}}{\Lambda_{2}}\right)^{1/\left(8\pi^{2}\right)}$
(6)

and it is straightforward to check that $g(\Lambda_{1})\to 0$ in the limit
$\Lambda_{1}\to\infty$, at fixed $g(\Lambda_{2})$, i.e. the coupling $g$ is
asymptotically free.

This remarkable result can be easily generalized to the case of different
$g_{n}$ if, at the scale $\Lambda_{1}$, the set of all $g_{n}$ is bounded by
$0\leq g_{n}\leq g$, where $g$ is a real number. Then, the exponential series
in (6) is greater than the series (5) term by term, and it is straightforward
to conclude $g_{n}(\Lambda_{1})\to 0$ for each $n$, when
$\Lambda_{1}\to\infty$. Conversely, if the potential term in Eq. (1) contains
only a finite set of non-zero couplings $g_{n}\neq 0$, the theory is
renormalizable, but no asymptotically free coupling is found. zappalast

##  III Fermionic Fields

In case of fermionic degrees of freedom, the higher derivative part of the 3+1
dimensional action, with $z=3$, is

$S_{F}=\\!\int\\!d^{3}x\,dt\;\bar{\psi}\left[i\gamma^{0}\partial_{0}-\left(b_{1}+\,\frac{\partial_{j}\partial^{j}}{M^{2}}\right)(i\gamma^{i}\partial_{i})-m_{f}\right]\,\psi$
(7)

where the summed indices $j$ and $i$ refer to the 3 space coordinates and, as
for the scalar case, we normalize $S_{F}$ by setting $b_{3}=1$ in front of the
higher derivative term. To preserve rotational invariance, the higher
derivative part in brackets does not contain a $b_{2}$ term (analogous of the
$a_{2}$ term in (1) ), although the essential terms, linear and cubic in the
space derivatives, are included in (7). The linear term in the time derivative
sets the scaling dimension of the fermion field, $[\psi]_{s}=3/2$, which,
unlike the scalar field, turns out to be equal to its canonical dimension.

Consequently, if we assume a Yukawa-like interaction with the scalar
excitations, the renormalizable interaction sector of the action is

$S_{Y}=-\\!\int\\!d^{3}x\,dt\;\sum_{n=1}^{\infty}\,y_{n}\;\frac{\bar{\psi}\psi\;\phi^{n}}{n!\;M^{n-1}}\,,$
(8)

where all powers of the scalar field are taken into account because
$[\phi]_{s}=0$ and, therefore, all these terms are renormalizable.

The analysis of the degree of divergence of the diagrams generated by $S_{Y}$
yields

$D_{\Lambda}=6-3\,\sum_{n}Y_{n}-\frac{3}{2}E_{f}$ (9)

where $Y_{n}$ is the number of vertices proportional to $y_{n}$ and $E_{f}$ is
the number of external fermionic legs. Then $D_{\Lambda}$ is non-negative if
$E_{f}=2$ and $\sum_{n}Y_{n}=1$, or $E_{f}=0$ and $\sum_{n}Y_{n}=1$, or
$E_{f}=0$ and $\sum_{n}Y_{n}=2$. In the first case, when $E_{f}=2$, diagrams
have, besides the two external fermionic legs, all possible multiple scalar
tadpoles insisting on the same vertex. In the second case, with $E_{f}=0$,
diagrams display the fermionic tadpole (b) in Fig. 1, together with all
possible multiple scalar tadpoles insisting on the same vertex. In the third
case, again with $E_{f}=0$, diagrams are of the (c) kind in Fig. 1, where not
only each vertex can be dressed with all possible multiple scalar tadpoles as
before, but also an arbitrary number of internal scalar propagators can link
the two vertices. In all three cases, an arbitrary number of scalar external
legs coming out of each vertex is allowed. Then, in the diagrams (b) and (c)
the vertices indicated by thick black dots, stand for the full sum of scalar
tadpoles together with the arbitrary number of scalar external legs.

By adopting the same regularization procedure employed for scalar fields, the
fermionic tadpole (b) gives the negative result $-4m_{f}I_{1}$, and the most
divergent part of the diagram (c), with $r\neq 0$ internal scalar propagators,
is $-2I_{1}^{r+1}$ while the pure fermionic one loop diagram (c) without
scalar propagators gives $-4I_{1}$. Therefore, in all cases the divergences
amount to powers of the logarithm of the cutoff.

It is now possible to analyze the renormalization of the couplings and we
consider the Yukawa sector in the simple case in which all couplings are
equal: $y=y_{n}$. For our purposes, together with $a_{3}=b_{3}=1$, which is
justified by the negligible running of these couplings in the region above
$M$, we can also take $a_{1}=a_{2}=b_{1}=0$ without loss of generality. The
sum of all diagrams contributing to the renormalization of $y$, integrated as
before between the two scales $\Lambda_{1}$ and $\Lambda_{2}$, gives

$y(\Lambda_{2})=y(\Lambda_{1})\;\exp\left[\frac{I_{1}(\Lambda_{1},\Lambda_{2},M)}{2M^{2}}\right]\equiv
y(\Lambda_{1})\;{\cal E}$ (10)

and this leads to the same conclusions discussed after Eq. (6), i.e. the
coupling $y$ is asymptotically free.

The renormalization of the scalar sector is more involved, because of the
diagrams in (c). In fact, by keeping different $g_{n}$ for different $n$, and
by retaining the leading divergences only, we find the following
generalization of Eq. (5)

$g_{n}(\Lambda_{2})=\sum_{m=0}^{\infty}\frac{g_{{}_{n+2m}}(\Lambda_{1})}{m!}\;\left[\frac{I_{1}(\Lambda_{1},\Lambda_{2},M)}{2M^{2}}\right]^{m}-c_{n}\;y(\Lambda_{1})^{2}\;{\cal
E}^{4}I_{1}$ (11)

where $c_{n}>0$ is associated with the combinatorial weight of the diagrams in
(c) and its value changes with $n$. Then, also the contribution to $g_{n}$
changes with $n$. Therefore, if all $g_{n}$ were equal to a single coupling
$g$ as for the case of the Yukawa coupling $y$, it would be impossible to
renormalize the full structure in Eq. (11).

In order to determine the $\beta$-function of $g_{n}$, it is sufficient to
take the derivative with respect to $\Lambda_{2}$ in Eq. (11), by keeping
fixed the values of all couplings at $\Lambda_{1}$, and it is easy to realize
that the two terms in the right hand side (rhs) of (11) give respectively a
negative and positive contribution, so that the sign of the $\beta$-function
depends on the balance of these contributions and, consequently, on the
particular values of the couplings, so that it is not uniquely established.

However, the following heuristic argument can be put forward: the stability of
the scalar potential requires $g_{n}\geq 0$, at least for all $n$ larger than
some $\bar{n}$. Then, for $n\geq\bar{n}$ the first term in the rhs of (11)
must be larger than the second. In addition, we noticed that the coupling $y$
is asymptotically free and it grows when the momentum scale is lowered and, in
case the $g_{n}$ are not sufficiently large at the scale $\Lambda_{1}$ to
ensure a negative $\beta$-function, $g_{n}$ decreases when the momentum scale
runs from $\Lambda_{1}$ towards $\Lambda_{2}$. Then, if these two scales are
sufficiently distant, $g_{n}$ could become negative at the infrared scale
$\Lambda_{2}$, yielding an unstable potential. This is avoided if the values
of $g_{n}(\Lambda_{1})$ are sufficiently large to ensure a negative
$\beta$-function, at least for all $g_{n}$ with $n\geq\bar{n}$.

It must be noticed that, in the analysis of the fermionic sector we
deliberately did not include 2n-fermion vertices. Actually, it is easy to show
by means of the analysis of the degree of divergence of diagrams, that the
only vertex of this kind that is renormalizable, is the 4-fermion vertex.
However, this vertex can always be reduced, through Hubbard-Stratonovich
transformation to a Yukawa-like interaction and, therefore, it should not be
treated as a fundamental interaction. In any case, the 4-fermion interaction
was already studied, Dhar:2009dx  and we do not expect that its inclusion
could substantially modify our conclusions.

A comment on the scalar square mass $g_{2}\;M^{2}$ ( as defined in (2) ),
whose correction is given in (11) for $n=2$, is in order. It is evident that,
unless a very subtle cancellation between the first and second term in the rhs
of (11) occurs, we get a finite non-vanishing correction to the scalar mass
from the momentum region between $\Lambda_{1}$ and $\Lambda_{2}$, with
$\Lambda_{1}\to\infty$. The correction is, in any case, proportional to the
scale $M^{2}$. This is analogous to the result obtained for the scalar case in
zappalast , and even if, as expected, the fermionic contribution to the
correction occurs here with opposite sign with respect to the scalar, an exact
cancellation of the two would be quite unnatural.

##  IV Gauge Fields

The inclusion of gauge fields can be realized by defining the appropriate
covariant derivative, which can be constructed by starting from the
generalized derivative $\widehat{\partial}_{\mu}$ already introduced in Eq.
(7) (here for simplicity, the coefficients of the various derivative terms,
analogous of $b_{1},\,b_{3}$, are taken equal to 1)

$\widehat{\partial}_{0}=\partial_{0}\;,\;\;\;\;\;\;\;\;\widehat{\partial}_{i}=\left(1+\,\frac{\partial_{j}\partial^{j}}{M^{2}}\right)\partial_{i}\;,$
(12)

and by including the gauge field in the usual way :
$\widehat{D}_{\mu}=\widehat{\partial}_{\mu}-iqA_{\mu}$. Accordingly, the
generalized electromagnetic tensor is
$\widehat{F}_{\mu\nu}=\widehat{\partial}_{\mu}A_{\nu}-\widehat{\partial}_{\nu}A_{\mu}$.
However, it is straightforward to realize that $\widehat{D}_{\mu}$, which is
non-linear in the derivatives, does not behave properly under gauge
transformations, so that the corresponding higher derivative electrodynamics
produces gauge violating contributions to the amplitudes, proportional to
powers of $p^{2}/M^{2}$, where $p^{2}=p^{\mu}p_{\mu}$ is some square momentum
scale, representative of the amplitude.

A different higher derivative formulation, that is gauge invariant, was
already proposed long agohorava:ym  but, in this case, the higher order space
derivatives appear in the form $(D_{h}D_{j}F_{ik})^{2}$ (where $D_{i}$ are the
space components of the standard covariant derivative), which contains terms
analogous to those proportional to $w_{m,k}$ in Eq. (2) for the scalar case.
As discussed above, in order to avoid model dependent modifications of the
light cone, we required to reject these kind of terms and, consequently, this
formulation cannot be accepted.

If we admit that, together with Lorentz symmetry, also gauge invariance is
violated, at least at high energy, and only recovered at low energies below
$M$, we can make use of the definition in (12) and of the corresponding
$\widehat{D}_{\mu}$, to include the dynamics of the field $A_{\mu}$. Then,
independently of gauge symmetry, the effect of the dynamics of $A_{\mu}$ in
Eq. (11) is, from an effective point of view, equivalent to the presence of
additional scalar degrees of freedom, with the consequence of modifying only
some coefficients in (11) but not the overall structure of the equation, so
that the conclusions discussed above are not qualitatively altered.

##  V Conclusions

We analyzed a realization of a higher derivative theory in 3+1 dimensions with
anisotropy exponent $z=3$ and, in particular, we considered the possible
generalization to fermionic and vector gauge fields of the specific
formulation for scalar fields studied in zappalast , with the main intent of
constructing a consistent UV completion of these fields and, ultimately, of
the Standard Model. Therefore, we focused on the high energy region that
corresponds to energies much larger than the scale $M$, which sets the scale
of Lorentz violating effects.

In fact, as evident from Eqs. (1) and (7), all modifications to the standard
Lorentz invariant dispersion relations of scalars and fermions, correspond to
Lorentz violating additive corrections proportional to powers of the ratio
$(p^{2}/M^{2})$, which become relevant only for $p^{2}>>M^{2}$. In addition,
an acceptable crossover to the region below $M^{2}$, where the usual
renormalization group flow of the various couplings is recovered and Lorentz,
as well as gauge symmetry are effectively restored, has been realized and
discussed in zappalast . As already mentioned, astrophysical high energy gamma
ray observations put a strong limit on $M$, that, for models with corrections
proportional to $(p^{2}/M^{2})$ as in our case, is about $M\geq 10^{10}$
GeV,ellis:2003 ; chenhuang  a very large value if compared with the
characteristic scales of the Standard Model.

From the point of view of renormalizability, there is a great improvement due
to the higher derivative terms. In particular, the scalar self-interacting
sector has a very simple UV structure that indicates that the couplings
$g_{n}$ are asymptotically free, at least when, at the scale $M$, all
couplings are positive and the set of all $g_{n}$ is bounded from above; in
addition, if $g_{n}=g$ for all $n$, the series of significant diagrams, due to
their combinatorial weights, sums up to a simple exponential form.

The interacting system of scalars and fermions reproduces the same structure
of divergences in the Yukawa sector with asymptotic freedom of the
corresponding coupling, while for the scalar couplings the pattern of
divergences becomes more complicated. However, under the hypothesis of
stability of the potential, also scalar couplings turn out to be
asymptotically free. On the other hand, we are not able to include gauge
fields without introducing interaction operators that produce physically
unacceptable modifications of the light cones for different species of
particles. Therefore in this approach, gauge symmetry can only be recovered as
an emergent low energy symmetry (below $M$).

Then, from the phenomenological point of view, if the violation of both
Lorentz and gauge symmetry could be kept under control by the large value of
the lower bound on $M$, the price to pay is that the strong suppression of
divergences is limited to the UV region above $M$, while standard scaling
holds at scales smaller than $M$. In particular, the correction to the scalar
square mass, although finite, is proportional to $g_{4}M^{2}$, that is very
large if compared for instance to the Higgs square mass. In this sense, the
positive features of this kind of theory, seem to have a very weak impact on
the physics of the Standard Model at present energies.

## Acknowledgments

DZ is grateful to A. Bonanno, V. Branchina and R. Percacci for useful comments
and suggestions. This work has been carried out within the INFN project FLAG.

## References

  * (1) W. Thirring, Phys. Rev. 77, 570 (1950).
  * (2) A. Pais and G. E. Uhlenbeck, Phys. Rev. 79, 145 (1950).
  * (3) K. S. Stelle, Phys. Rev. D 16, 953 (1977).
  * (4) F. J. de Urries and J. Julve, J. Phys.  A 31, 6949 (1998). [hep-th/9802115].
  * (5) S. W. Hawking and T. Hertog, Phys. Rev.  D 65, 103515 (2002). [hep-th/0107088].
  * (6) J. Collins, A. Perez, D. Sudarsky, L. Urrutia, and H. Vucetich, Phys. Rev. Lett. 93, 191301 (2004). [arXiv:gr-qc/0403053].
  * (7) R. Hornreich, M. Luban, and S. Shtrikman, Phys.Rev.Lett. 35, 1678 (1975).
  * (8) R. M. Hornreich, Journal of Magnetism and Magnetic Materials 15, 387 (1980).
  * (9) W. Selke, Physics Reports 170, 213 (1988).
  * (10) H. Diehl, Acta Phys.Slov. 52, 271 (2002).
  * (11) H. Diehl and M. Shpot, J.Phys. A35, 6249 (2002). [cond-mat/0204267].
  * (12) J. Ellis, N. Mavromatos, D.V. Nanopoulos, A.S. Sakharov Astron. Astrophys. 402, 409 (2003). [arXiv:astro-ph/0210124].
  * (13) B. Chen, Q.G. Huang, Phys. Lett. B 683, 108 (2010).
  * (14) J. Ellis, R. Konoplich, N. E. Mavromatos, L. Nguyen, A. S. Sakharov, and E. K. Sarkisyan-Grinbaum, Phys. Rev. D 99, 083009 (2019). [arXiv:1807.00189].
  * (15) A. Bonanno and D. Zappala, Nucl. Phys.  B 893, 501 (2015). [arXiv:1412.7046].
  * (16) D. Zappala, Phys. Lett.  B 773, 213 (2017). [arXiv:1703.00791].
  * (17) D. Zappala, Phys. Rev. D 98, 085005 (2018). [arXiv:1806.00043].
  * (18) D. Zappala, Int. J. Geom. Meth. Mod. Phys. 17, 2050053 (2020). [arXiv:1912.03071].
  * (19) N. Defenu, A. Trombettoni, and D. Zappala, Nucl. Phys. B 964, 115295 (2021). [arXiv:2003.04909].
  * (20) D. Buccio and R. Percacci, Renormalization Group flows between Gaussian Fixed Points, [arXiv:2207.10596].
  * (21) R. P. Woodard, Lect. Notes Phys. 720, 403 (2007). [astro-ph/0601672].
  * (22) P. Horava, Phys. Rev. D 79, 084008 (2009). [arXiv:0901.3775].
  * (23) A.O. Barvinsky, D. Blas, M. Herrero-Valea, S.M. Sibiryakov, and C.F. Steinwachs, Phys. Rev.  D 93, 064022 (2016). [arXiv:1512.02250].
  * (24) A.O. Barvinsky, D. Blas, M. Herrero-Valea, S.M. Sibiryakov, and C.F. Steinwachs, Phys. Rev. Lett. 119, 211301 (2017). [arXiv:1706.06809].
  * (25) A.O. Barvinsky, D. Blas, M. Herrero-Valea, and S.M. Sibiryakov, Phys. Rev.  D 100, 026012 (2019). [arXiv:1905.03798].
  * (26) A.O. Barvinsky, A.V. Kurov, and S.M. Sibiryakov, Phys. Rev.  D 105, 044009 (2022). [arXiv:2110.14688].
  * (27) D. Anselmi and M. Halat, Phys. Rev. D 76,125011 (2007). [arXiv:0707.2480].
  * (28) R. Iengo, J. G. Russo, and M. Serone, JHEP 11, 020 (2009). [arXiv:0906.3477].
  * (29) A. Dhar, G. Mandal, and S. R. Wadia, Phys. Rev. D 80, 105018 (2009). [arXiv:0905.2928].
  * (30) P. Horava, Phys. Lett. B 694, 172 (2011). [arXiv:0811.2217].
  * (31) M. Eune, W. Kim, and E. J. Son, Phys. Lett. B 703, 100 (2011). [arXiv:1105.5194].
  * (32) J. Alexandre, Int. J. Mod. Phys. A 26, 4523 (2011). [arXiv:1109.5629].
  * (33) K. Kikuchi, Prog. Theor. Phys. 127, 409 (2012). [arXiv:1111.6075].
  * (34) W. Chao, Commun. Theor. Phys. 65, 743 (2016). [arXiv:0911.4709].
  * (35) A. R. Solomon and M. Trodden, JCAP 02, 031 (2018). [arXiv:1709.09695].
  * (36) D. Zappala, Eur. Phys. J. C 82, 341 (2022). [arXiv:2111.08385].


# Cavity-modified electronic interactions in molecular ensembles under
vibrational strong coupling: Combined insights from cavity Born-Oppenheimer
perturbation theory and ab initio wave function theory

Eric W. Fischer<EMAIL_ADDRESS>Institut für Chemie, Humboldt-
Universität zu Berlin, Brook-Taylor-Straße 2, D-12489, Berlin, Germany

###### Abstract

Resonant vibrational strong coupling (VSC) between molecular vibrations and
quantized field modes of low-frequency optical cavities constitutes the
conceptual cornerstone of vibro-polaritonic chemistry. In this work, we
theoretically investigate complementary nonresonant interactions between
electrons and low-frequency cavity modes with a focus on cavity-induced
modifications of electronic interactions in molecular ensembles under VSC.
Methodologically, we combine cavity Born-Oppenheimer perturbation theory (CBO-
PT) [Fischer, Saalfrank, JCTC 19, 7215 (2023)] with ab initio wave function
approaches, specifically, CBO Hartree-Fock and coupled cluster theories. In a
first step, we derive cavity-induced energy corrections to intra- and
intermolecular electronic energies based on second-order CBO-PT. We find local
modifications of molecular PES for selected isomerization reactions to be
small and dominantly captured by the first-order dipole-fluctuation
correction. Excellent agreement between CBO-PT and non-perturbative methods is
obtained, which indicates potentially negligible state relaxation effects
under VSC. Secondly, we examine cavity-modified intermolecular interactions
revealing non-trivial corrections of dipole-dipole, dipole-induced dipole and
van-der-Waals type character emerging at second-order of CBO-PT. Besides an
explicit coupling of interacting dimers to cavity modes, we find a cavity-
polarization dependent enhanced long-range component of dipole-induced dipole
and van-der-Waals interactions. We present an illustrative numerical analysis
of interaction potentials for selected molecular dimer models based on the
CBO-CCSD approach, which provides a non-perturbative perspective on cavity-
modified intermolecular interactions under VSC.

##  I Introduction

The emerging field of vibro-polaritonic chemistry relies on resonant
vibrational strong coupling (VSC) between infrared active molecular
vibrational modes and confined light modes of an infrared Fabry-Pérot
cavity.[1, 2, 3] The VSC regime is characterized by light-matter hybrid states
known as vibrational polaritons[4], which have been experimentally reported to
profoundly impact vibrational spectroscopy[5, 6, 7] and more intriguingly
ground state chemistry[8, 9]. In this work, we theoretically investigate
complementary nonresonant interactions between low-frequency cavity modes and
electrons under VSC with a focus on cavity-induced modifications of electronic
interactions.

A theoretical framework consistently describing both resonant VSC and
nonresonant electron-photon interactions is provided by the cavity Born-
Oppenheimer (CBO) approach.[10, 11, 12] The interacting light-matter many-body
system is here decomposed into two coupled subsystems containing high energy
(“fast”) electrons and low energy (“slow”) nuclei and cavity modes. The
electronic subsystem is described by an CBO electronic Hamiltonian, which
accounts for nonresonant electron-photon interactions, while the related
nuclear-cavity problem fully accounts for VSC.[11] An adiabatic ground state
theory is obtained via the CBO approximation, which straightforwardly
generalizes ideas from molecular quantum mechanics to the strong light-matter
coupling scenario.[11, 12]

A valid CBO approximation indicates weak electron-photon interactions, which
we recently exploited to introduced a perturbative approach to the CBO
electronic ground state problem named cavity Born-Oppenheimer perturbation
theory (CBO-PT).[12, 13] We employ CBO-PT here to investigate cavity-modified
electronic interactions in dipole-coupled molecular ensembles under VSC. This
work extends recent related results obtained for the electronic strong
coupling (ESC) regime and its low-frequency limit.[14, 15] Methodologically,
we combine CBO-PT with ab initio wave function theory, specifically CBO
Hartree-Fock (CBO-HF)[16] and CBO coupled cluster theory (CBO-CC)[17], which
provide non-perturbative numerical approaches to the CBO electronic ground
state problem.

We first address cavity-induced intramolecular modifications of electronic
properties under VSC based on an intramolecular cavity reaction potential
(CRP)[18]. The CRP minimizes energy in the cavity subspace and therefore
provides access to cavity-induced corrections of stationary points on
molecular potential energy surfaces. From a combined analysis employing CBO-PT
and CBO wave function theory, we find intramolecular cavity-induced energy
corrections to be dominantly captured by the first-order dipole fluctuation
correction, which is small in magnitude relative to classical activation
energies. Moreover, we find negligibly small state relaxation effects under
VSC at both mean-field and correlated levels of wave function theory for the
studied examples.

In a second step, we investigate cavity-induced modifications of
intermolecular interactions and derive the corresponding CRP for an
interacting molecular dimer under VSC. From CBO-PT, we find cavity-induced
modifications of dipole-dipole, dipole-induced dipole and van-der-Waals
interactions with an enhanced long-range character partially found already for
the ESC scenario[14, 15]. We obtain non-trivial corrections to intermolecular
interactions only at second-order of CBO-PT in contrast to the locally
dominant first-order dipole fluctuation term. We provide a complementary
numerical analysis of interaction potentials for selected polar and nonpolar
molecular dimers on a CBO-CCSD level of theory. In contrast to local
properties discussed in the monomer example, we observe cavity-induced
modifications of intermolecular interactions not to be necessarily small. A
related relevance in the VSC regime is given by a potential yet unclear
connection to collective VSC effects.

The paper is structured as follows. In Sec.II, we introduce the CBO framework
besides the CBO-PT approach and provide up to second-order corrected intra-
and intermolecular potentials under VSC. In Sec.III, we introduce the idea of
a cavity reaction potential (CRP) followed by a brief overview of non-
perturbative CBO-HF and CBO-CCSD approaches in Sec.IV. In Secs.V andVI, we
provide detailed discussions of intramolecular and intermolecular interactions
under VSC by combining CBO-PT and ab initio methods. Sec.VII summarizes and
concludes this work.

##  II Perturbation Theory

###  II.1 The CBO Ground State Problem

We consider a light-matter hybrid system composed of molecules strongly
coupled to effective transverse field modes of a low-frequency optical cavity.
The corresponding Pauli-Fierz Hamiltonian in length gauge representation and
long-wavelength approximation reads[11]

$\displaystyle\hat{H}$
$\displaystyle=\hat{T}_{c}+\hat{T}_{n}+\hat{H}_{ec}\quad,$ (1)

with cavity and nuclear kinetic energy operators (KEO), $\hat{T}_{c}$ and
$\hat{T}_{n}$. The CBO electronic Hamiltonian, $\hat{H}_{ec}$, describes the
electronic subsystem in presence of both nuclei and the quantized cavity
radiation field, and enters the corresponding time-independent Schrödinger
equation (TISE)[11]

$\displaystyle\hat{H}_{ec}\ket{\Psi^{(ec)}_{\nu}(\underline{R},\underline{x})}$
$\displaystyle=E^{(ec)}_{\nu}(\underline{R},\underline{x})\ket{\Psi^{(ec)}_{\nu}(\underline{R},\underline{x})}\quad.$
(2)

The latter provides access to cavity potential energy surfaces (cPES),
$E^{(ec)}_{\nu}(\underline{R},\underline{x})$, and corresponding adiabatic
electronic states, $\ket{\Psi^{(ec)}_{\nu}(\underline{R},\underline{x})}$,
which depend on both nuclear and cavity displacement coordinates,
$\underline{R}$ and $\underline{x}$, respectively. In the following, we
restrict the discussion to the ground state cPES, $E^{(ec)}_{0}$, which
provides a potential for the strongly coupled nuclear-cavity subsystem[11, 12]

$\left(\hat{T}_{n}+\hat{T}_{c}+E^{(ec)}_{0}(\underline{R},\underline{x})\right)\chi^{(nc)}_{i0}(\underline{R},\underline{x})\\\
+\sum_{\mu>0}\hat{\mathcal{K}}^{(nc)}_{0\mu}\,\chi^{(nc)}_{i\mu}(\underline{R},\underline{x})=\varepsilon_{i}\,\chi^{(nc)}_{i0}(\underline{R},\underline{x})\quad,$
(3)

with (ro)vibro-polaritonic states,
$\chi^{(nc)}_{i0}(\underline{R},\underline{x})$, energies, $\varepsilon_{i}$,
and non-adiabatic coupling elements, $\hat{\mathcal{K}}^{(nc)}_{0\mu}$. In
analogy to the Born-Oppenheimer framework in molecular quantum mechanics,
TISEs (2) and (3) describe high energy (“fast”) and low energy (“slow”)
degrees of freedom. Eq.(2) relates to the dressed electronic subsystem subject
to nonresonant interactions with low-frequency cavity modes, which will be the
focus of this work, whereas Eq.(3) captures resonant VSC between cavity modes
and rovibrational degrees of freedom. We now assume that the ground state cPES
is energetically well separated from the adiabatic excited state manifold

$\displaystyle\Delta^{(ec)}_{1}\gg\Delta^{(nc)}_{1}\quad,$ (4)

with $\Delta^{(ec)}_{1}=E^{(ec)}_{1}-E^{(ec)}_{0}$ and
$\Delta^{(nc)}_{1}=\varepsilon_{1}-\varepsilon_{0}$, which motivates the
cavity Born-Oppenheimer approximation[11, 12]

$\displaystyle\hat{\mathcal{K}}^{(nc)}_{0\mu}\propto(\Delta^{(ec)}_{\mu})^{-1}\approx
0\quad,\quad\mu\geq 1\quad.$ (5)

An appealing consequence of the CBO approximation relates to the observation
that interactions between electrons and low-frequency cavity modes are weak.
This is in close analogy to the molecular BO approximation, which renders
respective vibronic interactions weak. Accordingly, we can exploit a
perturbative approach to investigate cavity-induced modifications of
electronic interactions in molecular ensembles under VSC, specifically, the
recently formulated CBO perturbation theory (CBO-PT)[12].

###  II.2 CBO Perturbation Theory for Molecular Ensembles

####  II.2.1 The CBO Electronic Hamiltonian for Dipole-Coupled Molecular
Ensembles

In the CBO approximation, the ground state problem for an ensemble of $M$
dipole-coupled molecules interacting with the quantized field of a low-
frequency optical cavity is described by an ensemble CBO electronic
Hamiltonian

$\displaystyle\hat{H}_{ec}$
$\displaystyle=\hat{H}^{(e)}_{\mathrm{ens}}+V_{c}+\hat{V}_{\mathrm{ens}}+\hat{W}_{\mathrm{ens}}\quad.$
(6)

The first term simply constitutes the electronic Hamiltonian of the non-
interacting molecular ensemble

$\displaystyle\hat{H}^{(e)}_{\mathrm{ens}}$
$\displaystyle=\sum^{M}_{a}\hat{H}^{(e)}_{a}\quad,$ (7)

with single-molecule electronic contributions

$\displaystyle\hat{H}^{(e)}_{a}$
$\displaystyle=\hat{T}^{(e)}_{a}+\hat{V}^{(ee)}_{a}+\hat{V}^{(en)}_{a}+\hat{V}^{(nn)}_{a}\quad,$
(8)

determined by the electronic KEO, $\hat{T}^{(e)}_{a}$, besides Coulomb
potentials for electron-electron, $\hat{V}^{(ee)}_{a}$, electron-nuclei,
$\hat{V}^{(en)}_{a}$, and inter-nuclear interaction, $\hat{V}^{(nn)}_{a}$. The
second term in Eq.(6) is the harmonic cavity mode potential

$\displaystyle V_{c}(\underline{x})$
$\displaystyle=\dfrac{1}{2}\sum^{2N_{c}}_{\lambda,k}\omega^{2}_{k}\,x^{2}_{\lambda
k}\quad,$ (9)

for $2N_{c}$ cavity modes characterized by harmonic frequencies, $\omega_{k}$,
and cavity displacement coordinates, $x_{\lambda k}$, with mode index, $k$,
and polarization index, $\lambda$. Every cavity mode, $k$, is doubly
degenerate with respect to two orthogonal polarization directions, $\lambda$.

The third term in Eq.(6) resembles intermolecular interactions

$\displaystyle\hat{V}_{\mathrm{ens}}$
$\displaystyle=\dfrac{1}{2}\sum^{M}_{a,b\neq a}\hat{V}_{ab}\quad,$ (10)

which we approximate by a dipole-dipole interaction [15]

$\displaystyle\hat{V}_{ab}$
$\displaystyle=-\sum_{\kappa,\kappa^{\prime}}\hat{d}^{a}_{\kappa}\mathcal{T}^{\kappa\kappa^{\prime}}_{ab}\,\hat{d}^{b}_{\kappa^{\prime}}\quad,$
(11)

with Cartesian indices, $\kappa,\kappa^{\prime}$, Cartesian dipole operator
components, $\hat{d}^{a}_{\kappa}$ and $\hat{d}^{b}_{\kappa^{\prime}}$, as
well as the intermolecular interaction tensor

$\displaystyle\mathcal{T}^{\kappa\kappa^{\prime}}_{ab}$
$\displaystyle=-\dfrac{1}{R^{3}_{ab}}\left(\delta_{\kappa\kappa^{\prime}}-\dfrac{3\,R^{\kappa}_{ab}R^{\kappa^{\prime}}_{ab}}{R^{2}_{ab}}\right)\quad.$
(12)

The latter is characterized by the intermolecular distance,
$R_{ab}=|\underline{R}_{ab}|$, and a direction vector, $\underline{R}_{ab}$,
with Cartesian components, $R^{\kappa}_{ab}$ and $R^{\kappa^{\prime}}_{ab}$.

The last term in Eq.(6) constitutes the ensemble light-matter interaction

$\displaystyle\hat{W}_{\mathrm{ens}}$
$\displaystyle=\sum^{M}_{a}\hat{W}_{a}+\sum^{M}_{a,b\neq a}\hat{W}_{ab}\quad,$
(13)

with one- and two-body contributions

$\displaystyle\hat{W}_{a}$
$\displaystyle=-g_{0}\sum^{2N_{c}}_{\lambda,k}\omega_{k}\,\hat{d}^{a}_{\lambda
k}\,x_{\lambda
k}+\dfrac{g^{2}_{0}}{2}\sum^{2N_{c}}_{\lambda,k}\hat{d}^{a}_{\lambda
k}\hat{d}^{a}_{\lambda k}\quad,\vspace{0.2cm}$ (14)
$\displaystyle\hat{W}_{ab}$
$\displaystyle=\dfrac{g^{2}_{0}}{2}\sum^{2N_{c}}_{\lambda,k}\hat{d}^{a}_{\lambda
k}\hat{d}^{b}_{\lambda k}\quad.$ (15)

The one-body term, $\hat{W}_{a}$, contains the light-matter interaction linear
in the light-matter interaction constant,
$g_{0}=\frac{1}{\sqrt{\epsilon_{0}V_{\mathrm{cav}}}}$, with cavity volume,
$V_{\mathrm{cav}}$, and the permittivity, $\epsilon_{0}$, which is determined
by a polarization-projected dipole operator, $\hat{d}^{a}_{\lambda
k}=\underline{e}_{\lambda k}\cdot\underline{\hat{d}}^{(en)}_{a}$, for a cavity
mode polarization vector, $\underline{e}_{\lambda k}$, and molecular dipole
operator, $\underline{\hat{d}}^{(en)}_{a}$. Further, the second term of
$\hat{W}_{a}$ is the intramolecular dipole self-energy quadratic in $g_{0}$,
whereas $\hat{W}_{ab}$ resembles an intermolecular contribution of the dipole
self-energy, which augments the molecular interaction, $\hat{V}_{ab}$, in
Eq.(11).

####  II.2.2 Ensemble CBO-PT

In CBO-PT, the CBO electronic Hamiltonian in Eq.(6) is partitioned as[12]

$\displaystyle\hat{H}_{ec}$
$\displaystyle=\hat{H}_{0}+\lambda\,\hat{V}_{1}\quad,$ (16)

with perturbation parameter $\lambda$, zeroth-order Hamiltonian,
$\hat{H}_{0}$, and perturbation, $\hat{V}_{1}$, given by

$\displaystyle\hat{H}_{0}$
$\displaystyle=\hat{H}^{(e)}_{\mathrm{ens}}+V_{c}\quad,\vspace{0.2cm}$ (17)
$\displaystyle\hat{V}_{1}$
$\displaystyle=\hat{V}_{\mathrm{ens}}+\hat{W}_{\mathrm{ens}}\quad.$ (18)

The low-frequency cavity field is here treated as perturbation of the
electronic subsystem. Moreover, $\hat{V}_{1}$ contains both the ensemble
light-matter interaction contributions, $\hat{W}_{a}$ and $\hat{W}_{ab}$, and
the molecular dipole-dipole interaction, $\hat{V}_{ab}$, which allows for
electron-photon correlation corrections on both intra- and intermolecular
electronic interactions.[14]

The zeroth-order Hamiltonian, $\hat{H}_{0}$, is separable in the electronic
subspace and the harmonic cavity potential, $V_{c}$, constitutes just an
additive constant energy shift, such that zeroth-order reference states are
direct products of single-molecule adiabatic states

$\displaystyle\\{\ket{\Psi^{(e)}_{0}},\ket{\Psi^{(e)}_{\mu_{a}}},\ket{\Psi^{(e)}_{\mu_{a}}\Psi^{(e)}_{\mu_{b}}}\\}\quad,$
(19)

with non-interacting ensemble adiabatic ground state, $\ket{\Psi^{(e)}_{0}}$,
monomer excited states, $\ket{\Psi^{(e)}_{\mu_{a}}}$ (only the
$a^{\mathrm{th}}$-molecule is excited), and dimer excited states,
$\ket{\Psi^{(e)}_{\mu_{a}}\Psi^{(e)}_{\mu_{b}}}$ (both $a^{\mathrm{th}}$ and
$b^{\mathrm{th}}$-molecules are excited), respectively. We write CBO-PT($n$)
at $n^{\mathrm{th}}$-order of CBO-PT.

###  II.3 Ensemble CBO-PT Energies

The second-order CBO-PT ground state energy,
$E^{(ec)}_{0}\approx\mathcal{E}^{(2)}_{0}$, for a dipole-coupled molecular
ensemble under VSC can be written as

$\displaystyle\mathcal{E}^{(2)}_{0}=\mathcal{E}^{(2)}_{\mathrm{intra}}+\mathcal{E}^{(2)}_{\mathrm{inter}}\quad.$
(20)

The first term contains local intramolecular energies

$\displaystyle\mathcal{E}^{(2)}_{\mathrm{intra}}$
$\displaystyle=\sum^{M}_{a}\left(E^{(e)}_{a}+E^{(1)}_{a}+E^{(2)}_{a}\right)+V_{c}\quad,$
(21)

and the second term resembles the non-local intermolecular contribution, which
can be written in analogy to a many-body expansion as

$\displaystyle\mathcal{E}^{(2)}_{\mathrm{inter}}$ $\displaystyle=\sum_{a\neq
b}\left(E^{(1)}_{ab}+E^{(2)}_{ab}+U^{(2)}_{ab}\right)+\sum_{a\neq b\neq
c}E^{(2)}_{abc}\quad.$ (22)

Here, two-body and three-body corrections, $E^{(n)}_{ab}$ and $E^{(2)}_{abc}$,
are related to interactions, where only excited states of one molecule
contribute. The pair-interaction energy, $U^{(2)}_{ab}$, is in contrast
determined by excited states from both molecules. In the following, we provide
explicit expressions for intra- and intermolecular second-order CBO-PT
energies.

####  II.3.1 Intramolecular CBO-PT Energies

The intramolecular energy, $\mathcal{E}^{(2)}_{\mathrm{intra}}$, collects
contributions from all $M$ molecules and their interactions with the cavity
field. Here, the first-order energy correction, $E^{(1)}_{a}$, contains the
light-matter interaction and the ground-state DSE term

$\displaystyle E^{(1)}_{a}$
$\displaystyle=-g_{0}\sum^{2N_{c}}_{\lambda,k}\omega_{k}\,d^{0_{a}}_{\lambda
k}\,x_{\lambda
k}+\dfrac{g^{2}_{0}}{2}\sum^{2N_{c}}_{\lambda,k}\mathcal{F}^{0_{a}}_{\lambda
k}\quad,$ (23)

with

$\displaystyle\mathcal{F}^{0_{a}}_{\lambda k}$
$\displaystyle=d^{0_{a}}_{\lambda k}d^{0_{a}}_{\lambda
k}+\tilde{\mathcal{F}}^{0_{a}}_{\lambda k}\quad.$ (24)

Here, $d^{0_{a}}_{\lambda k}$, is the polarization-projected ground state
dipole moment of the $a^{\mathrm{th}}$ molecule and
$\tilde{\mathcal{F}}^{0_{a}}_{\lambda k}$ is the corresponding dipole
fluctuation contribution

$\displaystyle\tilde{\mathcal{F}}^{0_{a}}_{\lambda k}$
$\displaystyle=\braket{\Psi^{(e)}_{0}}{\hat{d}^{(e)}_{\lambda
k}\hat{d}^{(e)}_{\lambda
k}}{\Psi^{(e)}_{0}}-\braket{\Psi^{(e)}_{0}}{\hat{d}^{(e)}_{\lambda
k}}{\Psi^{(e)}_{0}}^{2}\quad,$ (25)

which is purely electronic in character.[15] The dipole fluctuation term is
regularly neglected in a mean-field type approximation underlying effective
ground state Hamiltonians for vibro-polaritonic chemistry[12], which do not
capture cavity-induced changes of classical activation barrier but VSC related
distortions of cPES leading for example to barrier broadening[19, 18].

The second-order CBO-PT energy correction, $E^{(2)}_{a}$, decomposes into
three terms

$\displaystyle E^{(2)}_{a}$
$\displaystyle=-\dfrac{g^{2}_{0}}{2}\sum^{2N_{c}}_{\lambda,k}\sum^{2N_{c}}_{\lambda^{\prime},k^{\prime}}\omega_{k}\,x_{\lambda
k}\,\alpha^{0_{a}}_{\lambda\lambda^{\prime}}\,x_{\lambda^{\prime}k^{\prime}}\,\omega_{k^{\prime}}\vspace{0.2cm}$
(26) $\displaystyle\hskip
11.38092pt+\dfrac{g^{3}_{0}}{2}N_{c}\sum^{N_{c}}_{k}\sum^{2}_{\lambda,\lambda^{\prime}}\omega_{k}\,x_{\lambda
k}\,\mathcal{A}^{0_{a}}_{\lambda\lambda^{\prime}}\vspace{0.2cm}$ (27)
$\displaystyle\hskip
11.38092pt-\dfrac{g^{4}_{0}}{8}N^{2}_{c}\sum^{2}_{\lambda,\lambda^{\prime}}\mathcal{B}^{0_{a}}_{\lambda\lambda^{\prime}}\quad,$
(28)

which resembles a bare cavity correction, a light-matter interaction
correction and a bare molecular correction. The first term, determined by
polarization-projected static polarizability tensor elements,
$\alpha^{0_{a}}_{\lambda\lambda^{\prime}}$, constitutes a matter-induced
screening of harmonic cavity frequencies ($k=k^{\prime}$) besides a matter-
mediated inter-cavity mode coupling ($k\neq k^{\prime}$).[12] The interaction
and DSE corrections in Eqs.(27) and (28) contain linear response functions,
$\mathcal{A}^{0_{a}}_{\lambda\lambda^{\prime}}$ and
$\mathcal{B}^{0_{a}}_{\lambda\lambda^{\prime}}$, which are explicitly given in
Appendix .1.

####  II.3.2 Intermolecular CBO-PT Energies

Turning to the inter-molecular CBO-PT(2) energy in Eq.(22),
$\mathcal{E}^{(2)}_{\mathrm{inter}}$, the first-order contribution

$\displaystyle E^{(1)}_{ab}$ $\displaystyle=\dfrac{1}{2}V_{ab}+W_{ab}$ (29)

contains both a distance-dependent dipole-dipole interaction, $V_{ab}$,
besides a distance-independent DSE-induced bimolecular interaction, $W_{ab}$,
which explicitly read

$\displaystyle V_{ab}$
$\displaystyle=-\sum_{\kappa,\kappa^{\prime}}d^{0_{a}}_{\kappa}\mathcal{T}^{\kappa\kappa^{\prime}}_{ab}\,d^{0_{b}}_{\kappa^{\prime}}\quad,\vspace{0.2cm}$
(30) $\displaystyle W_{ab}$
$\displaystyle=\dfrac{g^{2}_{0}}{2}\sum^{2N_{c}}_{\lambda,k}d^{0_{a}}_{\lambda
k}d^{0_{b}}_{\lambda k}\quad.$ (31)

The second-order two-body contribution, $E^{(2)}_{ab}$, contains two distinct
terms

$\displaystyle E^{(2)}_{ab}$
$\displaystyle=\tilde{E}^{(2)}_{ab}+\tilde{E}^{(2)}_{abb}\quad,$ (32)

where, $\tilde{E}^{(2)}_{ab}$, decomposes into

$\displaystyle\tilde{E}^{(2)}_{ab}$
$\displaystyle=-g_{0}\sum^{2N_{c}}_{\lambda,k}\omega_{k}\,x_{\lambda
k}\,\sum_{\kappa,\kappa^{\prime}}\alpha^{0_{a}}_{\lambda\kappa}\,\mathcal{T}^{\kappa\kappa^{\prime}}_{ab}d^{0_{b}}_{\kappa^{\prime}}\vspace{0.2cm}$
(33) $\displaystyle\hskip
11.38092pt+\dfrac{g^{2}_{0}}{2}N_{c}\sum^{2}_{\lambda}\sum_{\kappa,\kappa^{\prime}}\mathcal{A}^{0_{a}}_{\lambda\kappa}\mathcal{T}^{\kappa\kappa^{\prime}}_{ab}d^{0_{b}}_{\kappa^{\prime}}\vspace{0.2cm}$
(34) $\displaystyle\hskip
11.38092pt+g^{3}_{0}\,N_{c}\sum^{2N_{c}}_{\lambda,k}\sum^{2}_{\lambda^{\prime}}\omega_{k}\,x_{\lambda
k}\,\alpha^{0_{a}}_{\lambda\lambda^{\prime}}\,d^{0_{b}}_{\lambda^{\prime}}\vspace{0.2cm}$
(35) $\displaystyle\hskip
11.38092pt-\dfrac{g^{4}_{0}}{2}N^{2}_{c}\sum^{2}_{\lambda,\lambda^{\prime}}\mathcal{A}^{0_{a}}_{\lambda\lambda^{\prime}}d^{0_{b}}_{\lambda^{\prime}}\quad.$
(36)

Here, the terms linear and quadratic in $g_{0}$ dependent via
$\mathcal{T}^{\kappa\kappa^{\prime}}_{ab}$ on the intermolecular distance,
$R_{ab}$, in contrast to distance-independent contributions cubic and quartic
in $g_{0}$. Further, terms scaling as $g_{0}$ and $g^{3}_{0}$ explicitly
couple interacting molecules to the cavity field via $x_{\lambda k}$. The two-
body three-index correction, $\tilde{E}^{(2)}_{abb}$, in Eq.(32) follows as
special case of the three-body term ($b=c$)

$\displaystyle E^{(2)}_{abc}$
$\displaystyle=-\dfrac{1}{2}\sum_{\kappa,\kappa^{\prime}}\sum_{\gamma,\gamma^{\prime}}\mathcal{T}^{\kappa\kappa^{\prime}}_{ab}\,\mathcal{T}^{\gamma\gamma^{\prime}}_{ac}d^{0_{b}}_{\kappa^{\prime}}\,\alpha^{0_{a}}_{\kappa\gamma}\,d^{0_{c}}_{\gamma^{\prime}}\vspace{0.2cm}$
(37) $\displaystyle\hskip
11.38092pt+g^{2}_{0}\,N_{c}\sum^{2}_{\lambda}\sum_{\kappa,\kappa^{\prime}}\mathcal{T}^{\kappa\kappa^{\prime}}_{ab}d^{0_{b}}_{\kappa^{\prime}}\alpha^{0_{a}}_{\kappa\lambda}\,d^{0_{c}}_{\lambda}\vspace{0.2cm}$
(38) $\displaystyle\hskip
11.38092pt-\dfrac{g^{4}_{0}}{2}N^{2}_{c}\sum^{2}_{\lambda,\lambda^{\prime}}d^{0_{b}}_{\lambda}\alpha^{0_{a}}_{\lambda\lambda^{\prime}}d^{0_{c}}_{\lambda^{\prime}}\quad,$
(39)

and can be interpreted as cavity-modified dipole-induced dipole interaction
(cf. Sec.VI). Finally, the second-order pair interaction correction reads

$\displaystyle U^{(2)}_{ab}$
$\displaystyle=U^{ab}_{\mathrm{vdW}}+U^{ab}_{\mathrm{cvdW}}+U^{ab}_{\mathrm{cav}}\quad,$
(40)

and contains three individual terms

$\displaystyle U^{ab}_{\mathrm{vdW}}$
$\displaystyle=-\dfrac{1}{2}\sum_{\kappa,\kappa^{\prime}}\sum_{\gamma,\gamma^{\prime}}\mathcal{T}^{\kappa\kappa^{\prime}}_{ab}\mathcal{T}^{\gamma\gamma^{\prime}}_{ab}\alpha^{ab}_{\kappa\gamma\kappa^{\prime}\gamma^{\prime}}\quad,\vspace{0.2cm}$
(41) $\displaystyle U^{ab}_{\mathrm{cvdW}}$
$\displaystyle=g^{2}_{0}\sum^{2}_{\lambda}\sum_{\kappa,\kappa^{\prime}}\mathcal{T}^{\kappa\kappa^{\prime}}_{ab}\alpha^{ab}_{\kappa\lambda\kappa^{\prime}\lambda}\quad,\vspace{0.2cm}$
(42) $\displaystyle U^{ab}_{\mathrm{cav}}$
$\displaystyle=-\dfrac{g^{4}_{0}}{2}N^{2}_{c}\sum^{2}_{\lambda,\lambda^{\prime}}\alpha^{ab}_{\lambda\lambda^{\prime}\lambda\lambda^{\prime}}\quad.$
(43)

The first term resembles the common van-der-Waals (vdW) interaction,
$U^{ab}_{\mathrm{vdW}}$, which is augmented by a cavity-induced long-range
vdW-type interaction, $U^{ab}_{\mathrm{cvdW}}$, and a purely cavity-related
distance independent term $U^{ab}_{\mathrm{cav}}$. Identical corrections to
the van-der-Waals interaction have been derived recently for intermolecular
interactions in the ESC regime.[14, 15]

With intra- and intermolecular second-order CBO-PT energies at hand, we will
now discuss cavity-induced modifications of electronic properties in
molecules. For this purpose, we first introduce the notion of a cavity
reaction potential.

##  III Cavity Reaction Potentials

Chemical reactions are understood as changes in nuclear configuration space,
which are regularly discussed in terms of stationary points on (cavity)
potential energy surfaces. In order to exploit this perspective in vibro-
polaritonic settings, we introduce the cavity reaction potential (CRP)[18] as
first term of the Taylor expansion

$\displaystyle E^{(ec)}_{0}(\underline{R},\underline{x})$
$\displaystyle=\mathcal{V}^{(ec)}_{\mathrm{crp}}(\underline{R})+\dfrac{1}{2}\Delta\underline{x}\,\underline{\underline{H}}^{(ec)}_{CC}(\underline{R})\Delta\underline{x}\quad,$
(44)

which expands CBO ground state cPES around a cavity reference coordinate,
$\underline{x}^{(ec)}_{0}$, with

$\displaystyle\mathcal{V}^{(ec)}_{\mathrm{crp}}(\underline{R})$
$\displaystyle=E^{(ec)}_{0}(\underline{R},\underline{x}^{(ec)}_{0})\quad.$
(45)

The second term in Eq.(44) contains the cavity-component of the CBO Hessian,
$\underline{\underline{H}}^{(ec)}_{CC}$[20], which we recently derived for
CBO-PT(2)[13], and displacements,
$\Delta\underline{x}=\underline{x}-\underline{x}^{(ec)}_{0}$. The CRP provides
a static perspective on cavity-altered reaction mechanisms, i.e., it fully
captures changes in nuclear configuration but no excitations in the cavity
mode subspace relevant for dynamical processes.[21, 18, 22, 23, 24, 25, 26,
27] The minimizing cavity coordinate, $\underline{x}^{(ec)}_{0}$, defining the
CRP satisfies two constraints[16]

$\displaystyle\dfrac{\partial}{\partial x_{\lambda
k}}E^{(ec)}_{0}=0\quad\forall\lambda,k\,\Leftrightarrow\,\braket{\Psi^{(ec)}_{0}}{\underline{\hat{E}}_{\perp}}{\Psi^{(ec)}_{0}}=0\quad,$
(46)

and is explicitly given by

$\displaystyle x^{(ec)}_{0}$
$\displaystyle=\dfrac{g_{0}}{\omega_{k}}\braket{\Psi^{(ec)}_{0}}{\hat{d}_{\lambda}}{\Psi^{(ec)}_{0}}\quad.$
(47)

The first condition in Eq.(46) minimizes the energy with respect to variations
along cavity displacement coordinates and the second one ensures a “non-
radiating” ground state: It connects to the physical condition of a vanishing
ground state expectation value of the transverse electrical component of the
cavity field, $\underline{\hat{E}}_{\perp}$.[28, 16]

In Sec.V, we will discuss the CRP of a monomer with a focus on cavity-induced
modifications of the intramolecular energy given by

$\displaystyle\Delta^{(ec)}_{\mathrm{crp}}(\underline{R};\lambda)$
$\displaystyle=\mathcal{V}^{(ec)}_{\mathrm{crp}}(\underline{R};\lambda)-E^{(e)}_{0}(\underline{R})\quad,$
(48)

which vanishes in the non-interacting limit, $g_{0}\to 0$. Specifically, we
will discuss different approximations to $\Delta^{(ec)}_{\mathrm{crp}}$ based
on either analytical CBO-PT results or non-perturbative numerical solutions of
the CBO electronic ground state problem obtained from ab initio wave function
methods to be introduced the next section IV. In addition, we discuss the CRP
of an interacting dimer in Sec.VI, which allows us to examine cavity-induced
modifications of intermolecular interactions.

##  IV Ab initio CBO Wave Function Theory

As non-perturbative reference for CBO-PT, we consider numerical solutions of
the ground state CBO electronic TISE (2) at two different levels of ab initio
wave function theory: CBO-Hartree Fock theory[16] and CBO coupled cluster
theory with singles and doubles (CBO-CCSD)[17]. Both methods have been
implemented via the PySCF software package[33, 34].

###  IV.1 CRP CBO Hartree-Fock Theory

Inspired by the concept of cavity reaction potentials, we introduce a simple
but beneficial reformulation of CBO-HF theory, which self-consistently
includes the minimizing cavity coordinate Eq.(47) in the CBO-Fockian. This
approach leads to an a priori origin-independent total energy and a cavity-
coordinate free form of the CBO-HF equations.

We consider a single Slater-Determinant, $\Psi^{(ec)}_{\mathrm{cbohf}}$, and
minimize the energy expectation value

$\displaystyle E^{(ec)}_{0}(\underline{R},\underline{x})$
$\displaystyle=\underset{\underline{x},\underline{\kappa}}{\mathrm{min}}\braket{\Psi^{(ec)}_{\mathrm{cbohf}}(\underline{\kappa})}{\hat{H}_{ec}(\underline{x})}{\Psi^{(ec)}_{\mathrm{cbohf}}(\underline{\kappa})}\,,$
(49)

simultaneously with respect to both orbital rotation parameters,
$\underline{\kappa}$, and cavity-displacement coordinates, $\underline{x}$.
From Eq.(49), we obtain a CRP-type CBO Fockian with matrix elements

$\displaystyle\tilde{f}^{pq}_{\lambda}$
$\displaystyle=\tilde{h}^{pq}_{\lambda}+\tilde{v}^{pq}_{\lambda}\quad,$ (50)

where we restrict the discussion to a single cavity mode with polarization,
$\lambda$, while generalization to multi-mode systems as discussed in Sec.II
being straightforward. Core and mean-field potential contributions for a
closed-shell system are given by

$\displaystyle\tilde{h}^{pq}_{\lambda}$
$\displaystyle=h_{pq}+\dfrac{g^{2}_{0}}{2}O^{pq}_{\lambda}\quad,\vspace{0.2cm}$
(51) $\displaystyle\tilde{v}^{pq}_{\lambda}$
$\displaystyle=J_{pq}-\dfrac{1}{2}\tilde{K}^{pq}_{\lambda}\quad,$ (52)

where the bare electronic core contribution, $h_{pq}$, is augmented by the
1-particle contribution of the DSE, $O^{pq}_{\lambda}$ (cf. Appendix .2).
Moreover, the self-consistent inclusion of Eq.(47) leads to a mean-field
potential determined by the bare Coulomb contribution, $J_{pq}$, and a DSE-
corrected exchange term

$\displaystyle\tilde{K}^{pq}_{\lambda}=\sum_{rs}D_{rs}\left(g_{prsq}-g^{2}_{0}\,d^{pr}_{\lambda}d^{sq}_{\lambda}\right)\quad.$
(53)

The total CRP CBO-RHF ground state energy reads

$\displaystyle\mathcal{V}^{\mathrm{cbohf}}_{\mathrm{crp}}(\underline{R})$
$\displaystyle=E^{(e)}_{\mathrm{cbohf}}+\dfrac{g^{2}_{0}}{2}\tilde{\mathcal{F}}^{\mathrm{cbohf}}_{\lambda}+V_{nn}\quad,$
(54)

which is independent of both cavity mode frequencies and cavity displacement
coordinates. Here, $E^{(e)}_{\mathrm{cbohf}}$ is the RHF energy evaluated with
respect to CBO-RHF molecular orbitals, which account for DSE-induced
relaxation effects, and $\tilde{\mathcal{F}}^{\mathrm{cbohf}}_{\lambda}$
resembles the CBO-RHF dipole-fluctuation correction (cf. Appendix .2)

$\displaystyle\tilde{\mathcal{F}}^{\mathrm{cbohf}}_{\lambda}$
$\displaystyle=2\sum_{i}O^{ii}_{\lambda}-2\sum_{ij}d^{ij}_{\lambda}d^{ji}_{\lambda}\quad.$
(55)

This expression is formally equivalent to Eq.(.2.11) but evaluated from
relaxed CBO-RHF molecular orbitals in contrast to bare molecular orbitals.
Moreover, from Eq.(54) we obtain the corresponding mean-field energy
correction, $\Delta^{\mathrm{cbohf}}_{\mathrm{crp}}$.

###  IV.2 CRP CBO Coupled Cluster Theory

We consider a CRP-type CBO coupled cluster approach with singles and doubles
corrections (CBO-CCSD)[17] to account for electron correlation under VSC as
determined by the wave function ansatz

$\displaystyle\ket{\Psi^{(ec)}_{0}}$
$\displaystyle=\exp\left(\hat{T}_{1}+\hat{T}_{2}\right)\ket{\Psi^{(ec)}_{\mathrm{cbohf}}}\quad,$
(56)

with common singles and doubles excitation operators, $\hat{T}_{1}$ and
$\hat{T}_{2}$. The only difference to molecular CCSD is the CRP CBO-HF
reference state, $\Psi^{(ec)}_{\mathrm{cbohf}}$. We employ the CBO-CCSD
approach to numerically calculate the non-perturbative cavity reaction
potential

$\displaystyle\mathcal{V}^{(ec)}_{\mathrm{crp}}(\underline{R})$
$\displaystyle=\braket{\Psi^{(ec)}_{0}}{\hat{H}^{0}_{ec}}{\Psi^{(ec)}_{0}}\quad,\vspace{0.2cm}$
(57)
$\displaystyle=\braket{\Psi^{(ec)}_{0}}{\hat{H}_{e}}{\Psi^{(ec)}_{0}}+\dfrac{g^{2}_{0}}{2}\mathcal{F}^{(ec)}_{\lambda}+V_{nn}\quad,$
(58)

where the CBO adiabatic ground state, $\Psi^{(ec)}_{0}$, is here approximated
by Eq.(56). The second term in Eq.(58) corresponds to the dipole-fluctuation
contribution

$\displaystyle\mathcal{F}^{(ec)}_{\lambda}$
$\displaystyle=\braket{\Psi^{(ec)}_{0}}{\hat{d}_{\lambda}\hat{d}_{\lambda}}{\Psi^{(ec)}_{0}}-\braket{\Psi^{(ec)}_{0}}{\hat{d}_{\lambda}}{\Psi^{(ec)}_{0}}^{2}\quad,$
(59)

evaluated with respect to the CBO-CCSD adiabatic ground state and the CRP-type
CBO electronic Hamiltonian is obtained with Eq.(47) as

$\displaystyle\hat{H}^{0}_{ec}$
$\displaystyle=\sum_{pq}\tilde{h}^{0}_{pq}\hat{E}_{pq}+\dfrac{1}{2}\sum_{pqrs}\tilde{g}_{pqrs}\hat{e}_{pqrs}+\tilde{V}_{0}\quad.$
(60)

Indices $p,q,r,s$ run over the molecular orbital basis and the constant energy
shift is given by

$\displaystyle\tilde{V}_{0}$
$\displaystyle=V_{nn}+\dfrac{g^{2}_{0}}{2}\braket{\Psi^{(ec)}_{0}}{\hat{d}_{\lambda}}{\Psi^{(ec)}_{0}}^{2}\quad.$
(61)

As pointed out in Ref.[17], $\hat{H}^{0}_{ec}$ allows for a direct
implementation of (CRP) CBO-CCSD methods by augmenting 1- and 2-electron
integrals as

$\displaystyle\tilde{h}^{0}_{pq}$
$\displaystyle=\hat{h}_{pq}+\dfrac{g^{2}_{0}}{2}O^{pq}_{\lambda}-g^{2}_{0}\,d^{pq}_{\lambda}\braket{\Psi^{(ec)}_{0}}{\hat{d}_{\lambda}}{\Psi^{(ec)}_{0}}\quad,$
(62)

and

$\displaystyle\tilde{g}_{pqrs}$
$\displaystyle=g_{pqrs}-g^{2}_{0}\,d^{pq}_{\lambda}d^{rs}_{\lambda}\quad,$
(63)

where $\tilde{h}^{0}_{pq}$ depends on the ground state, which requires a self-
consistent solution of the problem. Practically, the dipole expectation value
in $\tilde{h}^{0}_{pq}$ will be first approximated via CRP CBO-HF to generate
a correlated solution, which is subsequently used to update the Hamiltonian
towards a CBO-CC solution.

##  V Cavity-Modified Intramolecular Interactions

We comparatively discuss cavity-induced modifications of intramolecular
electronic interactions by analysing VSC-related energy changes of molecular
PES based on both CBO-PT and non-perturbative CBO-RHF/CBO-CCSD approaches in
their CRP-type formulation.

Figure 1: Torsion PES $E^{(e)}_{0}(\theta)$ (TPSSh/cc-pVTZ) in black of
1,2-DCE (a,b) and 1,2-DFE(c,d) with cavity-induced energy corrections
evaluated via CBO-PT(1)-RHF/cc-pVDZ and CBO-RHF/cc-pVDZ (a,c) besides CBO-
PT(1)-CCSD/cc-pVDZ and CBO-CCSD/cc-pVDZ (b,d) at light-matter coupling
$g_{0}=0.03\sqrt{E_{h}}/e\,a_{0}$ for cavity-mode polarizations $\lambda=x$
(red), $\lambda=y$ (green) and $\lambda=z$ (blue).

###  V.1 Unimolecular Cavity Reaction Potential

We first derive the unimolecular CRP at second-order of CBO-PT. We consider a
single molecule coupled to a single cavity mode with coordinate $x_{c}$ and
polarization, $\lambda$. The corresponding CBO-PT(2) cPES reads

$\displaystyle\mathcal{E}^{(2)}_{0}$
$\displaystyle=E^{(e)}_{0}+V_{c}+E^{(1)}_{0}+E^{(2)}_{0}\quad,$ (64)

where first- and second-order CBO-PT energy corrections are explicitly given
by

$\displaystyle E^{(1)}_{0}$
$\displaystyle=-g_{0}\,\omega_{c}\,d^{0}_{\lambda}x_{c}+\dfrac{g^{2}_{0}}{2}d^{0}_{\lambda}d^{0}_{\lambda}+\dfrac{g^{2}_{0}}{2}\tilde{\mathcal{F}}^{0}_{\lambda}\quad,$
(65)

and

$\displaystyle E^{(2)}_{0}$
$\displaystyle=-\dfrac{g^{2}_{0}}{2}\omega^{2}_{c}\,\alpha^{0}_{\lambda\lambda}\,x^{2}_{c}+\dfrac{g^{3}_{0}}{2}\omega_{c}\,x_{c}\,\mathcal{A}^{0}_{\lambda\lambda}-\dfrac{g^{4}_{0}}{8}\mathcal{B}^{0}_{\lambda\lambda}\quad.$
(66)

The second-order minimizing cavity coordinate, $x^{(2)}_{0}$, for the
unimolecular problem respecting both conditions in Eq.(46) is obtained as

$\displaystyle x^{(2)}_{0}$
$\displaystyle=\dfrac{g_{0}\,d^{0}_{\lambda}}{\omega_{c}\left(1-g^{2}_{0}\,\alpha^{0}_{\lambda\lambda}\right)}-\dfrac{g^{3}_{0}\,\alpha^{0}_{\lambda\lambda}d^{0}_{\lambda}}{2\omega_{c}\left(1-g^{2}_{0}\,\alpha^{0}_{\lambda\lambda}\right)}\quad,$
(67)

where the non-radiating ground state condition imposes a mean-field
approximation of the interaction-DSE cross term in Eq.(66) scaling as
$g^{3}_{0}$

$\displaystyle\mathcal{A}^{0}_{\lambda\lambda}$
$\displaystyle\approx\alpha^{0}_{\lambda\lambda}d^{0}_{\lambda}\quad.$ (68)

The unimolecular CBO-PT(2) cavity reaction potential is then obtained from
Eqs.(64) and (67) as

$\displaystyle\mathcal{V}^{(2)}_{\mathrm{crp}}(\underline{R};\lambda)=E^{(e)}_{0}(\underline{R})+\dfrac{g^{2}_{0}}{2}\tilde{\mathcal{F}}^{0}_{\lambda}(\underline{R})-\dfrac{g^{4}_{0}}{8}\mathcal{B}^{0}_{\lambda\lambda}(\underline{R})\quad,$
(69)

where we made the parametric dependence of $\mathcal{V}^{(2)}_{\mathrm{crp}}$
on the cavity-mode polarization, $\lambda$, explicit. The first term on the
right-hand side is simply the molecular ground state PES, $E^{(e)}_{0}$, which
is in leading-order corrected by the dipole-fluctuation term,
$\tilde{\mathcal{F}}^{0}_{\lambda}$. The negative third term emerges at second
order of CBO-PT and is determined by the linear response function,
$\mathcal{B}^{0}_{\lambda\lambda}$ (cf. Appendix .1).

###  V.2 Unimolecular Isomerization under VSC

We will now analyse cavity reaction potentials for the isomerization of
1,2-Dichloroethane (1,2-DCE) and its less polarizable fluorinated equivalent
1,2-Difluoroethane (1,2-DFE) under VSC with $g_{0}=0.03\sqrt{E_{h}}/e\,a_{0}$.
In both cases, we have calculated chemically accurate molecular PES,
$E^{(e)}_{0}(\theta)$, as function of the torsion angle $\theta$ (cf. Fig.1)
via TPSSh/cc-pVTZ as implemented in ORCA 5.0.4[29, 30].

Specifically, we compare cavity-induced energy corrections obtained from CBO-
PT and non-perturbative wave function approaches. To this end, we rewrite the
non-perturbative energy correction in Eq.(48) as

$\displaystyle\Delta^{(ec)}_{\mathrm{crp}}(\underline{R};\lambda)$
$\displaystyle=\dfrac{g^{2}_{0}}{2}\mathcal{\tilde{F}}^{0}_{\lambda}(\underline{R})+\Delta^{(ec)}_{\mathrm{rel}}(\underline{R};\lambda)\quad,$
(70)

such that the right-hand side contains simply the CBO-PT(1) energy correction
given in Eq.(25) and a second term

$\displaystyle\Delta^{(ec)}_{\mathrm{rel}}(\underline{R};\lambda)$
$\displaystyle=\mathcal{V}^{(ec)}_{\mathrm{crp}}(\underline{R};\lambda)-\mathcal{V}^{(1)}_{\mathrm{crp}}(\underline{R};\lambda)\quad,$
(71)

which captures an energy correction due to cavity-induced relaxation effects
of adiabatic states. Eq.(70) indicates that cavity-induced intramolecular
corrections to molecular PES are formally captured already at first-order of
CBO-PT with higher-order corrections being implicitly contained in
$\Delta^{(ec)}_{\mathrm{rel}}$.

We now compare both CBO-RHF/CBO-PT(1)-RHF and CBO-CCSD/CBO-PT(1)-CCSD levels
of theory, where we restrict the discussion to the first-order dipole-
fluctuation correction. For both 1,2-DCE and 1,2-DFE (cf. Fig.1), we find a
cavity-induced energy shift of the torsion PES whose magnitude and dependence
of cavity mode polarization is slightly more pronounced for the more
polarizable 1,2-DCE. In contrast to previous results[31], which did not take
into account the dipole self-energy term[32, 28], we do not observe
significant changes in classical activation barriers. Polarization-dependence
of the energy corrections is small

At both levels of theory, we find an excellent agreement between CBO-PT(1) and
ab initio wave function theory (cf. Fig.1a,b vs. Fig.1c,d), which indicates
that cavity-induced intramolecular energy corrections are already captured at
first-order of CBO-PT. In addition, this observation points at negligible
state relaxation effects for the present examples, i.e.,
$\Delta^{(ec)}_{\mathrm{rel}}\ll\frac{g^{2}_{0}}{2}\mathcal{\tilde{F}}^{0}_{\lambda}$
in Eq.(70). This statement translates into small cavity-induced changes in MO
coefficients for CBO-HF and cluster amplitudes for CBO-CCSD approaches. Thus
local cavity-induced modifications of electronic interactions seem to be
captured by CBO-PT(1), which is not sufficient for the related Hessian, whose
leading correction appears only in CBO-PT(2).[12]

##  VI Cavity-Modified Intermolecular Interactions

We now illustrate that second-order CBO-PT corrections are essential to
capture non-trivial cavity-induced modifications of intermolecular
interactions, which might be relevant for understanding collective strong
coupling effects under VSC.

Figure 2: Cavity-modified interaction potentials for (a) a H2 dimer at light-
matter coupling $g_{0}=0.03\sqrt{E_{h}}/e\,a_{0}$ and (b) a H2O dimer at
light-matter coupling $g_{0}=0.05\sqrt{E_{h}}/e\,a_{0}$ evaluated via CBO-
CCSD/aug-cc-pVDZ for cavity-mode polarizations $\lambda=x$ (red), $\lambda=y$
(green) and $\lambda=z$ (blue). Interaction potentials are obtained by
subtracting energies at $50$ Å.

###  VI.1 Interacting Dimer under VSC

We start by deriving the CRP for an interacting dimer composed of molecules A
and B at distance $R$ under VSC with a single cavity mode. The corresponding
CBO-PT(2) cPES reads

$\displaystyle\mathcal{E}^{(2)}_{\mathrm{dimer}}$
$\displaystyle=\tilde{\mathcal{E}}^{(2)}_{A}+\tilde{\mathcal{E}}^{(2)}_{B}+V_{c}+W^{(2)}_{AB}\quad,$
(72)

with unimolecular contributions, $\tilde{\mathcal{E}}^{(2)}_{A}$ and
$\tilde{\mathcal{E}}^{(2)}_{B}$, in Eq.(64) excluding the cavity potential,
$V_{c}$. The intermolecular interaction

$\displaystyle W^{(2)}_{AB}$
$\displaystyle=V^{(2)}_{AB}+\mathcal{E}^{(2)}_{AB}+U^{(2)}_{AB}\quad,$ (73)

contains three distinct contributions, which we consider along a fixed
Cartesian axis, $\kappa$, such that the intermolecular interaction tensor in
Eq.(12) turns into $\mathcal{T}^{\kappa\kappa}_{AB}=2R^{-3}$. The first term,
$V^{(2)}_{AB}$, relates to a cavity-altered dipole-dipole interaction

$\displaystyle V^{(2)}_{AB}(R)$
$\displaystyle=-\dfrac{D_{3}}{R^{3}}-g_{0}\dfrac{D^{\prime}_{3\lambda}\,\omega_{c}}{R^{3}}x_{c}+g^{2}_{0}\dfrac{D^{\prime\prime}_{3\lambda}}{R^{3}}\quad,$
(74)

where the second and third term follow from the leading-order contributions of
Eq.(36) with coefficients

$\displaystyle D_{3}$
$\displaystyle=2\,d^{0_{A}}_{\kappa}d^{0_{B}}_{\kappa}\quad,\vspace{0.2cm}$
(75) $\displaystyle D^{\prime}_{3\lambda}$
$\displaystyle=2\left(\alpha^{0_{A}}_{\lambda\kappa}d^{0_{B}}_{\kappa}+\alpha^{0_{B}}_{\lambda\kappa}d^{0_{A}}_{\kappa}\right)\quad,\vspace{0.2cm}$
(76) $\displaystyle D^{\prime\prime}_{3\lambda}$
$\displaystyle=\alpha^{0_{A}}_{\lambda\kappa}d^{0_{A}}_{\kappa}d^{0_{B}}_{\kappa}+\alpha^{0_{B}}_{\lambda\kappa}d^{0_{B}}_{\kappa}d^{0_{A}}_{\kappa}\quad.$
(77)

Here, we assume the mean-field approximation in Eq.(.1.5) for
$D^{\prime\prime}_{3}$. The first correction ($D^{\prime}_{3}$) is linear in
the light-matter coupling constant, $g_{0}$, and constitutes a cavity-mediated
dipole-polarizability type interaction. This term directly couples the
molecular dimer to the cavity field via the coordinate, $x_{c}$, which will
directly influence the dynamics of the cavity subsystem due to energy transfer
from the molecules. The second correction ($D^{\prime\prime}_{3}$) resembles a
higher-order interaction between two polar molecules, which additionally
involves polarizability effects.

The second interaction term in Eq.(73) is distance-independent and exclusively
contains cavity-induced contributions

$\displaystyle\mathcal{E}^{(2)}_{AB}$
$\displaystyle=g^{2}_{0}\,\tilde{D}_{3\lambda}+g^{3}_{0}\,\tilde{D}^{\prime}_{3\lambda}\,\omega_{c}\,x_{c}-\dfrac{g^{4}_{0}}{2}\tilde{D}^{\prime\prime}_{3\lambda}\quad,$
(78)

with coefficients

$\displaystyle\tilde{D}_{3\lambda}$
$\displaystyle=d^{0_{A}}_{\lambda}d^{0_{B}}_{\lambda}\quad,\vspace{0.2cm}$
(79) $\displaystyle\tilde{D}^{\prime}_{3\lambda}$
$\displaystyle=\alpha^{0_{A}}_{\lambda\lambda}d^{0_{B}}_{\lambda}+\alpha^{0_{B}}_{\lambda\lambda}d^{0_{A}}_{\lambda}\quad,\vspace{0.2cm}$
(80) $\displaystyle\tilde{D}^{\prime\prime}_{3\lambda}$
$\displaystyle=\alpha^{0_{A}}_{\lambda\lambda}d^{0_{A}}_{\lambda}d^{0_{B}}_{\lambda}+\alpha^{0_{B}}_{\lambda\lambda}d^{0_{B}}_{\lambda}d^{0_{A}}_{\lambda}\quad,$
(81)

which are closely related to their counterparts in Eq.(74). The first term
($\tilde{D}_{3}$) constitutes the first-order DSE-induced intermolecular
interaction, whereas the second and third term scaling as $g^{3}_{0}$ and
$g^{4}_{0}$ stem from the second-order correction in Eq.(36). All coefficients
entering Eq.(78) are fully determined by the cavity-mode polarizations in
contrast to the Cartesian axis determining coefficients in Eq.(74).

The third interaction term in Eq.(73) can be decomposed into two distance-
dependent contributions

$\displaystyle U^{(2)}_{AB}(R)$
$\displaystyle=\tilde{U}^{(2)}_{AB}(R)+\tilde{u}^{(2)}_{AB}(R)\quad.$ (82)

The first term relates to a cavity-modified van-der-Waals interaction with
long-range character

$\displaystyle\tilde{U}^{(2)}_{AB}(R)$
$\displaystyle=-\dfrac{C_{6}}{R^{6}}+g^{2}_{0}\dfrac{C_{3\lambda}}{R^{3}}-g^{4}_{0}\dfrac{C_{0\lambda}}{R^{0}}\quad,$
(83)

with coefficients

$\displaystyle C_{6}$
$\displaystyle=2\alpha^{AB}_{\kappa\kappa\kappa\kappa}\quad,\vspace{0.2cm}$
(84) $\displaystyle C_{3\lambda}$
$\displaystyle=2\alpha^{AB}_{\kappa\lambda\kappa\lambda}\quad,\vspace{0.2cm}$
(85) $\displaystyle C_{0\lambda}$
$\displaystyle=\dfrac{\alpha^{AB}_{\lambda\lambda\lambda\lambda}}{2}\quad.$
(86)

which has been previously derived by Koch and coworkers in an electronic
strong coupling scenario[14, 15]. Here, a long-range correction scaling as
$g^{2}_{0}\,R^{-3}$ enters besides a weak distance-independent contribution
scaling as $g^{4}_{0}$. The latter can be traced back to the dipole
approximation.

The second term in Eq.(82), $\tilde{u}^{(2)}_{AB}$, constitutes a cavity-
modified dipole-induced dipole interaction

$\displaystyle u^{(2)}_{AB}(R)$
$\displaystyle=-\dfrac{\tilde{C}_{6}}{R^{6}}+g^{2}_{0}\dfrac{\tilde{C}_{3\lambda}}{R^{3}}-g^{4}_{0}\dfrac{\tilde{C}_{0\lambda}}{R^{0}}\quad,$
(87)

which is determined by coefficients

$\displaystyle\tilde{C}_{6}$
$\displaystyle=2\left(\alpha^{0_{A}}_{\kappa\kappa}d^{0_{B}}_{\kappa}d^{0_{B}}_{\kappa}+d^{0_{A}}_{\kappa}d^{0_{A}}_{\kappa}\alpha^{0_{B}}_{\kappa\kappa}\right)\quad,\vspace{0.2cm}$
(88) $\displaystyle\tilde{C}_{3\lambda}$
$\displaystyle=2\left(\alpha^{0_{A}}_{\kappa\lambda}d^{0_{B}}_{\kappa}d^{0_{B}}_{\lambda}+d^{0_{A}}_{\kappa}d^{0_{A}}_{\lambda}\alpha^{0_{B}}_{\kappa\lambda}\right)\quad,\vspace{0.2cm}$
(89) $\displaystyle\tilde{C}_{0\lambda}$
$\displaystyle=\dfrac{1}{2}\left(\alpha^{0_{A}}_{\lambda\lambda}d^{0_{B}}_{\lambda}d^{0_{B}}_{\lambda}+d^{0_{A}}_{\lambda}d^{0_{A}}_{\lambda}\alpha^{0_{B}}_{\lambda\lambda}\right)\quad,$
(90)

and results from the two-body contributions in $E^{(2)}_{abb}$ (cf. Eqs.(37)
to (39)). At this point, we like to note that there is only a single first-
order contribution to the cavity-modified intermolecular interaction energy,
which is the distance-independent DSE-interaction term in Eq.(78) with
coefficient $\tilde{D}_{3\lambda}$. Thus, in contrast to intramolecular
corrections, second-order CBO-PT is crucial for the intermolecular interaction
scenario.

###  VI.2 Bimolecular Cavity Reaction Potential

In analogy to the unimolecular scenario in Sec.V.1, we shall now discuss the
cavity reaction potential for the interacting dimer with a focus on cavity-
modified intermolecular interactions. The minimizing cavity coordinate
obtained for $\mathcal{E}^{(2)}_{\mathrm{dimer}}$ in Eq.(72) is given by

$\displaystyle x^{(2)}_{0}$
$\displaystyle=x^{(2)}_{\mathrm{loc}}+x^{(2)}_{\mathrm{int}}\quad,$ (91)

with an intramolecular term

$\displaystyle x^{(2)}_{\mathrm{loc}}$
$\displaystyle=\dfrac{g_{0}}{\omega_{c}}\left(d^{0_{A}}_{\lambda}+d^{0_{B}}_{\lambda}\right)+\dfrac{g^{3}_{0}}{2\omega_{c}}\left(\alpha^{0_{A}}_{\lambda\lambda}d^{0_{A}}_{\lambda}+\alpha^{0_{B}}_{\lambda\lambda}d^{0_{B}}_{\lambda}\right)\quad,$
(92)

and an intermolecular distance-dependent contribution

$\displaystyle x^{(2)}_{\mathrm{int}}$
$\displaystyle=\dfrac{g_{0}\,D^{\prime}_{3}}{\omega_{c}}\dfrac{1}{R^{3}}+\dfrac{g^{3}_{0}\,D^{\prime}_{3}\left(\alpha^{0_{A}}_{\lambda\lambda}+\alpha^{0_{B}}_{\lambda\lambda}\right)}{\omega_{c}}\dfrac{1}{R^{3}}\quad.$
(93)

Those two expressions provide access to a cavity reaction potential accurate
up to terms scaling as $g^{4}_{0}$ in line with the second-order cPES (cf.
Eq.(72)). In the dilute-gas limit, the interaction contribution in Eq.(93)
vanishes

$\displaystyle x^{(2)}_{\mathrm{int}}\to 0\quad,\quad R\to\infty\quad,$ (94)

which leads with Eq.(72) to the corresponding bimolecular CRP

$\displaystyle\mathcal{V}^{(2)}_{\mathrm{crp}}$
$\displaystyle=\sum_{i=A,B}\left(E^{(e)}_{0_{i}}+\dfrac{g^{2}_{0}}{2}\tilde{\mathcal{F}}^{0_{i}}_{\lambda}-\dfrac{g^{4}_{0}}{8}\mathcal{B}^{0_{i}}_{\lambda\lambda}\right)+\mathcal{W}^{(2)}_{AB}\quad.$
(95)

Here, the summation runs simply over individual contributions of Monomer A and
B in line with the unimolecular result in Eq.(69). In addition, we find here
an interaction term

$\displaystyle\mathcal{W}^{(2)}_{AB}$
$\displaystyle=-g^{4}_{0}\dfrac{C_{0}+\tilde{C}_{0}}{R^{0}}\quad,$ (96)

which contains both the distance-independent van-der-Waals-type correction in
agreement with Ref.[15] and the related dipole-induced dipole correction in
Eq.(87).

In presence of distance-dependent intermolecular interactions,
$\mathcal{W}^{(2)}_{AB}$ turns into

$\displaystyle\mathcal{W}^{(2)}_{AB}(R)$
$\displaystyle=U^{(2)}_{AB}(R)+\mathcal{V}^{(2)}_{AB}(R)\quad,$ (97)

where the first term is given by Eq.(82) and the second term reads

$\mathcal{V}^{(2)}_{AB}(R)=-\dfrac{D_{3}}{R^{3}}-g^{2}_{0}\left(\dfrac{V_{3\lambda}}{R^{3}}+\dfrac{V_{6\lambda}}{2R^{6}}\right)\\\
-\dfrac{g^{4}_{0}}{2}\left(\dfrac{V^{\prime}_{3\lambda}}{R^{3}}+\dfrac{V^{\prime}_{6\lambda}}{R^{6}}\right)\quad,$
(98)

with coefficients

$\displaystyle V_{3\lambda}$
$\displaystyle=D^{\prime}_{3\lambda}\left(d^{0_{A}}_{\lambda}+d^{0_{B}}_{\lambda}\right)-D^{\prime\prime}_{3\lambda}\quad,\vspace{0.2cm}$
(99) $\displaystyle V^{\prime}_{3\lambda}$
$\displaystyle=D^{\prime}_{3\lambda}\left(\alpha^{0_{A}}_{\lambda\lambda}d^{0_{A}}_{\lambda}+\alpha^{0_{B}}_{\lambda\lambda}d^{0_{B}}_{\lambda}\right)\quad,\vspace{0.2cm}$
(100) $\displaystyle V_{6\lambda}$
$\displaystyle=D^{\prime\,2}_{3\lambda}\quad,\vspace{0.2cm}$ (101)
$\displaystyle V^{\prime}_{6\lambda}$
$\displaystyle=D^{\prime\,2}_{3\lambda}\left(\alpha^{0_{A}}_{\lambda\lambda}+\alpha^{0_{B}}_{\lambda\lambda}\right)\quad.$
(102)

The leading-order correction ($g^{2}_{0}$) in $\mathcal{V}^{(2)}_{AB}$
contains interactions scaling with distance as both $R^{-3}$ and $R^{-6}$,
which is in contrast to the second-order cPES, where $R^{-6}$-contributions
are independent of light-matter coupling. We moreover recall that
$\mathcal{W}^{(2)}_{AB}(R)$ describes cavity-modified intermolecular
interactions under the constraint that the total energy is minimized with
respect to cavity modes, i.e., static modifications of intermolecular
interaction potentials.

###  VI.3 Interaction Potentials under VSC

In Fig.(2), we present numerical results obtained via CBO-CCSD/aug-cc-pVDZ for
interaction potentials of H2 and H2O dimers under VSC. Those systems have
already been studied in earlier work on intermolecular interactions under
electronic strong coupling.[17, 35] In agreement with Refs.[17, 35], we
observe a cavity-induced destabilization of the interaction potential for
cavity mode polarizations parallel to the dissociation direction and a
stabilizing effect for orthogonal polarization. Cavity-induced modifications
of van-der-Waals interactions in the H2 dimer are significantly more
pronounced than the dipole-dipole dominated H2O interaction due to their
inherently weaker character. We note, that the H2O dimer has been investigated
at a higher coupling strength due to this reason to gain insight into the
qualitative nature of the cavity-modified interaction potential.

We close this section by a brief discussion of our results in the context of
findings obtained in the ESC regime[35, 14, 15]. The caviy-modified van-der-
Waals interaction is traced back to a perturbative second-order correction
emerging from dipole self-energy and the molecular dipole-dipole interaction.
Since the cavity modes are not explicitly contributing to this interaction, it
is equivalent for ab initio cavity QED and CBO frameworks employed to study
ESC and VSC regimes, respectively. This is no longer the case for interaction
contributions involving an explicit coupling to the cavity subsystem, which
can be understood as follows: In ab initio cavity QED as formulated for
electronic strong coupling problems[35], cavity modes are grouped with
electrons as fast degrees of freedom. This leads to a mixed fermionic-bosonic
TISE fully accounting for excitations in both electronic and cavity mode
subspaces. This is not the case in the CBO framework, where Eq.(2) is purely
electronic in character and excitations of the low-frequency cavity subsystem
are captured by Eq.(3). A connection between the CBO and ab initio cavity QED
approaches has been recently proposed by taking the low-frequency limit of the
latter.[15] We argue that this does not lead to formally equivalent
descriptions, which becomes apparent when one consistently takes the low-
frequency limit before invoking the adiabatic approximation. In this case, the
mixed electron-cavity TISE in ab initio cQED leads to electron-polariton
states with a fundamental excitation energy similar to the cavity mode
frequency as recently illustrate in a numerical study[12]. Therefore, an
adiabatic approximation might no longer be a valid since nonadiabatic coupling
elements can become significant. In addition, dynamical processes in the low-
frequency limit of ab initio cavity QED are by construction non-adiabatic in
nature since cavity mode excitations are accounted for in the electron-
polariton problem. We therefore conclude that the CBO framework is well suited
for the low-frequency cavity VSC regime and ab initio cavity QED successfully
addresses the high-frequency cavity ESC regime as already argued in Ref.[17]

##  VII Conclusion

We theoretically discussed the role of nonresonant electron-photon
interactions in dipole-coupled molecular ensembles under vibrational strong
coupling (VSC) with a focus on modifications of intra- and intermolecular
electronic interactions. Methodologically, we combined cavity Born-Oppenheimer
perturbation theory (CBO-PT) with ab initio wave function theory,
specifically, CBO Hartree-Fock (CBO-HF) and coupled cluster (CBO-CCSD)
theories, which allowed us to obtain analytical and numerical results for
cavity-modified ground state energies. Specifically, we discussed intra- and
intermolecular effects from the perspective of cavity reaction potentials
(CRP), which minimize the energy in the cavity subspace and satisfy an
approximate non-radiating ground state condition in CBO-PT. We explicitly
derived the CRP for a molecular monomer, which we applied to study cavity-
induced modifications of isomerization potentials under VSC, besides the CRP
for an interacting molecular dimer, which allowed us to examine VSC-modified
intermolecular interactions.

We find local electronic properties to be dominantly corrected by the first-
order CBO-PT dipole-fluctuation term, where we obtained excellent agreement
with non-perturbative mean-field (CBO-HF) and correlated (CBO-CCSD) results.
Higher-order perturbative corrections are identified as cavity-induced state
relaxation effects, i.e., changes in molecular orbital coefficients or cluster
amplitudes in the present study. For investigated torsion potentials of
chlorinated and fluorinated ethane derivatives, we do not find substantial
cavity-induced modifications of classical activation barriers and identify
state relaxation effects due to VSC to be negligible. We based our non-
perturbative study on a simple reformulation of the CBO-HF approach in spirit
of CRPs, which self-consistently accounts for the constrains determining the
minimizing cavity coordinate. This approach renders the related CBO-HF
equations a priori cavity-coordinate free and the resulting total energy
origin-independent.

In contrast to small first-order modifications of local electronic properties,
we obtain cavity-altered dipole-dipole, dipole-induced dipole and van-der-
Waals interactions at second-order of CBO-PT. Dipole-dipole interactions are
found to contain an explicit coupling to the cavity subsystem, which indicates
potentially altered energy transfer processes of an interacting molecular
dimer under VSC. Further, we obtained dipole-induced dipole and van-der-Waals
interactions with enhanced long-range character scaling as $R^{-3}$ and
$R^{0}$, i.e., the latter being distance-independent under the dipole
approximation. We complemented our analytical CBO-PT results via a numerical
non-perturbative analysis of selected dimer interaction potentials based on
CBO-CCSD. Van-der-Waals interactions are here observed to be altered most
easily under VSC due to their comparatively weak character, while all
interaction types show a dependence on cavity-mode polarization but are
independent of the cavity mode frequency.

The presented analysis naturally rises the question of how cavity-modified
intermolecular interactions in the VSC regime might be related to the
experimentally relevant collective VSC. Furthermore, small local modifications
of electronic interactions supports the hypothesis of a dynamical mechanism
underlying microscopic process in vibro-polaritonic chemistry. An integrated
approach accounting for both aspects in the spirit of very recent related
investigations[36, 37] might be a promising next step.

## Acknowledgements

E.W. Fischer acknowledges valuable discussions with Prof. Michael Roemelt and
funding by the Deutsche Forschungsgemeinschaft (DFG, German Research
Foundation) through DFG project 536826332.

## Appendix

###  .1 Details on CBO-PT(2) Corrections

Linear response functions in Eqs.(27) and (28) are explicitly given by

$\displaystyle\mathcal{A}^{0}_{\lambda\lambda^{\prime}}$
$\displaystyle=2\sum_{\mu\neq
0}\dfrac{\braket{\Psi_{0}}{\hat{d}_{\lambda}\hat{d}_{\lambda}}{\Psi_{\mu}}\braket{\Psi_{\mu}}{\hat{d}_{\lambda^{\prime}}}{\Psi_{0}}}{|\Delta^{(e)}_{\mu}|}\quad,\vspace{0.2cm}$
(.1.1) $\displaystyle\mathcal{B}^{0}_{\lambda\lambda^{\prime}}$
$\displaystyle=2\sum_{\mu\neq
0}\dfrac{\braket{\Psi_{0}}{\hat{d}_{\lambda}\hat{d}_{\lambda}}{,\Psi_{\mu}}\braket{\Psi_{\mu}}{\hat{d}_{\lambda^{\prime}}\hat{d}_{\lambda^{\prime}}}{\Psi_{0}}}{|\Delta^{(e)}_{\mu}|}\quad,$
(.1.2)

which can equivalently be written as

$\displaystyle\mathcal{A}^{0}_{\lambda\lambda^{\prime}}$
$\displaystyle=\alpha^{0}_{\lambda\lambda^{\prime}}d^{0}_{\lambda^{\prime}}+\tilde{\mathcal{A}}^{0_{a}}_{\lambda\lambda^{\prime}}\quad,\vspace{0.2cm}$
(.1.3) $\displaystyle\mathcal{B}^{0}_{\lambda\lambda^{\prime}}$
$\displaystyle=d^{0}_{\lambda}\alpha^{0}_{\lambda\lambda^{\prime}}d^{0}_{\lambda^{\prime}}+\tilde{\mathcal{B}}^{0_{a}}_{\lambda\lambda^{\prime}}\quad,$
(.1.4)

with fluctuation terms,
$\tilde{\mathcal{A}}^{0_{a}}_{\lambda\lambda^{\prime}}$ and
$\tilde{\mathcal{B}}^{0_{a}}_{\lambda\lambda^{\prime}}$. The latter
partitioning follows the same idea as Eq.(24), with fluctuation terms
resembling contributions from excited adiabatic states in analogy to the
dipole fluctuation term, $\tilde{\mathcal{F}}^{0}_{\lambda k}$. Mean-field
approximations are obtained by neglecting the fluctuation terms

$\displaystyle\mathcal{A}^{0}_{\lambda\lambda^{\prime}}$
$\displaystyle\approx\alpha^{0}_{\lambda\lambda^{\prime}}d^{0}_{\lambda^{\prime}}\,\Leftrightarrow\tilde{\mathcal{A}}^{0}_{\lambda\lambda^{\prime}}\approx
0\vspace{0.2cm}$ (.1.5)
$\displaystyle\mathcal{B}^{0}_{\lambda\lambda^{\prime}}$ $\displaystyle\approx
d^{0}_{\lambda}\alpha^{0}_{\lambda\lambda^{\prime}}d^{0}_{\lambda^{\prime}}\,\Leftrightarrow\tilde{\mathcal{B}}^{0}_{\lambda\lambda^{\prime}}\approx
0\quad.$ (.1.6)

###  .2 Ab initio Dipole-Fluctuation Correction

The dipole-fluctuation correction, $\tilde{\mathcal{F}}^{0}_{\lambda}$,
contains ground state expectation values of 1- and 2-electron terms. In in
first quantization, we have

$\braket{0^{(e)}}{\hat{d}^{(e)}_{\lambda}\hat{d}^{(e)}_{\lambda}}{0^{(e)}}=\braket{0^{(e)}}{\sum^{N_{e}}_{i}r^{2}_{i\lambda}}{0^{(e)}}\\\
+\braket{0^{(e)}}{\sum^{N_{e}}_{i\neq
j}r_{i\lambda}r_{j\lambda}}{0^{(e)}}\quad,$ (.2.1)

with polarization-projected electronic coordinates,
$r_{i\lambda}=\underline{e}_{\lambda}\cdot\underline{r}_{i}$. In second
quantization, the expectation value turns into

$\braket{0^{(e)}}{\hat{d}^{(e)}_{\lambda}\hat{d}^{(e)}_{\lambda}}{0^{(e)}}=\sum_{pq}O^{pq}_{\lambda}D_{pq}\\\
+\sum_{pqrs}d^{pq}_{\lambda}d^{rs}_{\lambda}d_{pqrs}\quad,$ (.2.2)

where $p,q,r,s$ run over the molecular orbital basis and 1-electron integrals
read

$\displaystyle
d^{pq}_{\lambda}=-\braket{\psi_{p}}{r_{1\lambda}}{\psi_{q}}\quad,\quad
O^{pq}_{\lambda}=\braket{\psi_{p}}{r^{2}_{1\lambda}}{\psi_{q}}\quad,$ (.2.3)

such that

$\displaystyle\tilde{\mathcal{F}}^{0}_{\lambda}$
$\displaystyle=\sum_{pq}O^{pq}_{\lambda}D_{pq}+\sum_{pqrs}d^{pq}_{\lambda}d^{rs}_{\lambda}\left(d_{pqrs}-D_{pq}D_{rs}\right)\,.$
(.2.4)

In case of restricted Hartree-Fock, 1- and 2-particle reduced density
matrices, $D_{pq}$ and $d_{pqrs}$, simplify as

$\displaystyle D_{pq}$ $\displaystyle=2\delta_{pq}\quad,\vspace{0.2cm}$ (.2.5)
$\displaystyle d_{pqrs}$
$\displaystyle=D_{pq}D_{rs}-\dfrac{1}{2}D_{ps}D_{rq}\quad,$ (.2.6)

which leads to

$\displaystyle\sum_{pq}O^{pq}_{\lambda}D_{pq}$
$\displaystyle=2\sum_{p}O^{pp}_{\lambda}\quad,\vspace{0.2cm}$ (.2.7)
$\displaystyle\sum_{pqrs}d^{pq}_{\lambda}d^{rs}_{\lambda}d_{pqrs}$
$\displaystyle=4\sum_{pr}d^{pp}_{\lambda}d^{rr}_{\lambda}-2\sum_{pr}d^{pr}_{\lambda}d^{rp}_{\lambda}\quad,\vspace{0.2cm}$
(.2.8) $\displaystyle\sum_{pq}d^{pq}_{\lambda}D_{pq}$
$\displaystyle=2\sum_{p}d^{pp}_{\lambda}\quad.$ (.2.9)

Thus, the dipole fluctuation correction turns into

$\tilde{\mathcal{F}}^{0}_{\lambda}=2\sum_{p}O^{pp}_{\lambda}+4\sum_{pr}d^{pp}_{\lambda}d^{rr}_{\lambda}-2\sum_{pr}d^{pr}_{\lambda}d^{rp}_{\lambda}\\\
-4\sum_{pr}d^{pp}_{\lambda}d^{rr}_{\lambda}\quad,$ (.2.10)

where the forth term stemming from the product of dipole expectation values
cancels the second term, which finally gives

$\displaystyle\tilde{\mathcal{F}}^{\mathrm{rhf}}_{\lambda}$
$\displaystyle=2\sum_{p}O^{pp}_{\lambda}-2\sum_{pr}d^{pr}_{\lambda}d^{rp}_{\lambda}\quad,$
(.2.11)

in MO representation.

## References

  * [1] Hirai, K.; Hutchison, J.A.; Uji-i, H. Recent progress in vibropolaritonic chemistry. ChemPlusChem 85, 1981 (2020).
  * [2] Nagarajan, K.; Thomas, A.; Ebbesen, T.W. Chemistry under vibrational strong coupling. J. Am. Chem. Soc. 143, 16877 (2021).
  * [3] Dunkelberger, A.D.; Simpkins, B.S.; Vurgaftman, I.; Owrutsky, J.C. Vibration-cavity polariton chemistry and dynamics. Ann. Rev. Phys. Chem. 73, 429 (2022).
  * [4] Ebbesen, T.W. Hybrid light-matter states in a molecular and material science perspective. Acc. Chem. Res. 49, 2403 (2016).
  * [5] Shalabney, A.; George, J.; Hutchison, J.A.; Pupillo, G.; Genet, C.; Ebbesen, T.W. Coherent coupling of molecular resonators with a microcavity mode. Nat. Commun. 6, 5981 (2015).
  * [6] George, J.; Shalabney, A.; Hutchison, J.A.; Genet, C.; Ebbesen, T.W. Liquid-phase vibrational strong coupling. J. Phys. Chem. Lett. 6, 1027 (2015).
  * [7] Long, J.P.; Simpkins, B.S. Coherent coupling between a molecular vibration and Fabry-Perot optical cavity to give hybridized states in the strong coupling limit. ACS Photonics 2, 130 (2015).
  * [8] Thomas, A.; George, J.; Shalabney, A.; Dryzhakov, M.; Varma, S.J.; Moran, J.; Chervy, T.; Zhong, X.; Devaux, E.; Genet, C.; Hutchison, J.A.; Ebbesen, T.W. Ground-state chemical reactivity under vibrational coupling to the vacuum electromagnetic field. Angew. Chem. Int. Ed. 55, 11462 (2016).
  * [9] Thomas, A.; Lethuillier-Karl, L.; Nagarajan, K.; Vergauwe, R.M.A.; George, J.; Chervy, T.; Shalabney, A.; Devaux, E.; Genet, C.; Moran, J.; Ebbesen, T.W. Tilting a ground-state reactivity landscape by vibrational strong coupling. Science 363, 615 (2019).
  * [10] Flick, J.; Ruggenthaler, M.; Appel, H.; Rubio, A. Atoms and molecules in cavities, from weak to strong coupling in quantum-electrodynamics (QED) chemistry. PNAS 114, 3026 (2017).
  * [11] Flick, J.; Appel, H.; Ruggenthaler, M.; Rubio, A. Cavity Born-Oppenheimer approximation for correlated electron-nuclear-photon systems. J. Chem. Theory Comput. 13, 1616, (2017).
  * [12] Fischer, E.W.; Saalfrank, P. Beyond Cavity Born-Oppenheimer: On Nonadiabatic Coupling and Effective Ground State Hamiltonians in Vibro-Polaritonic Chemistry. J. Chem. Theory Comput. 19, 7215 (2023).
  * [13] Fischer, E.W.; Syska, J.A.; Saalfrank, P. A Quantum Chemistry Approach to Linear Vibro-Polaritonic Infrared Spectra with Perturbative Electron-Photon Correlation. J. Phys. Chem. Lett. 15, 2262 (2024).
  * [14] Philbin, J.P.; Haugland, T.S.; Ghosh, T.K.; Ronca, E.; Chen, M.; Narang, P.; Koch, H. Molecular van der Waals Fluids in Cavity Quantum Electrodynamics. J. Phys. Chem. Lett. 14, 8988 (2023).
  * [15] Haugland, T.S.; Philbin, J.P.; Ghosh, T.K.; Chen, M.; Koch, H.; Narang, P. Understanding the polaritonic ground state in cavity quantum electrodynamics. arXiv preprint, arXiv:2307.14822 (2023).
  * [16] Schnappinger, T.; Sidler, D.; Ruggenthaler, M.; Rubio, A., Kowalewski, M. Cavity Born-Oppenheimer Hartree-Fock Ansatz: Light-Matter Properties of Strongly Coupled Molecular Ensembles. J. Phys. Chem. Lett. 14, 8024 (2023).
  * [17] Angelico, S.; Haugland, T.S.; Ronca, E.; Koch, H. Coupled cluster cavity Born-Oppenheimer approximation for electronic strong coupling. J. Chem. Phys. 159, 214112 (2023).
  * [18] Fischer, E.W.; Anders, J.; Saalfrank, P. Cavity-altered thermal isomerization rates and dynamical resonant localization in vibro-polaritonic chemistry. J. Chem. Phys. 156, 154305 (2022).
  * [19] Fischer, E.W.; Saalfrank, P. Ground state properties and infrared spectra of anharmonic vibrational polaritons of small molecules in cavities. J. Chem. Phys. 154, 104311 (2021).
  * [20] Bonini, J.; Flick, J. Ab Initio Linear-Response Approach to Vibro-Polaritons in the Cavity Born-Oppenheimer Approximation. J. Chem. Theory Comput. 18, 2764 (2022).
  * [21] Li, X.; , A.; Huo, P. Cavity frequency-dependent theory for vibrational polariton chemistry. Nat. Commun. 12, 1315 (2021).
  * [22] Sun, J.; Vendrell, O. Suppression and enhancement of thermal chemical rates in a cavity. J. Phys. Chem. Lett. 13, 4441 (2022).
  * [23] Schäfer, C.; Flick, J., Ronca, E.; Narang, P.; Rubio, A. Shining light on the microscopic resonant mechanism responsible for cavity-mediated chemical reactivity. Nat. Commun. 13, 7817 (2022).
  * [24] Fischer, E.W.; Saalfrank, P. Cavity-catalyzed hydrogen transfer dynamics in an entangled molecular ensemble under vibrational strong coupling. Phys. Chem. Chem. Phys. (2023).
  * [25] Lindoy, L.P.; Mandal, A.; Reichman, D.R. Quantum dynamical effects of vibrational strong coupling in chemical reactivity. Nat. Commun. 14, 2733 (2023).
  * [26] Ying, W.; H. Pengfei. Resonance theory and quantum dynamics simulations of vibrational polariton chemistry. J. Chem. Phys. 159, 084104 (2023).
  * [27] Yaling, K.; Richardson, J.O. Insights into the mechanisms of optical cavity-modified ground-state chemical reactions. J. Chem. Phys. 160, 224704 (2024).
  * [28] Schäfer, C.; Ruggenthaler, M.; Rokaj, V.; Rubio, A. Relevance of the quadratic diamagnetic and self-polarization terms in cavity quantum electrodynamics. ACS Photonics 7, 975 (2020).
  * [29] Neese, F. The ORCA program system. WIREs Comput. Mol. Sci. 2, 73 (2012).
  * [30] Neese, F. Software update: The ORCA program system-Version 5.0. WIREs Comput. Mol. Sci. 5, e1606 (2022).
  * [31] Galego, J.; Climent, C; Garcia-Vidal, F.J.; Feist, J. Cavity Casimir-Polder forces and their effects in ground-state chemical reactivity.  Phys. Rev. X 9, 021057 (2019).
  * [32] Rokaj, V.; Welakuh, D.M.; Ruggenthaler, M.; Rubio, A. Light-matter interaction in the long-wavelength limit: no ground-state without dipole self-energy. J. Phys. B: At. Mol. Opt. Phys. 51, 034005 (2018).
  * [33] Q. Sun, T. C. Berkelbach, N. S. Blunt, G. H. Booth, S. Guo, Z. Li, J. Liu, J. McClain, S. Sharma, S. Wouters, and G. K.-L. Chan. PySCF: the Python-based simulations of chemistry framework. WIREs Comput. Mol. Sci. 8, e1340 (2018).
  * [34] Q. Sun, X. Zhang, S. Banerjee, P. Bao, M. Barbry, N. S. Blunt, N. A. Bogdanov, G. H. Booth, J. Chen, Z.-H. Cui, J. J. Eriksen, Y. Gao, S. Guo, J. Hermann, M. R. Hermes, K. Koh, P. Koval, S. Lehtola, Z. Li, J. Liu, N. Mardirossian, J. D. McClain, M. Motta, B. Mussard, H. Q. Pham, A. Pulkin, W. Purwanto, P. J. Robinson, E. Ronca, E. R. Sayfutyarova, M. Scheurer, H. F. Schurkus, J. E. T. Smith, C. Sun, S.-N. Sun, S. Upadhyay, L. K. Wagner, X. Wang, A. White, J. Daniel Whitfield, M. J. Williamson, S. Wouters, J. Yang, J. M. Yu, T. Zhu, T. C. Berkelbach, S. Sharma, A. Yu. Sokolov, and G. K.-L. Chan. Recent developments in the PySCF program package. J. Chem. Phys. 153, 024109 (2020).
  * [35] Haugland, T.; Schäfer, C.; Ronca, E.; Rubio, A.; Koch, H. Intermolecular interactions in optical cavities: An ab initio QED study. J. Chem. Phys. 154, 094113 (2021).
  * [36] Sidler, D.; Schnappinger, T.; Obzhirov, A.; Ruggenthaler, R.; Kowalewski, M.; Rubio, A. Unraveling a Cavity-Induced Molecular Polarization Mechanism from Collective Vibrational Strong Coupling. J. Phys. Chem. Lett. 15, 5208 (2023).
  * [37] Castagnola, M.; Haugland, T.S.; Ronca, E.; Koch, H.; Schäfer, C. Collective strong coupling modifies aggregation and solvation. J. Phys. Chem. Lett. 15, 1428 (2024).