Predicting the Post‐Hartree‒Fock Electron Correlation Energy of Complex Systems with the Information‐Theoretic Approach

1. Introduction

Electron correlation energy lies at the heart of quantum chemistry. [1,2] However, the computational cost of high-level post-Hartree–Fock methods skyrockets with system size. In this context, there is a pressing need for alternative lower-scaling cost-efficient methods across broad classes of systems. In recent years, the information-theoretic approach (ITA) [3,4,5,6] has emerged as a promising framework for understanding and predicting the electron correlation energy from the perspective of information theory. By treating the electron density as a continuous probability distribution, ITA introduces a set of descriptors—such as Shannon entropy [7] and Fisher information [8]—that encode global and local features of the electron density distribution. These quantities are inherently basis-agnostic and physically interpretable, providing a new lens through which quantum chemical problems can be approached.

In continuation with our previous work by employing the simple density-based ITA quantities to appreciate response properties [9,10,11,12,13] (such as molecular polarizability and NMR chemical shielding constant) and energetics of elongated hydrogen chains, [14] in this work, we aim to predict the post-Hartree–Fock (see Figure 1) electron correlation energies of various molecular clusters and linear or quasi-linear organic polymers with increasing cluster size and polymer length. These systems including 24 octane isomers (see Figure 2); [11,15] (ii) polymeric structures (see Figure 3), polyyne, polyene, all-trans-polymethineimine, and acene; [11] (iii) molecular clusters (see Figure 4), such as metallic Be_nand Mg_n, [16,17] covalent S_n, [18,19] hydrogen-bonded protonated water clusters H⁺(H₂O)_n, [20] and dispersion-bound carbon dioxide (CO₂)_n, [21] and benzene clusters (C₆H₆)_n. [22] We construct strong linear relationships between the low-cost Hartree–Fock [23] ITAs and the electron correlation energies from post-Hartree–Fock methods, such as MP2 or RI-MP2, [24,25] CCSD, [26] and CCSD(T). [27] It is noteworthy to mention that MP2 is mainly used here only as a proof-of-concept; Hartree–Fock can be simply replaced with any approximate functionals of density functional theory (DFT) [28,29].

By examining trends across increasing cluster size and polymer length, we assess the transferability, scalability, and physical insights provided by ITA features in capturing electron correlation. Our findings highlight not only the feasibility of ITA-driven correlation energy prediction, but also reveal key descriptors that most strongly govern correlation effects in extended systems. These results suggest that ITA may serve as a promising direction for developing efficient, interpretable, and physically grounded models in quantum machine learning and electronic structure theory.

2. Results

To validate the accuracy of the LR(ITA) method, we choose a total of 24 octane isomers as shown in Figure 2. MP2, CCSD, and CCSD(T) are used to generate the electron correlation energies and ITA quantities are obtained at the Hartree‒Fock level at the same basis set 6-311++G(d,p) level. More details can be found in the Supplementary Information (SI, Table S1). Table 1 shows the linear relationships and RMSDs between the LR(ITA)-predicted and calculated electron correlation energies. For S_S, I_F, and S_GBP, the RMSDs are < 2.0 mH, indicating that LR(ITA) should be accurate enough to predict the electron correlation energies. Because CCSD and CCSD(T) are too computationally-intensive and intractable, only MP2 is used hereafter as proof-of-concept.

In Table 2, Table 3, Table 4 and Table 5, we have collected the linear correlation coefficients (R² = 1.000) and RMSDs (root mean squared deviations) between the calculated correlation energies at the MP2/6-311++G(d,p) level and those predicted based on the ITA quantities at the HF/6-311++G(d,p) level for polyyne, polyene, all-trans-polymethineimine, and acene, respectively. More details can be found in Tables S2–S5. It is clearly showcased that R² is close to 1 for most ITA quantities. More strikingly, based on the linear regression (LR) equations of ITA quantities, the predicted electron correlations deviate from the calculated ones only by ~1.5 mH for polyyne, ~3.0 mH for polyene, and < 4.0 mH for all-trans-polymethineimine. For acene, the RMSDs are reasonably satisfactory by ~ 10 ‒ 11 mH. These results collectively reveal that ITA quantities are indeed good descriptors of electron correlations for those linear or quasi-linear polymeric systems with delocalized electronic structures. For more challenging acenes, a single ITA quantity fails to capture sufficient amount of information of more delocalized electronic structures.

