Energy decompositions according to physical space partitioning schemes: Treatments of the density cumulant

This article is a continuation of our previous paper on schemes of energy decompositions of molecular systems in the real space [D. R. Alcoba et al., J. Chem. Phys. 122, 074102 (2005)] now using correlated state functions. We study, according to physical arguments, the appropriate management of the density cumulant arising from the second-order reduced density matrix at correlated level, whose contributions can be assigned to one-center or to two-center terms in the energy partitioning. Our treatments are applied within two physical space partitioning schemes: the Bader partitioning into atomic basins and the fuzzy atom procedure. The results obtained in selected molecules are analyzed and discussed in detail.


I. INTRODUCTION
The decomposition of the electronic energy of N-electron systems into contributions which can be associated with one, two, or more atom components has received significant attention from many authors recently. 1-12 This interest arises from the possibility of determining satisfactorily the strength of the chemical bondings between the atoms which compose a given molecule, atomic group, cluster, etc. The information provided by this type of studies is of paramount importance for chemists since it allows one to know the interactions within a system. As is well known, this kind of methodology has been tackled by different approaches. On the one hand, some methods perform the analyses of the one-center, two-center, etc., energy terms in the Hilbert space of basis functions, each of which is assigned to one of the atoms. 1-5 Alternatively, other methods use schemes based on the partitioning of the three-dimensional physical space into atomic domains. [5][6][7][8][9][10][11] This division of the physical space has usually been accomplished by means of Bader's atoms in molecules 13 ͑AIM͒ theory or by the fuzzy atom procedure. [14][15][16] In Ref. 10, we have described a general algorithm to decompose the energy of a molecule into one-and twocenter contributions, according to a partitioning of the real space. The procedure has been applied within the framework of the AIM theory and within the fuzzy atom approach. Moreover, our algorithm can be applied to any level of theory, that is, at Hartree-Fock ͓self-consistent field ͑SCF͔͒ level or at correlated one. However, we limited our calculations to the Hartree-Fock level, in order to establish an appropriate comparison with the results arising from procedures reported by other authors. The aim of this paper is to analyze the influence of many-body effects on the partitioning scheme itself and to explore all the capabilities of our method. For theoretical and practical reasons, the elements of the second-order reduced density matrix required for the energy expression are now formulated in terms of the elements of the first-order reduced density matrix and the elements of the density cumulant matrix. [17][18][19] Hence, the application of this formulation at correlated level rises the question of how to manage the matrix elements of the density cumulant within the partitioning of the energy into onecenter and two-center terms.
The paper is organized as follows. Section II describes the energy partitioning into one-center and two-center terms by means of two different treatments of the matrix elements of the density cumulant. In a first procedure, these elements are assigned to the two-center terms while in another treatment, according to physical reasons, 20 they are managed as possessing one-center character. Both treatments are developed within the AIM theory and within the fuzzy atom approach. Section III reports the computational aspects and the numerical determinations performed in selected molecules, which provide the comparison with previous results at SCF level and the corresponding discussion. Finally, in Sec. IV we report the remarks and conclusions of this work. a͒

II. TREATMENTS OF THE DENSITY CUMULANT
We will express the energy of an N-electron molecule with clamped nuclei as

͑1͒
in which i , j , k , l ,... are orbital functions of an orthonormal basis set, A , B ,... are the nuclei of the molecule, R AB the distance between those nuclei, and Z A , Z B ,... the corresponding nuclear charges. T j i mean the matrix elements of the kinetic energy operator, ͑ − 1 2 ٌ 2 ͒ . A V j i are the matrix elements corresponding to the ͑−Z A / r A ͒ operator and B jl ik = ͗ik ͉ jl͘ the standard two-electron integrals in the ͑12͉12͒ convention. 1 D j i and 2 D jl ik denote the matrix elements of the spin-free first-order and second-order reduced density matrices, 21 respectively.
For practical purposes this equation will be rewritten as and the Kronecker deltas will be formulated in a general way as in which ͗i ͉ j͘ are the standard overlap integrals ͑where the integration is performed over the total space͒. These integrals are partitioned according to the S j i ͑A͒ quantities related with the atomic domains of the centers A.
Within the AIM scheme, we will formulate the quantities S j i ͑A͒ as where ͗i ͉ j͘ A are the overlap integrals over the Bader atomic domain associated with the nucleus A. 13 Alternatively, in the fuzzy atom approach 14-16 the quantities S j i ͑A͒ will be formulated as in which a non-negative continuous weight function w A ͑r͒ is introduced for each nucleus A. These weight functions measure the degree in which a given point of space r is considered to belong to atom A, fulfilling the conditions The introduction of the Kronecker deltas in Eq. ͑2͒ according to Eq. ͑3͒ allows us to express the energy as a sum of one-center and two-center terms in which Equations ͑8͒-͑10͒ constitute a compact formulation for the partitioning of the electronic energy of an N-electron molecule in the three-dimensional physical space, which is expressed by means of one-and two-center terms ͑E A and E AB , respectively͒. From a mathematical point of view, both approaches, Bader and fuzzy atom ones, only differ in the formulation of the quantities S j i ͑A͒. Further manipulations of Eq. ͑2͒ can lead to three-and four-center energy components which have no clear chemical significance. 2,3 Hence, we propose the partitioning scheme into one-and two-center contributions of Eq. ͑8͒.
The elements of the spin-free second-order reduced density matrix are formulated as 4,18,19,22,23 where ⌳ jl ik are the elements of the spin-free density cumulant matrix. All these matrix elements are zero at the SCF level ͑for singlet states͒. Consequently, the use of correlated functions implies the appearance of new contributions in the formulation of the partitioning of the molecular energy whose treatment requires to be studied.
The substitution of Eq. ͑11͒ into Eqs. ͑9͒ and ͑10͒ leads one to The partitioning of the energy expressed by Eqs. ͑12͒ and ͑13͒ will be denominated two-center model because the elements of the density cumulant ⌳ are distributed between the one-center and two-center energy components. Alternatively, as a consequence of the physical meaning of the cumulant describing effectively unpaired electrons, which do not contribute to bondings, 20 its matrix elements ⌳ jl ik should be assigned to the one-center components. From this point of view, the elements ⌳ jl ik will only contribute to one-center energies; this treatment of the density cumulant will be called one-center model. In this case, the energy partitioning leads one to the formulas As has been mentioned in the Introduction, the aim of this paper is to study both treatments of the density cumulant. In Sec. III, we report results arising from both models.

