11.5.1 Theoretical Background

The total energy of the system consists of the interaction energy of the effective fragments ($E^{\mathrm{ef}-\mathrm{ef}}$) and the energy of the ab initio (i.e., QM) region in the field of the fragments. The former includes electrostatics, polarization, dispersion and exchange-repulsion contributions (the charge-transfer term, which might be important for description of the ionic and highly polar species, is omitted in the current implementation):

$$E^{\mathrm{ef}\text{-}\mathrm{ef}}=E_{\mathrm{elec}}+E_{\mathrm{pol}}+E_{\mathrm{disp}}+E_{\mathrm{ex}\text{-}\mathrm{rep}}.$$

(11.49)

The QM-EF interactions are computed as follows. The electrostatics and polarization parts of the EFP potential contribute to the quantum Hamiltonian via one-electron terms,

$$H_{pq}^{\prime }=H_{pq}+⟨p|{\widehat{V}}^{\mathrm{elec}}+{\widehat{V}}^{\mathrm{pol}}|q⟩$$

(11.50)

whereas dispersion and exchange-repulsion QM-EF interactions are treated as additive corrections to the total energy.

The electrostatic component of the EFP energy accounts for Coulomb interactions. In molecular systems with hydrogen bonds or polar molecules, this is the leading contribution to the total intermolecular interaction energy.¹⁰⁹ An accurate representation of the electrostatic potential is achieved by using multipole expansion (obtained from the Stone’s distributed multipole analysis) around atomic centers and bond midpoints (i.e., the points with high electronic density) and truncating this expansion at octupoles.^{924, 923, 208, 318} The fragment-fragment electrostatic interactions consist of charge-charge, charge-dipole, charge-quadrupole, charge-octupole, dipole-dipole, dipole-quadrupole, and quadrupole-quadrupole terms, as well as terms describing interactions of electronic multipoles with the nuclei and nuclear repulsion energy.

Electrostatic interaction between an effective fragment and the QM part is described by perturbation ${\widehat{V}}^{\mathrm{elec}}$ of the ab initio Hamiltonian (see Eq.  (11.50)). The perturbation enters the one-electron part of the Hamiltonian as a sum of contributions from the expansion points of the effective fragments. Contribution from each expansion point consists of four terms originating from the electrostatic potential of the corresponding multipole (charge, dipole, quadrupole, and octupole).

The multipole representation of the electrostatic density of a fragment breaks down when the fragments are too close. The multipole interactions become too repulsive due to significant overlap of the electronic densities and the charge-penetration effect. The magnitude of the charge-penetration effect is usually around 15% of the total electrostatic energy in polar systems, however, it can be as large as 200% in systems with weak electrostatic interactions.⁸⁸⁸ To account for the charge-penetration effect, the simple exponential damping of the charge-charge term is used.^{266, 888} The charge-charge screened energy between the expansion points $k$ and $l$ is given by the following expression, where ${\alpha }_k$ and ${\alpha }_l$ are the damping parameters associated with the corresponding expansion points:

	$E_{kl}^{\mathrm{ch}\text{-}\mathrm{ch}}$	$=$	$\left[1-(1+{\alpha }_k{R}_{kl}/2)e^{-{\alpha }_k{R}_{kl}}\right]q^k{q}^l/R_{kl},\text{if }{\alpha }_k={\alpha }_l$		(11.51)
	or	$=$	$\left(1-{\displaystyle \frac{{\alpha }_l^2}{{\alpha }_l^2-{\alpha }_k^2}}{e}^{-{\alpha }_k{R}_{kl}}-{\displaystyle \frac{{\alpha }_k^2}{{\alpha }_k^2-{\alpha }_l^2}}{e}^{-{\alpha }_l{R}_{kl}}\right)q^k{q}^l/R_{kl},\text{if }{\alpha }_k\ne {\alpha }_l$		(11.52)

Damping parameters are included in the potential of each fragment, but QM-EFP electrostatic interactions are currently calculated without damping corrections.

Alternatively, one can obtain the short-range charge-penetration energy using the spherical Gaussian overlap (SGO) approximation:⁸⁸⁹