Shown in Table 6, Table 7 and Table 8 are the results of the linear correlation coefficients (R²) and RMSDs (root mean squared deviations) between the calculated correlation energies at the MP2/6-311++G(d,p) level and those predicted based on the ITA quantities at the HF/6-311++G(d,p) level for neutral metallic Be_n, Mg_n, and covalent S_n systems, respectively. More details can be found in Tables S6–S11. One can see that strong correlations exist (R² > 0.990) between ITA quantities and MP2 correlation energies, indicating that they are extensive in nature. However, the predicted electron correlation energies deviate much from the calculated ones by ~28 ‒ 37 mH for Be_n, ~17 ‒ 33 mH for Mg_n, and ~26 ‒ 42 mH for S_n, respectively. These results collectively showcase that for 3-dimensional metallic clusters, Be_n and Mg_n, and covalent S_n, a single ITA quantity fails to quantitively capture enough amount information of electron energies of complex systems.

Shown in Table 9 are the results of the linear correlation coefficients (R²), the corresponding regression coefficients, and RMSDs (root mean squared deviations) between the calculated correlation energies at the MP2/6-311++G(d,p) level and those predicted based on the ITA quantities at the HF/6-311++G(d,p) level for hydrogen-bonded protonated water clusters. One can see that strong correlations exist (R² = 1.000) between (8 out of 11) ITA quantities and MP2 correlation energies, indicating that they are extensive in nature. The RMSDs range from 2.1 ($E_2$ and $E_3$) to 9.3 ($G_3$) mH, indicating that ITA quantities are good descriptors of the post-Hartree‒Fock electron correlation energies of hydrogen-bonded systems.

Finally, we will switch our gear to two dispersion-bound clusters, (CO₂)_n and (C₆H₆)_n. Table 10 gives the strong correlations (R² = 1.000) and RMSDs between the RI-MP2 correlation energies and Hartree‒Fock ITA quantities at the same basis set 6-311++G(d,p) for (CO₂)_n(n = 4 ‒ 40). More details can be found in Table S12. The RMSDs vary from 6.3 ($E_2$ and $E_3$) to 10.8 ($G_3$) to 14.6 ($S_{\mathrm{S}}$) mH. For (C₆H₆)_n (n = 4 ‒ 14) clusters, we have calculated the linear correlations (R² = 1.000) and RMSDs between the MP2/6-311++G(d,p) electron correlation energies and HF/6-311++G(d,p) ITA quantities, as collected in Table 11. More details can be found in Tables S13 and S14. The RMSDs range from 2.8 ($G_3$) to 6.9 ($E_3$) to 10.7 ($S_{\mathrm{S}}$) mH. The RMSD results collectively suggest (8 out of 11) ITA quantities are reasonably good descriptors of the post-Hartree‒Fock electron correlation energies of dispersion-bound clusters.

To further verify the accuracy of the LR(ITA) method, we employ some relatively larger (C₆H₆)_n (n = 15 ‒ 30) clusters to this end. Plus, conventional MP2/6-311++G(d,p) calculations are too computationally-intensive, we employ GEBF [30,31,32,33] to obtain the MP2-level electron correlation energies as reference. Finally, as the linear regression based on the ITA quantity $G_3$ has the least RMSD value, we choose LR(G₃) to make predictions of electron correlation energies of benzene clusters. More details can be found in Tables 15 and 16. Figure 5 shows a comparison of the LR(G₃)-predicted and GEBF-calculated MP2 electron correlation energies for benzene clusters. The RMSD is 8.6 mH, indicative that the LR(ITA) method has a comparable performance to the linear-scaling GEBF method. In addition, we have found that when subsystem wavefunctions (thus electron density and ITA quantities) are used to obtain the subsystem electron correlation energies, the final total electron correlation energies deviate from GEBF by 40.0 mH in terms of RMSD as shown in Table S16. One possible for this observation may come from the error accumulation, rather than error cancellation.

3. Discussion

