The accuracy of molecular bond lengths computed by multireference electronic structure methods
Introduction
Extensive comparisons of computed equilibrium bond lengths for 20 molecules using various single-reference methods, ranging in computational level from SCF to CCSD(T), have been reported by Helgaker et al. [1] and by Bak et al. [2]. In the present paper, we compare these results with bond lengths computed for the same molecules using multireference methods based on MCSCF reference functions of a simple GVB-RCI direct-product form [3], [4], [5]. Systematic sequences of atomic orbital basis sets are used to extrapolate to the infinite basis set limits for all of the methods.
This study explores several general trends. First, the computed MCSCF Re values tend to be longer than the experimental results. The MCSCF wave functions used in this study include only valence correlation effects. The extra-valence dynamical correlation is described using multireference CI methods, and these more flexible wave functions tend to produce a shortening of the computed bond lengths, relative to the MCSCF values, for the molecules in this study. Second, we observe that more flexible orbital basis sets also tend to result in shorter computed bond lengths for the multireference methods. Thus, we observe a general trend of systematic convergence from long computed bond lengths toward the experimental values as the orbital basis sets are improved and as more flexible wave function expansions are used. This is in contrast to the SCF computed bond lengths [1] which are systematically shorter than the experimental values, and to the corresponding single-reference methods which often do not show smooth convergence to the experimental values as the basis sets are improved.
We attribute these bond length trends to the ionic contamination of SCF wave functions. As is well known in the simple example of H2, the SCF wave function dissociates to an equal mixture of neutral and ionic products, introducing a −1/R contribution to the dissociation potential energy function. This spurious ionic contamination tends to make the computed supermolecule dissociation energy asymptote too high (relative to the neutral open-shell RHF fragments) and the −1/R contribution tends to make the computed equilibrium bond length too short. This same general principle also applies to the SCF description of single and multiple bonds in molecules. In molecules with multiple bonds, or in the simultaneous distortion of several bonds, each bond contributes its own separate ionic contamination to the SCF dissociation energy curve. Compared to a correctly dissociating MCSCF wave function, the ionic contamination results in a function Δ(R) = Escf(R) − Emcscf(R) that is rigorously nonnegative (due to the Ritz variational principle), and at long bond distances is characterized asymptotically as the monotonically increasing function Δ(R) ≈ Δ(∞) − α/R, in which the scalar α depends on the nuclear charges of the fragments and on the nature of the distorted bond(s). Near the equilibrium structure, Δ(R) does not display this simple asymptotic functional form, but it tends nonetheless to be increasing with R (i.e. dΔ/dR > 0), showing the remnants of this long-range asymptotic behavior. This feature results in computed SCF Re values that tend to be shorter than the MCSCF Re values, although due to various other features of the wave functions near the equilibrium structures (such as orbital hybridization differences), this is rather more of a general trend than a rigorous relation.
The molecules studied are all closed-shell singlet molecules at their equilibrium structures. While the SCF wave function, used as the zero-order reference function for the single-reference methods, suffers from ionic contamination, the correlated methods based on this reference function are usually expected to perform well for such systems near their equilibrium geometry. Nevertheless, traces of the ionic contamination of the SCF reference function persist at some levels of the single-reference correlated methods, becoming more pronounced as bonds are stretched, and thus affecting the calculated equilibrium bond lengths. The additional wave function flexibility of the correlated methods must eliminate this spurious ionic contamination in the SCF reference wave function (which would tend to increase the computed bond lengths toward the experimental values), while simultaneously describing the dynamical correlation energy (which would tend to decrease the computed bond lengths). Due to these counteracting effects, the single-reference methods do not generally show the same kind of systematic convergence with orbital basis sets and with increasing wave function flexibility as the multireference methods. An important objective of the multireference methods is to eliminate the incorrect dissociation behavior of the reference function, and thus improve the performance of the correlated methods.
Improved treatment of the electron correlation beyond the MCSCF description can consist of including additional configuration state functions (CSFs) in the wave function expansion. In the present work, these additional CSFs are generated by allowing single and double excitations from the occupied orbitals in the MCSCF reference into the set of external correlating orbitals. There are two energy functions of this form that are used in the present work, the MR-SDCI energy and the MR-AQCC energy. The MR-SDCI energy is the straightforward application of the Ritz variational principle in which the ground state energy is both an upper bound to the exact (full-CI) energy within the orbital basis, and it is variational with respect to the wave function expansion coefficients [6]. The MR-AQCC energy, which accounts for the most important size-extensivity effects [7], is variational in the wave function expansion coefficients, but it is not a rigorous upper bound to the full-CI energy. The variational nature of both of these methods allows the Hellmann–Feynman theorem to be exploited, allowing for very efficient computation of analytic energy gradients [8], [9], [10], [11].
The purpose of this work is to extend the single-reference-based geometry optimization results of Refs. [1], [2] to include MCSCF, MR-SDCI, and MR-AQCC results. These data have their own practical scope, completely independent of the analysis in terms of ionic contamination. Nonetheless, the above general trends regarding dynamical correlation and orbital basis set effects will be demonstrated, reasons for the exceptions will be discussed, and the accuracy of the computed geometries will be compared quantitatively to those of the single-reference methods.
Section snippets
Methodology
The first issue to be addressed is the specific choice of MCSCF and MR-SDCI wave function expansions. It is desirable to choose a general wave function expansion form and then use that same form for all of the molecules. One of the simplest forms is the GVB-RCI wave function expansion [3], [4], [5]. This is a direct-product wave function expansion in which pairs of electrons are assigned either to a single orbital with double occupancy, or to a pair of active orbitals whose individual
Calculations
The MCSCF, MR-SDCI, and MR-AQCC results for each of the molecules in the sample set will now be discussed. These molecules are presented in the order of increasing complexity of the MCSCF expansion space. A direct-product wave function expansion is defined by disjoint subsets of invariant orbital subspaces [8], [3]. A feature of such expansions is that the orbitals within a subspace may be transformed arbitrarily among themselves without changing the wave function. The MCSCF expansion space for
Overall statistics
The optimized bond lengths are determined for the MCSCF, the MR-SDCI, and the MR-AQCC wave functions for each of the selected molecules and for each of the three orbital basis sets as discussed in the previous section. The computed bond lengths Re,i are compared to the experimental values, , and the mean errors , the standard deviation of these errors , the absolute errors ∣Δi∣, and the mean absolute errors , are determined for the
Timings
A molecular structure optimization consists of a sequence of single point calculations, each of which involves the computation of an energy and of the analytic gradient of the energy with respect to the 3 Natom cartesian coordinates of the molecule. The computation of the MCSCF energy consists of an AO integral computation step and a wave function optimization step in which both the molecular orbitals and the CSF expansion coefficients are optimized. The associated analytic energy gradient step
Conclusions
Although the multireference methods used in this study would be expected to produce accurate and reliable global potential energy surfaces, they are lacking with respect to the better single-reference methods in an important aspect – they do not account for any correlation effects that involve three or more electrons in the external (extra-valence) orbital space. For this set of closed-shell molecules near their equilibrium structures, valence correlation effects would not be expected to be
Acknowledgements
The authors gratefully acknowledge communication of the full data set of bond lengths from T. Helgaker for the earlier single-reference calculations. R.S. acknowledges support from the Office of Basic Energy Sciences, Division of Chemical Sciences, Geosciences, and Biosciences, U.S. Department of Energy under Contract Number DE-AC02-06CH11357. H.L. acknowledges support by the Austrian Science Fund within the framework of the Special Research Program F16 (Advanced Light Sources). P.G.S.
References (63)
- et al.
Chem. Phys. Lett.
(1993)J. Chem. Phys.
(1995) - et al.
J. Chem. Phys.
(2001) - et al.
J. Chem. Phys.
(2000) - et al.
Molecular Spectra and Molecular Structure Constants of Diatomic Molecules
(1979) - et al.
J. Chem. Phys.
(2005) Can. J. Phys.
(1958)- et al.
J. Mol. Spectrosc.
(1978) Mol. Phys.
(1974)- et al.
Mol. Phys.
(1986) - et al.
J. Chem. Phys.
(1992)