A chemical Hamiltonian approach study of the basis set superposition error changes on electron densities and one-and two-center energy components

The basis set superposition error-corrected first-order electron densities of several hydrogen bonded complexes of increasing molecular size have been obtained with the Hartree–Fock and density-functional theory versions of the chemical Hamiltonian approach ~CHA! methodology. A detailed analysis of the local basis set superposition error ~BSSE! effects has been carried out by comparing the uncorrected electron densities and energy components with the CHA ones. Topological analysis of the electron density through the atoms in molecules theory is used in order to obtain a quantitative measure of the BSSE effects in terms of the characterization of the critical points of the electron density. Density difference isocontour maps are also depicted in order to show the local electron density redistributions induced by the BSSE-correction. We show that the effects of the BSSE are common for all the complexes studied, namely water dimer, formic acid dimer and uracil–water complex. The formic acid dimer and uracil–water density difference maps at frozen geometry reveal that the effects of the BSSE do not extend significantly beyond the atoms involved in the interaction and their first neighbors. The main redistribution effects are not strictly localized on the intermolecular region and mostly take place in the valence shells of the heavy atoms directly involved in the intermolecular interaction. These trends are also confirmed by means of an energy decomposition analysis performed at the Hartree–Fock level of theory with the recently proposed chemical energy component analysis ~CECA! method. In agreement to previous results, we found that inclusion of diffuse functions is of utmost importance in order to minimize the magnitude of the BSSE. However, both the electron density difference maps and the CECA analysis confirm that the local effects of the BSSE are very different when diffuse functions are present in the calculation. © 2002 American Institute of Physics. @DOI: 10.1063/1.1463439#