To accurately and efficiently predict the post-Hartree‒Fock electron correlation energy at a relatively low cost is a hot area in the community of quantum chemistry. Starting from Hartree‒Fock molecular orbitals, there exist two typical methods. One is to calculate the local correlation energy, whose early development is due to Pulay and Sæbø; [34,35,36] the other is to predict the correlation energy with the aid of deep learning (DL). [37,38,39,40,41,42,43,44,45,46] Our proposed LR(ITA) method is a special favor of DL. Suffice to note that an inherent drawback of local correlation methods is to perform orbital localization. This problem is also encountered by the DL-driven method. For our LR(ITA) method, only the molecular orbitals (thus electron density) are required without any manipulation. Very recently, we have showcased the good accuracy of LR(ITA) and its variant DL(ITA). With LR(ITA), one can even predict the FCI-level electron correlation with the DMRG (density matrix renormalization group) algorithm as a solver for elongated hydrogen chain, [14] and the RMSD is only a few mH. Moreover, with DL(ITA) where a total 11 ITA quantities are used as input, [13] we have predicted the DLPNO-MP2 electron correlation energy for a database of > 90 K real organic molecules and the RMSD is about 6.8 mH. In addition, LR(ITA) is not limited to any post-Hartree‒Fock electronic structure methods; MP2 is used here as a proof-of-concept. Thus, we have showcased that LR(ITA) is designed with architectural and conceptual simplicity and is numerically shown to be a good protocol to predict the electron correlation energies of various systems.

Admittedly, using LR(ITA) to accurately and efficiently predict the electron correlation energy is still in its infancy. In the near future, we will implement a new concept of “ITL-DL Loop”. The physics behind is simple: low-tier (such as semiempirical PM7 or even promolecular) electron densities are used as input for ITA quantities, and DL is introduced to obtain high-tier (such as DFT) electron densities. Based on the newly generated electron densities, ITA quantities are obtained and used as input for another either classical or quantum DL model to predict the electron correlation energies of electrons of physicochemical properties of molecules. Work along this line is in progress and the results will be presented elsewhere.

4. Materials and Methods

4.1. Information-Theoretic Approach Quantities

Shannon entropy $S_{\mathrm{S}}$ [7] and Fisher information $I_{\mathrm{F}}$ [8] are two foundational quantities in information theory. They are defined as Equations (1) and (2), respectively.

$$S_{\mathrm{S}}=-\int \rho \left(r\right)\ln \rho \left(r\right)dr$$

(1)

$$I_{\mathrm{F}}=\int {\displaystyle \frac{{\left|\nabla \rho \left(r\right)\right|}^2}{\rho \left(r\right)}}dr$$

(2)

where $\rho \left(r\right)$ is the electron density and $\nabla \rho \left(r\right)$is the density gradient. Physically, $S_{\mathrm{S}}$ characterizes the spatial delocalization of the electron density, while $I_{\mathrm{F}}$ reflects its sharpness or localization. Of note, $S_{\mathrm{S}}$ and $I_{\mathrm{F}}$ are not mutually exclusive and but always intercorrelated.

Beyond the total electron density, additional quantities such as kinetic-energy density can be incorporated into the formulation of information-theoretic approaches (ITA). Utilizing both electron density and kinetic-energy density, Ghosh, Berkowitz, and Parr introduced an entropy functional known as ($S_{\mathrm{GBP}}$), [47]

$$S_{\mathrm{GBP}}=-\int {\displaystyle \frac{3}{2}}k\rho \left(r\right)\left[c+\ln {\displaystyle \frac{t\left(r;\rho \right)}{t_{\mathrm{TF}}\left(r;\rho \right)}}\right]dr$$

(3)

where t(r; ρ) and t_TF(r; ρ) represent the non-interacting and Thomas–Fermi (TF) kinetic energy density, respectively. The constants are defined as follows: k is the Boltzmann constant, c = (5/3) +ln(4πc_K/3), and c_K = (3/10)(3π²)^2/3]. The non-interacting kinetic energy density $t\left(r;\rho \right)$ integrates to give the total kinetic energy $T_{\mathrm{S}}$,

$$\int t\left(r;\rho \right)dr=T_{\mathrm{S}}$$

(4)

It can be computed from the canonical orbital densities as,

$$t\left(r;\rho \right)={\displaystyle \sum }_i{\displaystyle \frac{1}{8}}{\displaystyle \frac{\nabla {\rho }_i·\nabla {\rho }_i}{{\rho }_i}}-{\displaystyle \frac{1}{8}}{\nabla }^2\rho$$