$$E_{kl}^{\mathrm{pen}}=-2{\left(\frac{1}{-2ln|S_{kl}|}\right)}^{\frac{1}{2}}\frac{S_{kl}^2}{R_{kl}}$$

(11.53)

where $S_{kl}$ is the overlap integral between localized MOs $k$ and $l$, calculated for the exchange-repulsion term, Eq. (11.65). This charge-penetration energy is calculated and printed separately from the rest of the electrostatic energy. Using overlap-based damping generally results in a more balanced description of intermolecular interactions and is recommended.

Polarization accounts for the intramolecular charge redistribution in response to external electric field. It is the major component of many-body interactions responsible for cooperative molecular behavior. EFP employs distributed polarizabilities placed at the centers of valence LMOs. Unlike the isotropic total molecular polarizability tensor, the distributed polarizability tensors are anisotropic.

The polarization energy of a system consisting of an ab initio and effective fragment regions is computed as²⁰⁸

$$E^{\mathrm{pol}}=-\frac{1}{2}\sum _k{\mu }^k(F^{\mathrm{mult},k}+F^{\mathrm{ai},\mathrm{nuc},k})+\frac{1}{2}\sum _k{\overline{\mu }}^k{F}^{\mathrm{ai},\mathrm{elec},k}$$

(11.54)

where ${\mu }^k$ and ${\overline{\mu }}^k$ are the induced dipole and the conjugated induced dipole at the distributed point $k$; $F^{\mathrm{mult},k}$ is the external field due to static multipoles and nuclei of other fragments, and $F^{\mathrm{ai},\mathrm{elec},k}$ and $F^{\mathrm{ai},\mathrm{nuc},k}$ are the fields due to the electronic density and nuclei of the ab initio part, respectively.

The induced dipoles at each polarizability point $k$ are computed as

$${\mu }^k={\alpha }^k{F}^{\mathrm{total},k}$$

(11.55)

where ${\alpha }^k$ is the distributed polarizability tensor at $k$. The total field $F^{\mathrm{total},k}$ comprises from the static field and the field due to other induced dipoles, $F_k^{\mathrm{ind}}$, as well as the field due to nuclei and electronic density of the ab initio region:

$$F^{\mathrm{ai},\mathrm{total},k}=F^{\mathrm{mult},k}+F^{\mathrm{ind},k}+F^{\mathrm{ai},\mathrm{elec},k}+F^{\mathrm{ai},\mathrm{nuc},k}$$

(11.56)

As follows from the above equation, the induced dipoles on a particular fragment depend on the values of the induced dipoles of all other fragments. Moreover, the induced dipoles on the effective fragments depend on the ab initio electronic density, which, in turn, is affected by the field created by these induced dipoles through a one electron contribution to the Hamiltonian:

$${\widehat{V}}^{\mathrm{pol}}=-\frac{1}{2}\sum _k\sum _a^{x,y,z}\frac{({\mu }_a^k+{\overline{\mu }}_a^k)a}{R^3}$$

(11.57)

where $R$ and $a$ are the distance and its Cartesian components between an electron and the polarizability point $k$. In sum, the total polarization energy is computed self-consistently using a two level iterative procedure. The objectives of the higher and lower levels are to converge the wave function and induced dipoles for a given fixed wave function, respectively. In the absence of the ab initio region, the induced dipoles of the EF system are iterated until self-consistent with each other.

Self-consistent treatment of polarization accounts for many-body interaction effects. Polarization energy between EFP fragments is augmented by gaussian-like damping functions with default parameter $\alpha =\beta =0.6$, applied to electric field $F$:⁸⁸⁹

$$F=F_0{f}^{\mathrm{damp}}$$

(11.58)

$$f^{\mathrm{damp}}=1.0-exp(-\sqrt{\alpha \beta }{r}^2)(1+\sqrt{\alpha \beta }{r}^2)$$

(11.59)

Dispersion provides a leading contribution to van der Waals and $\pi$-stacking interactions. The dispersion interaction is expressed as the inverse $R$ dependence:

$$E^{\mathrm{disp}}=\sum _n{C}_6{R}^{-6}$$

(11.60)