INTRODUCTION
Accurate ab initio or density-functional theory ͑DFT͒ calculations of weakly bound molecular complexes must take into account the so-called basis set superposition error ͑BSSE͒. 1 BSSE is caused by the imbalance in the description of the monomers forming the complex and the same monomers isolated.Briefly, in the whole complex calculation, the intramolecular operators associated to each molecule or fragment can be expanded to some extent in the basis functions of the other molecule.In contrast, for the description of the isolated monomers, each molecule can use only its own basis functions.
It is well known that BSSE leads to an overestimation of the interaction energy between the monomers and a shortening of the intramolecular distance. 2BSSE is caused by the use of finite basis sets to expand the molecular orbitals in ab initio and Kohn-Sham DFT calculations.Therefore, in principle, the use of large basis sets should diminish the magnitude of the error.Eventually, in the complete basis set ͑CBS͒ limit, the BSSE should be zero.Nevertheless, in practice, for weakly bound molecular complexes, the magnitude of the BSSE can be of the same order of the interaction energy.
There are several methods available for the correction of the BSSE.The counterpoise method ͑CP͒ 3 has been widely used in the last decades in order to obtain BSSE-free energies at any level of theory.The CP method has also been generalized and applied successfully to the correction of molecular geometries, vibrational frequencies and, in general, any molecular property that depends on derivatives of the energy. 4The CP method shows the desired asymptotic behavior in the CBS limit.Furthermore, the use of CPcorrected values allow for a meaningful extrapolation of supermolecular properties 5 to the CBS limit.However, CP does not provide a corrected wave function.There have been attempts to obtain CP-corrected electron densities, but only the so-called interaction density, that is, the difference between the complex electron density and that of the monomers at a given geometry, can be determined.Since there is no CPcorrected Hamiltonian, no CP-corrected molecular orbitals are available for the molecular complex and hence no molecular electron density.
An alternative to the CP method is the chemical Hamiltonian approach ͑CHA͒. 6The CHA method eliminates the unphysical terms in the Hamiltonian that cause BSSE delocalizations between the fragments within the molecular complex.The CHA methodology has been developed and successfully applied at several levels of theory, namely Hartree-Fock ͓CHA-SCF [7][8][9] and CHA/F 10 ͔, DFT ͑CHA/DFT 11 ͒, and second-order Møller-Plesset perturbation theory, ͑CHA-MP2͒. 12Several studies have found that CP and CHA lead to very similar corrections to both the interaction energies and molecular geometries 13,14 for several van der Waals or hydrogen bonded complexes.The main advantage of CHA with respect to the CP method is that CHA corrects the BSSE a priori, yielding a BSSE-free wave function and, therefore, also a BSSE-free electron density, (r).
There is still another a priori BSSE-correction method, namely the self-consistent field for molecular interactions ͑SCF-MI͒ [15][16][17] that has recently been used to obtain BSSEfree electron densities of hydrogen bonded complexes. 15owever, it has been pointed out that the SCF-MI eliminates some true interaction terms and thus overestimates the BSSE correction. 18The consequence is that, as opposed to CP and CHA, SCF-MI does not exhibit the correct asymptotic behavior near to the CBS limit. 15ecently, we have carried out an analysis of the effects of the BSSE on the electron density of the hydrogen fluoride dimer, by comparing electron densities obtained with conventional Hartree-Fock ͑HF͒ and DFT methods, using BLYP and B3LYP ͑Becke three-parameter Lee-Yang-Parr͒ functionals, to the corresponding BSSE-free densities obtained with the CHA/F and CHA/DFT methods, respectively, using several basis sets of different size. 19We have formally considered the BSSE correction as a perturbation to the system.In this way, we have been able to analyze the effects of BSSE correction both in the nuclear geometries and electron densities, according to Scheme 1.
X//Y in the scheme above refers to a single-point calculation using method X at the molecular geometry optimized with method Y. SCF and CHA are used to denote conventional and BSSE-corrected calculations, respectively, at any level of theory.Thus, SCF refers to conventional HF or DFT calculations, while CHA denotes the corresponding BSSEcorrected CHA/F and CHA/DFT calculations, respectively.According to this notation, SCF//SCF in Scheme 1 corresponds to the conventional or unperturbed calculation, that is, a HF or DFT calculation with no BSSE correction.With respect to this reference calculation, the correction of the BSSE with no nuclear relaxation induces a rearrangement of the electron density.This is denoted as CHA//SCF.Finally, geometry optimization with the CHA method leads to the fully corrected situation, CHA//CHA.The cycle can be closed with an SCF//CHA calculation, meaning an uncorrected calculation at the molecular geometry obtained with the CHA method.One should take into account that the BSSE correction is not a perturbation itself.However, this scheme allows to analyze the effects of removing the BSSE on the electron density with and without nuclear relaxation.In order to quantify the effect of the BSSE-correction on the electron density several tools can be used.First, according to the theory of atoms in molecules ͑AIM͒, 20 the interaction between two molecules can be characterized by analyzing the properties of the critical points in the electron density associated to the intermolecular interaction.This will provide us with quantitative comparisons of the uncorrected and the BSSE-corrected electron densities, as well as both the nuclear and electronic relaxation contributions separately.Additionally, when no nuclear relaxation is allowed, one can depict density difference maps between uncorrected and corrected densities to analyze visually the effects of BSSE correction.Finally, the quantum molecular similarity 21 theory can provide us with quantitative measures of the global similarity of two electron density distributions.The determination of the similarity indexes between the uncorrected and the BSSE-corrected densities of the molecular complexes studied here as well as the hydrogen fluoride dimer have been recently carried out in our laboratory. 22he geometry of the hydrogen fluoride dimer has been found to be very dependent on the level of theory and basis set.However, some systematic trends have been found for the effects of BSSE correction on the charge density at a frozen geometry.SCF//SCF-CHA//SCF density differences reveal that the BSSE correction generally depletes the electron density from the intermolecular region.However, the stronger effects take place in the valence shells of the two heavy atoms.For each F atom, the BSSE correction overestimates and underestimates the electron density in two regions which are approximately perpendicular to each other.For the donor F, the BSSE correction removes the electron density along the intramolecular F-H axis while, for the acceptor F, the intermolecular F¯H axis corresponds to a zone of density overestimation.These trends are common to most of the combinations of method and basis set, except when diffuse functions are used.Inclusion of diffuse functions leads to an overall decrease of the BSSE by one order of magnitude, and hence to smaller corrections of the electron density, as revealed also by means of a quantum molecular similarity study. 22Moreover, when diffuse functions are used, difference maps reveals that the BSSE correction actually induces to an increase of the charge density in the intermolecular region, and the patterns of electron density redistribution around the heavy atoms are lost.
The present paper extends this kind of analysis to other molecular complexes, namely, the water dimer, the formic acid dimer and the water-uracil complex.Its first objective is to assess whether the trends found for the BSSE correction with the CHA method for the hydrogen fluoride dimer are common to other systems.A second goal is the study of larger systems, specially the water-uracil complex, that will allow to analyze whether the effects of BSSE are restricted Nuclear and electronic relaxation paths.
strictly to the atoms directly involved in the intermolecular interaction or are extended to other regions of the molecule.As a third purpose, the density analysis is complemented by using the recently proposed chemical energy component analysis ͑CECA͒. 23The CECA uses a projected integral expansion technique to approximately decompose the molecular energy into one-and two-center terms.This decomposition allows to analyze the local BSSE effects in terms of atomic ͑one-center͒ and interatomic ͑two-center͒ energy components.