(5)

while the Thomas–Fermi expression is given by,

$$t_{\mathrm{TF}}\left(r;\rho \right)=c_{\mathrm{K}}{\rho }^{5/3}\left(r\right)$$

(6)

It is important to note that kinetic-energy density may take different forms depending on context. [48,49,50,51,52,53,54,55] Nonetheless, S_GBP satisfies the maximum-entropy principle from a rigorous mathematical viewpoint. [47]

Expanding further, several ITA descriptors have been proposed to characterize chemical reactivity. Within the framework of conceptual density functional theory (CDFT), [56,57,58,59,60] one such example is relative Rényi entropy [60] of order n

$$R_n^r={\displaystyle \frac{1}{n-1}}\ln \left[\int {\displaystyle \frac{{\rho }^n\left(r\right)}{{\rho }_0^{n-1}\left(r\right)}}dr\right]$$

(7)

Another related measure is information gain ($I_{\mathrm{G}}$) [61], also called Kullback−Leibler divergence or relative Shannon entropy,

$$I_{\mathrm{G}}=\int \rho \left(r\right)\ln {\displaystyle \frac{\rho \left(r\right)}{{\rho }_0\left(r\right)}}dr$$

(8)

In both expressions (7) and (8),${\rho }_0\left(r\right)$ is a reference-state density, and both ${\rho }_0\left(r\right)$ and $\rho \left(r\right)$ are normalized.

More recently, [62] one of the present authors introduced three ITA descriptors G₁, G₂, and G₃, applicable at both atomic and molecular levels:

$$G_1={\displaystyle \sum }_{\mathrm{A}}\int {\nabla }^2{\rho }_{\mathrm{A}}\left(r\right){\displaystyle \frac{{\rho }_{\mathrm{A}}\left(r\right)}{{\rho }_{\mathrm{A}}^0\left(r\right)}}dr$$

(9)

$$G_2={\displaystyle \sum }_{\mathrm{A}}\int {\rho }_{\mathrm{A}}\left(r\right)\left[{\displaystyle \frac{{\nabla }^2{\rho }_{\mathrm{A}}\left(r\right)}{{\rho }_{\mathrm{A}}\left(r\right)}}-{\displaystyle \frac{{\nabla }^2{\rho }_{\mathrm{A}}^0\left(r\right)}{{\rho }_{\mathrm{A}}^0\left(r\right)}}\right]dr$$

(10)

$$G_3={\displaystyle \sum }_{\mathrm{A}}\int {\rho }_{\mathrm{A}}\left(r\right){\left[\nabla \ln {\displaystyle \frac{{\rho }_{\mathrm{A}}\left(r\right)}{{\rho }_{\mathrm{A}}^0\left(r\right)}}\right]}^{2}dr$$

(11)

Finally, to partition electron density into atomic contributions within a molecule, the Hirshfeld stockholder approach [63,64] is frequently adopted. It is defined as:

$${\rho }_{\mathrm{A}}\left(r\right)={\omega }_{\mathrm{A}}\left(r\right)\rho \left(r\right)={\displaystyle \frac{{\rho }_{\mathrm{A}}^0\left(r\right)\left(r-R_{\mathrm{A}}\right)}{{{\displaystyle \sum }}_{\mathrm{B}}{\rho }_{\mathrm{B}}^0\left(r-R_{\mathrm{B}}\right)}}\rho \left(r\right)$$

(12)

Here,${\rho }_{\mathrm{A}}\left(r\right)$ is the atomic (Hirshfeld) density, ${\omega }_{\mathrm{A}}\left(r\right)$ is the weight or “sharing” function, ${\rho }_{\mathrm{B}}^0\left(r-R_{\mathrm{B}}\right)$ represents the reference (typically spherically averaged) atomic density centered at $R_{\mathrm{B}}$. The denominator is known as the promolecular density. The stockholder method naturally aligns with ITA due to its information-theoretic foundation. Alternative partitioning schemes include Becke’s fuzzy atom method [65] and Bader’s atoms-in-molecules (AIM) approach based on zero-flux surfaces [66]. A summary of our recent work in this direction is available in Refs [67,68,69].