where coefficients $C_6$ are derived from the frequency-dependent distributed polarizabilities with expansion points located at the LMO centroids, i.e., at the same centers as the polarization expansion points. The higher-order dispersion terms (induced dipole-induced quadrupole, induced quadrupole/induced quadrupole, etc.) are approximated as $1/3$ of the $C_6$ term.⁶

For small distances between effective fragments, dispersion interactions are corrected for charge penetration and electronic density overlap effect either with the Tang-Toennies damping formula⁹⁴⁸ with parameter $b=1.5$,

$$C_6^{kl}\to \left(1-e^{-bR}\sum _{k=0}^6\frac{{(bR)}^k}{k!}\right)C_6^{kl},$$

(11.61)

or else using interfragment overlap (so-called overlap-based damping):⁸⁸⁹

$$C_6^{kl}\to \left(1-S_{kl}^2\left(1-2\log |S_{kl}|+2{\log }^2|S_{kl}|\right)\right)C_6^{kl}$$

(11.62)

QM-EFP dispersion interactions are currently disabled.

Exchange-repulsion originates from the Pauli exclusion principle, which states that the wave function of two identical fermions must be anti-symmetric. In traditional classical force fields, exchange-repulsion is introduced as a positive (repulsive) term, e.g., $R^{-12}$ in the Lennard-Jones potential. In contrast, EFP uses a wave function-based formalism to account for this inherently quantum effect. Exchange-repulsion is the only non-classical component of EFP and the only one that is repulsive.

The exchange-repulsion interaction is derived as an expansion in the intermolecular overlap, truncated at the quadratic term,^{430, 431} which requires that each effective fragment carries a basis set that is used to calculate overlap and kinetic one-electron integrals for each interacting pair of fragments. The exchange-repulsion contribution from each pair of localized orbitals $i$ and $j$ belonging to fragments $A$ and $B$, respectively, is:

	$E_{ij}^{\mathrm{exch}}$	$=$	$-4\sqrt{{\displaystyle \frac{-2\ln |S_{ij}|}{\pi }}}{\displaystyle \frac{S_{ij}^2}{R_{ij}}}$		(11.65)
			$-2S_{ij}\left({\displaystyle \sum _{k\in A}}{F}_{ik}^A{S}_{kj}+{\displaystyle \sum _{l\in B}}{F}_{jl}^B{S}_{il}-2T_{ij}\right)$
			$+2S_{ij}^2\left({\displaystyle \sum _{J\in B}}-Z_J{R}_{iJ}^{-1}+2{\displaystyle \sum _{l\in B}}{R}_{il}^{-1}+{\displaystyle \sum _{I\in A}}-Z_I{R}_{Ij}^{-1}+2{\displaystyle \sum _{k\in A}}{R}_{kj}^{-1}-R_{ij}^{-1}\right)$

where $i$, $j$, $k$ and $l$ are the LMOs, $I$ and $J$ are the nuclei, $S$ and $T$ are the intermolecular overlap and kinetic energy integrals, and $F$ is the Fock matrix element.

The expression for the $E_{ij}^{\mathrm{exch}}$ involves overlap and kinetic energy integrals between pairs of localized orbitals. In addition, since Eq. (11.65) is derived within an infinite basis set approximation, it requires a reasonably large basis set to be accurate [6-31+G* is considered to be the smallest acceptable basis set, 6-311++G(3df,2p) is recommended]. These factors make exchange-repulsion the most computationally expensive part of the EFP energy calculations of moderately sized systems.

Large systems require additional considerations. Since total exchange-repulsion energy is given by a sum of terms in Eq. (11.65) over all the fragment pairs, its computational cost formally scales as $𝒪(N^2)$ with the number of effective fragments $N$. However, exchange-repulsion is a short-range interaction; the overlap and kinetic energy integrals decay exponentially with the inter-fragment distance. Therefore, by employing a distance-based screening, the number of overlap and kinetic energy integrals scales as $𝒪(N)$. Consequently, for large systems exchange-repulsion may become less computationally expensive than the long-range components of EFP (such as Coulomb interactions).

The QM-EFP exchange-repulsion energy is currently disabled.