The CHAÕF and CHAÕDFT methods
A detailed explanation of the CHA methodology is available in a recent review. 6Briefly, the CHA is based on the removal from the Hamiltonian of the intermolecular contributions that cause the BSSE, while keeping all the physically true interaction terms.Thus, the CHA ensures that the description of each monomer or fragment isolated and within the complex is consistent.The separation of the intramolecular, pure intermolecular and BSSE contributions in the Hamiltonian is carried out by means of a mixed secondquantization formalism.The resulting Hamiltonian is not Hermitian, so the application of the variational principle to derive the SCF equations is not trivial.However, by making use of the Brillouin theorem it is still possible to derive a set of SCF equations and obtain the corresponding CHA-SCF canonical orbitals. 7The CHA/F 10 method is a variation of CHA-SCF.In this case ͑CHA/F͒, the modifications are performed to the standard Fock equations, so the one-and twoelectron integrals remain unchanged, with respect to the original SCF method.The CHA/F equations can be considered as an approximation to the CHA-SCF ones; in practice, it has been shown that CHA-SCF and CHA/F results are practically equivalent.The main advantage of using the CHA/F scheme is that the same philosophy can be used in the context of Kohn-Sham DFT.In this case, the equivalent to the CHA/F method is CHA/DFT.Only CHA/F and CHA/ DFT were used throughout the present paper.
As the CHA Hamiltonian is not Hermitian, its eigenfunctions are not orthonormalized.This does not represent a problem in order to obtain the electron density since Slater determinants built from orthogonal or nonorthogonal orbital sets differ only by a physically irrelevant phase factor.However, it is more convenient to orthogonalize the occupied orbitals in order to construct the first-order density matrix like in a conventional calculation.Therefore, at each iteration, the canonic CHA orbitals are orthonormalized in order to build up the density matrix for the next iteration.After self-consistency, the CHA orbitals obtained are also orthonormalized and the final density matrix is obtained.
Finally, it must be taken into account that the evaluation of analytic gradients of the CHA energy is not straightforward.Actually, the corresponding equations have been developed but not yet implemented. 24Therefore, the location of stationary points on the BSSE-corrected potential energy surface ͑PES͒ was performed with numerical gradient methods. 13e chemical energy component analysis "CECA… With semiempirical quantum methods, the molecular energy can be expressed in terms of one-and two-center contributions only.This property allows for a straightforward decomposition of the molecular energy into atomic ͑one-center͒ and interatomic ͑two-center͒ contributions.In contrast, calculation of molecular energies with ab initio methods involves also three-and four-center terms.These terms may contribute significantly to the total energy and cannot be generally ignored; however, it is difficult to attach a direct chemical significance to those multicenter terms.
The basic idea in the CECA is to use a projective integral expansion scheme that allows to express approximately each multicenter term as a summation of one-and two-center terms.The theory behind the CECA is tightly related to the CHA ͑see Ref. 23 for a detailed description͒.The main advantage of this approach compared to the classical bond order analysis is that bonding interactions with formally the same multiplicity ͑single, double, triple bonds͒ can be now clearly differentiated in terms of energetic ͑static͒ contributions to the overall energy of the system.Also, it allows to distinguish between bonding and antibonding interactions.However, one must take into account that this decomposition is exact only for diatomic molecules, where no multicenter contributions are possible.In the general case, the sum of all the energy contributions does not match exactly the total molecular energy.Mayer states that the CECA reproduces HF energies of small molecules with a precision between 10 and 40 mHartrees ͑Ϸ6 -25 kcal•mol Ϫ1 ͒. 23 This is a relatively small error, compared to the total molecular energy.Furthermore, one should take into account that this error arises from the summation of all the one-and two-center energy contributions.One can expect that the individual components have much smaller errors.Therefore, the accuracy of the CECA decomposition scheme should be sufficient in most cases.
A second problem of the CECA is related to the fact that the energy decomposition is performed in the Hilbert space.That is, each atomic orbital or basis function is assigned to a single atom and the partition of each energetic term is carried out by projecting over the subspace spawned by these atomic orbitals.Practically, the basis functions are supposed to belong to the atom in which they are centered.Then, one-and two-center contributions are related to atomic and diatomic components.However, the results of this kind of decomposition schemes can be very dependent on the basis set used for the calculation.Particularly, these analysis may lose significance when diffuse functions are included in the basis set.For instance, it is well known that Mulliken charges, also based on the formal partition of the atomic orbital's space, have little chemical meaning when diffuse basis functions are used.This problem can be solved by performing a similar decomposition in the Euclidean space, for instance in the frame of the AIM theory. 25It can be shown that a partitioning of that kind, that can be connected to the CECA one by a simple mapping of the integrals, 26 can provide an exact decomposition of the HF molecular energy.Several test calculations revealed some differences between the CECA and the AIM-based energy partitioning. 25Unfortunately, the huge computational cost that implies the evaluation of double in-tegrations over two AIM atomic basins prevents the use of this methodology for the present work, and only CECA approximated results will be shown.

Computational details
Molecular wave functions, energies and first-order electron densities were calculated using the GAUSSIAN 94 package, 27 using conventional HF and DFT methods.The CHA/F and CHA/DFT calculations were done using a modified version of the GAUSSIAN 92 package. 28The B3LYP functional was used for the DFT and CHA/DFT calculations, henceforth CHA/B3LYP.Six different standard basis sets were used for the water dimer: 6-31G, 6-31G(d), 6-31G(d,p), 6-31ϩϩG(d, p), 6-311G(d, p) and 6-311ϩ ϩG(3d f ,2pd).For this system, molecular geometries were completely optimized for each combination of computational method and basis set in a previous study. 13For all the basis sets, additional single-point CHA/F and CHA/B3LYP calculations were carried out at the optimized HF and B3LYP geometries, respectively.For the formic acid dimer and uracil-water complexes, only the HF method and two different basis sets for each complex were used.Moreover, stationary points for these complexes were located only on the conventional PES.Therefore, SCF//SCF and CHA//SCF results are available for the three complexes studied and, additionally, CHA//CHA results also for ( For all the calculations, the critical points on (r) were located and characterized using the AIMPAC package. 29Special attention was paid to intermolecular bond critical points ͑bcp͒ and ring critical points ͑rcp͒.Difference maps correspond to differences between SCF//SCF and CHA//SCF electron densities.For the water dimer, the maps were performed on the symmetry plane of the complex.For the formic acid dimer and uracil-water complex, they were calculated in the molecular plane containing all the heavy atoms.Positive values in the maps ͑solid lines͒ correspond to zones where the uncorrected calculation overestimates the electron density, whereas negative values ͑dashed lines͒ correspond to zones where the electron density is underestimated, with respect to the corresponding CHA electron density.The energy decomposition analysis was carried out for the three complexes at the HF level of theory using the program APOST. 30 The current implementation of the program does not allow for f and higher angular momentum basis functions so the results for the 6-311ϩϩG(3d f ,2pd) basis set could not be obtained.Instead, the values for the 6-311ϩϩG(d,p) are reported.