4.2. An Outline of GEBF

In the generalized energy-based fragmentation (GEBF) method, [30,31,32,33] the total energy of a large system—such as a macromolecule or molecular aggregate—is expressed as a linear combination of the energies of smaller embedded subsystems, as given in Equation (13):

$$E_{\mathrm{tot}}={\displaystyle \sum }_m{C}_m{E}_m-\left[\left({\displaystyle \sum }_m{C}_m\right)-1\right]{\displaystyle \sum }_{\mathrm{A}}{\displaystyle \sum }_{\mathrm{B}>\mathrm{A}}{\displaystyle \frac{Q_{\mathrm{A}}{Q}_{\mathrm{B}}}{R_{\mathrm{AB}}}}$$

(13)

Here, E_m and C_m stand for the total energy and the coefficient of the mth subsystem, respectively. Q_A, is the atomic charge on atom A. R_AB is the interatomic distance between atoms A and B.

The general procedure for performing GEBF calculations involves several steps. Employing a molecular cluster of benzene (C₆H₆) as illustrated in Figure 4f, each benzene molecule is treated as a fragment. Primitive subsystems are then constructed centered at each fragment, defined by a distance threshold (ζ). These primitive subsystems are assigned coefficients C_m = +1. Due to the spatial overlap among primitive subsystems, smaller derivative subsystems are generated. The coefficients of these derivative subsystems are determined automatically using the principle of inclusion and exclusion, ensuring proper energy accounting. Another parameter, γ_max, representing the maximum number of fragments allowed in a subsystem, is introduced to control subsystem size.

All quantum chemical calculations for the subsystems are carried out using the GEBF method as implemented in the LSQC (low scaling quantum chemistry) package. [70] In this work, the two key GEBF parameters, (ζ, γ_max) are set to be (4.0, 6).

4.3. Computational Details

A total of 24 of octane isomers, metallic clusters Be_n (n =3 to 25), Mg_n (n = 3 to 20 and 28), (CO₂)_n (n = 4 to 40), organic clusters of (C₆H₆)_n (n = 4 to 30), covalent S_n (n =2 to 18); polymeric structures (see Figure 2) of polyyne, polyene, all-trans-polymethineimine, and acene, were taken from our previous publication. For the protonated clusters [(H₂O)_n(H₃O)]⁺ , they were taken from Ref 20. For cluster sizes n = 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, and 20, there are 74, 79, 113, 119, 108, 140, 121, 138, 114, 125, and 143 structures, respectively.

Molecular wavefunctions for all the systems were obtained at the HF/6-311++G(d,p) level. The Multiwfn 3.8 [71,72] program was utilized to calculate all ITA quantities by using the Gaussian 16 checkpoint or wavefunction file as the input. The stockholder Hirshfeld partition scheme of atoms in molecules was employed when atomic contributions were concerned. The reference-state density was the neutral atom calculated at the restricted open-shell ROHF/6-311++G(d,p) level. CCSD and CCSD(T) calculations for octane isomers were performed with the Gaussian 16 [73] package. For RI-MP2 calculations, Hartree‒Fock (HF) orbitals from the Gaussian 16 calculations were then transformed into the ORCA [74] format by using the MOKIT [75] program. The frozen core formalism [76] was used throughout this work, unless otherwise stated.

5. Conclusions

To summarize, in this work, we have applied the information-theoretic approach (ITA) quantities to appreciate the post-Hartree‒Fock (such as MP2 or RI-MP2) correlation energies for various molecular clusters and polymeric systems with both localized and delocalized electronic structures. We have found that for linear or quasi-linear polymeric systems, such as polyyne and polyene, the predicted results based on the Hartree‒Fock ITA quantities, are in excellent agreement with the calculated MP2 correlation energies. For other systems, such as hydrogen-bonded protonated water clusters and dispersion-bound carbon dioxide and benzene clusters, satisfactory results can be obtained with the LR(ITA) protocol. For metallic Be_n and Mg_n, as well as covalent S_n, one can still obtain reasonable results. In addition, for relatively larger benzene clusters, we compare the LR(ITA) results with those from the GEBF method and similar accuracy is observed. Our results collectively showcase that LR(ITA) is a promising method as a cost-efficient tool in predict the electron correlation energy.