III. COMPUTATIONAL ASPECTS AND RESULTS
The numerical determinations performed in this work have required the calculation of the one-electron integrals T j i and A V j i as well as that of the two-electron integrals B jl ik . All of them were computed with the basis sets 6-31G using a modified version of GAMESS program. 24 The overlap integrals ͗i ͉ j͘ A ͑calculated over the Bader domains͒ were evaluated with the PROAIM code, 25 while the integrals ͗i͉w A ͑r͉͒j͘ were obtained with the modified code cited in Ref. 26 which follows a Becke integration scheme 27 based on the weight functions w A ͑r͒ originally proposed by this author. The atomic radii used for determining the w A ͑r͒ weights were identical to those employed in our previous work. 10 The PSI3 program package 28 provided the first-order and second-order reduced density matrices derived from N-electron state functions calculated at the level of configuration interaction ͑CI͒ with single and double excitations ͑CISD͒. All the calculations were implemented with the experimental geometries in order to compare more directly the resulting energies with the experimental ones arising from dissociation measures. 29 To test the behavior of the above described models we have performed energy partitionings in simple systems ͑H 2 ͒ and hydrocarbons with nonpolar bonds ͑CH 4 , C 2 H 6 , C 2 H 4 , and C 2 H 2 ͒ as well as in hydrides with growing polarity bonds ͑NH 3 , H 2 O, and HF͒. Tables I and II report the results within the schemes AIM and fuzzy atoms, respectively. Both tables collect results arising from the SCF approximation and from the two-center and one-center models of treatment of the density cumulant ͑calculated within the CISD approxi-mation͒. In the SCF case, 10 all the elements of the density cumulant matrix are zero and consequently, no distinction between both models can be made. However, we have included the SCF results ͑obtained with the experimental ge-ometries͒ for the sake of comparison with the results at correlated level. As is well known, the two-center energy components ͑E AB ͒ characterize the strengths of the bonds between the individual atoms within the molecule but they cannot be related immediately to the dissociation energies.  Table I ͑partitioning of the physical space according to the AIM picture͒ shows that the bonding energy numerical values obtained with the one-center model for the density cumulant are closer to the experimental ones than those obtained with the two-center model. This behavior is observed in all the values related with classical chemical bonds. Moreover, these one-center model results are also closer to the SCF ones than their two-center counterparts and the differences ͉E AB SCF − E AB one-center ͉ lies in the interval ͑10 −4 , 3 ϫ 10 −2 ͒ a.u. The lower limit corresponds to the CC bond in C 2 H 6 molecule, while the higher one corresponds to the bonds OH and HF in the systems H 2 O and HF, respectively. These results seem to indicate that the differences between both bonding energy values grow up from pure covalent bonds to polar ones passing through multiple bonds in the C 2 H 4 and C 2 H 2 molecules. Consequently, such differences, although small, tend to increase in systems in which the effectively unpaired electron density concentrates in atomic sites. 30,31 The results reported in Table II ͑partitioning of the physical space according to the fuzzy atom picture͒ show a similar trend, although the values E CC ͑in the C 2 H 6 , C 2 H 4 , and C 2 H 2 systems͒ turn out to be higher ͑in all the cases͒ than in the AIM picture. As in the previous table, the results derived from the one-center model are also closer to the SCF ones than those arising from the two-center model although the behavior of the differences ͉E AB SCF − E AB one-center ͉ is less clear. This may be interpreted as a consequence of the bias of the surfaces defining the atoms due to the semiempirical character of the fuzzy atom methodology choosing the same atomic radii for the elements in different compounds.
The physical meaning of the density cumulant and the results described in this paper confirm our previous suggestions on the role played by the density cumulant. In fact, in Ref. 20 we pointed out the nonpairing character of the density cumulant in the description of bond orders. Similar arguments can be exhibited from the fact that a tensorial con-taction of the ⌳ matrix 19 leads to the effectively unpaired electron matrix 32,33 which as is well known describes unpaired electrons assigned to determined domains.

IV. CONCLUDING REMARKS
In this work, we have proposed and examined two possible treatments for the density cumulant in the studies of energy decompositions according to partitionings of the physical space at correlated level. In the two-center model the elements of the density cumulant matrix are assigned to one-center and to two-center energy components, while in the one-center model those elements are assigned only to one-center components. The numerical determinations, performed within the AIM and fuzzy atom schemes of partitioning of the physical space, point out our one-center model as the most appropriate treatment, mainly in the AIM scheme. These results agree with the physical meaning of the density cumulant and confirm our previous predictions on the management of this device. We are currently studying the treatment of the correlation effects within energy partitionings carried out in the Hilbert space of atomic orbitals and the results will be discussed in a forthcoming paper. 34