RESULTS
This section presents the electron density analysis for the three systems: the water dimer, the formic acid dimer and the uracil-water complex.For each case, the effects of the BSSE correction were assessed: ͑i͒ By comparing the properties of the intermolecular bcp's and rcp's with and without BSSE correction, ͑ii͒ by analyzing ͑SCF//SCF-CHA//SCF͒ density difference maps.The energy decomposition analysis for the three complexes is presented last.

Water dimer
All the calculations refer to the trans-linear water complex, with C s symmetry and a single H bond between the two water molecules ͑see Scheme 2͒.This structure is preserved in all the calculations; however, the level of theory, basis set and BSSE correction have a significant impact on the molecular geometry.A detailed discussion about the effects of the method of calculation, basis set size and BSSE on the geometry of this complex is available in Ref. 13.In general, the rO-O distance is larger at the SCF level of theory, compared to B3LYP.At both levels, increasing the basis set size or correcting the BSSE leads also to a lengthening of the rO-O distance.As for the ␣ and ␤ angles that define the mutual orientation of the two water monomers ͑see Scheme 2͒, it is generally necessary to employ relatively large basis sets or BSSE correction to obtain reliable values, specially for the B3LYP calculations.In general, as the size of the basis set is increased, the corrected and uncorrected geometrical parameters converge to values near to the experimental ones. 31 all levels of theory, the two water molecules are connected through a O-H¯O bond.The properties of the corresponding bcp reflect the characteristics of the water-water interaction ͑see Tables I and II, for the HF and B3LYP results, respectively͒.According to the most accurate calculation, B3LYP/6-311ϩϩG(3d f ,2pd), the values of the charge density, bcp (r), and its Laplacian, ٌ 2 bcp (r), at the bcp are 0.024 and 0.076, characteristic of a strong O-H¯O hydrogen bond.At the HF/6-311ϩϩG(3d f ,2pd) level, bcp (r) and ٌ 2 bcp (r) are significantly smaller: 0.016 and 0.066, respectively.In general, the correction of the BSSE without allowing nuclear relaxation ͑CHA//SCF͒ does not change significantly the bcp (r) and ٌ 2 bcp (r) values.Indeed, corrections in bcp (r) and ٌ 2 bcp (r), defined as SCF// SCF-CHA//SCF, are usually smaller than 10 Ϫ3 and 10 Ϫ2 , respectively.The location of the bcp does not change appreciably upon BSSE correction either.Moreover, the sign of the correction depends on the basis set.For instance, with the 6-31G basis set, the BSSE correction to bcp (r) is positive, leading to higher bcp (r) values, both at the HF and B3LYP levels.Addition of polarization functions ͓6-31G(d) and 6-31G(d,p) basis sets͔ makes the correction negative.In contrast, with respect to the 6-31G(d,p) basis set, addition of more valence functions ͓6-311G(d, p)͔ or diffuse functions ͓6-31ϩϩG(d, p)͔ leads to positive corrections again.Finally, for the largest basis set, 6-311Gϩϩ(3d f ,2pd), the correction is positive.Note that, in general, the use of larger basis sets does not lead to smaller corrections in bcp (r) and ٌ 2 bcp (r).When nuclear relaxation is allowed ͑CHA//CHA calcu-lations͒, bcp (r) and ٌ 2 bcp (r) values always decrease, at the same time that the ellipticity and the distance of the bcp to the O of the acceptor molecule also increase slightly.These trends are common for the HF and B3LYP calculations, and for all the basis sets; however, the differences between the corrected and uncorrected calculations become progressively smaller as the basis set size is increased.In general, these trends, combined with the behavior of the rO-O distance discussed above, reveal that both the increase in basis set size and the BSSE correction work in the direction of weakening the H bond interaction.
Figures 1 and 2 collect SCF//SCF-CHA//SCF density difference maps at the HF and B3LYP levels of theory, respectively, for all the basis sets used in the study.The posi-tions of atomic nuclei and intermolecular bcp's ͑only those corresponding to the SCF//SCF electron density are shown͒ are denoted with circles and stars, respectively.All the difference maps corresponding to HF or B3LYP calculations using basis sets without diffuse functions ͓Figs.1͑a͒-1͑c͒, 1͑e͒, 2͑a͒-2͑c͒, and 2͑e͔͒ are similar.In general, the intermolecular bcp is located close to the zero isodensity contour.For the calculations with diffuse functions and also with the 6-31G and 6-311G(d,p) basis sets, ͓Figs.1͑a͒, 1͑d͒-1͑f͒, 2͑a͒, and 2͑d͒-2͑f͔͒, the intermolecular region is quite flat and exhibits negative values.Thus, depending on the combination of level of theory and basis set, the intermolecular bcp may be found in positive or negative zones of the SCF// SCF-CHA//SCF density difference maps.Actually, according to the density difference maps, the main effects of the BSSE correction in the electron density take place in the  tor molecule, the effect of BSSE correction is just the opposite.Indeed, the subtle density differences found in the intermolecular region may actually be just a consequence of the redistributions taking place around the heavy atoms.Finally, density difference maps corresponding to calculations with diffuse functions, Figs.1͑d͒, 1͑f͒, 2͑d͒, and 2͑f͒ exhibit also maximal density differences around the heavy atoms, but not the polarization patterns characteristic of the maps in Figs.

Formic acid dimer and uracil-water complex
Results for the hydrogen fluoride 19 and water dimers suggest that inclusion of diffuse functions in the basis set is the main factor influencing the magnitude of the BSSE, while the level of theory and inclusion of more valence or polarization functions has a minor impact.Therefore, only the HF method and two different basis sets were used for each system, namely, 6-31G(d, p) and 6-31ϩϩG(d,p) for the formic acid dimer, and 6-31G(d) and 6-31ϩG(d) for the uracil-water complex.Tables III and IV gather the properties of the intermolecular critical points for the formic acid dimer and water-uracil complex, respectively, while Figs. 3  and 4 depict the SCF//SCF-CHA//SCF density difference maps for the two complexes.
The optimized structures belong to the C 2h and C s symmetries, for the formic acid dimer and uracil-water complex, respectively ͑see Scheme 3͒.In the case of the uracil-water system, the structure reported corresponds in fact to one of several minima which are very close in energy.In all cases, there are two H-bonds, which are equivalent in the formic acid dimer.Accordingly, two intermolecular bcp's exist, together with a ring critical point ͑rcp͒.Since no nuclear relax-ation is allowed, the differences between corrected and uncorrected bcp (r) and ٌ 2 bcp (r) values are again quite small ͑see Tables III and IV͒.For the formic acid dimer, BSSE correction increases bcp (r) and ٌ 2 bcp (r), for both basis sets.The same trend is found for the two bcp's in uracilwater with the 6-31ϩG(d) basis set.However, for the 6-31G(d) basis set, BSSE correction decreases the bcp (r) values and increases the ٌ 2 bcp (r) ones.As for the rcp's, the BSSE correction decreases rcp (r) for the two systems, when using the smaller basis sets.However, when diffuse functions are considered, rcp (r) shows no variation ͑for the formic acid dimer͒ or increases slightly ͑for uracil-water͒.ٌ 2 rcp (r) increases for all the calculations, except uracilwater with the 6-31ϩG(d) basis set.In general, the changes induced by BSSE correction in the rcp's properties are even smaller than the changes in the bcp's.
Formic acid dimer ͑a͒ and Uracil-water ͑b͒ complex.Scheme 3. Density difference maps give further insight on the local effects of BSSE on the electron densities.The formic acid dimer exhibits the main trends found for the water dimer ͑and previously for the hydrogen fluoride dimer 19 ͒.Thus, for the 6-31G(d,p) basis set, Fig. 3͑a͒, there is a narrow inter-molecular region where the BSSE correction decreases the electron density.This region includes the rcp but not the two bcp's.The main density redistribution effects take place in the valence shells of the heavy atoms.In particular, the O atoms in the hydroxyl and carbonyl moieties exhibit the den-FIG.3. Formic acid dimer SCF//SCF-CHA//SCF density difference isocontour maps.͑a͒ 6-31G(d,p), ͑b͒ 6-31 ϩϩG(d,p).Isodensity contours at Ϯ2.10 Ϫ4 , Ϯ4.10 Ϫ4 , Ϯ8.10 Ϫ4 , etc. sity redistribution patterns directed along the bonding axes characteristic of H-donor and acceptor atoms, respectively.The C-H moiety, which is the only one that does not participate directly in the water-water interaction, exhibits also some density redistribution due to the polarization of the neighboring atoms.Indeed, it appears that the BSSE correc-tion underestimates the electron density on the C atom.
When diffuse functions are added ͓see Fig. 3͑b͔͒, SCF// SCF-CHA//SCF density differences in the intermolecular region become slightly negative.In this case, all the intermolecular cp's fall into this negative zone.The strong redistribution patterns associated to the H-donor and acceptor atoms FIG. 4. Uracil-water complex SCF// SCF-CHA//SCF density difference isocontour maps.͑a͒ 6-311G(d), ͑b͒ 6-31ϩG(d).Isodensity contours at Ϯ2.10 Ϫ4 , Ϯ4.10 Ϫ4 , Ϯ8.10 Ϫ4 , etc. in Fig. 3͑a͒ are not found in this case.It appears that the negative region in the intermolecular zone is followed by alternating positive and negative regions at each side, with some positive regions focused strictly on the heavy nuclei.Nevertheless, the density difference decreases dramatically when including the diffuse functions.The maximum density differences observed with the 6-31G(d, p) basis set were Ϫ0.0190 and 0.0076 a.u., whereas for the 6-31ϩϩG(d,p) these values decrease to Ϫ0.0012 and 0.0015 a.u., respectively.
The uracil-water complex has some features that may add interesting insights.First of all, it is a relatively large system, which allows to study the scope of the BSSE effects on molecular electron densities.Second, the O in the water moiety acts as H-donor and acceptor at the same time.Figure 4͑a͒ corresponds to the SCF//SCF-CHA//SCF map with the 6-31G(d) basis set.In this case, a positive and a negative region are found in the intermolecular region.All the intermolecular cp's fall into the negative one.In the uracil molecule, the O 7 and the N 2 atoms exhibit the directional density redistribution patterns characteristic of H-acceptor and donor systems, respectively.Thus, the BSSE correction underestimates the electron density along the N 2 -H 8 bond and in the middle of the O 7 -H 13 intermolecular H bond.The O atom in the water molecule combines both features: The BSSE correction underestimates the density along the O 14 -H 13 intermolecular bond, but there is an increase in the direction of the intramolecular O 14 -H 8 bond.The C 1 atom, which is bonded to an acceptor and to a donor atom, exhibits minor density redistributions, similar to the ones in the C atom in (HCOOH) 2 .Finally, the effect of the BSSE in the rest of atoms is practically negligible, except for the carbonyl O 9 atom.
The 6-31ϩG(d,p) difference map ͓Figure 4͑b͔͒ presents a relatively large intermolecular region with negative values, which encloses all the intermolecular cp's.Significant density redistribution takes place only around the atoms directly involved in the H-bond interactions.As usual, the highly directional density redistribution patterns around heavy atoms found in the 6-31G(d) difference map are lacking in the 6-31ϩG(d) one.Atoms not involved in the intermolecular interaction do not exhibit appreciable density redistributions, except for O 9 .Again, the maximum density differences are about one order of magnitude larger when no diffuse functions are included.

Energy decomposition analysis
Tables V-VIII gather the results of the CECA decomposition for all the HF calculations.Table V collects the oneand two-center energy components obtained for the water dimer with the 6-31G(d,p) and 6-31ϩϩG(d,p) basis sets.Atomic energies are always negative ͑stabilizing͒, as well as the interaction between bonded atom pairs.Some terms like the O-O and H-H interactions are repulsive, which agrees with the chemical intuition.Direct comparison of the energy components obtained with different basis sets is not very convenient since the total molecular energy can be very different.Hence, since we are interested in the analysis of the effects of the BSSE in the energy, only selected differences between the SCF//SCF and CHA//SCF values for each basis set are discussed.
We use E BSSE (A) and E BSSE (A,B) to denote the effect of the BSSE correction in one-and two-center energy components involving atom A and the atomic pair A, B, respectively.Negative E BSSE (A) and E BSSE (A,B) values correspond to energy components that are too stabilizing in the uncorrected calculation, because of the BSSE.That is, the given one or two-center component is more stabilizing ͑less destabilizing͒ for the SCF//SCF than for the CHA//SCF calculation.Inversely, positive values correspond to energy components that are lower in energy for the CHA than for the SCF calculations.Tables VI-VIII also list the total energy difference for each monomer as well as the correction to the static interaction energy, computed by summing up all the corresponding CECA one-and two-center terms.
Note that the sum of all the CECA intermolecular energy components must be clearly distinguished from the conventional stabilization and interaction energies, the former reported in Tables I and II.In the supermolecular approach, the interaction energy is defined as the difference between the energy of the complex and the energies of the monomers at the complex's geometry.The stabilization energy holds for the global stabilization of a complex with respect to the isolated ͑noninteracting͒ fragments.Hence, both the interaction and the stabilization energies take into account the electronic relaxation, as the wave function of the monomers is computed to obtain the corresponding energies.In the case of the stabilization energy, the nuclear relaxation of the monomers is also taken into account.In contrast, the static interaction energy account only for local energetic interactions extracted uniquely from the complex's wave function, and is easily obtained as the summation of all the CECA energy components associated to intermolecular two-center interactions.In a similar way, the static monomer energies can also be obtained by collecting all the one-and two-center CECA components involving the atoms of the given monomer.The summation of all BSSE corrections to each ͑static͒ monomer energy and to the ͑static͒ intermolecular component yields the total correction to the complex energy.The overall BSSE correction calculated as the difference between the uncorrected ͑SCF//SCF͒ and corrected ͑CHA//SCF͒ energies is also reported.Comparison of these values gives a measure of the accuracy of the CECA partition in each case.
For the 6-31G basis sets, the BSSE correction is manifested mainly in the energy components related to H 5 , which is the H participating in the intermolecular bond.The principal stabilizing contribution comes from the one-center component in H 5 , while the major destabilizing contributions correspond to two-center components involving H 5 and other atoms.Thus, E BSSE (H 5 ) is ϩ10.0 kcal•mol Ϫ1 , while E BSSE (O 1 ,H 5 ) and E BSSE (O 4 ,H 5 ) are Ϫ9.4 and Ϫ7.0 kcal •mol Ϫ1 , respectively.However, these trends are not general for all the calculations.For all the basis sets without diffuse functions, ͑a͒, ͑b͒, ͑c͒, and ͑e͒, E BSSE (O 1 ,H 5 ) is ϳϪ10 kcal•mol Ϫ1 .E BSSE (H 5 ) and E BSSE (O 4 ,H 5 ) also contribute to the BSSE, but to a small extent, compared to the 6-31G results.E BSSE (H 5 ) is always positive, but E BSSE (O 4 ,H 5 ) can be positive or negative, depending on the basis set.In some cases, other components exhibit also significant BSSE.In general, for all these basis sets, the BSSE correction stabilizes the two water monomers, especially the donor one, but makes the intermolecular component less stable.The overall effect of the BSSE in the static interaction TABLE VII.CECA analysis of the formic acid dimer.Given values represent energetic differences between SCF and CHA at the SCF ͑uncorrected͒ geometries in kcal/mol.See Scheme 3͑a͒ for the selected one-and two-center energy differences.⌬E A-D , and ⌬E int are the static differences on the formic acid moiety and the interaction energy, respectively, computed from the CECA one-and two-center terms.BSSE and BSSE C give the exact and the CECA approximated difference between the CHA and the SCF energies.
Basis   Ϫ1 , respectively.In terms of molecular and intermolecular components, the BSSE-correction contribution to the intermolecular term is always favorable ͑ϳϩ9 kcal•mol Ϫ1 in both cases͒.For the 6-31ϩϩG(d,p) and 6-311ϩϩG(d, p) basis sets, the overall BSSEcorrection contribution is negative for both the donor and acceptor molecules.Altogether, the effect of BSSE correction on the molecular static interaction energy is always destabilizing, but the sign of the contributing terms is reversed, compared to the calculation with no diffuse functions.
Table V also lists also the BSSE calculated ͑i͒ as the difference between SCF//SCF and CHA//SCF energies, and ͑ii͒ as the summation of the BSSE in each one-and twocenter energy component.The difference between the two values can be used to estimate the accuracy of the CECA.In general, the differences are significant, taking into account that the BSSE is generally small.For the basis sets with no diffuse functions, the difference is always less than 0.8 kcal •mol Ϫ1 , and the CECA always underestimates the magnitude of the BSSE.On the contrary, for the 6-31ϩ ϩG(d,p) and 6-311ϩϩG(d, p) basis set, the CECA overestimates the magnitude of the BSSE by ϳ1.5 and 3.5 kcal •mol Ϫ1 , respectively.Similar conclusions can be drawn when comparing SCF//SCF and CHA//CHA energies, so the results are not reported.
The results for the formic acid dimer are presented in Table VII The overall picture is the same in both calculations: The BSSE correction destabilizes the two-center components related to the H-bond interactions, as well as the one-center components in the acceptor atoms, but stabilizes the H atoms participating in the intermolecular bond.Altogether, the BSSE correction stabilizes each formic acid monomer but decreases the attractive intermolecular energy component.Hence, the overall contribution of the BSSE correction to the molecular interaction is destabilizing for both basis sets.For this complex, the error in the CECA analysis is quite small, compared to (H 2 O) 2 .
The results of the CECA analysis for the uracil-water complex are collected in

DISCUSSION
It is interesting to remark that the main effects of the BSSE correction on the electron density of the water dimer are very similar to those found previously for the hydrogen fluoride dimer.The patterns of electron redistribution caused by the removal of the BSSE at frozen geometries for (HF) 2 and (H 2 O) 2 are very similar.Indeed, for calculations without diffuse functions, the main feature of the difference maps is the redistribution of electron density in the valence shells of the heavy atoms in both cases.Moreover, similar trends are found for the 6-31G(d,p) and 6-31G(d) calculations on the formic acid dimer and uracil-water complexes, respectively.Furthermore, addition of diffuse functions leads to similar effects for all the systems analyzed: An overall decrease of the differences between corrected and uncorrected densities, negative differences in the intermolecular region, and lack of the highly directional density redistribution patterns in heavy donor and acceptor atoms that are observed with smaller basis sets.
In fact, some of the differences between the SCF//SCF and CHA//SCF electron densities appear to be at odds with simple chemical intuition.For instance, it might be expected that the BSSE correction should weaken the intermolecular interaction and therefore lead to a decrease of the electron density in the intermolecular region.Actually, in many cases, the BSSE correction works in the opposite direction, leading to an accumulation of electron density in the intermolecular region.Moreover, the BSSE correction also decreases the electron density in the intramolecular bonds of the donor moieties.In general, it should be taken into account that the CHA//SCF results used in the difference maps do not correspond to stationary points on the BSSE-corrected surface.It is well-known that geometry relaxation is necessary for fully correcting the BSSE.In fact, when nuclear relaxation is taken into account, there is always a depletion of the electron density in the intermolecular region, as reflected in the properties of the intermolecular bcp's.
The study of larger systems, like the formic acid dimer, and especially the uracil-water complex, reveals that the effects of BSSE on the electron density are generally restricted to the intermolecular region and especially to the atoms directly involved in the intermolecular interaction and their first-neighbors.
The CECA decomposition scheme has been found to be a valuable tool for analyzing the effects of BSSE correction in terms of atomic and interatomic contributions.However, one has to be aware that the CECA decomposition is not exact.Therefore, the applicability of this method to analyze the subtle effects of the BSSE correction on the molecular energy depends on the accuracy of the approximation.In general, for the calculations reported in this paper, the accuracy of the decomposition, calculated as the difference to the true BSSE correction, is good or acceptable when basis sets without diffuse functions are used.In these cases, the results of the CECA analysis are in agreement with chemical intuition: The BSSE correction generally stabilizes the purely intramolecular energies of the two molecules forming the complex, but it dwindles the intermolecular energy component.The final result is that BSSE correction always leads to less attractive interaction energies.In general, it is worth to note that, although the BSSE in the total molecular energy is usually small, the individual atomic or interatomic contributions can be quite large.
When diffuse functions are taken into account, the results of the analysis are just the opposite.That is, in general, small energy destabilization results from the combination of a large destabilization of the intramolecular energies and a stabilization of the intermolecular term upon BSSE correction.However, the validity of the CECA analysis in these cases is questionable for two reasons.First, the accuracy of the CECA decomposition is very low ͑the BSSE is overestimated by several kcal•mol Ϫ1 , except for the formic acid dimer͒.Second, the identification of the one-and two-center components with atomic and interatomic contributions is doubtful when diffuse functions are involved.Nevertheless, the CECA results in these particular cases seem to agree with the density difference maps in the sense that the differences tend to be smaller when diffuse functions are included but they are more delocalized.Indeed, the H atoms involved in the intermolecular H-bonds have similar BSSE effect on the two-center components with the acceptor atom and the neighbors of this atom.In other words, since the diffuse functions are so spread in space and hardly assigned to a given nuclei, the BSSE is not only energetically localized in the bonds but also to the same extent in nonbonded interactions.This effect is in agreement with the observations of corresponding isocontour density difference maps.

CONCLUSIONS
We have carried out a detailed comparison of the local effects of the BSSE on the electron densities and energy components of three representative H-bonded complexes.These results complement previous studies of the effects of the BSSE on the geometries, energies and electron densities of these and other complexes.In general, we have found that the effects of the BSSE are common for all the complexes studied.The elimination of the BSSE by means of the CHA always leads to lower interaction energies.When nuclear relaxation is taken into account, the BSSE correction also leads to larger interatomic distances and a decrease of the electron density at all the intermolecular critical points on the electron density.Density difference maps at frozen geometry reveal that the effects of the BSSE are not limited to the intermolecular region.Rather, the main redistribution effects take place in the valence shells of the heavy atoms directly involved in the intermolecular interaction.For the larger complex, uracil-water, the effects of the BSSE do not extend significantly beyond the atoms involved in the interaction and their first neighbors.
These trends are confirmed by means of an energy decomposition analysis performed with the CECA method.In general, the BSSE effects on the energy are centered also in the components involving the atoms participating in intermolecular interactions.In general, two-center terms related to intermolecular components and one-center terms centered on the H bonding account for a large part of the BSSE.However, other components can also make non-negligible contributions to the total BSSE.
The BSSE is inherent to the expansion of the molecular wave function in terms of basis functions.Therefore, the size and characteristics of the basis set is one of the main factor influencing the magnitude of the BSSE.In agreement to previous results, we found that inclusion of diffuse functions is of utmost importance in order to minimize the magnitude of the BSSE.Moreover, the present study confirms that the origin of the BSSE in terms of one-and two-center contributions is very different depending on the inclusion of diffuse functions.In general, the energy decomposition analysis reveals that the small BSSE correction is the result of a near cancellation between larger errors in the intramolecular and intermolecular components.When no diffuse functions are considered, the BSSE correction is destabilizing for the intermolecular component and stabilizing for the intramolecular components.When diffuse functions are used, exactly the opposite is found.One should take into account that the CECA analysis has probably only a qualitative value when diffuse functions are considered.Anyway, since density difference maps do also reveal systematic differences between calculations with and without diffuse functions, the CECA analysis may well be meaningful as a complement for the understanding of intermolecular interactions.

TABLE I .
Geometrical parameters ͑Å and degrees͒, stabilization energies ͑kcal/mol͒ and topological parameters of the intermolecular critical points of the electron density for the (H 2 O) 2 calculated in five different basis sets at SCF and CHA/F levels of theory.

TABLE II .
Geometrical parameters ͑Å and degrees͒, stabilization energies ͑kcal/mol͒ and topological parameters of the intermolecular critical points of the electron density for the (H 2 O) 2 calculated in five different basis sets at B3LYP and CHA/B3LYP levels of theory.

TABLE III .
Electron density and its Laplacian at the intermolecular critical points located on the first-order electron density of the formic acid dimer, using two different basis sets.SCF and CHA/F values ͑in parenthesis͒ are reported.All the calculations have been performed at the geometry optimized on the SCF PES.
a Reference 22.

TABLE IV .
Electron density and its Laplacian at the intermolecular critical points located on the first-order electron density of the uracil-water complex ͓see Scheme 3͑b͔͒, using two different basis sets.SCF and CHA/F values ͑in parenthesis͒ are reported.All the calculations have been performed at the geometry optimized on the SCF PES.
a Reference 22.

TABLE V .
SCF one-and two-center energy components for the (H 2 O) 2 complex for the 6-31G(d,p) and 6-31ϩϩG(d,p) ͑lower triangle, in italics͒ in a.u..The values in parentheses correspond to the H 2 -H 3 diatomic term.

TABLE VI .
CECA analysis of the (H 2 O) 2 at the SCF level of theory for several basis sets.Given values represent energetic differences between SCF//SCF and CHA//SCF calculations in kcal/mol.E BSSE (O 1 -H 5 ), (O 4 -H 5 ) and (H 5 ), hold for the BSSE contribution on selected two-and one-center interactions ͑see Scheme 2͒.⌬E D , ⌬E A , and ⌬E int are the static BSSE contributions on donor, acceptor and interaction energies, respectively, computed from the CECA one-and two-center terms.Last two columns give the exact and the CECA approximated BSSE.
Table VIII.The main trends are very similar to those found for the formic acid dimer.Thus, for both basis sets, the BSSE correction introduces a large destabilizing contribution to the two-center components related to the H bonds, reflected in the large negative values of E BSSE (O 7 ,H 13 ) and E BSSE (O 14 ,H 8 ), while E BSSE (H 8 ) and E BSSE (H 13 ) are positive.The main difference between the