DEVELOPMENT AND APPLICATIONS OF COUPLED-CLUSTER METHODS AND POTENTIAL ENERGY SURFACE EXTRAPOLATION SCHEMES By Jesse J. Lutz A DISSERTATION Submitted to Michigan State University in partial fulfillment of the requirements for the degree of DOCTOR OF PHILOSOPHY Chemistry 2011 ABSTRACT DEVELOPMENT AND APPLICATIONS OF COUPLED-CLUSTER METHODS AND POTENTIAL ENERGY SURFACE EXTRAPOLATION SCHEMES By Jesse J. Lutz The generation of highly accurate potential energy surfaces (PESs) for reactive processes represents a difficult challenge for modern electronic structure theory. Since chemical reactions often involve breaking and forming bonds or intermediate and transition state species, one must employ a methodology that provides a balanced and highly accurate description of varying levels of electronic degeneracy, but that is also practical enough to be applied to a wide range of chemical problems. Using small to medium sized systems, we examine the performance of two classes of coupled-cluster (CC) methods which are capable of accounting for the diverse electron correlation effects encountered in the majority of groundand excited-state PES considerations. The first class of methods are the size-extensive completely renormalized (CR) CC approaches for ground-states and their equation of motion (EOM) CC extensions for excited-states, in which noniterative corrections due to higherorder excitations are added to the energies obtained with the standard CC and EOMCC approximations, such as CCSD (CC with singles and doubles) or EOMCCSD (EOMCC with singles and doubles), respectively. In particular, we focus on the left-eigenstate CR-CC(2,3) and CR-EOMCC(2,3) methods, in which a noniterative correction due to triple excitations is added to the CCSD or EOMCCSD energy, respectively, and, when necessary, a noniterative correction for quadruple excitations is also included via the CR-CC(2,3)+Q approach. A new variant of the CR-EOMCC(2,3) method, abbreviated as δ-CR-EOMCC(2,3), that can pro- vide a size-intensive treatment of excitation energies, is discussed as well. The second class of methods considered here is the active-space variants of the electron-attached (EA) and ionized (IP) EOMCC theories, which utilize the idea of applying a linear electron-attaching or ionizing operator to the correlated, ground-state CC wave function of an N -electron closed-shell system in order to generate the ground and excited states of the related (N ± 1)electron radical species of interest. These approaches use a physically motivated set of active orbitals to a priori select the dominant higher-order correlation effects to be included in the calculation, which significantly reduces the costs of the high-level EA- and IP-EOMCC approximations needed for obtaining accurate results for open-shell species without sacrificing accuracy. We have also developed a general extrapolation strategy for reducing the cost of generating PESs with correlated electronic structure methods using the concept of correlation energy scaling. Benchmark studies were performed to demonstrate typical accuracies for two types of PES extrapolation schemes, namely, the single-level PES extrapolation schemes, in which the essential quantity, the correlation energy scaling factor, is generated using only the quantum chemistry method of interest, and the dual-level PES extrapolation schemes, where lower-order approaches are used to estimate the correlation energy scaling factor corresponding to the method of interest. Unifying features of these PES extrapolation techniques are discussed, including the role of pivot geometries and base wave functions, and PES extrapolation to the complete basis set limit is examined as well. Finally, the most essential details of the new open-shell EOMCCSD and EA- and IP-EOMCC computer codes for the GAMESS software package, developed as part of this thesis research, are described. Copyright by JESSE J. LUTZ 2011 This dissertation is dedicated to my loving and ever supportive parents who taught me the most important things I could ever learn. v ACKNOWLEDGMENT I owe a great deal of gratitude to my Ph.D. advisor, Professor Piotr Piecuch for countless opportunities to contribute important research which never would have been possible otherwise. He and his group were nurturing and supportive in many ways, particularly regarding research, but never limited to it. I would also like to thank the rest of my guidance committee, namely my second reader, Professor Katharine L. C. Hunt, as well as Professor Robert I. Cukier and Professor James E. Jackson for lending their time and support, and for all their comments and advice regarding research and life in general. It is also important to mention my great appreciation for the efforts of Professor Marta Wloch and Dr. Jeffery R. Gour towards my professional development. As a result of their able guidance, I progressed much more quickly then I would have without them, and I am deeply grateful for their help and collaboration in many of the projects described in this dissertation. I would also like to thank other members of the Piecuch group, both past and present, for the help and companionship they provided, including Dr. Maricris Lodriguito, Dr. Wei Li, Dr. Jun Shen, and Mr. Jared Hansen. Finally, I would like to acknowledge Michigan State University and the Department of Chemistry for several fellowships that supported part of my research. I would also like to acknowledge the MSU High Performance Computing Center for providing resources which were utilized for part of this research. Most of the research discussed in this thesis was supported by the U. S. Department of Energy, Office of Science, Office of Basic Energy Sciences, through grants awarded to my advisor and research fellowships to me, and I am appreciative for this support as well. vi TABLE OF CONTENTS List of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . List of Figures ix . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . xiv 1 Introduction 1 2 Project Objectives 17 3 Applications of Coupled-Cluster and Equation-of-Motion Coupled-Cluster Methods 19 3.1 Theory . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19 3.1.1 Single-Reference Coupled-Cluster Theory for Ground States . . . . . 20 3.1.2 Equation-of-Motion Coupled-Cluster Theory for Electronically Excited, Electron-Attached, and Ionized States . . . . . . . . . . . . . . . . . 25 3.1.2.1 Excited States . . . . . . . . . . . . . . . . . . . . . . . . . 26 3.1.2.2 Electron-Attached and Ionized States . . . . . . . . . . . . . 29 3.1.3 The Method of Moments of Coupled-Cluster Equations . . . . . . . . 32 3.1.3.1 Completely Renormalized Coupled-Cluster and Equation-ofMotion Coupled-Cluster Approaches . . . . . . . . . . . . . 35 3.1.3.2 Biorthogonal Formulation of the MMCC Equations and the CR-CC(2,3) and CR-EOMCC(2,3) Methods . . . . . . . . . 43 3.1.3.3 A Size-Intensive Variant of CR-EOMCC(2,3): The δ-CREOMCC(2,3) Method . . . . . . . . . . . . . . . . . . . . . 49 3.1.4 The Active-Space Coupled-Cluster and Equation-of-Motion CoupledCluster Approaches . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51 3.2 Applications . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 56 3.2.1 The DBH24 Benchmark Database for Thermochemical Kinetics . . . 57 3.2.2 Addition Reactions of Ethylene and Acetylene to Ozone . . . . . . . 67 3.2.3 Mechanism of the Isomerization of Bicyclobutane to Butadiene . . . . 75 3.2.4 Excited-State Potential Energy Surfaces for the Dissociation of Water 81 3.2.5 Excitation Energies and Hydrogen-Bonding-Induced Spectral Shifts in Complexes of cis-7-Hydroxyquinoline . . . . . . . . . . . . . . . . . . 91 3.2.6 Geometries and Adiabatic Excitation Energies of CNC, C2 N, N3 , and NCO . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 100 4 Potential Energy Surface Extrapolation Schemes 129 4.1 Motivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 129 vii 4.2 4.3 4.4 Theory . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 134 Single-Level Potential Energy Surface Extrapolation Schemes . . . . . . . . . 137 4.3.1 Potential Energy Surface Extrapolation to Larger Basis Sets . . . . . 137 4.3.2 Potential Energy Surface Extrapolation to the Complete Basis Set Limit138 4.3.3 The Role of Pivot Geometries and Base Wave Functions in Single-Level PES Extrapolations . . . . . . . . . . . . . . . . . . . . . . . . . . . . 140 4.3.4 Application to the Isomerization of Bicyclobutane to Butadiene . . . 141 Dual-Level Potential Energy Surface Extrapolation Schemes . . . . . . . . . 153 4.4.1 Using Lower-Order Methods to Produce Correlation Energy Scaling Factors in High-Level Calculations . . . . . . . . . . . . . . . . . . . 154 4.4.2 Application to Single Bond-Breaking Potential Energy Curves . . . . 156 4.4.3 Application to the Isomerization of Bicyclobutane to Butadiene . . . 163 4.4.4 The Role of Pivot Geometries and Base Wave Functions in Dual-Level PES Extrapolations . . . . . . . . . . . . . . . . . . . . . . . . . . . . 166 5 Development of Computer Codes for the GAMESS software Package 5.1 Key Details of Efficient Computer Implementations . . . . . . . . . . . . . . 5.1.1 Standard Equation-of-Motion Coupled-Cluster Theory with Singles and Doubles for Open-Shell Systems . . . . . . . . . . . . . . . . . . 5.1.2 Electron-Attached and Ionized Equation-of-Motion Coupled-Cluster Theories . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 174 175 6 Summary and Concluding Remarks 187 References . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 191 viii 176 184 LIST OF TABLES 3.1 Representative barrier heights database DBH24 taken from Ref. [204]. . . . 59 3.2 Mean signed error (MSE) and mean unsigned error (MUE) of coupled cluster methods calculated with all electrons correlated compared against the DBH24 benchmark database (in kcal/mol). . . . . . . . . . . . . . . . . . . . . . . . 61 Mean signed error (MSE) and mean unsigned error (MUE) of coupled cluster methods calculated with frozen core approximation compared against the DBH24 benchmark database (in kcal/mol). . . . . . . . . . . . . . . . . . . 62 Errors resulting from CR-CC(2,3) calculations of forward and reverse barrier = = = = heights, Vf and Vr respectively, reported as ǫ(Vf ) / ǫ(Vr )) at varying basis set levels compared against DBH24 benchmark values (in kcal/mol). . . . . 66 Mean signed errors (MSE) and mean unsigned errors (MUE) resulting from CR-CC(2,3) calculations of all DBH24 forward and reverse barrier heights, = = = = Vf and Vr respectively, reported as ǫ(Vf ) / ǫ(Vr )) at varying basis set levels (in kcal/mol). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 68 3.6 Benchmark values for the C2 H2 and C2 H4 ozonolysis reaction pathways. . . 70 3.7 Energetics of stationary points relative to reactants, in kcal/mol, produced at various levels of coupled-cluster theory for the ozonolysis of ethylene and acetylene. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74 Relative enthalpies, in kcal/mol, with respect to reactant as calculated at various levels of theory for the conrotatory and disrotatory bicbut → tbut isomerization pathways. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 78 Ground-state X 1 A′ energies for the asymmetric single-bond breaking of H2 O at a series of the O-H bond lengths, R, in atomic units. The last two rows give the mean unsigned errors (MUE) and non-parallelity errors (NPE) relative to the full CI PES obtained in Ref. [271]. . . . . . . . . . . . . . . . . . . . . . 83 3.3 3.4 3.5 3.8 3.9 3.10 Same as Table (3.9) for the three lowest-lying 1 A′ states. ix . . . . . . . . . . 86 3.11 Same as Table (3.9) for the three lowest-lying 3 A′ states. . . . . . . . . . . 87 3.12 Same as Table (3.9) for the two lowest-lying 1 A′′ states. . . . . . . . . . . . 88 3.13 Same as Table (3.9) for the two lowest-lying 3 A′′ states. . . . . . . . . . . . 89 3.14 The basis-set dependence of the vertical excitation energies ωπ→π ∗ and the environment-induced shifts ∆ωπ→π ∗ (in cm−1 ) obtained with the EOMCCSD approach corresponding to the lowest π → π ∗ transition in the cis-7HQ chromophore and its complexes with the water and ammonia molecules. . . . . 95 3.15 The vertical excitation energies ωπ→π ∗ and the environment-induced shifts ∆ωπ→π ∗ (in cm−1 ) obtained with the EOMCCSD/6-31+G(d), EOMCCSD/6311+G(d), δ-CR-EOMCC(2, 3), A/6-31+G(d), and δ- CR-EOMCC(2, 3), D/631+G(d) approaches, and their composite EOMCC corresponding to the lowest π → π ∗ transition in the cis-7HQ chromophore and its various complexes. 97 3.16 Total and adiabatic excitation energies for the ground and low-lying excited states of CNC, as obtained with the different EA-EOMCC approaches using the DZP [4s2p1d] and cc-pVXZ (X = D, T, Q) basis sets and extrapolating to the CBS limit. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 123 3.17 Total and adiabatic excitation energies for the ground and low-lying excited states of C2 N, as obtained with the different EA-EOMCC approaches using the DZP [4s2p1d] and cc-pVXZ (X = D, T, Q) basis sets and extrapolating to the CBS limit. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 124 3.18 Comparison of the optimized equilibrium geometries for the low-lying states of CNC and C2 N, as obtained with the EA-EOMCC and SAC-CI-SDT-R/PS approaches using the DZP[4s2p1d] and cc-pVXZ (X = D, T, Q) basis sets. . 125 3.19 Total and adiabatic excitation energies for the ground and low-lying excited states of NCO, as obtained with the different IP EOMCC approaches using the DZP [4s2p1d] and cc-pVXZ (X = D, T, Q) basis sets and extrapolating to the CBS limit. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 126 3.20 Total and adiabatic excitation energies for the ground and low-lying excited states of N3 , as obtained with the different IP EOMCC approaches using the DZP [4s2p1d] and cc-pVXZ (X = D, T, Q) basis sets and extrapolating to the CBS limit. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 127 x 3.21 Comparison of the optimized equilibrium geometries for the low-lying states of N3 and NCO, as obtained with the IP EOMCC and SAC-CI-SDT-R/PS approaches using the DZP [4s2p1d] and cc-pVXZ (X = D, T, Q) basis sets. 128 4.1 Computer costs of typical ab initio wave function calculations at the aug-ccpVTZ basis set level, taken from Ref. [204] . . . . . . . . . . . . . . . . . . 133 4.2 The calculated CR-CC(2,3)/cc-pVDZ and CR-CC(2,3)/cc-pVTZ energies and the calculated and extrapolated CR-CC(2,3)/cc-pVQZ energies at the stationary points defining the conrotatory and disrotatory pathways characterizing the bicbut→t-but isomerization. The bicbut reactant defines the pivot geometry for the PES extrapolations. . . . . . . . . . . . . . . . . . . . . . . . 143 4.3 The calculated CR-CC(2,3)/cc-pVDZ and CR-CC(2,3)/cc-pVTZ energies and the calculated and extrapolated CR-CC(2,3)/cc-pVQZ energies at the stationary points defining the conrotatory and disrotatory pathways characterizing the bicbut→t-but isomerization. The con TS transition state defines the pivot geometry for the PES extrapolations. . . . . . . . . . . . . . . . . . . 144 4.4 The calculated CR-CC(2,3)/cc-pVDZ and CR-CC(2,3)/cc-pVTZ energies and the calculated and extrapolated CR-CC(2,3)/cc-pVQZ energies at the stationary points defining the conrotatory and disrotatory pathways characterizing the bicbut→t-but isomerization. The dis TS transition state defines the pivot geometry for the PES extrapolations. . . . . . . . . . . . . . . . . . . 145 4.5 The calculated CR-CC(2,3)/cc-pVDZ and CR-CC(2,3)/cc-pVTZ energies and the calculated and extrapolated CR-CC(2,3)/cc-pVQZ energies at the stationary points defining the conrotatory and disrotatory pathways characterizing the bicbut→t-but isomerization. The g-but intermediate defines the pivot geometry for the PES extrapolations. . . . . . . . . . . . . . . . . . . . . . 146 4.6 The calculated CR-CC(2,3)/cc-pVDZ and CR-CC(2,3)/cc-pVTZ energies and the calculated and extrapolated CR-CC(2,3)/cc-pVQZ energies at the stationary points defining the conrotatory and disrotatory pathways characterizing the bicbut→t-but isomerization. The gt TS transition state defines the pivot geometry for the PES extrapolations. . . . . . . . . . . . . . . . . . . 147 4.7 The calculated CR-CC(2,3)/cc-pVDZ and CR-CC(2,3)/cc-pVTZ energies and the calculated and extrapolated CR-CC(2,3)/cc-pVQZ energies at the stationary points defining the conrotatory and disrotatory pathways characterizing the bicbut→t-but isomerization. The t-but product defines the pivot geometry for the PES extrapolations. . . . . . . . . . . . . . . . . . . . . . . . 148 xi 4.8 The calculated CR-CC(2,3)/cc-pVDZ, CR-CC(2,3)/cc-pVTZ, and CR-CC(2,3)/ccpVQZ energies and the CBS values of the CR-CC(2,3) energies obtained using the point-wise extrapolations exploiting Eq. (3.98) and the CBS extrapolation procedure combining Eqs. (4.1)–(4.3) with Eq. (3.98) discussed in the text at the stationary points defining the conrotatory and disrotatory pathways characterizing the bicbut→t-but isomerization. The bicbut reactant was used to define the pivot geometry for the PES extrapolations based on Eqs. (4.1)–(4.3). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 152 4.9 A comparison of calculated and extrapolated CR-CC(2,3)/aug-cc-pVQZ energies of the H2 O molecule in which one of the two O-H bonds (R) is stretched, while keeping the other O-H bond at the equilibrium length and the H-O-H angle fixed at 104.5 ◦ . The equilibrium geometry defines the pivot geometry Re and RHF defines the base wave function for the PES extrapolations. In all post-RHF calculations, the lowest orbital, correlating with the 1s shell of the oxygen atom was kept frozen. . . . . . . . . . . . . . . . . . . . . . . . 160 4.10 A comparison of calculated and extrapolated CR-CC(2,3)/aug-cc-pVQZ energies for several internuclear separations RH-Cl of the HCl molecule. The equilibrium geometry defines the pivot geometry Re and RHF defines the base wave function for the PES extrapolations. In all post-RHF calculations, the lowest five orbitals, correlating with the 1s, 2s, and 2p shells of Cl, were kept frozen. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 161 4.11 A comparison of calculated and extrapolated CR-CC(2,3)/aug-ccpVQZ energies for several internuclear separations RF-F of the F2 mole- cule. The equilibrium geometry defines the pivot geometry Re and RHF defines the base wave function for the PES extrapolations. In all post-RHF calculations, the lowest two orbitals, correlating with the 1s shells of the F atoms, were kept frozen. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 162 4.12 A comparison of calculated and extrapolated CR-CC(2,3)/cc-pVQZ energies at the stationary points defining the conrotatory and disrotatory pathways characterizing the bicbut→t-but isomerization. The bicbut reactant defines the pivot geometry Re and RHF defines the base wave function for the PES extrapolations. In all post-RHF calculations, the lowest four orbitals, correlating with the 1s shells of the carbon atoms, were kept frozen. . . . . 165 xii 4.13 A comparison of calculated and extrapolated CR-CC(2,3)/aug-cc-pVQZ energies of the H2 O molecule, in which one of the two O-H bonds R is stretched, while keeping the other O-H bond at the equilibrium length and the H-O-H angle fixed at 104.5 ◦ . Each geometry serves as its own pivot geometry and RHF defines the base wave function for the PES extrapolations. In all postRHF calculations, the lowest orbital, correlating with the 1s orbital of the oxygen atom, was kept frozen. . . . . . . . . . . . . . . . . . . . . . . . . . 169 4.14 A comparison of calculated and extrapolated CR-CC(2,3)/aug-cc-pVQZ energies for several internuclear separations R of the HCl molecule. Each geometry serves as its own pivot geometry and RHF defines the base wave function for the PES extrapolations. In all post-RHF calculations, the lowest five orbitals, correlating with the 1s, 2s, and 2p shells of Cl, were kept frozen. . . . . . . 170 4.15 A comparison of calculated and extrapolated CR-CC(2,3)/aug-cc-pVQZ energies for several internuclear separations R of the F2 molecule. Each geometry serves as its own pivot geometry and RHF defines the base wave function for the PES extrapolations. In all post-RHF calculations, the lowest two orbitals, correlating with the 1s orbitals on the fluorine atoms, were kept frozen. . . 171 4.16 A comparison of calculated and extrapolated CR-CC(2,3)/aug-cc-pVQZ energies at the stationary points defining the conrotatory and disrotatory pathways characterizing the bicbut→t-but isomerization. Each geometry serves as its own pivot geometry and RHF defines the base wave function for the PES extrapolations. In all post-RHF calculations, the lowest four orbitals, correlating with the 1s orbitals of the carbon atoms, were kept frozen. . . . 172 4.17 A summary of the necessary calculations and the corresponding computer resources required to utilize different tiers of the PES extrapolation scheme based on Eqs. (4.1) – (4.3) to scale the bicyclobutane isomerization pathway from the CR-CC(2,3)/cc-pVTZ level of theory to the CR-CC(2,3)/cc-pVQZ level, along with the corresponding extrapolation errors. The bicbut structure is used to provide the pivot geometry. . . . . . . . . . . . . . . . . . . 173 5.1 Explicit algebraic expressions for the one- and two-body matrix elements of ¯ (CCSD) HN,open . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 177 xiii LIST OF FIGURES 3.1 3.2 Stationary points along the C2 H2 (top row) and C2 H4 (bottom row) ozonolysis reaction pathways. In each row, the structures from left to right represent the van der Waals minimum (vdW), transition state (TS), and the cycloadduct, respectively. The oxygen, carbon, and hydrogen atoms are represented by the red, yellow, and grey spheres, respectively. For interpretation of the references to color in this and all other figures, the reader is referred to the electronic version of this dissertation. . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69 The conrotatory and disrotatory pathways describing the isomerization of bicyclo[1.1.0]butane to trans-buta-1,3-diene, along with the enthalpy values relative to the reactant at all stationary points obtained in the explicit CRCC(2,3)/cc-pVQZ//CASSCF(10,10)/cc-pVDZ calculations (numbers in roman) and the enthalpy values obtained with the PES extrapolation procedure using the CR-CC(2,3)/cc-pVDZ, CR-CC(2,3)/cc-pVTZ, and RHF/cc-pVQZ energies at all stationary points resulting from the CASSCF(10,10)/cc-pVDZ optimizations, and the CR-CC(2,3)/cc-pVQZ correlation energy at the reactant geometry (numbers in bold italics). The available experimental enthalpy values are in parentheses. All enthalpies are in kcal/mol. . . . . . . . . . . . 77 3.3 Cuts of the PESs for a few ground- and excited-states of a single-bond stretching model of H2 O as obtained with (a) the CCSD/EOMCCSD and (b) the CR-CC(2,3)/CR-EOMCC(2,3) approaches, and the TZ basis set. Lines are used to represent the CC/EOMCC data, while the corresponding full CI values are represented by points (which are identically replicated in (a) and (b)). 92 3.4 The eight hydrogen-bonded complexes of the cis-7HQ molecule. . . . . . . . xiv 94 Chapter 1 Introduction Potential energy surfaces (PESs) play a central role in the theoretical description of molecular structures, properties, and reactivities, making them of great utility in many areas of chemical research, including spectroscopy and kinetics, investigations of reaction mechanisms, and the design of force fields for biological and materials science applications. Unfortunately, obtaining a molecular PES is typically very difficult due to the enormous mathematical complexity one faces when trying to solve the electronic Schr¨dinger equation. Development o efforts in quantum chemistry have focused on the design of computationally manageable, yet reliable, approximation methods for the generation of energies and properties that can be applicable to a wide range of molecular systems. In the interest of reducing mathematical complexity some simplifying approximations are introduced from the outset, including the neglect of relativistic effects, which are small in systems with light atoms [1], and the decoupling of the nuclear and electronic wave functions, accomplished under the Born-Oppenheimer approximation [2], which allows the electronic energy of a molecule to be expressed as a function of its geometry, thus providing the conceptual basis for a PES. Despite these and 1 other commonly employed simplifications, the electronic Schr¨dinger equation remains too o formidable to be solved exactly for any system with two or more electrons [3]. One of the principal challenges in the conventional determination of a PES, in which one is required to solve the electronic Schr¨dinger equation point-by-point for each fixed nuclear geometry o of interest, is the development, implementation, and benchmarking of new approximate ab initio electronic structure methods that can provide suitably accurate results for the wide variety of chemical species commonly encountered when studying reactive PESs, while requiring only modest computational resources. Thus, the first goal of this dissertation is to make a significant contribution toward the development and benchmarking of leading modern methods of electronic structure theory. It is difficult at the outset to determine what characteristics are important when developing new ab initio methods of electronic structure theory. One of the ultimate goals of quantum chemistry is to provide a predictive tool which can guide experimental efforts. In order to be considered a quantitatively predictive model, an ab initio method must be able to reliably produce numerical values for reaction energies which are accurate to within ≃ 1 kcal/mol of widely acknowledged benchmark data, a threshold often called “chemical accuracy”. It should also be clear what cases a method is and is not appropriate for, and, when necessary, it should be clear what to do to improve a poor result. A shortcoming of many ab initio methods is that they do not offer a balanced treatment of the short-range or “dynamical” and long-range or “non-dynamical” electronic interactions and therefore are unreliable for applications involving chemically reactive processes, since the relative importance of these competing effects can vary rapidly moving from one region to another on the corresponding PES. The goal is to develop a method which can give a balanced treatment 2 of both types of electronic interactions, while also having a straightforward and logical way of improving the result in the case that the predictions are shown to be inadequate. The starting point for single-reference (SR) electronic structure methods is the independent particle model (IPM), usually Hartree-Fock (HF), approximation [4–7], and the more practical discretized algebraic form of the resulting equations based on the linear combination of atomic orbitals (LCAO) self-consistent field (SCF) formalism, which is the origin of the basis set approximation defining ab initio models [8]. The HF method produces the best single determinant description of the electronic wave function of interest and although it is by now well established that the HF approximation yields over 99% of the total energy, its ability to provide a reliable description of chemical phenomena is very limited due to the need to describe relatively small energy differences in characterizing chemical processes. The inadequacies of HF have been shown time and again, but as an outstanding early example, Wahl’s study of the F2 molecule demonstrated that an SCF description predicts that F2 is unbound [9]. Despite the completely unphysical relative energetics often produced by the HF approximation, the LCAO SCF formalism is still a valuable tool as it provides the HF model, a critical element to our conceptual understanding of electronic structure and chemical reactivity. To retain the HF picture and be able to make predictive calculations, methods which improve upon the basic molecular orbital approximation must be developed. The focus of these so-called ’post-HF methods’ is to accurately recover the small amount of energy neglected at the HF level, a quantity which is known as the correlation energy. The most popular ab initio post-HF approaches can be grouped into the variational and perturbative classes. Variational approaches provide an upper bound to the exact energy, as in the configuration interaction (CI) theory [10–14], while perturbative approaches, such as the 3 many-body perturbation theory (MBPT) [15–22] and coupled cluster (CC) methods [23–27], have other advantages which will be discussed shortly. To understand how these approaches may be utilized to produce powerful new ab initio electronic structure approaches, their individual strengths and weaknesses must first be understood. The CI method is the most straightforward way of accounting for electronic correlation. The CI wave function is a linear combination of the reference (usually HF) determinant and Slater determinants obtained by exciting electrons from occupied to unoccupied orbitals. If this determinantal expansion is complete (full CI), the correlation energy and total energy of the molecular system will be exact within whatever basis set approximation was used used in the calculation. Due to the fact that the number of Slater determinants grows factorially with the size of the system, full CI calculations are impractical for systems of more than a few electrons. Approximations are commonly made in order to reduce this expense, for example, by truncating the CI expansion to include only the simplest classes of excitations, singles and doubles (SD), which is called the CISD method. A particularly appealing feature of the CI methods is that, by systematically adding classes of excitations, e.g., CI with singles, doubles, and triples, CI with singles, doubles, triples, and quadruples, and so on, one may obtain a series of increasingly accurate energies approaching the full CI value. Without this systematically improvable nature, it is difficult to claim a particular level of convergence in electronic structure calculations unless there is benchmark data available for comparison. Unfortunately, truncated CI methods lack size-extensivity, i.e., the correct linear scaling with the number of electrons, and the hierarchy converges slowly with the rank of excitation included. While the lack of size-extensivity can be approximately accounted for at the CISD level by adding Davidson corrections [28, 29], the slow convergence with level of excitation 4 to the full CI limit makes the CI method comparatively inefficient in the context of other modern electronic structure methods. Perturbation theory offers another systematically improvable hierarchy of methods for determining the correlation energy. Although these methods are not variational, they are strictly size-extensive at every level. In the SR-MBPT approach, the electronic Hamiltonian is divided into an unperturbed part, corresponding to a single-determinantal reference description, usually the HF determinant, and a perturbation part, which usually describes the electronic correlation. In this form of MBPT, introduced by Møller and Plesset (MP) [15], the perturbation corrections to the reference wave function and energy are then calculated using the Rayleigh-Schr¨dinger perturbation theory. Systematic improvements in the eleco tronic energy may be made by considering perturbative corrections of increasing orders, i.e., second-order MP (MP2), third-order MP (MP3), and so on. Additionally, comparisons of each order of perturbation with various CI ranks have made possible further division within levels of perturbation, allowing for the selection of certain classes of excitations. As an example, one may choose to calculate only certain components within the fourth-order of MP theory, e.g., MP4 with only contribution from doubles (MP4D), doubles and quadruples (MP4DQ), singles, doubles, and quadruples (MP4SDQ), or singles, doubles, triples, and quadruples (full MP4). This allows for certain computationally inexpensive higher-order corrections to be included without significantly increasing time requirements. In general, if the unperturbed HF determinant describing the molecular system is close to the exact wave function for that system, the convergence of the MBPT series is usually very rapid; however, when chemical bonds are stretched, the MBPT series becomes divergent. Thus, MBPT produces very good energies near the equilibrium geometries of molecules; however, 5 unlike CI, in regions of PESs where bonds are broken or formed MBPT energies become unphysical, making them inadequate for describing reactive PESs. The conventional SR coupled-cluster (CC) theory is currently considered the preeminent ab initio method for ground-state calculations of modest-sized molecular systems. Contrary to the CI approach, which is characterized by a linear expansion of configuration state functions, the CC approach uses an exponential ansatz for the wave function, inherently assures that truncated forms of CC theory remain size-extensive, and produces a much faster convergence to the full CI wave function. It can be shown using a perturbation theory analysis that the improved convergence of CC as compared with the same level of truncation in CI theory is due to higher-order excitations being folded in as products of lower-order excitations by the exponential form of the CC ansatz. At the same time, thanks to the use of diagram factorization techniques commonly employed in efficient computer implementations of CC methods, the computer costs of CC calculations are similar to those characterizing the CI approaches truncated at the same excitation levels. This is why CC methods can offer higher accuracy at relatively lower costs as compared with CI or MBPT methods and even though the energies produced are not variational, they are typically considerably more accurate than those produced by CI at the same level of truncation. Excited states may also be accessed in CC theory through the equation-of-motion (EOM) CC formalism [30–34] or its symmetry-adapted-cluster configuration-interaction (SAC-CI) [35–39] or linear-response CC analogs [40–44]. The most appealing electronic structure methods come as a result of combining the different abovementioned approaches. EOMCC, for example, is just a CIlike expansion starting from an approximate ground-state wave function obtained using CC theory. Approaches based on the CC and EOMCC formalisms which utilize perturbation 6 theory to obtain inexpensive corrections due to higher-order correlation effects also exist and will be discussed in more detail shortly. All of these appealing features of CC theory make it the most promising electronic structure approach of those considered so far, but there are still open problems yet to be completely addressed in CC theory. For instance, the most widely used and computationally practical CC approximation, the CC with singles and doubles (CCSD) method [45–48], fails completely when applied to describe a PES involving bond breaking. The CCSD method, which is based on an iterative procedure with central processing unit (CPU) steps scaling as n2 n4 , where no and nu are the numbers of occupied and unoccupied orbitals in the refero u ence model, respectively, or as N 6 , where N is the size of the system expressed as the sum of the exponents of the basis functions, neglects the important triply excited, quadruply excited, and other higher-order clusters needed to describe bond breaking. Unfortunately, the CC method with singles, doubles, and triples (CCSDT) [49, 50] and the CC approach with singles, doubles, triples, and quadruples (CCSDTQ) [51–54] which include these clusters are prohibitively expensive for anything but small molecules, as they require iterative steps that scale as n3 n5 (N 8 ), and n4 n6 (N 10 ), respectively. A parallel problem exists o u o u in the case of the EOMCC or response CC methods where the basic singles and doubles approximation (EOMCCSD) [31–33], which is characterized by iterative n2 n4 (N 6 ) scaling o u steps, fails to describe excited states dominated by two-electron and other many-electron transitions out of the ground state. Again, the EOMCC theory with singles, doubles, and triples (EOMCCSDT) [55–58] and the EOMCC approach with singles, doubles, triples, and quadruples (EOMCCSDTQ) [59], which can describe such states, do not offer a suitable alternative in the majority of applications because of their prohibitively expensive iterative N 8 7 and N 10 scaling steps, respectively. The massive success of the SRCC theories in molecular applications [60–66], particularly for nondegenerate ground-states of molecules near their equilibrium geometry and electronically excited states dominated by one-electron transitions, in addition to substantial progress in recent years in code parallelization [67–77], and in the development of various local correlation CC techniques (see, e.g., Refs. [78–83]), including, for example the cluster-in-molecule CC method developed by the Piecuch group [84–86], continues to stimulate parallel effort toward the extension of CC theory to handle quasidegerate states characterizing bond breaking and many-electron excitations with the same or close to the same level of computational effort as that required by CCSD. The most natural way to handle quasi-degenerate electronic states is to either turn to one of the variants of multi-reference (MR) perturbation theory, such as the popular MCQDPT2 [87,88] and CASPT2 [89] methods, which are designed to handle large nondynamical correlation effects and low-order dynamical correlation effects, the MRCI approaches, such as the popular MRCI(Q) approximation [90, 91], or the genuine MRCC methods, which, in analogy to the previously discussed SR analogs, offer a better treatment of the dynamical correlation effects as compared with the MRMBPT or MRCI methods. The genuine MRCC methods can be categorized into two types. The first is the hierarchy of the Fock-space or valence-universal methods [92, 93], in which a single valence-universal wave operator operates on the system of interest and its ions which are obtained by removing one, two, etc. active electrons from active orbitals. This is a convenient formalism when one is interested in ionization potentials and electron affinities, however, these methods are unfavorable when a wide range of geometries must be considered, as is the case in the generation of reactive PESs. The other category of genuine MRCC methods, the Hilbert-space or state-universal 8 (SU) approaches [94], employ the Jeziorski-Monkhorst wave function ansatz in conjunction with the multi-root Bloch wave-operator formalism. The SU-MRCC approaches have been shown to produce very accurate results for the ground and excited states of systems undergoing severe geometrical transformations (see, e.g., Refs. [95–100]), making them a very attractive choice within the context of PES generation. However, their routine use for such applications is complicated by many factors. First, a model space of reference states must be pre-defined by the user, which often contains many states irrelevent to a given problem, and then a truncation scheme must be chosen which is compatible with the model space. This requires expert-level decisions on a case-by-case basis for each molecular problem. Then, even when appropriate choices are made, intruder states may appear [95–97], which can severely complicate interpretation of results. The aforementioned ambiguous parameters become even more undesirable when considered in the context of practical implementation. The need to accomodate such general choices in combination with the massive number of cluster amplitudes which must be computed make writing general SU-MRCC computer programs excessively difficult. While widespread routine use of genuine MRCC methods seems unlikely in the near future, these methods have provided great insight into ways to improve existing SRCC formalisms and inspired activity toward MRCC approaches dealing with a single quantum state. We refer the reader to Ref. [101], and references therein for further detailed discussion. Due to the difficulties with implementing practical, user-friendly MRCC methods and the prohibitive computational expenses characterizing the SR CCSDT and CCSDTQ approaches and their EOM analogs, substantial research effort has been directed toward developing approximate approaches for including higher-order correlation effects within a SRCC formal9 ism. This led to both approximate iterative approaches, such as CCSDT-n [102–104], and CCSDTQ-1 [105], as well as the more popular non-iterative approaches, such as CCSD[T] [104, 106] and CCSD(T) [107], where triply excited T3 clusters are approximated using perturbative arguments which take the form of a relatively inexpensive non-iterative n3 n4 (N 7 ) o u scaling step in addition to the cost of the underlying CCSD calculation. Methods were also derived for cases when quadruple excitations should not be disregarded, where perturbative noniterative corrections for both triples and quadruples are added, e.g., the CCSD(TQf ) method [108], characterized by a n2 n5 (N 7 ) scaling step. These and other noniterative o u perturbative CC approaches have become extremely popular because they efficiently account for most of the important connected triple or triple and quadruple excitations in a user-friendly (“black-box”) fashion, while avoiding the steep iterative N 8 or N 10 scaling steps required by full CCSDT or CCSDTQ. Unfortunately, their applicability is limited to molecules near the equilibrium geometries, since the perturbative arguments used to derive the noniterative corrections of CCSD(T) and similar approaches fail when bonds are streched or broken. A parallel challenge is found within the development of the standard response CC and EOMCC approximations for excited states, in which the effects of triply or triply and quadruply excited configurations are estimated using arguments originating from MBPT (see, e.g., Refs. [109–114]). The perturbatively corrected EOMCC methods, ˜ such as, for example, EOMCCSD(T) [110], where the basic EOMCCSD approximation is corrected for triples via a noniterative n3 n4 or N 7 scaling correction, break down for excited o u states having larger contributions due to doubly excited configurations, which is particularly common for excited-state PESs along bond breaking coordinates (see, e.g., Ref. [115]). Quite a few approaches have been suggested in recent years which attempt to overcome 10 the failures of the conventional perturbative SR CC and EOMCC methods at larger internuclear separations and for excited states dominated by many-electron transitions, while avoiding the complexity of the genuine MRCC approaches. Examples which will be considered here include the externally corrected CC methods, such as the reduced MRCC (RMRCC) approaches [116–122], the completely renormalized (CR) CC methods [65,115,123–132], and the active-space CC theories [54–57, 133–149] (see Ref. [101] for a recent review). All of these methods are related in that they are concerned with improving the description of bond-breaking processes and other cases involving electronic quasidegeneracies, while relying on a SR-like formulation. Each of these methods will now be discussed, paying particular attention to the strengths of each approach and making mention of selected recent applications in the literature. The externally corrected SRCC methods represent an alternative to the perturbativelyderived corrections in which the CCSD equations corrected for terms containing the triply (T3 ) and quadruply (T4 ) excited clusters are solved after replacing T3 and T4 amplitudes by their values obtained in the cluster analysis of some non-CC wave function which exhibits good behavior at large internuclear separations, such as the projected unrestricted HF [150, 151], valence bond [152–154], multiconfigurational SCF or complete active-space SCF (CASSCF) [155–157], or MRCI wave functions [116–121]. Although all of these methods help in bond breaking situations, the latter, MRCI-corrected CCSD approach, referred to as RMR-CCSD, and its RMR-CCSD(T) extension [121, 122], have shown the most substantial improvements in the CCSD results. Unfortunately, the generation of the MRCISD wave function is very expensive when compared to the CCSD approach and MRCI is not sizeextensive, which limits the applicability of the RMR-CC approaches to smaller systems. One 11 can partly address these deficiencies by turning to the so-called (N ,M )-CCSD methods [158]. This notation implies that an M -reference general-model-space (GMS) SU-MRCCSD calculation [99,100] is corrected for higher-order clusters by drawing the relevent information from M pertinent wave functions [159] produced by an N -reference MRCISD calculation. Extensive benchmarking and testing have shown the RMR-CCSD method and its aforementioned extension with a perturbative correction for triples, RMR-CCSD(T) [121, 122, 160–172], as well as the (N ,M )-CCSD methods [158, 173–181] to be quite accurate in practice, although all of these approaches are very complex and require a lot of expertise. The ground-state CR-CC [65, 123–131] and excited-state CR-EOMCC [65, 115, 123, 126, 129, 132] methods represent another class of approaches which were developed with the intention of removing the pervasive failings of the conventional CC/EOMCC perturbative methods. These approaches are based on the more general formalism of the method of moments of CC equations (MMCC) [65, 123–125, 127–129, 131, 182, 183]. In analogy to the conventional perturbative methods, the CR-CC and CR-EOMCC methods allow one to calculate noniterative state-specific energy corrections corresponding to selected higher-order excitations which are added to the energies obtained from conventional CC/EOMCC calculations, such as CCSD or EOMCCSD. These CR-CC and CR-EOMCC corrections are based on the asymmetric energy expressions and resulting moment expansions which form the underlying framework for all MMCC methods. The original CR-CC approaches, such as CR-CCSD(T) [124,125] suffered from small errors due to a lack of strict size-extensivity [65], but these issues were addressed in the more recent left-eigenstate CR-CC approaches, including CR-CC(2,3), which are based on yet another form of the of the moment expansion of the full CI energy that defines the biorthogonal MMCC formalism [127–132]. The main 12 advantage of the CR-CC, CR-EOMCC, and other MMCC approaches is that the resulting ground- and excited-state energies are not dependent on any choice of active orbitals or other subjective parameters that one has to choose in MR calculations. Another advantage is the computational expense of the CR-CC and CR-EOMCC methods, which is on the order of the conventional perturbative SRCC methods. For example, the noniterative triples correction in the CR-CC(2,3) approach scales as n3 n4 or N 7 , in analogy to the scalo u ing of CCSD(T). The CR-CC(2,3) method has already been proven to be very accurate and robust, particularly in applications involving single bond breaking [127–130,184–188], mechanistic studies involving biradicals [127, 128, 130, 131, 189–193], and singlet-triplet gaps in biradical/magnetic systems [130, 131, 192]. Similar successes have been reported for the CREOMCC(2,3) method [129, 132] and its CR-EOMCCSD(T) predecessor [65, 115, 123, 126], where noniterative, n3 n4 or N 7 scaling corrections for triples are added to the underlyo u ing EOMCCSD energies. For example, the CR-EOMCC(2,3) method has recently been shown to accurately reproduce adiabatic excitation energies for various closed- and openshell molecules which are believed to be dominated by two-electron transitions out of the ground-state [132]. The third and final effort toward constructing a SR formalism capable of handling stronger non-dynamical correlations which will be considered here are the active-space CC and EOMCC methods [54, 133–149] (see Ref. [101] for a recent review). By specifically targeting the higher-order cluster and excitation amplitudes which become large in such situations by assigning a small subset of active orbitals defining these excitations, highly accurate results may be obtained while avoiding much of the expense of the higher-order parent methods. For example, the CCSDt and CCSDtq methods are based on the idea of 13 selecting T3 and T4 clusters within the CCSDT and CCSDTQ systems of equations using active orbitals [54, 133–146]. When properly implemented, the most expensive CPU steps of 2 2 the full CCSDt and CCSDtq approaches scale as No Nu n2 n4 and No Nu n2 n4 , respectively, o u o u where No and Nu are the numbers of active occupied and unoccupied orbitals, respectively. All active-space approaches have a few distinct advantages over competing methods. As an example, take again the CCSDt and CCSDtq approaches. These methods recover the exact results of their parent CCSDT and CCSDTQ approaches, respectively, in the limit that all orbitals are assigned as active. They are also systematically improvable, approaching the CCSDT and CCSDTQ limits as the number of active orbitals is increased. Given an appropriate selection of a usually small number of active orbitals, the active-space CC results are typically virtually perfect when compared to the values produced by the parent methods. Another advantage is the fact that all active-space CC methods are characterized by relatively low computational scalings, which are small prefactors times the n2 n4 steps o u of the CCSD type. Naturally, following the development of the ground-state CCSDt and CCSDtq active-space approaches, the EOMCCSDt and EOMCCSDtq methods were developed, with the first implementation of EOMCCSDt and the proposal for all such EOMCCbased methods occurring in the Piecuch group [55–58]. In addition to the EOMCCSDt approach, another class of the active-space EOMCC methods was also developed by the Piecuch group, namely the active-space electron attached (EA) and ionized (IP) EOMCC approaches [147–149], which may be used to generate ground and excited states of valence open-shell systems out of a related closed-shell ground-state, as in radical species. Generally, the EA- and IP-EOMCC methods (see Refs. [66,147–149] and references therein for information) have the distinct advantage over the traditional open-shell EOMCC approaches based 14 on restricted open-shell HF (ROHF) or unrestricted HF (UHF) references in the fact that they automatically generate orthogonally spin-adapted wave functions (a difficult thing to accomplish in CC theory). For comparison, traditional ROHF- or UHF-based open-shell CC/EOMCC methods introduce spin-contamination to the resulting wave functions. Spincontamination can have non-negligible effects on the energy of the states generated, but, more importantly, it prohibits the designation of the spin-symmetry of a given state, causing the identification of a particular state to become most inconvenient. This is an important issue, particularly when the excited states of open-shell systems are examined, which will be returned to periodically throughout this dissertation. While promising electronic structure approaches have been and continue to be developed, many of which provide a highly accurate and balanced description of chemical species typically encountered while scanning molecular PESs, reliable methods are still prohibitively expensive when one faces the calculation of the hundreds or thousands of points which are typically needed to sample these PESs. Electronic structure methods could certainly benefit from an auxiliary approach for predicting points on the PES based on inexpensive calculations with either less expensive quantum chemistry approaches or smaller basis sets, which would ameliorate the staggering computational expense of generating hundreds or thousands of points using high levels of electronic structure theory. To address this problem, an ab initio extrapolation scheme has been proposed that predicts the PES corresponding to expensive high-level calculations from the results of a series of comparatively inexpensive lower-level calculations using the concept of correlation energy scaling. The PES extrapolation scheme of this type was originally suggested in Ref. [184] and was designed such that it can be used in conjunction with any typical post-HF or post-CASSCF electronic structure method. 15 Thus, the second principal goal of this dissertation is to develop and perform benchmark studies on the PES extrapolation schemes based on correlation energy scaling ideas, which allow one to generate chemically reactive molecular PESs with very little or almost no information from high-level calculations while retaining the desired accuracy that high-level ab initio electronic structure methods provide. Finally, after new methods are developed and shown to be useful in benchmark studies, general-purpose computer programs should be written and distributed based on the successful theories to give scientists and engineers around the world a tool for interpreting, or even predicting, experimental results. The third and final goal in this dissertation is to outline the key details for the computer implementations of the open-shell EOMCCSD and EAand IP-EOMCC methods which which were written for the GAMESS electronic structure software package [194], a freely available suite of computer codes with tens of thousands of registered users. The three goals of this dissertation, summarized in Chapter 2, are addressed in this thesis, with a chapter being devoted to each. Thus, Chapter 3 is devoted to the development and application of new ground and excited-state CC/EOMCC methods for chemically reactive systems including the CR-CC/EOMCC methods and the active-space EA- and IP-EOMCC approaches. Chapter 4 begins with a discussion motivating the need for auxiliary methods to aid in the generation of molecular PESs and moves toward a presentation of the theory and various applications which help demonstrate various ways the PES extrapolation schemes based on correlation energy scaling may be used. Chapter 5 covers the development of computer codes for the GAMESS software package in detail, presenting programmable equations for a few of the theories implemented in this work as well as a brief discussion of how they are solved in practice. 16 Chapter 2 Project Objectives The main objectives of this work are: A. Performing benchmark calculations using the new generations of CR-CC approaches, including barrier heights of hydrogen transfer, heavy-atom transfer, nucleophilic substitution, and unimolecular and association reactions, and PESs for addition and isomerization reactions involving species with varying degrees of electronic degeneracy in order to demonstrate what levels of theory are appropriate in different situations. B. Developing and performing benchmark applications for the new generations of CREOMCC approaches, including the calculations of vertical excitation energies and environment-induced spectral shifts of organic chromophores and the calculation of excited-state PESs along bond-breaking channels. C. Performing benchmark applications for the EA- and IP-EOMCC methods including geometry optimizations and the calculation of adiabatic excitation energies of small open-shell molecules. 17 D. Developing the PES extrapolation schemes based on correlation energy scaling and performing benchmark applications to demonstrate the full range of capabilities offered and typical accuracies which should be expected in practice. E. Outlining the key details of computer implementations of the ROHF-based EOMCCSD and RHF- or ROHF-based EA- and IP-EOMCC programs recently developed for GAMESS. 18 Chapter 3 Applications of Coupled-Cluster and Equation-of-Motion Coupled-Cluster Methods 3.1 Theory As explained in the Introduction, the SR CC and EOMCC methods are the preeminent methods for the determination of electronic energies and properties in chemistry. The majority of this dissertation is concerned with the new generations of the CC and EOMCC methods which are useful in situations where the conventional CC and EOMCC approximations fail. We begin by reviewing the conventional CC and EOMCC theories in Sects. (3.1.1) and (3.1.2), respectively. The CR-CC and CR-EOMCC approaches and the active-space CC and EOMCC methods, with special emphasis on their EA- and IP-EOMCC extensions are discussed afterwards, in Sects. (3.1.3) and (3.1.4), respectively. 19 3.1.1 Single-Reference Coupled-Cluster Theory for Ground States In the SRCC theory, the ground-state wave function |Ψ0 of an N -electron system is expressed using the exponential ansatz, |Ψ0 = eT |Φ , (3.1) where T is the cluster operator and |Φ is an IPM reference configuration, e.g., the HF determinant (throughout this thesis, the RHF or ROHF determinant). Typically, we truncate the many–body expansion of T at a conveniently chosen excitation level mT , to obtain an approximate T , i.e., T ≃ T (A) , hoping that one can reach the desired accuracies with mT << N . The truncated cluster operator T (A) defining the approximate CC method A is given by T (A) = mT Tn , (3.2) n=1 with i ...i Tn = i1 <··· 1, we solve the eigenvalue problem given by Eq. (3.18) by diagonalizing the similarity-transformed a ...a ¯ Hamiltonian H (A) , Eq. (3.6), in a space spanned by the excited determinants |Φi 1...inn , 1 (A) with n = 1, ..., mR , corresponding to the many-body excitation operators included in Rµ . In general, in order for Eq. (3.18) to hold and to obtain a size-intensive description [44, 196] of vertical excitation energies, mR should not exceed mT [34], but, as already mentioned, one typically chooses mR = mT . For example, in the EOMCCSD method (note that the EE-EOMCC methods are often abbreviated as simply EOMCC), where Rµ is approximated as (CCSD) Rµ = Rµ,0 + Rµ,1 + Rµ,2 , (3.20) with i rµ,a aa ai (3.21) ij Rµ,1 = (3.22) i,a and rµ,ab aa ab aj ai , Rµ,2 = i mA disregarded in the 1 conventional CC calculations a ...a ¯ i1 ...in M0,a ...an (mA ) = Φi 1...inn |H (A) |Φ , 1 1 (3.35) where mA = mT is the maximum level of excitation included in the CC calculation being corrected. The excited-state moments needed to correct the EOMCC energies resulting from truncating T and Rµ at the mA -body components (so that mT = mR = mA ) are projections a ...a of the EOMCC equations on the excited determinants |Φi 1...inn , 1 i1 ...in a ...a ¯ (A) (A) Mµ,a1 ...an (mA ) = Φi 1...inn |Hopen Rµ,open |Φ . 1 (3.36) These moments are central to the ground- and excited-state MMCC theory and will be shown to be crucial quantities for evaluating the desired CR-CC and CR-EOMCC energy corrections in the following sections. (A) Several ways of expressing the δµ i ...i 1 n corrections in terms of moments Mµ,a1 ...an (mA ) have been proposed to date [65,123–125,127–129,131,182,183]. The original and historically oldest formula, obtained in [124] for the ground states and [182] for excited states has the 33 following form: (A) δµ (A) ≡ Eµ − Eµ N n Ψµ |Cn−k (mA ) Mµ,k (mA )|Φ / = n=mA +1 k=mA +1 (A) Ψµ |Rµ eT (A) |Φ . (3.37) Here, Cn−k (mA ) = (eT (A) )n−k (3.38) (A) are the (n − k)-body components of the exponential wave operator eT , defining the CC method A, |Ψµ is the full CI ground- (µ = 0) or excited- (µ > 0) state, and i ...i Mµ,k (mA ) = i1 <··· mA , to determine the noniterative energy correction δµ , Eq. (3.33). The Cn−k (mA ) terms are very easy to calculate. The zero-body term, C0 (mA ), equals 1; the 2 one-body term, C1 (mA ), equals T1 ; the two-body term, C2 (mA ), equals T2 + 1 T1 if mA ≥ 2; 2 3 3 the three-body term C3 (mA ) equals T1 T2 + 1 T1 if mA = 2 and T3 + T1 T2 + 1 T1 if mA ≥ 3, 6 6 i ...i 1 n etc. The computation of moments Mµ,a1 ...an (mA ) for the most interesting cases of correct- ing the CCSD or EOMCCSD energies (mA = 2) is straightforward too, particularly if we limit ourselves to the corrections due to triples (k = 3) or quadruples (k = 4). As an example, if one is interested in recovering the exact ground-state energy E0 through 34 (CCSD) the addition of the full correction δ0 (CCSD) to the CCSD energy Eµ (where mA = 2), i1 ...ik one has to consider the generalized moments of the CCSD equations M0,a ...a (2) with 1 k k > 2. After a quick diagrammatic analysis, we can show that this seemingly long expansion i1 ...ik contains moments M0,a ...a (2) with k ≤ 6 only, since the electronic Hamiltonian contains 1 k (CCSD) only up to two-body interactions. Thus, all terms which make up the correction δ0 approaching the exact energy contain relatively few moments, namely, a ...a ¯ i1 ...ik M0,a ...a (2) = Φi 1...i k |H (CCSD) |Φ , 1 k 1 k k = 3 − 6. (3.40) The projections of the CCSD equations on higher-than-hextuply excited configurations do not have to be calculated, since for Hamiltonians containing up to two-body interactions i1 ...ik the generalized moments M0,a ...a (2) with k > 6 vanish. Similar simplifications occur 1 k in the case of correcting the excited-state EOMCCSD energies, where the only moments i1 ...ik M0,a ...a (2) that matter are those with k = 3 − 8. Although this is a considerable reduction 1 k of the computer effort, it is usually not computationally feasible to calculate up to six-body or higher moments to obtain a given MMCC correction. The CR-CC and CR-EOMCC methods discussed in Sects. (3.1.3.1)-(3.1.3.3) address this issue by focusing on the approximate i ...i 1 k corrections due to triples and quadruples which use moments Mµ,a1 ...ak (2) with k = 3 and 4 only and simplify Eq. (3.37) accordingly. 3.1.3.1 Completely Renormalized Coupled-Cluster and Equation-of-Motion CoupledCluster Approaches In general, the CR-CC and CR-EOMCC approaches are obtained by approximating |Ψµ in Eq. (3.37) by a quasi-perturbative form that brings information about the desired cor35 (A) relation effects we want Eµ to be corrected for. In particular, the CR-CCSD(T) [124], CR-CCSD(TQ) [124, 125] and CR-EOMCCSD(T) [115] methods, which are of interest in this dissertation, are obtained when the wave functions |Ψµ in Eq. (3.37) are approximated by low-order MBPT-like expressions. In the CR-CCSD(T) method which corrects (CCSD) the ground-state CCSD energy E0 for triply excited clusters, |Ψ0 in Eq. (3.37), where mA = 2, is replaced by the following second-order-type, MBPT(2)[SDT]-like expression CCSD(T) |Ψ0 [2] = (1 + T1 + T2 + T3 + Z3 )|Φ , (3.41) where T1 and T2 are the singly and doubly excited clusters obtained in the CCSD calculations, the [2] (3) T3 |Φ = R0 (VN T2 )C |Φ (3.42) term is an approximation of the connected triples (T3 ) contribution, which is correct through second order, and (3) Z3 |Φ = R0 VN T1 |Φ (3.43) is the disconnected triples correction, which distinguishes the CCSD(T) approach from its (3) CCSD[T] predecessor. Conventional MBPT notation is used here, in which R0 designates the three-body component of the MBPT reduced resolvent and VN is the two-body part of HN . In the CR-CCSD(TQ) method (note that here and elsewhere the so-called CRCCSD(TQ),b variant is implied and that a discussion of the other variants, which can be found, for example, in [65, 123], will be omitted for the sake of brevity), which corrects the ground-state CCSD energy for the combined effect of triples and quadruples, |Ψ0 in Eq. (3.37), where mA = 2, is replaced by the following second-order-type, MBPT(2)[SDTQ]-like 36 expression CCSD(TQ) CCSD(T) |Ψ0 CCSD(T) where |Ψ0 = |Ψ0 1 2 + 2 T2 |Φ , (3.44) is given by Eq. (3.41). Once these approximations have been established, the following compact formulas for the CR-CCSD(T) and CR-CCSD(TQ) energies can be written: (CR-CCSD(T)) E0 (CCSD) = E0 + N CR(T) /D(T) (3.45) and (CR-CCSD(TQ)) E0 (CCSD) = E0 + N CR(TQ) /D(TQ) , (3.46) where the N CR(T) and N CR(TQ) numerators are defined as [2] N CR(T) = Φ|(T3 )† M0,3 (2)|Φ + Φ|(Z3 )† M0,3 (2)|Φ (3.47) † 1 N CR(TQ) = N CR(T) + 2 Φ|(T2 )2 [T1 M0,3 (2) + M0,4 (2)]|Φ , (3.48) and CCSD(T) and the D(T) and D(TQ) denominators, representing the overlaps between the |Ψ0 CCSD(TQ) and |Ψ0 wave functions, Eqs. (3.41), and (3.44), respectively, with the CCSD 37 ground state, as in Eq. (3.37), are calculated as CCSD(T) T1 +T2 |e |Φ D(T) ≡ Ψ0 † † 1 2 = 1 + Φ|T1 T1 |Φ + Φ|T2 (T2 + 2 T1 )|Φ [2] 1 3 + Φ|(T3 )† (T1 T2 + 6 T1 )|Φ † 1 3 + Φ|Z3 (T1 T2 + 6 T1 )|Φ (3.49) and D(TQ) ≡ CCSD(TQ) T1 +T2 |e |Φ Ψ0 † 1 2 1 2 1 4 = D(T) + 1 Φ|(T2 )2 ( 2 T2 + 2 T1 T2 + 24 T1 )|Φ . 2 (3.50) The quantities M3 (2) and M4 (2) in Eqs. (3.47) and (3.48) are expressed in terms of the triply and quadruply excited moments of the CCSD equations, easily calculated as in Eqs. (3.39) and (3.40), with k set at 3 and 4, respectively. Specifically, ijk M0,abc (2) |Φabc , ijk M0,3 (2)|Φ = (3.51) i 0) CR-EOMCC(2,3) energies, Eµ (2, 3), can be calculated in the following manner: (CCSD) Eµ (2, 3) = Eµ + δµ (2, 3), (3.68) where ijk ℓabc Mµ,abc (2), µ,ijk δµ (2, 3) = Φ|Lµ,3 Mµ,3 (2)|Φ = (3.69) i 0) equations defined by Eqs. (3.52) and (3.61), respectively. Since the triply excited moments of the CCSD/EOMCCSD equations are already welldefined, the determination of the three-body amplitudes ℓabc entering Eq. (3.69) becomes µ,ijk the primary focus. Since the exact values of these amplitudes may not be obtained without solving the full CI problem for the exact bra state Ψµ |, a method for the approximate determination of the ℓabc amplitudes has to be proposed. Following Refs. [127, 128], the µ,ijk derivation of the ℓabc amplitudes used in practical CR-CC(2,3) and CR-EOMCC(2,3) calµ,ijk culations begins by defining the approximate form of the deexcitation operator Lµ which 44 parameterizes the full CI state Ψµ | as (CCSD) Lµ ≈ L µ + Lµ,3 , (3.70) where Lµ,3 is the three-body component of Lµ of interest, which is expressed in the usual way as ℓabc ai aj ak ac ab aa , µ,ijk Lµ,3 = (3.71) i 0 case, the EOMCCSD excitation energies ωµ , and in the : ijk abc ℓabc = Nµ,ijk /Dµ,abc , µ,ijk (3.76) ijk abc where the numerator Nµ,ijk and denominator Dµ,abc are defined as follows: abc Nµ,ijk = = (CCSD) ¯ (CCSD) abc H |Φijk Φ|Lµ (CCSD) ¯ Φ|[(Lµ,1 H2 (CCSD) ¯ +(Lµ,2 H2 )C ]|Φabc , ijk 46 (CCSD) ¯ )DC + (Lµ,2 H1 )DC (3.77) and (CCSD) ¯ − Φabc |H (CCSD) |Φabc ijk ijk (CCSD) ijk − Dµ,abc = Eµ = ωµ 3 (CCSD) ¯ Φabc |Hn ijk |Φabc . ijk (3.78) n=1 ijk Note that the denominator Dµ,abc used here is the same as that of the CR-EOM-CCSD(T) approach (see Eq. 3.62). The CR-CC(2,3) (µ = 0) and CR-EOMCC(2,3) (µ > 0) approaches ijk are obtained by substituting Eqs. (3.52) and (3.61) for Mµ,abc (2) and Eq. (3.76) for ℓabc , µ,ijk ijk abc where Nµ,ijk and Dµ,abc are given by Eqs. (3.77) and (3.78), respectively, into the triples correction formula, Eq. (3.69), which is subsequently added to the CCSD/EOMCCSD energy (CCSD) Eµ to obtain the total energy Eµ (2, 3), as in Eq. (3.68). In both the CR-CC(2,3) and CR-EOMCC(2,3) theories, the exact treatment of the ijk ijk Epstein-Nesbet-like denominator Dµ,abc , as in Eq. (3.78), where no terms in Dµ,abc are neglected, characterizes the most complete variant of the CR-CC(2,3) and CR-EOMCC(2,3) approaches designated as CR-CC(2,3),D or CR-EOMCC(2,3),D, respectively. By neglectijk ing selected terms in Eq. (3.78) for Dµ,abc , we obtain approximate CR-CC(2,3) and CREOMCC(2,3) schemes. Let us focus on CR-EOMCC(2,3) for this discussion, which contains CR-CC(2,3) as a special case corresponding to µ = 0. Variant C of the CR-EOMCC(2,3) theory, designated as the CR-EOMCC(2,3),C approach, is obtained by ignoring the three-body ¯ (CCSD) |Φabc term, while keeping the ¯ component of H (CCSD) in Eq. (3.78), i.e., the Φabc |H3 ijk ijk ijk ¯ contributions to Dµ,abc from the one- and two-body components of H (CCSD) intact. The CR-EOMCC(2,3),B approach is obtained by ignoring the two- and three-body components ¯ (CCSD) |Φabc in ¯ of H (CCSD) in Eq. (3.78), leaving only the one-body contribution Φabc |H1 ijk ijk 47 ijk Dµ,abc . Finally, variant A of the CR-EOMCC(2,3) approach is obtained by replacing the ijk Epstein-Nesbet-like denominator Dµ,abc , Eq. (3.78), by the Møller-Plesset-like denominator (CCSD) for triple excitations, ωµ − (ǫa + ǫb + ǫc − ǫi − ǫj − ǫk ), where ǫp ’s are the single-particle energies associated with spin-orbitals p (diagonal elements of the Fock matrix). An analogous discussion may be made for the ground-state case resulting in the A, B, C, and D variants of CR-CC(2,3). In analogy to the CR-CCSD(T) approach discussed in the previous subsection, one can extend the CR-CC(2,3) method to higher-than-triple excitations as in, for example, the CR-CC(2,4) scheme [127, 128, 166]. The CR-CC(2,4) approach, when implemented fully, combines the N 6 -type steps of CCSD with the N 7 (n3 n4 ) steps of CR-CC(2,3) needed to o u determine the triples correction, and the N 9 (n4 n5 ) steps needed to calculate the analogous o u correction due to quadruples. In order to address this CPU-time increase, Piecuch et al. proposed the so-called CR-CC(2,3)+Q method, also tested in this thesis, where one calculates the ground-state energy as follows [186] E CR-CC(2,3)+Q = E CR-CC(2,3) + E CR-CCSD(TQ) − E CR-CCSD(T) , (3.79) i.e., one adds the quadruples correction extracted from the CR-CCSD(TQ) calculations to the CR-CC(2,3) energy. This has an advantage over CR-CC(2,4) in the fact that the CPUtime costs of the CR-CCSD(TQ) calculations in the quadruples correction part scale as N 7 (n2 n5 ) with the system size N , as opposed to the N 9 steps of CR-CC(2,4). As a result, o u the CR-CC(2,3)+Q approach is almost as affordable as the CR-CC(2,3) method itself, while bringing information about connected quadruply excited clusters that become important in multiple bond breaking situations [186, 187]. 48 The CR-CC(2,3) approach is capable of breaking bonds and, unlike its CR-CCSD(T) predecessor, is size-extensive. Unfortunately, the δµ (2, 3) corrections to the EOMCCSD energies, defining CR-EOMCCSD(2,3), violate the property of size-intensivity of the EOMCC excitation energies [44, 115, 196]. Although this violation is often unimportant, it is useful to consider the possibility of restoring size intensivity in CR-EOMCC(2,3). This aspect is discussed next. 3.1.3.3 A Size-Intensive Variant of CR-EOMCC(2,3): The δ-CR- EOMCC(2,3) Method As shown in Refs. [132, 200], although the ground-state variants of CR-CC(2,3) are sizeextensive, their excited-state CR-EOMCC(2,3) analogs do not satisfy the property of sizeintensivity [44,126,196], i.e., the vertical excitation energy of a non-interacting system A+B, in which fragment A is excited, resulting from the CR-EOMCC(2,3) calculations, is not the same as that obtained for the isolated system A. The lack of size intensivity of the CREOMCC(2,3) excitation energies can be traced back to the presence of the size-extensive contribution [200, 201] ijk (rµ,0 ℓabc − ℓabc ) M0,abc (2) 0,ijk µ,ijk βµ = (3.80) ij>k,a>b>C t4 = I>J>k>l,a>b>C>D where I and J are summed over only active occupied orbitals and C and D are summed over only active unoccupied orbitals. The CCSDtq system of equations for the relevant ti , a ij Ijk tab , tabC , and tIJkl amplitudes has the form of CCSDTQ equations in which T3 and T4 are abCD replaced by t3 and t4 , respectively. We obtain Φa |CCSD + (HN t3 )C |Φ = 0, i (3.92) Φab |CCSD + [HN (t3 + T1 t3 + t4 )]C |Φ = 0, ij (3.93) 2 ΦabC |CCSD + [HN (t3 + T1 t3 + t4 + T1 t4 + T2 t3 + 1 T1 t3 )]C |Φ = 0, Ijk 2 53 (3.94) 1 2 ΦabCD | CCSD + [HN (t3 + T1 t3 + t4 + T1 t4 + T2 t3 + 2 T1 t3 IJkl 1 2 1 3 +T2 t4 + 1 t3 2 + 2 T1 t4 + T1 T2 t3 + 6 T1 t3 )]C |Φ = 0 , 2 (3.95) where t3 and t4 are defined by Eqs. (3.90) and (3.91), respectively. The CCSDt method could be obtained from the CCSDtq equations described here by setting t4 = 0 and solving the system of equations given by Eqs. (3.92) – (3.94) only. The active-space methods of the CCSDt and CCSDtq types were shown to be effective for excited-state theories as well, for example, in Refs. [55–57], where the EOMCCSDt approach, an active-space variant of EOMCCSDT was reported for the first time. While it is rather straightforward to generalize the active-space methods to particle conserving EOMCC theories, such as EOMCCSDT and EOMCCSDTQ, in this dissertation the focus is on the active-space variants of particle non-conserving EOMCC theories, in particular, the EA-EOMCCSD(3p-2h) and IP-EOMCCSD(3h-2p) methods. The full inclusion of the (N +1) Rµ,3p-2h and Rµ,3h-2p components of the Rµ (N −1) and Rµ operators in the EA- and IP- EOMCC calculations needed to obtain an accurate description of electronic excitations in radicals comes at a high price, increasing the N 5 -like no n4 and n2 n3 operations defining the u o u iterative diagonalization steps of the base EA-EOMCCSD(2p-1h) and IP-EOMCCSD(2h-1p) schemes to the N 7 -like n2 n5 and n3 n4 steps, respectively. One way to retain the accuracy of o u o u the higher-order EA- and IP-EOMCC schemes with the 3p-2h and 3h-2p excitations, while avoiding this steep computer cost increase, is to use the active-space variants of the EA/IPEOMCC methods with higher-than 2p-1h/2h-1p excitations described in Refs. [147–149]. In analogy to the ground-state active-space CC approaches described above, in the activespace EA- and IP-EOMCC methods one divides the available orbitals of the N -electron 54 reference system into core, active occupied, active unoccupied, and virtual categories, and (N +1) uses active orbitals to define the electron attaching and ionizing operators Rµ (N −1) Rµ and , respectively. In particular, the active-space EA-EOMCCSD(3p-2h){Nu } approach using Nu active unoccupied orbitals is obtained by replacing the 3p-2h component Rµ,3p-2h (N +1) of the electron attaching operator Rµ , Eq. (3.28), by jk rAbc aA ab ac ak aj . rµ,3p-2h = (3.96) j>k,Aj>k,b 2Re region, with CCSD appearing to be the only lower-order method which may be used successfully at all internuclear separations. This is especially evident for the F2 curve where results based on correlation energy scaling factors obtained with MP2, MP3, MP4D, MP4DQ, and MP4SDQ diverge badly at 4Re from the true CR-CC(2,3)/aug-cc-pVQZ curve, and while the use of CCSD instead leads to reasonable behavior. In fact, the maximum reported extrapolation errors are again found to be very close at all internuclear distances when the CCSD and CR-CC(2,3) correlation energy scaling factors are used to scale the CRCC(2,3)/aug-cc-pVTZ HCl and F2 PESs. For HCl, they are 1.891 and 1.664 millihartree, (CR-CC(2,3)) respectively, each corresponding to a recovery of about 90% of ∆E4,3 (CCSD) extrapolation errors found using χ4,3 (CR-CC(2,3)) ∆E4,3 (R). For F2 , (R) show a recovery of 93–99% of the corresponding (CR-CC(2,3)) (R) values, which may be compared to 95–100% obtained with χ4,3 159 . Table 4.9: A comparison of calculated and extrapolated CR-CC(2,3)/aug-cc-pVQZ energies of the H2 O molecule in which one of the two O-H bonds (R) is stretched, while keeping the other O-H bond at the equilibrium length and the H-O-H angle fixed at 104.5 ◦ . The equilibrium geometry defines the pivot geometry Re and RHF defines the base wave function for the PES extrapolations. In all post-RHF calculations, the lowest orbital, correlating with the 1s shell of the oxygen atom was kept frozen. (B) ǫ4 (χ4,3 )d R/Re a 0.75 0.90 1.00 1.10 1.25 1.50 2.00 3.00 4.00 E4 b ∆E4,3 c MP2 MP3 MP4D MP4DQ MP4SDQ CCSD CR-CC(2, 3) -76.259133 -76.351620 -76.363051 -76.355931 -76.329620 -76.276456 -76.199626 -76.162662 -76.161060 -22.676 -21.605 -21.059 -20.641 -20.278 -19.952 -19.168 -18.479 -18.303 -2.993 -2.759 -2.742 -2.723 -2.624 -2.639 -3.613 -4.086 -2.232 -0.763 -0.573 -0.576 -0.580 -0.509 -0.513 -1.217 -0.601 1.965 -0.611 -0.425 -0.428 -0.430 -0.349 -0.316 -0.788 0.729 3.688 -0.257 -0.046 -0.033 -0.019 0.080 0.130 -0.359 0.774 2.953 -0.073 0.121 0.131 0.145 0.253 0.338 -0.071 1.450 4.798 -0.186 0.031 0.054 0.085 0.225 0.394 0.294 0.912 1.124 -0.140 0.012 0.000 0.009 0.133 0.314 0.244 0.744 0.760 a The equilibrium value of R used here is R = 0.95785 ˚. b The calculated CR-CC(2,3)/aug-cc-pVQZ total energies in A e c Differences, in millihartree, between the actual CR-CC(2,3)/aug-cc-pVQZ and CR-CC(2,3)/aug-cc-pVTZ energies. hartree. d Differences, in millihartree, between the calculated and extrapolated CR-CC(2,3)/aug-cc-pVQZ energies, where the latter (B) energies were generated by applying the correlation energy scaling factors χ4,3 (R) obtained with B = MP2, MP3, MP4D, MP4DQ, MP4SDQ, CCSD, and CR-CC(2,3). The choice of B = CR-CC(2,3) is equivalent to the single-level PES extrapolation scheme. 160 Table 4.10: tions RH-Cl function for shells of Cl, A comparison of calculated and extrapolated CR-CC(2,3)/aug-cc-pVQZ energies for several internuclear separaof the HCl molecule. The equilibrium geometry defines the pivot geometry Re and RHF defines the base wave the PES extrapolations. In all post-RHF calculations, the lowest five orbitals, correlating with the 1s, 2s, and 2p were kept frozen. (B) ǫ4 (χ4,3 )d R/Re a 0.75 0.90 1.00 1.10 1.25 1.50 2.00 3.00 4.00 E4 b ∆E4,3 c MP2 MP3 MP4D MP4DQ MP4SDQ CCSD CR-CC(2, 3) -460.245772 -460.351216 -460.364178 -460.356856 -460.329720 -460.277589 -460.212229 -460.192428 -460.192532 -22.053 -20.979 -20.704 -20.550 -20.327 -19.823 -18.747 -18.345 -17.809 -4.017 -3.055 -2.610 -2.330 -2.178 -2.429 -3.755 -2.889 -0.558 -1.415 -0.537 -0.127 0.128 0.261 0.065 -0.828 0.849 2.958 -1.254 -0.387 0.017 0.268 0.407 0.276 -0.293 2.050 3.966 -0.834 0.033 0.439 0.692 0.834 0.697 0.055 1.675 1.670 -0.854 0.022 0.437 0.703 0.868 0.793 0.456 4.602 9.747 -0.896 -0.003 0.427 0.713 0.926 0.991 0.966 1.891 1.573 -1.273 -0.410 0.000 0.270 0.469 0.580 0.730 1.664 1.049 a The equilibrium value of R used here is R = 1.27455 ˚. b The calculated CR-CC(2,3)/aug-cc-pVQZ total energies in A e hartree. c Differences, in millihartree, between the actual CR-CC(2,3)/aug-cc-pVQZ and CR-CC(2,3)/aug-cc-pVTZ energies. d Differences, in millihartree, between the calculated and extrapolated CR-CC(2,3)/aug-cc-pVQZ energies, where the latter (B) energies were generated by applying scaling factors χ4,3 (R) obtained with B = MP2, MP3, MP4D, MP4DQ, MP4SDQ, CCSD, and CR-CC(2,3). The choice of B=CR-CC(2,3) is equivalent to the single-level extrapolation scheme. 161 Table 4.11: A comparison of calculated and extrapolated CR-CC(2,3)/aug-ccpVQZ energies for several internuclear separations RF-F of the F2 mole- cule. The equilibrium geometry defines the pivot geometry Re and RHF defines the base wave function for the PES extrapolations. In all post-RHF calculations, the lowest two orbitals, correlating with the 1s shells of the F atoms, were kept frozen. (B) ǫ4 (χ4,3 )d R/Re a 0.75 0.90 1.00 1.10 1.25 1.50 2.00 3.00 4.00 E4 b ∆E4,3 c MP2 MP3 MP4D MP4DQ MP4SDQ CCSD CR-CC(2, 3) -199.196897 -199.349657 -199.364732 -199.356794 -199.333764 -199.307123 -199.296727 -199.297648 -199.297985 -56.537 -53.009 -51.743 -51.172 -50.464 -49.741 -49.300 -49.368 -49.422 -5.098 -3.991 -3.835 -3.871 -4.209 -4.370 -1.232 6.748 12.149 -2.516 -1.752 -1.737 -1.793 -2.010 -1.655 3.454 13.398 18.901 -1.915 -1.132 -1.028 -0.929 -0.775 0.587 8.189 18.543 22.880 -1.208 -0.444 -0.378 -0.349 -0.389 0.339 5.241 0.429 30.269 -0.689 0.065 0.191 0.284 0.313 1.119 6.107 4.979 31.933 -1.120 -0.197 0.036 0.236 0.396 0.944 2.453 3.239 3.585 -1.011 -0.258 0.000 0.227 0.455 0.965 2.146 2.460 2.673 a The equilibrium value of R used here is R = 0.988351 ˚. b The calculated CR-CC(2,3)/aug-cc-pVQZ total energies in A e hartree. c Differences, in millihartree, between the actual CR-CC(2,3)/aug-cc-pVQZ and CR-CC(2,3)/aug-cc-pVTZ energies. d Differences, in millihartree, between the calculated and extrapolated CR-CC(2,3)/aug-cc-pVQZ energies, where the latter (B) energies were generated by applying the correlation energy scaling factors χ4,3 (R) obtained with B = MP2, MP3, MP4D, MP4DQ, MP4SDQ, CCSD, and CR-CC(2,3). The choice of B=CR-CC(2,3) is equivalent to the single-level PES extrapolation scheme. 162 In summary, MP2 consistently fails to predict the proper information about the CRCC(2,3) correlation energy scaling, while MP3 typically does much better. This leads to the conclusion that at least third-order correlation energy effects must be included for the proper description of the correlation energy scaling in the near-equilibrium PES regions. In regions further away from equilibrium, the third and even the partial fourth-order perturbative corrections still cannot produce the proper correlation energy scaling factors that would be compatible with those of CR-CC(2,3), which is a consequence of the inability of MBPT to describe bond-breaking, but one can use the CCSD approach instead, which is qualitatively correct at larger internuclear separations in single-bond breaking situations, providing a reasonable estimate of the CR-CC(2,3) correlation energy scaling with the basis set. 4.4.3 Application to the Isomerization of Bicyclobutane to Butadiene The results of dual-level PES extrapolations on the reaction profiles for the isomerization of bicyclobutane to butadiene are given in Table (4.12). The format of this table is similar to that of Tables (4.9)–(4.11), except that here the geometries of interest are the stationary points along the conrotatory and disrotatory pathways shown in Figure (3.2) and the cc-pVXZ rather than aug-cc-pVXZ basis sets are employed throughout. From the results presented in Table (4.12) it can be seen that the MP3 and various MP4 methods may be successfully used to probe the CR-CC(2,3) correlation energy scaling for the energetically favored conrotatory reaction profile, which consists entirely of species with correlation energy dominated by lower-order excitations, but in every case the MBPT methods produce a significantly larger error for the highly biradical dis TS geometry. On the other hand, when 163 (CCSD) χ4,3 (R) is used to extrapolate the CR-CC(2,3)/cc-pVQZ PES, the reported extrapo- lation errors remain within a millihartree of the explictly calculated CR-CC(2,3)/cc-pVQZ energy values at every stationary point, rivaling the sub-millihartree accuracies found using (CR-CC(2,3)) the single-level PES extrapolation scheme where χ4,3 (R) is employed instead. This is another case where a quasi-degeneracy, in this case resulting from the biradical nature of the dis TS configuration, inhibits the MBPT methods from producing the correct correlation energy scaling information. It is clear, of the methods considered here, that only the CCSD approach can offer a correlation energy scaling factor compatible with CR-CC(2,3), although the extrapolation errors obtained with MP4SDQ, which are 1.950 millihartree or ∼ 1.5 kcal/mol for the strongly biradical dis TS structure (located over 60 kcal/mol above the reactant) and less than 1 millihartree for the remaining structures, are excellent as well. Also, once again, MP2 correlation energies do not contain sufficient information to model the CR-CC(2,3) correlation energy scaling at any geometry. 164 Table 4.12: A comparison of calculated and extrapolated CR-CC(2,3)/cc-pVQZ energies at the stationary points defining the conrotatory and disrotatory pathways characterizing the bicbut→t-but isomerization. The bicbut reactant defines the pivot geometry Re and RHF defines the base wave function for the PES extrapolations. In all post-RHF calculations, the lowest four orbitals, correlating with the 1s shells of the carbon atoms, were kept frozen. (B) ǫ4 (χ4,3 )c Structure bicbut con TS dis TS g-but gt TS t-but E4 a ∆E4,3 b MP2 MP3 MP4D MP4DQ MP4SDQ CCSD CR-CC(2, 3) -155.695497 -155.625392 -155.582481 -155.733236 -155.728073 -155.738043 -44.149 -43.727 -42.447 -43.909 -43.764 -43.963 -9.521 -10.169 -13.436 -9.730 -9.815 -9.614 -0.860 -1.578 -4.643 -0.879 -0.934 -0.786 -0.397 -0.993 -3.890 -0.325 -0.374 -0.234 0.845 0.292 -2.671 0.899 0.844 0.972 0.760 0.397 -1.950 0.906 0.837 0.987 0.582 0.466 0.096 0.927 0.876 1.008 0.000 0.245 0.631 0.482 0.410 0.563 a The calculated CR-CC(2,3)/cc-pVQZ total energies in hartree. b Differences, in millihartree, between the actual CRCC(2,3)/cc-pVQZ and CR-CC(2,3)/cc-pVTZ energies. c Differences, in millihartree, between the calculated and extrapolated CR-CC(2,3)/cc-pVQZ energies, where the latter energies were generated by applying the correlation energy scaling factors (B) χ4,3 (R) obtained with B = MP2, MP3, MP4D, MP4DQ, MP4SDQ, CCSD, and CR-CC(2,3). The choice of B = CR-CC(2,3) is equivalent to the single-level extrapolation scheme. 165 4.4.4 The Role of Pivot Geometries and Base Wave Functions in Dual-Level PES Extrapolations In Tables (4.9)–(4.12) it was clear that regardless of the correlation energy scaling factor used in the PES extrapolation scheme based on Eqs. (4.1) – (4.3), the extrapolation errors are largest in regions where the nature of the electron correlation effects differ most from those found at the pivot geometry. To see how much the reported errors could be reduced by employing additional pivot geometries, the same sets of PES extrapolations were considered, but this time with all nuclear configurations treated as pivot geometries. In addition to introducing a new extrapolation method which should yield improved accuracies, this approach provides a direct measure of the error introduced when lower-order correlation energy scaling factors are employed to predict the results of higher-order calculations with a larger basis set, since any error due to the earlier assumption of the approximate transferability of the scaling factor from one geometry to another is eliminated. Additional calculations required to perform these extrapolations, when compared to those required for the extrapolations of the previous section (done with a single pivot geometry), consist of the remaining calculations required to obtain each PES using the aug-cc-pVQZ basis for the bond-breaking curves or the cc-pVQZ basis for the bicbut→t-but isomerization pathways using the lowerlevel methodology B. This is still relatively inexpensive, since we never have to perform the high-level CR-CC(2,3) calculation with the largest basis set employed at any geometry. A comparison of Tables (4.9) and (4.13), (4.10) and (4.14), (4.11) and (4.15), and (4.12) and (4.16) demonstrates the full extent of the accuracy which may be gained by using additional pivot geometries for each of the systems considered in this study. The most interesting detail to note in these comparisons is that when the same method is used to 166 predict the correlation energy scaling factor, increasing the number of pivot geometries does not necessarily result in an overall improvement in the quality of the extrapolated surface. Extrapolation errors are, for the most part, improved when the number of pivot geometries is increased, especially when MP4SDQ and CCSD are used to generate χ, but there are no significant benefits from switching from the previously discussed single-point approach to its multi-point analog. As stated before, the many options inherently included in the PES extrapolation scheme based on Eqs. (4.1) – (4.3) allow one to tailor it to make it more accurate or more affordable, as required by a given application. To demonstrate clearly and concisely the different levels of accuracy and savings in the computer effort which may be obtained using different tiers of the PES extrapolation scheme, the results of four different PES extrapolations examined here, using the bicbut→t-but isomerization as an example, are collected in Table (4.17). We recall that the goal of each extrapolation in this case is to predict the CR-CC(2,3)/cc-pVQZ PES (as represented by six structures on the corresponding conrotatory and disrotatory pathways) from the results of lower-level calculations. The most efficient PES extrapolation considered for this table, where a lower-order MP4SDQ scaling factor is used to extrapolate the CR-CC(2,3)/cc-pVTZ PES to the level of CR-CC(2,3)/cc-pVQZ calculations using one pivot geometry (bicbut) and RHF base energies, requires only 5% of the CPU time of the conventional method and the mean unsigned error (MUE) is slightly below 1 millihartree. The MUE is reduced by 0.314 millihartree when the correlation energy scaling factor obtained with MP4SDQ is replaced by that produced by CCSD, but the required CPU time is also increased to 10% of that required by the conventional CR-CC(2,3)/cc-pVQZ calculations. The results obtained with the original single-level PES extrapolation procedure are also 167 presented, where the time required to perform these calculations is 22% of that required conventionally, but the mean unsigned error falls to below 0.4 millihartree. Finally, when (CR-CC(2,3)) the base energy is obtained in CCSD calculations and χ4,3 (R) is used to scale the remaining correlation energy, an MUE of 0.170 is attained, but now the computational time savings amounts only to about one third of the conventional procedure. It can be seen in Table (4.17) that regardless of the PES extrapolation approach used, the MUE remains below 1 millihartree, which is a relatively insignificant loss in accuracy compared to the conventional, point-wise CR-CC(2,3)/cc-pVQZ calculations. By using the dual-level PES extrapolation scheme to extrapolate CR-CC(2,3)/cc-pVQZ energies, we have reduced the time required to accurately construct a PES for this problematic polyatomic isomerization by more than an order of magnitude. These computational savings would only grow larger if another system were considered which contained a greater number of electrons or a larger number of points on the PES were considered. 168 Table 4.13: A comparison of calculated and extrapolated CR-CC(2,3)/aug-cc-pVQZ energies of the H2 O molecule, in which one of the two O-H bonds R is stretched, while keeping the other O-H bond at the equilibrium length and the H-O-H angle fixed at 104.5 ◦ . Each geometry serves as its own pivot geometry and RHF defines the base wave function for the PES extrapolations. In all post-RHF calculations, the lowest orbital, correlating with the 1s orbital of the oxygen atom, was kept frozen. (B) ǫ4 (χ4,3 )d R/Re a 0.75 0.90 1.00 1.10 1.25 1.50 2.00 3.00 4.00 E4 b ∆E4,3 c MP2 MP3 MP4D MP4DQ MP4SDQ CCSD -76.259133 -76.351620 -76.363051 -76.355931 -76.329620 -76.276456 -76.199626 -76.162662 -76.161060 -22.676 -21.605 -21.059 -20.641 -20.278 -19.952 -19.168 -18.479 -18.303 -2.766 -2.737 -2.742 -2.769 -2.835 -3.072 -4.107 -5.201 -3.619 -0.504 -0.534 -0.576 -0.636 -0.733 -0.951 -1.651 -1.501 1.054 -0.360 -0.391 -0.428 -0.475 -0.538 -0.668 -1.073 0.009 3.143 -0.014 -0.016 -0.033 -0.060 -0.096 -0.190 -0.591 0.207 2.604 0.129 0.138 0.131 0.114 0.097 0.033 -0.301 0.856 4.400 0.059 0.063 0.054 0.041 0.034 0.021 -0.080 0.101 0.164 a The equilibrium value of R used here is R = 0.95785 ˚. b The calculated CR-CC(2,3)/aug-cc-pVQZ total energies in A e c Differences, in millihartree, between the actual CR-CC(2,3)/aug-cc-pVQZ and CR-CC(2,3)/aug-cc-pVTZ energies. hartree. d Differences, in millihartree, between the calculated and extrapolated CR-CC(2,3)/aug-cc-pVQZ energies, with the latter (B) energies generated by applying the correlation energy scaling factors χ4,3 (R) obtained with B = MP2, MP3, MP4D, MP4DQ, MP4SDQ, and CCSD. 169 Table 4.14: A comparison of calculated and extrapolated CR-CC(2,3)/aug-cc-pVQZ energies for several internuclear separations R of the HCl molecule. Each geometry serves as its own pivot geometry and RHF defines the base wave function for the PES extrapolations. In all post-RHF calculations, the lowest five orbitals, correlating with the 1s, 2s, and 2p shells of Cl, were kept frozen. (B) ǫ4 (χ4,3 )d R/Re a 0.75 0.90 1.00 1.10 1.25 1.50 2.00 3.00 4.00 E4 b ∆E4,3 c MP2 MP3 MP4D MP4DQ MP4SDQ CCSD -460.245772 -460.351216 -460.364178 -460.356856 -460.329720 -460.277589 -460.212229 -460.192428 -460.192532 -22.053 -20.979 -20.704 -20.550 -20.327 -19.823 -18.747 -18.345 -17.809 -2.707 -2.625 -2.610 -2.680 -2.885 -3.370 -5.151 -5.225 -2.748 -0.044 -0.087 -0.127 -0.220 -0.404 -0.739 -1.906 -0.655 2.016 0.123 0.056 0.017 -0.062 -0.204 -0.426 -1.204 0.927 3.643 0.516 0.466 0.439 0.373 0.249 0.046 -0.736 0.985 1.988 0.491 0.453 0.437 0.384 0.284 0.146 -0.265 4.375 11.510 0.447 0.424 0.427 0.401 0.362 0.386 0.336 0.933 0.533 a The equilibrium value of R used here is R = 1.27455 ˚. b The calcuated CR-CC(2,3)/aug-cc-pVQZ total energies in A e hartree. c Differences, in millihartree, between the actual CR-CC(2,3)/aug-cc-pVQZ and CR-CC(2,3)/aug-cc-pVTZ energies. d Differences, in millihartree, between the calculated and extrapolated CR-CC(2,3)/aug-cc-pVQZ energies, with the latter (B) energies generated by applying the correlation energy scaling factors χ4,3 (R) obtained with B = MP2, MP3, MP4D, MP4DQ, MP4SDQ, and CCSD. 170 Table 4.15: A comparison of calculated and extrapolated CR-CC(2,3)/aug-cc-pVQZ energies for several internuclear separations R of the F2 molecule. Each geometry serves as its own pivot geometry and RHF defines the base wave function for the PES extrapolations. In all post-RHF calculations, the lowest two orbitals, correlating with the 1s orbitals on the fluorine atoms, were kept frozen. (B) ǫ4 (χ4,3 )d R/Re a 0.75 0.90 1.00 1.10 1.25 1.50 2.00 3.00 4.00 E4 b ∆E4,3 c MP2 MP3 MP4D MP4DQ MP4SDQ CCSD -199.196897 -199.349657 -199.364732 -199.356794 -199.333764 -199.307122 -199.296727 -199.297648 -199.297985 -56.537 -53.009 -51.743 -51.172 -50.464 -49.741 -49.300 -49.368 -49.422 -3.415 -3.570 -3.835 -4.116 -4.758 -5.496 -3.680 4.354 9.504 -1.068 -1.389 -1.737 -2.030 -2.539 -2.676 1.008 11.655 17.002 -0.764 -0.877 -1.028 -1.063 -1.062 -0.019 6.381 18.778 22.892 -0.128 -0.215 -0.378 -0.466 -0.646 -0.237 3.741 -4.021 30.905 0.310 0.290 0.191 0.176 0.107 0.664 4.751 2.203 32.328 0.190 0.149 0.036 0.009 -0.115 -0.119 0.250 0.723 0.876 a The equilibrium value of R used here is R = 0.988351 ˚. b The calculated CR-CC(2,3)/aug-cc-pVQZ total energies in A e hartree. c Differences, in millihartree, between the actual CR-CC(2,3)/aug-cc-pVQZ and CR-CC(2,3)/aug-cc-pVTZ energies. d Differences, in millihartree, between the calculated and extrapolated CR-CC(2,3)/aug-cc-pVQZ energies, with the latter (B) energies generated by applying scaling factors χ4,3 (R) obtained with B = MP2, MP3, MP4D, MP4DQ, MP4SDQ, and CCSD. 171 Table 4.16: A comparison of calculated and extrapolated CR-CC(2,3)/aug-cc-pVQZ energies at the stationary points defining the conrotatory and disrotatory pathways characterizing the bicbut→t-but isomerization. Each geometry serves as its own pivot geometry and RHF defines the base wave function for the PES extrapolations. In all post-RHF calculations, the lowest four orbitals, correlating with the 1s orbitals of the carbon atoms, were kept frozen. (B) ǫ4 (χ4,3 )c Structure bicbut con TS dis TS gbut gt TS tbut E4 a ∆E4,3 b MP2 MP3 MP4D MP4DQ MP4SDQ CCSD -155.695497 -155.625392 -155.582481 -155.733236 -155.728073 -155.738043 -44.149 -43.727 -42.447 -43.909 -43.764 -43.963 -9.521 -10.455 -14.074 -10.245 -10.366 -10.220 -0.860 -1.614 -4.825 -1.009 -1.085 -1.007 -0.397 -1.036 -4.164 -0.473 -0.548 -0.472 0.845 0.322 -2.824 0.781 0.702 0.776 0.760 0.373 -2.131 0.728 0.652 0.729 0.582 0.295 -0.529 0.699 0.641 0.703 a The calculated CR-CC(2,3)/cc-pVQZ total energies in hartree. b Differences, in millihartree, between the actual CRCC(2,3)/cc-pVQZ and CR-CC(2,3)/cc-pVTZ energies. c Differences, in millihartree, between the calculated and extrapolated (B) CR-CC(2,3)/cc-pVQZ energies, with the latter energies generated by applying the correlation energy scaling factors χ4,3 (R) obtained with B = MP2, MP3, MP4D, MP4DQ, MP4SDQ, and CCSD. 172 Table 4.17: A summary of the necessary calculations and the corresponding computer resources required to utilize different tiers of the PES extrapolation scheme based on Eqs. (4.1) – (4.3) to scale the bicyclobutane isomerization pathway from the CR-CC(2,3)/cc-pVTZ level of theory to the CR-CC(2,3)/cc-pVQZ level, along with the corresponding extrapolation errors. The bicbut structure is used to provide the pivot geometry. Base Energya RHF RHF RHF CCSD Correlation Energy cc-pVQZ Calculations Required c CPU Time Mean Unsigned Scaling Factorb RHF CCSD CR-CC(2,3) (t/tconv )d Error (millihartree)e MP4SDQ CCSD CR-CC(2, 3) CR-CC(2, 3) 6 6 6 6 0 1 1 6 0 0 1 1 0.05 0.10 0.22 0.63 0.973 0.659 0.389 0.170 6 6 6 1 — Conventional Calculationf a The method used to generate the base energy. b The method B used to generate the correlation energy scaling factor χ(B) (R) 4,3 in Eq. (4.1). c The number of cc-pVQZ basis set calculations which must be performed using a given base wave function and a given correlation energy scaling factor to extrapolate the CR-CC(2,3)/cc-pVQZ PES. d The CPU time needed to perform the necessary calculations for each PES extrapolation type relative to the time needed to generate the true CR-CC(2,3)/ccpVQZ PES. e The mean unsigned error representing an extrapolated CR-CC(2,3)/cc-pVQZ reaction pathway generated using the designated base energy and correlation energy scaling factor. f Characteristics of the true PES calculation, in which each stationary point on the conrotatory and disrotatory pathways of the bicbut→t-but isomerization is calculated at teh CR-CC(2,3)/cc-pVQZ level. 173 Chapter 5 Development of Computer Codes for the GAMESS software Package A large portion of the current doctoral effort was devoted to writing EOMCC computer codes for the GAMESS software package. In this section, we discuss the highly efficient GAMESS implementations of the open-shell EOMCCSD and IP-EOMCCSD(2h-1p) methods developed as part of this thesis project, based on theory discussed in Sect. (3.1.2) and the corresponding factorized equations, in terms of recursively generated intermediates that lead to the vectorized computer codes through the use of fast matrix multiplication rountines from the BLAS library. The open-shell EOMCCSD and IP-EOMCCSD(2h-1p) codes were interfaced with previously existing ROHF and RHF/ROHF integral routines, respectively, available in the GAMESS software package [194], as well as the CC programs and routines for the generation of matrix elements of the similarity-transformed Hamiltonian of CCSD originally developed for GAMESS by the Piecuch group at Michigan State University. In Sect. (5.1), we begin our discussion of the implementation of these programs, with specific 174 details outlined for the open-shell EOMCCSD and IP-EOMCCSD(2h-1p) methods in Sects. (5.1.1) and (5.1.2), respectively. 5.1 Key Details of Efficient Computer Implementations Both the open-shell EOMCCSD and IP-EOMCCSD(2h-1p) codes must begin by solving the usual CCSD equations for the ground-state of the N -electron reference system in order to ij obtain the singly and doubly excited cluster amplitudes, ti and tab , respectively. In both a cases, this is done using the general ROHF-based CCSD codes included in GAMESS that work for closed- and open-shells, developed by the Piecuch group and described in [130]. ij Following the CCSD calculation, the converged ti and tab amplitudes are used to contruct a the one- and two-body matrix elements of the CCSD similarity-transformed Hamiltonian ¯ ¯ (CCSD) ¯ q HN,open , hp and hrs , respectively, which define the one- and two-body components of pq ¯ (CCSD) HN,open within the second quantized formalism, ¯p ¯ (CCSD) = hq ap aq , H1 (5.1) ¯ ¯ (CCSD) = 1 hrs N [ap aq as ar ], H2 4 pq (5.2) and respectively. N [. . .] is the normal product of the operators between the brackets and the Einstein summation convention over repeated upper and lower indices is assumed throughout. The explicit equations defining these matrix elements in terms of the matrix elements q rs of the Hamiltonian in the normal-ordered form fp = p|f |q and vpq = pq|v|rs − pq|v|sr , ij and CCSD cluster amplitudes ti and tab , are given in Table (5.1). Once these common a 175 initial steps are complete, the appropriate expressions for solving the EOMCCSD and IPEOMCCSD(2h-1p) eigenvalue problems have to be constructed. These steps, which are specific to the open-shell EOMCCSD and IP-EOMCCSD(2h-1p) codes considered here, are outlined in Sects. (5.1.1) and (5.1.2), respectively. Once the suitable EOMCCSD and IPEOMCCSD(2h-1p) equations are constructed, another common feature in our implementation of the EOMCCSD and IP-EOMCCSD(2h-1p) approaches is that we rely on the HiraoNakatsuji generalization [303] of the Davidson diagonalization algorithm [304] to solve the resulting non-Hermitian eigenvalue problems, Eqs. (3.24) and (3.25) in Sect. (3.1.2.1). 5.1.1 Standard Equation-of-Motion Coupled-Cluster Theory with Singles and Doubles for Open-Shell Systems ¯ (CCSD) The left-hand sides of the EOMCCSD equations are calculated by projecting [HN,open α (Rµ,1 + Rµ,2 )]C |Φ onto the subspace of all singly- and doubly-excited determinants, Φaα |, i aβ aβ bβ aα bβ aα α Φi |, Φiα jbα |, Φi j |, and Φiα j |, to obtain the following expressions: β β β β α ¯ (CCSD) Φaα |[HN,open (Rµ,1 + Rµ,2 )]C |Φ i i mα i ¯ ¯ ¯ = heα reα − hiαα raα + heαα raα mα aα α m m α eα e i e m β m ¯ β i ¯i ¯α β +hmβ raα mα + haα eαα reαα + haα mβ reβ eβ αm β i f e f m i 1¯ α 1¯ − 2 hmααα ra α nα + 2 haα nα reα nα n α α α fα α fα iα fβ mα nβ eα fβ iα nβ ¯ ¯ −hmα nβ ra f + haα nβ re f , α β α β 176 (5.3) Table 5.1: Explicit algebraic expressions for the one- and two-body matrix elements of ¯ ¯ (CCSD) ¯ q HN,open (hp and hrs , respectively) taken from Refs. [126, 132]. pq Intermediate ¯ ha i j ¯ hi ¯ hb a ¯ hbc ai ¯ ka hij ¯ hcd ab ¯ hkl ij ¯ jb hia ¯ hic ab ¯ jk hia ′b Ia b Ia ′jb Iia Expressiona ae fia + vim tm e ef mj j je ¯ j fi + vim tm + 1 vmi tef + he te e 2 i b − hb tm Ia ¯ m a bc bc vai − vmi tm a ka ea vij + vij tk e cd cd cd ¯ vab + 1 vmn tmn − hcd tm + vbm tm am b a 2 ab ef ke kl ¯ vij + 1 vij tkl − hle tk + vij tl e ij e 2 ef ′jb jm ¯ jb I − v eb tea − h tm im ia im a ic + v ec ti − hic tm + I ′ic tm − hc tim + ¯ ¯ vab ma b m ab ab e mb a ce tim − v ce tim + 1 hic tnm ¯ ¯ hbm ae am be 2 nm ab jk ke j ¯ jk ¯ je via + hmi tm − via te + A jk him tkm a ae ¯ e tjk + I ′je tk − 1 v ef tjk +hi ea 2 ai ef ia e b + v be tm fa am e ′b − 1 v eb tmn Ia 2 mn ea jb eb j via + via te a Summation over repeated upper and lower indices is assumed. f q = p|f |q and v rs = p pq pq|v|rs − pq|v|sr are the one- and two-body matrix elements of the Hamiltonian in the ij normal-ordered form (one- and two-electron integrals), and the ti and tab are the singly a and doubly excited cluster amplitudes defining the ground-state CCSD wave function of the N -electron reference system. The antisymmetrizer Apq = 1 − (pq) operator is also used, where (pq) is the transposition of indices p and q. 177 a β ¯ (CCSD) Φi |[HN,open (Rµ,1 + Rµ,2 )]C |Φ β e i i m e i m β ¯ β β ¯β ¯ β β β = haβ reβ − hmβ raβ + hmβ raβ eβ i m i e m i e β α β ¯ ¯β β ¯β α m +heαα raβ eα + haβ mβ reβ + haβ mα reαα m iβ fβ mβ nβ eβ fβ iβ nβ 1¯ ¯ − 1 hmβ nβ ra f + 2 haβ nβ re f 2 β β β β i f e f m n i n β α ¯ β α β α ¯β α −hmβ nα ra f + haβ nα re f , β α β α α α ¯ (CCSD) Φaα jbα |[HN,open (Rµ,1 + Rµ,2 )]C |Φ i (5.4) m j i j ¯ ¯ = −Aiα jα hiαα ra α α + Aaα bα heα reα bα m aα b α α α α mα ¯e j i ¯ iα j −Aaα bα hm α raα + Aiα jα haα bα reα α α α α bα e f i j i j m 1¯ α 1¯ + 2 haα bα reα fα + 2 hmαα α ra α nα n α α α α α bα m j ¯ iα e −Aiα jα Aab hmααα re α α a b α α i e mβ jα ¯α β +Aiα jα Aab haα mβ re b β α mβ i j e m ¯ iα e j ¯αα β −hm α α reαα + ha b m reβ α α β α aα bα m ¯ iα j f − 1 Aaα bα hm α α ra α nα 2 α bα nα α fα iα jα fβ mα nβ ¯ −Aaα bα hm b n ra f α α β α β e j f α α α α α e j f i n ¯ αα β α β −Aiα jα ha b n re f , α α β α β i 1 ¯ + 2 Aiα jα haα bα nα reα nα f 178 (5.5) a b β β ¯ (CCSD) Φi j |[HN,open (Rµ,1 + Rµ,2 )]C |Φ β β iβ mβ jβ eβ iβ jβ mβ eβ jβ iβ ¯ ¯ = −Aiβ jβ hmβ ra b + Aa b haβ re b β β β β β β iβ jβ ¯ ¯ −Aa b hm b raβ + Aiβ jβ ha b reβ β β β β β β eβ fβ iβ jβ iβ jβ mβ nβ 1¯ ¯ + 1 ha b re f + 2 hmβ nβ ra b 2 β β β β β β iβ eβ mβ jβ iβ eα mα jβ ¯ −Aa b Aiβ jβ hmβ aβ re b β β β β ¯ +Aa b Aiβ jβ haβ mα re b α β β β i e j m i j e β ¯β β α m ¯β β β −hm a b reβ + ha b m reαα β β α β β β iβ jβ fβ mβ nβ 1 ¯ − 2 Aa b hm b n ra f β β β β β β β iβ jβ fα mβ nα ¯ −Aa b hm b n ra f β β β β α β α eβ jβ fβ iβ nβ 1 ¯ + 2 Aiβ jβ ha b n re f β β β β β eβ jβ fα iβ nα ¯ −Aiβ jβ ha b n re f , β β α β α 179 (5.6) and a b α β ¯ (CCSD) Φiα j |[HN,open (Rµ,1 + Rµ,2 )]C |Φ β j m j i j i n α β α β ¯ β α β ¯ ¯ = −hiαα ra b − hnβ ra b + heα re b aα α m α β α β β i j i j i j f n ¯ α β mα ¯ α β β ¯ β αβ +hb ra f − hm b raα − haα nβ rb α β α β β β e j i f j ¯ αβ i ¯α β β +ha b reα + ha b rf α α β α β β eα fβ iα jβ iα jβ mα nβ e j i n i f m j i f n j e j ¯ ¯ +ha b re f + hmα nβ ra b α β α β α β α β ¯ αβ α β ¯α β −haα nβ re b − hm b ra f α β α β α β ¯α β β β ¯ αβ m +haα nβ rf b + hm b reαα iα aα α β β β i e j i j f n ¯α αβ m ¯αβ β β −hm a b reαα − ha n b rf α β β α β β β m j f j i n α β ¯ iα e ¯ β β α β −hmααα re b − hb n ra f a α β β β α β i j f i j e n ¯αβ β β ¯αβ α m +ha b m reαα + ha b n rf α β β β α β α i j f i j f m n iα jβ eβ mβ nβ α β 1¯ α β α m ¯αβ β − 2 hm b n ra α nα − hm b n ra f α β β α β α β α α fα iα jβ eα mα nβ 1¯ ¯ −haα nβ mα re b − 2 haα nβ mβ re b α β β β e j f e j f i n 1¯ α β α i ¯ αβ β α β + 2 ha b n reα nα + ha b n re f α β α α fα α β β α β iα fβ eα jβ mα ¯ +ha b m rf e α β α β α iα fβ eβ jβ mβ 1¯ + 2 ha b m rf e , α β β β β 180 (5.7) where the antisymmetrizer Apq = 1 − (pq) is used, with (pq) designating the transposition of indices p and q, and i j f eα fα iα jβ ¯αβ α hm b n = vmα nα te b , α β α β α i j f eα fβ iα jβ ¯αβ β hm b n = vmα nβ te b , α β β α β i j e e f i j α β α β ¯αβ α haα nβ mα = vmα nβ ta f , α i j e ¯αβ β haα nβ mβ e j f ¯ αβ α ha b n α β α eα jβ fβ ¯ ha b n α β β β eβ fβ iα jβ = vmβ nβ ta f , α β eα fα mα jβ = vmα nα ta b , α β eα fβ mα jβ ¯α β α ha b m i f e eα fβ iα nβ i f e eβ fβ iα nβ = vmα nβ ta b , α β = vmα nβ ta b , α β α α β and ¯α β β ha b m = vmβ nβ ta b , α β β α β (5.8) where the µ subscript was dropped from the r amplitudes for clarity. By substituting the three-body components of the similarity transformed Hamiltonian of CCSD given in Eq. (5.8) and factorizing the resulting equations, the open-shell EOMCCSD equations projected on doubly excited determinants, Eqs. (5.5), (5.6), and (5.7), may be rewritten in the following 181 way: aα α ¯ (CCSD) Φiα jbα |[HN,open (Rµ,1 + Rµ,2 )]C |Φ m j i j 1¯ ¯ = Aiα jα Aaα bα (− 1 hiαα ra α α + 2 heα reα bα aα α α 2 m α bα mα 1 ¯e j i ¯ iα j − 1 hm α raα + 2 haα bα reα 2 α bα α α α e f i j i j m 1¯ α ¯ + 1 haα bα reα fα + 8 hmαα α ra α nα n 8 α α α α α bα i e mβ jα m j ¯α β ¯ e −hiαααα re α α + haα mβ re b m a b β α α α i j m j 1 1 − 2 χmα ta α α + 2 χeα teα bα ) aα α α iα α bα a b β β ¯ (CCSD) Φi j |[HN,open (Rµ,1 + Rµ,2 )]C |Φ β β iβ mβ jβ (5.9) eβ iβ jβ 1¯ ¯ = Aa b Aiβ jβ (− 2 hmβ ra b + 1 haβ re b 2 β β β β β β iβ jβ mβ eβ jβ iβ 1¯ 1¯ − 2 hm b raβ + 2 ha b reβ β β β β eβ fβ iβ jβ iβ jβ mβ nβ 1¯ ¯ + 8 ha b re f + 1 hmβ nβ ra b 8 β β β β β β i e m j iβ mβ jβ i e m j α β β β ¯β β ¯β α −hmβ aβ re b + haβ mα re b α β β β eβ iβ jβ 1 1 − 2 χmβ ta b + 2 χaβ te b ) β β β β 182 (5.10) and a b α β ¯ (CCSD) Φiα j |[HN,open (Rµ,1 + Rµ,2 )]C |Φ β m j j i j i n α β α β ¯ β α β ¯ ¯ = −hiαα ra b − hnβ ra b + heα re b aα α m α β α β β i j i j i j f n ¯ α β mα ¯ α β β ¯ β αβ +hb ra f − hm b raα − haα nβ rb α β α β β β e j i f j ¯ αβ i ¯α β β +ha b reα + ha b rf α α β α β β eα fβ iα jβ iα jβ mα nβ e j i n i f m j i f n j e j ¯ ¯ +ha b re f + hmα nβ ra b α β α β α β α β ¯ αβ α β ¯α β −haα nβ re b − hm b ra f α β α β α β ¯α β β β ¯ αβ m +haα nβ rf b + hm b reαα iα aα α β β β m j f j i n α β ¯ iα e ¯ β β α β −hmααα re b − hb n ra f a α β β β α β mα jβ iα jβ −χiαα ta b + χeα te b m aα α α β β jβ iα nβ fβ iα jβ −χnβ ta b + χb ta f α β β α β (5.11) where iα fβ i f nβ e f α n α n χiαα = −vnα mα rf α + vmα nβ rf − vnα mα rf α tiα m α α α eα α β eα fβ nβ eα fβ e f iα nβ α α i +vmα nβ rf tiα + 1 vmα nα reα nα + vmα nβ re f , eα 2 α β α fα β eα fβ nβ f e (5.12) e f n α α n χeα = −vaα nα rf α + vaα nβ rf + vnα mα rf α tmα aα α α α α aα β eα fβ nβ e f eα fβ mα nβ 1 α α m −vmα nβ rf tmα − 2 vmα nα ra α nα − vmα nβ ra f , a α β α fα β α 183 (5.13) iβ iβ fβ iβ fα nβ eβ fβ nβ iβ n χmβ = −vnβ mβ rf + vmβ nα rf α − vnβ mβ rf teβ α β β eβ fα iβ eβ fβ iβ nβ eβ fα iβ nα n 1 +vmβ nα rf α teβ + 2 vmβ nβ re f + vmβ nα re f , α β β β α (5.14) and eβ fβ eβ nβ eβ fα eβ fβ nβ mβ n χaβ = −vaβ nβ rf + vaβ nα rf α + vnβ mβ rf taβ α β β eβ fα mβ eβ fβ mβ nβ eβ fα mβ nα n −vmβ nα rf α taβ − 1 vmβ nβ ra f − vmβ nα ra f , 2 α β α β β (5.15) which is the final form of the open-shell EOMCCSD equations used in the efficient vectorized GAMESS code. Once the singly and doubly excited amplitudes defining the EOMCCSD (CCSD) ij i excitation operator, rµ,a and rµ,ab , respectively, and the vertical excitation energy ωµ have been determined by solving Eqs. (3.24) (3.25), rµ,0 is calculated a posteriori from the following expression: (CCSD) ¯ (CCSD) (CCSD) + R(CCSD) )]C |Φ /ωµ rµ,0 = Φ|[HN,open (Rµ,1 . µ,2 5.1.2 (5.16) Electron-Attached and Ionized Equation-of-Motion CoupledCluster Theories The key difference between the open-shell EOMCC theory and EA- or IP-EOMCC is the sector of the Fock space the similarity-transformed Hamiltonian of CCSD is diagonalized within. As an example, in the IP-EOMCC approaches we diagonalize the similarity-transformed Hamiltonian obtained in calculations for an N -electron reference system in the sector of the Fock space corresponding to (N − 1) electrons. Thus, the left-hand sides of the IP184 EOMCCSD(2h-1p) eigenvalue problem, (CCSD) (2h-1p) ¯ (HN,open Rµ (2h-1p) (2h-1p) Rµ |Φ )C |Φ = ωµ (5.17) bβ are obtained by projecting Eq. (5.17) onto all Φiβ |, Φi j |, and Φi bjα | determinants. In β α β β this way the following equations are obtained: ¯ (CCSD) (2h-1p) Φiβ |(HN,open RN,open )C |Φ i iβ mα e i m ¯β m ¯ ¯ β β β = −hmβ r β + heαα reα + hmβ reβ m iβ fα mβ nα ¯ −hmβ nα rf α bβ ¯ (CCSD) (2h-1p) Φi j |(HN,open RN,open )C |Φ β β − jβ iβ nβ ¯ = −Aiβ jβ hnβ rb β 1 hiβ fβ r mβ nβ , ¯ 2 mβ nβ fβ (5.18) fβ iβ jβ ¯ + hb rf β β i j i m e j ¯β β m ¯ β α β α −hm b r β + Aiβ jβ hm b reα β α β β eβ jβ iβ mβ ¯ −Aiβ jβ hb m reβ β β eβ fβ i j mβ nβ ¯β β + 1 hmβ nβ reβ 2 iβ jβ mβ nβ − 1 vmβ nβ te b rf 2 β β β eβ fα iβ jβ mβ nα −vmβ nα te b rf β β α 185 , (5.19) and ¯ (CCSD) (2h-1p) Φi bjα |(HN,open RN,open )C |Φ α β i j iβ nα ¯f β α + hb α rf α α α i j iβ mβ jα ¯β α m ¯ − hm b r β −hmβ rb α β α j ¯ = −hnα rb α eβ bα iβ mβ ¯ α − hb αm reα α α iβ eα mβ jα ¯ + hmβ nα rb ¯ −hm jα reβ β ¯ −hm b reα β α eβ fβ e j iβ mα iβ jα mβ nα α iβ jα mβ nβ 1 − 2 vmβ nβ te b rf β α β eβ fα iβ jα mβ nα −vmβ nα te b rf β α α . (5.20) As in the case of EOMCCSD, the IP-EOMCCSD(2h-1p) equations are solved using the Hirao-Nakatsuji algorithm [303]. 186 Chapter 6 Summary and Concluding Remarks In this dissertation, we addressed the problem of generating highly accurate potential energy surfaces (PESs) for reactive processes by introducing and demonstrating the performance of electronic structure methodologies that can provide a balanced description of chemical species with varying levels of electronic degeneracy, but are also practical enough to be applied to a wide range of chemical problems, as well as extrapolation techniques which facilitate the generation of PESs corresponding to high-level electronic structure calculations in a much more efficient manner than that offered by conventional and laborious point-wise computations. In particular, we examined the performance of two classes of coupled-cluster (CC) methods which are capable of accounting for the diverse electron correlation effects encountered in the majority of ground- and excited-state PES considerations. The first class of methods consisted of the size-extensive completely renormalized (CR) CC approaches for ground states and their equation-of-motion (EOM) CC extensions for excited states, in which noniterative corrections due to higher-order correlation effects are added to the energies obtained with the standard CC and EOMCC approximations, such as CCSD or EOMCCSD, respectively. 187 We showed that the left-eigenstate CR-CC(2,3) and CR-EOMCC(2,3) methods that belong to this category offer excellent performance for a diverse range of applications, including a benchmark database of barrier heights for thermochemical kinetics, a pair of bimolecular association mechanisms involving the ozone molecule, competing intramolecular reaction mechanisms describing the isomerization of bicyclobutane to butadiene, and the ground- and excited-state PES cuts for the water molecule. When necessary, corrections for quadruple excitations were also included via the CR-CC(2,3)+Q method which usually improved the performance of the CR-CC(2,3) methods from chemical to sub-chemical accuracies for many of the studied systems. A new variant of the CR-EOMCC(2,3) method was also presented and discussed, namely, the δ-CR-EOMCC(2,3) approach that can provide a size-intensive treatment of excitation energies. This method was applied to describe excitation energies and hydrogen-bonding-induced spectral shifts in complexes of 7-Hydroxyquinoline with considerable success, helping to explain problems with time-dependent density functional theory. The second class of methods considered here were the active-space variants of the electron attached (EA) and ionized (IP) EOMCC theories. The EA- and IP-EOMCC approaches were shown to be an excellent alternative to open-shell CC and EOMCC methods and their perturbative extensions for describing open-shell molecular systems, providing spin-adapted results while their active-space variants proved to be extremely efficient, significantly reducing the costs of the high-level parent EA- and IP-EOMCC approximations without sacrificing accuracy. We also developed a general strategy for reducing the cost of generating PESs with correlated electronic structure methods via the concept of correlation energy scaling. In order to demonstrate typical accuracies one may expect when using the two types of PES extrapolation schemes presented here, namely, the single-level and dual-level schemes, a number 188 of benchmark applications were presented, such as the previously mentioned bicyclobutane isomerization and single-bond breaking potential energy curves for the H2 O, HCl, and F2 molecules. The single-level extrapolation schemes were shown to reproduce PESs obtained in laborious high-level point-by-point computations to within fractions of a millihartree in most cases, even when used to extrapolate the PES to the CBS-limit. Meanwhile, the duallevel PES extrapolation schemes were shown to be capable of producing similar accuracies at a tiny fraction of the computational cost of their single-level analogs. The insensitivity of the results to the choice of pivot geometry and improvements in accuracy available when a higher-order base wave function is chosen were also demonstrated. Finally, the development of new open-shell EOMCCSD and IP-EOMCCSD(2h-1p) computer codes for the GAMESS software package, along with the corresponding programmable equations, was discussed. 189 REFERENCES 190 REFERENCES [1] P. A. M. Dirac, P. R. Soc. London-ContA. 123, 714 (1929). [2] M. Born, R. Oppenheimer, Annalen der Physik 389, 457 (1927). [3] H. M. James, A. S. Coolidge, J. Chem. Phys. 1, 825 (1933). [4] D. R. Hartree, Proc. Cam. Phil. Soc. 24, 89 (1928). [5] D. R. Hartree, Proc. Cam. Phil. Soc. 24, 111 (1928). [6] D. R. Hartree, Proc. Cam. Phil. Soc. 24, 426 (1928). [7] V. Fock, Z. Physik. 61, 126 (1930). [8] C. C. J. Roothaan, Rev. Mod. Phys. 23, 69 (1951). [9] A. C. Wahl, J. Chem. Phys. 41, 2600 (1964). [10] S. F. Boys, Proc. Roy. Soc. A201, 125 (1950). [11] P. O. L¨wdin, Phys. Rev. 97, 1474 (1955). o [12] P. O. L¨wdin, Phys. Rev. 97, 1490 (1955). o [13] P. O. L¨wdin, Phys. Rev. 97, 1509 (1955). o [14] J. A. Pople, J. S. Binkley, R. Seeger, Int. J. Quantum Chem. Symp. 10, 1, (1976). 191 [15] C. Møller, M. S. Plesset, Phys. Rev. 46, 618 (1934). [16] K. A. Brueckner, Phys. Rev. 97, 1353 (1955). [17] K. A. Brueckner, Phys. Rev. 100, 36 (1955). [18] J. Goldstone, Proc. Roy. Soc. A239, 267 (1957). [19] J. Hubbard, Proc. Roy. Soc. A240, 539 (1957). [20] J. Hubbard, Proc. Roy. Soc. A243, 336 (1958). [21] J. Hubbard, Proc. Roy. Soc. A244, 199 (1958). [22] N. M. Hugenholtz, Physica 23, 481 (1957). [23] F. Coester, Nucl. Phys. 7, 421 (1958). [24] F. Coester, H. K¨mmel, Nucl. Phys. 17, 477 (1960). u ˇ ıˇ [25] J. C´zek, J. Chem. Phys. 45, 4256 (1966). ˇ ıˇ [26] J. C´zek, Adv. Chem. Phys. 14, 35 (1969). ˇ ıˇ [27] J. C´zek, J. Paldus, Int. J. Quantum Chem. 5, 359 (1971). [28] S. R. Langhoff, E. R. Davidson, Int. J. Quantum Chem. 8, 61 (1974). [29] E. R. Davidson, In The World of Quantum Chemistry, edited by R. Daudel, B. Pullmsn, pages 17–30 (Reidel, Dordrecht, 1974). [30] K. Emrich, Nucl. Phys. A 351, 379 (1981). [31] J. Geertsen, M. Rittby, and R. J. Bartlett, Chem. Phys. Lett. 164, 57 (1989). [32] D. C. Comeau and R. J. Bartlett, Chem. Phys. Lett. 207, 414 (1993). [33] J. F. Stanton and R. J. Bartlett, J. Chem. Phys. 98, 7029 (1993). 192 [34] P. Piecuch and R. J. Bartlett, Adv. Quantum Chem. 34, 295 (1999). [35] H. Nakatsuji and K. Hirao, Chem. Phys. Lett. 47, 569 (1977). [36] H. Nakatsuji and K. Hirao, J. Chem. Phys. 68, 4279 (1978). [37] H. Nakatsuji, Chem. Phys. Lett. 59, 362 (1978). [38] H. Nakatsuji, In Computational Chemistry: Reviews of Current Trends, edited by J. Leszczy´ski, vol. 2, pages 62–124 (World Scientific, Singapore, 1997), and references n therein. [39] H. Nakatsuji, Bull. Chem. Soc. Jpn. 78, 1705 (2005), and references therein. [40] H. Monkhorst, Int. J. Quantum Chem. Symp. 11, 421 (1977). [41] E. Dalgaard and H. Monkhorst, Phys. Rev. A 28, 1217 (1983). [42] M. Takahashi and J. Paldus, J. Chem. Phys. 85, 1486 (1986). [43] H. Koch and P. Jørgensen, J. Chem. Phys. 93, 3333 (1990). [44] H. Koch, H. J. A. Jensen, P. Jørgensen, and T. Helgaker, J. Chem. Phys. 93, 3345 (1990). [45] G. D. Purvis III, R. J. Bartlett, J. Chem. Phys. 76, 1910 (1982). [46] J. M. Cullen and M. C. Zerner, J. Chem. Phys. 77, 4088 (1982). [47] G. E. Scuseria, A. C. Scheiner, T. J. Lee, J. E. Rice, and H. F. Schaefer III, J. Chem. Phys. 86, 2991 (1987). [48] P. Piecuch and J. Paldus, Int. J. Quantum Chem. 36, 429 (1989). [49] J. Noga and R. J. Bartlett, J. Chem. Phys. 86, 7041 (1987); 89, 3401 (1988) [Erratum]. [50] G. E. Scuseria and H. F. Schaefer, III, Chem. Phys. Lett. 152, 382 (1988). [51] N. Oliphant and L. Adamowicz, J. Chem. Phys. 95, 6645 (1991). 193 [52] S. A. Kucharski and R. J. Bartlett, Theor. Chim. Acta 80, 387 (1991). [53] S. A. Kucharski and R. J. Bartlett, J. Chem. Phys. 97, 4282 (1992). [54] P. Piecuch and L. Adamowicz, J. Chem. Phys. 100, 5792 (1994). [55] K. Kowalski and P. Piecuch, J. Chem. Phys. 113, 8490 (2000). [56] K. Kowalski and P. Piecuch, J. Chem. Phys. 115, 643 (2001). [57] K. Kowalski and P. Piecuch, Chem. Phys. Lett. 347, 237 (2001). [58] S. A. Kucharski, M. Wloch, M. Musial, and R. J. Bartlett, J. Chem. Phys. 115, 8263 (2001). [59] S. Hirata, J. Chem. Phys. 121, 51 (2004). [60] J. Paldus, In Methods in computational molecular physics, edited by S. Wilson and G. H. F. Diercksen, NATO ASI series B: physics, vol 293, page 99, (Plenum, New York, 1992). [61] R. J. Bartlett, In Advanced series in physical chemistry. Modern electronic structure theory, part I, edited by D. R. Yarkony, page 1047, (World Scientific, Singapore, 1995). [62] J. Paldus and X. Li, Adv. Chem. Phys. 110, 1 (1999). [63] T. D. Crawford and H. F. Schaefer III, Rev. Comp. Chem. 14, 33 (2000). [64] J. Gauss In Encyclopedia of Computational Chemistry, edited by P. v. R. Schleyer, N. L. Allinger, T. Clark, J. Gasteiger, P. A. Kollman, H. F. Schaefer III, and P. R. Schreiner, vol. 1, page 615 (Wiley, Chichester, 1998). [65] P. Piecuch, K. Kowalski, I. S. O. Pimienta, and M. J. McGuire, Int. Rev. Phys. Chem. 21 527 (2002). [66] M. Musial and R. J. Bartlett, Rev. Mod. Phys. 79, 291 (2007). [67] A. P. Rendell, T. J. Lee, and A. Komornicki, Chem. Phys. Lett. 178, 462 (1991). 194 [68] A. P. Rendell, T. J. Lee, and R. Lindh, Chem. Phys. Lett. 194, 84 (1992). [69] A. P. Rendell, M. F. Guest, and R. A. Kendall, J. Comput. Chem. 14, 1429 (1993). [70] R. Kobayashi and A. P. Rendell, Chem. Phys. Lett. 265, 1 (1997). [71] P. Piecuch and J. I. Landman, Parallel Comp. 26, 913 (2000). [72] S. Hirata, J. Phys. Chem. A 107, 9887 (2003). [73] G. Baumgartner, A. Auer, D. E. Bernholdt, A. Bibireata, V. Choppella, D. Cociorva, X. Gao, R. J. Harrison, S. Hirata, S. Krishnamoorthy, S. Krishnan, C. Lam, Q. Lu, M. Nooijen, R. M. Pitzer, J. Ramanujam, P. Sadayappan, and A. Sibiryakov, Proc. IEEE 93, 276 (2005). [74] P. Piecuch, S. Hirata, K. Kowalski, P.-D. Fan, and T. L. Windus, Int. J. Quantum Chem. 106, 79 (2006). [75] T. Janowski, A. R. Ford, and P. Pulay, J. Chem. Theory Comput. 3, 1368 (2007). [76] R. M. Olson, L. J. Bentz, R. A. Kendall, M. W. Schmidt, and M. S. Gordon, J. Chem. Theory Comput. 3, 1312 (2007). [77] M. E. Harding, T. Metzroth, J. Gauss, and A. A. Auer, J. Chem. Theory Comput. 4, 64 (2008). [78] M. Sch¨tz and H.-J. Werner, J. Chem. Phys. 114, 661 (2001). u [79] M. Sch¨tz, J. Chem. Phys. 114, 113 9986 (2000). u [80] M. Sch¨tz and H.-J. Werner, Chem. Phys. Lett. 318, 370 (2000). u [81] M. Sch¨tz, J. Chem. Phys. 116, 8772 (2002). u [82] J. E. Subotnik and M. Head-Gordon, J. Chem. Phys. 123, 064108 (2005). [83] J. E. Subotnik, A. Sodt, and M. Head-Gordon, J. Chem. Phys. 125 074116 (2006). [84] W. Li, J.R. Gour, P. Piecuch, and S. Li, J. Chem. Phys. 131, 114109 (2009). 195 [85] W. Li and P. Piecuch, J. Phys. Chem. A 114, 6721 (2010). [86] W. Li and P. Piecuch, J. Phys. Chem. A 114, 8644 (2010). [87] H. J. Nakano, Chem. Phys. 99, 7983 (1993). [88] H. J. Nakano, Chem. Phys. Lett. 207, 372 (1993). [89] K. Andersson, P. A. Malmqvist, and B. O. Roos, J. Chem. Phys. 96, 1218 (1992). [90] H.-J. Werner and P. J. Knowles, J. Chem. Phys. 89, 5803 (1988). [91] P. J. Knowles and H.-J. Werner, Chem. Phys. Lett. 145, 514 (1988). [92] I. Lindgren and D. Mukherjee, Phys. Rep. 151, 93 (1987). [93] D. Mukherjee and S. Pal, Adv. Quantum Chem. 20, 291 (1989). [94] B. Jeziorski and H.J. Monkhorst, Phys. Rev. A 24, 1668 (1981). [95] J. Paldus, P. Piecuch, L. Pylypow, and B. Jeziorski, Phys. Rev. A 47, 2738 (1993). [96] P. Piecuch, R. Tobola, and J. Paldus, Chem. Phys. Lett. 210, 243 (1993). [97] P. Piecuch and J. Paldus, Phys. Rev. A 49, 3479 (1994). [98] P. Piecuch and K. Kowalski, Int. J. Mol. Sci. 3, 676 (2002). [99] X. Li and J. Paldus, Int. J. Quantum Chem. 110, 2734 (2010). [100] X. Li and J. Paldus, J. Phys. Chem. A 114, 8591 (2010). [101] P. Piecuch, Mol. Phys. 108, 2987 (2010). [102] Y.S. Lee and R.J. Bartlett, J. Chem. Phys, 80, 4371 (1984). [103] Y.S. Lee, S.A. Kucharski, and R.J. Bartlett, J. Chem. Phys., 81, 5906 (1984); 82, 5761 (1985) [Erratum]. 196 [104] P. Piecuch and J. Paldus, Theor. Chim. Acta., 78, 65 (1990). [105] S.A. Kucharski and R.J. Bartlett, Chem. Phys. Lett., 158, 550 (1989). [106] M. Urban, J. Noga, S. J. Cole, and R. J. Bartlett, J. Chem. Phys. 83, 4041 (1985). [107] K. Raghavachari, G. W. Trucks, J. A. Pople, and M. Head-Gordon, Chem. Phys. Lett. 157, 479 (1989). [108] S. A. Kucharski and R. J. Bartlett, J. Chem. Phys., 108, 9221 (1998). [109] J. D. Watts and R. J. Bartlett, Chem. Phys. Lett. 233, 81 (1995). [110] J. D. Watts and R. J. Bartlett, Chem. Phys. Lett. 258, 581 (1996). [111] H. Koch, O. Christiansen, P. Jørgensen, and J. Olsen, Chem. Phys. Lett. 244, 75 (1995). [112] O. Christiansen, H. Koch, and P. Jørgensen, J. Chem. Phys. 103, 7429 (1995). [113] O. Christiansen, H. Koch, and P. Jørgensen, J. Chem. Phys. 105, 1451 (1996). [114] O. Christiansen, H. Koch, P. Jørgensen, and J. Olsen, Chem. Phys. Lett. 256, 185 (1996). [115] K. Kowalski and P. Piecuch, J. Chem. Phys. 120, 1715 (2004). [116] X. Li and J. Paldus, J. Chem. Phys. 107, 6257 (1997). [117] X. Li and J. Paldus, J. Chem. Phys. 108, 637 (1998). [118] X. Li and J. Paldus, J. Chem. Phys. 113, 9966 (2000). [119] X. Li and J. Paldus, Mol. Phys. 98, 1185 (2000). [120] X. Li and J. Paldus, Int. J. Quant. Chem. 80, 743 (2000). [121] X. Li and J. Paldus, J. Chem. Phys. 128, 144118 (2008). 197 [122] X. Li and J. Paldus, J. Chem. Phys. 128, 144119 (2009). [123] P. Piecuch, K. Kowalski, I. S. O. Pimienta, P. D. Fan, M. Lodriguito, M. J. McGuire, S. A. Kucharski, T. Kus, and M. Musial, Theor. Chem. Acc. 112, 349 (2004). [124] K. Kowalski and P. Piecuch, J. Chem. Phys. 113, 18 (2000). [125] K. Kowalski and P. Piecuch, J. Chem. Phys. 113, 5644 (2000). [126] M. Wloch, J. R. Gour, K. Kowalski, and P. Piecuch J. Chem. Phys. 122, 214107 (2005). [127] P. Piecuch and M. Wloch, J. Chem. Phys. 123, 224105 (2005). [128] P. Piecuch, M. Wloch, J. R. Gour, and A. Kinal, Chem. Phys. Lett. 418, 463 (2005). [129] M. Wloch, M. D. Lodriguito, P. Piecuch, and J. R. Gour, Mol. Phys. 104, 2149 (2006). [130] M. Wloch, J. R. Gour, and P. Piecuch, J. Phys. Chem. A 111, 11359 (2007). [131] P. Piecuch, J. R. Gour, and M. Wloch, Int. J. Quantum Chem. 108, 2128 (2008). [132] P. Piecuch, J. R. Gour, and M. Wloch, Int. J. Quantum Chem. 109, 3268 (2009). [133] N. Oliphant and L. Adamowicz, J. Chem. Phys. 94, 1229 (1991). [134] N. Oliphant and L. Adamowicz, J. Chem. Phys. 96, 3739 (1992). [135] N. Oliphant and L. Adamowicz, Int. Rev. Phys. Chem. 12, 339 (1993). [136] P. Piecuch, N. Oliphant, and L. Adamowicz, J. Chem. Phys. 99, 1875 (1993). [137] P. Piecuch and L. Adamowicz, Chem. Phys. Lett. 221, 121 (1993). [138] P. Piecuch and L. Adamowicz, J. Chem. Phys. 102, 898 (1995). [139] K. B. Ghose, P. Piecuch, and L. Adamowicz, J. Chem. Phys. 103, 9331 (1995). [140] K. B. Ghose and L. Adamowicz, J. Chem. Phys. 103, 9324 (1995). 198 [141] V. Alexandrov, P. Piecuch, and L. Adamowicz, J. Chem. Phys. 102, 3301 (1995). [142] K. B. Ghose, P. Piecuch, S. Pal, and L. Adamowicz, J. Chem. Phys. 104, 6582 (1996). [143] L. Adamowicz, P. Piecuch, and K. B. Ghose, Mol. Phys. 94, 225 (1998). [144] P. Piecuch, S. A. Kucharski, and R. J. Bartlett, J. Chem. Phys. 110, 6103 (1999). ˇ [145] P. Piecuch, S. A. Kucharski, and V. Spirko, J. Chem. Phys. 111, 6679 (1999). [146] K. Kowalski and P. Piecuch, Chem. Phys. Lett. 344, 165 (2001). [147] J. R. Gour, P. Piecuch, and M. Wloch, J. Chem. Phys. 123, 134113 (2005). [148] J. R. Gour, P. Piecuch, and M. Wloch, Int. J. Quantum Chem. 106, 2854 (2006). [149] J. R. Gour and P. Piecuch, J. Chem. Phys. 125, 234107 (2006). [150] P. Piecuch, R. Tobola, and J. Paldus, Phys. Rev. A 54, 1210 (1996). ˇ ıˇ [151] J. Paldus, J. C´zek, and M. Takahashi, Phys. Rev. A 30, 2193 (1884). [152] J. Paldus and J. Planelles, Theor. Chim. Acta 89, 13 (1994). [153] J. Paldus , and X. Li, Theor. Chim. Acta 89, 33 (1994). [154] J. Paldus , and X. Li, Theor. Chim. Acta 89, 59 (1994). [155] L. Stolarczyk, Chem. Phys. Lett. 217, 1 (1994). [156] G. Peris, J. Planelles, F. Rajadall, and J. Paldus, J. Chem. Phys. 62, 137 (1997). [157] X. Li, G. Peris, J. Planelles, F. Rajadall, and J. Paldus, J. Chem. Phys. 107, 90 (1997). [158] X. Li and J. Paldus, J. Chem. Phys. 119, 5334 (2003) [159] J. Paldus and X. Li, J. Chem. Phys. 118, 6769 (2003). 199 [160] X. Li and J. Paldus, J. Chem. Phys. 124, 034112 (2006). [161] X. Li and J. Paldus, J. Chem. Phys. 125, 064107 (2006). [162] X. Li and J. Paldus, Chem. Phys. Lett. 431, 179 (2006). [163] J. Paldus and X. Li, Collect. Czech. Chem. Commun. 72, 100 (2007). [164] X. Li and J. Paldus, J. Chem. Phys. 126, 224304 (2007). [165] X. Li and J. Paldus, J. Phys. Chem. A 111, 11189 (2007). [166] X. Li, J.R. Gour, J. Paldus, and P. Piecuch, Chem. Phys. Lett. 461, 321 (2008). [167] X. Li and J. Paldus, J. Theor. Comp. Chem. 7, 805 (2008). [168] X. Li and J. Paldus, J. Chem. Phys. 131, 114103 (2009). [169] X. Li and J. Paldus, Int. J. Quantum Chem. 109, 3305 (2009). [170] X. Li and J. Paldus, J. Chem. Phys. 132, 114103 (2010). [171] X. Li and J. Paldus, J. Chem. Phys. 134, 074301 (2011). [172] V. Spirko, X. Li, and J. Paldus, Collect. Czech. Chem. Commun. 76, 327 (2011). [173] X. Li and J. Paldus, J. Chem. Phys. 119, 5320 (2003). [174] X. Li and J. Paldus, J. Chem. Phys. 119, 5346 (2003). [175] X. Li and J. Paldus, J. Chem. Phys. 120, 5890 (2004). [176] X. Li and J. Paldus, Mol. Phys. 104, 661 (2006). [177] X. Li, Collect. Czech. Chem. Commun. 70, 755 (2005). [178] X. Li and J. Paldus, J. Chem. Phys. 124, 034112 (2006). 200 [179] X. Li and J. Paldus, J. Phys. Chem. A 114, 8591 (2010). [180] X. Li and J. Paldus, Int. J. Quantum Chem. 110, 2734 (2010). [181] X. Li and J. Paldus, J. Chem. Phys. 133, 024102 (2010). [182] K. Kowalski and P. Piecuch, J. Chem. Phys., 115, 2966 (2001). [183] K. Kowalski and P. Piecuch, J. Chem. Phys., 116, 7411 (2002). [184] A. J. C. Varandas and P. Piecuch, Chem. Phys. Lett. 430, 448 (2006). [185] Y. Ge, M. S. Gordon, and P. Piecuch, J. Chem. Phys. 2007, 127, 174106 (2007). [186] P. Piecuch, M. Wloch, and A. J. C. Varandas, Theor. Chem. Acc. 120, 59 (2008). [187] Y. Z. Song, A. Kinal, P. J. S. B. Caridade, A. J. C. Varandas, and P. Piecuch, J. Mol. Struct.: THEOCHEM. 859, 22 (2008). [188] Y. Ge, M. S. Gordon, P. Piecuch, M. Wloch, and J. R. Gour, J. Phys. Chem. A 112, 11873 (2008). [189] A. Kinal and P. Piecuch, J. Phys. Chem. A 111, 734 (2007). [190] C. J. Cramer, M. Wloch, P. Piecuch, C. Puzzarini, and L. Gagliardi, J. Phys. Chem. A 110, 1991 (2006); 111, 4871 (2007) [Erratum]. [191] C. J. Cramer, A. Kinal, M. Wloch, P. Piecuch, and L. Gagliardi, J. Phys. Chem. A 110, 11557 (2006); 111, 4871 (2007) [Erratum]. [192] C. J. Cramer, J. R. Gour, A. Kinal, M. Wloch, P. Piecuch, A. R. M. Shahi, and Gagliardi, L. J. Phys. Chem. A 112, 3754 (2008). [193] J. J. Lutz and P. Piecuch, J. Chem. Phys. 128, 154116 (2008). [194] M. W. Schmidt, K. K. Baldridge, J. A. Boatz, S. T. Elbert, M. S. Gordon, J. H. Jensen, S. Koseki, N. Matsunaga, K. A. Nguyen, S. J. Su, T. L. Windus, M. Dupuis, and J. A. Montgomery, J. Comput. Chem. 14, 1347 (1993). 201 [195] M. S. Gordon and M. W. Schmidt, In: Theory and Applications of Computational Chemistry: The First Forty Years, edited by C. E. Dykstra, G. Frenking, K. S. Kim, and G. E. Scuseria, pages 1167–1190, (Elsevier, Amsterdam, 2005). [196] L. Meissner and R. J. Bartlett J. Chem. Phys. 102, 7490 (1995). [197] E. A. Salter, G. W. Trucks, and R. J. Bartlett, J. Chem. Phys. 90, 1752 (1989). [198] P. S. Epstein, Phys. Rev. 28, 695 (1934). [199] R. D. Nesbet, Proc. Roy. Soc. (London) A230, 312 (1955). [200] T. Shiozaki, K. Hirao, S. Hirata, J. Chem. Phys. 126, 244106 (2007). [201] G. Fradelos, J. J. Lutz, T. A. Wesolowski, P. Piecuch, M. Wloch, J. Chem. Theory Comput. 7, 1647 (2011). [202] B. J. Lynch and D. G. Truhlar, J. Phys. Chem. A 107, 3898 (2003). [203] Y. Zhao, N. Gonzalez-Garcia, and D. G. Truhlar, J. Phys. Chem. A 109, 2012 (2005). [204] J. Zheng, Y. Zhao, and D. G. Truhlar, J. Chem. Theory Comput. 3, 569 (2007). [205] J. M. L. Martin and G. de Oliverira, J. Chem. Phys. 111, 1843 (1999). [206] J. Zheng, J. R. Gour, J. J. Lutz, M. Wloch, P. Piecuch, D. G. Truhlar, J. Chem. Phys. 128, 044108 (2008). [207] B. J. Lynch, Y. Zhao, and D. G. Truhlar, J. Phys. Chem. A 107, 1384 (2003). [208] T. H. Dunning, J. Chem. Phys. 90, 1007 (1989). [209] D. E. Woon and T. H. Dunning, J. Chem. Phys. 98, 1358 (1993). [210] R. A. Kendall, T. H. Dunning, and R. J. Harrison, J. Chem. Phys. 96, 6796 (1992). [211] T. H. Dunning, K. A. Peterson, and A. K. Wilson, J. Chem. Phys. 114, 9244 (2001). [212] D. E. Woon and T. H. Dunning, J. Chem. Phys. 103, 4572 (1995). 202 [213] K. A. Peterson and T. H. Dunning, J. Chem. Phys. 117, 10548 (2002). [214] L. A. Curtiss, K. Raghavachari, C. Redfern, V. Rassolov, and J. A. Pople, J. Chem. Phys. 109, 7764 (1998). [215] P. R. Taylor, In Lecture Notes in Quantum Chemistry, edited by B. O. Roos, page 406, (Springer-Verlag, Berlin, 1992). [216] A. Halkier, T. Helgaker, P. Jorgensen, W. Klopper, and J. Olsen, Chem. Phys. Lett. 302, 437 (1999). [217] T. Helgaker, W. Klopper, H. Koch, and J. Noga, J. Chem. Phys. 106, 9639 (1997). [218] W. D. Laidig and H. F. Schaefer, J. Chem. Phys. 74, 3411 (1981). [219] P. Borowski, K. Andersson, P.A. Malmqvist and B. O. Roos, J. Chem. Phys. 97, 5568 (1992). [220] M. L. Leininger and H. F. Schaefer, J. Chem. Phys. 107, 9059 (1997). [221] X. Li and J. Paldus, J. Chem. Phys. 110, 2884 (1999). [222] S. P. Walch and W. A. Goddard III, J. Am. Chem. Soc. 97, 5319 (1975). [223] S. D. Kahn, W. J. Hehre, and J. A. Pople, J. Am. Chem. Soc. 109, 1871 (1987). [224] P. C. Hiberty and C. J. Leforestier, J. Am. Chem. Soc. 100, 2012 (1978). [225] S. E. Wheeler, D. G. Ess, and K. N. Houk, J. Phys. Chem. A 112, 1798 (2008). [226] W. B. DeMore and C. L. Lin, J. Org. Chem. 38, 985 (1973). [227] W. B. DeMore, Int. J. Chem. Kinet. 1, 209 (1969). [228] D. J. Miller, T. E. Nemo, and L. A. Hull, J. Org. Chem. 40, 2675 (1975). [229] S. Jackson and L. A. Hull, J. Org. Chem. 41, 3340 (1976). [230] K. Griesbaum and Y. Dong, J. Prakt. Chem. 339, 575 (1997). 203 [231] O. Horie and G. K. Moortgat, Acc. Chem. Res. 31, 381 (1998). [232] J. M. Anglada, R. Crehuet, and J. M. Bofill, Chem. Eur. J. 5, 1809 (1999). [233] L. C. Li, P. Deng, M. H. Xu, and N. B. Wong, Int. J. Quantum Chem. 98, 309 (2004). [234] W. T. Chan and I. P. Hamilton, J. Chem. Phys. 118, 1688 (2003). [235] I. Ljubic and A. Sabljic J. Phys. Chem. A 106, 4745 (2002). [236] J. Z. Gillies, C. W. Gillies, F. J. Lovas, K. Matsumura, R. D. Suenram, E. Kraka, and D. Cremer, J. Am. Chem. Soc. 113, 6408 (1991). [237] M. Olzmann, E. Kraka, D. Cremer, R. Gutbrod, and S. Anderson, J. Phys. Chem. A 101, 9421 (1997). [238] D. Cremer, R. Crehuet, J. Anglada, J. Am. Chem. Soc. 123, 6127 (2001). [239] D. Cremer, E. Kraka, R. Crehuet, J. Anglada, and J. Grafenstein, Chem. Phys. Lett. 347, 268 (2001). [240] J. Z. Gillies, C. W. Gillies, F. J. Lovas, K. Matsumura, R. D. Suenram, E. Kraka, and D. Cremer, J. Am. Chem. Soc. 113, 2412 (1991). [241] W. T. Chan, C. Weng, and J. D. Goddard, J. Phys. Chem. A 111, 4792 (2007). [242] Y. Zhao, O. Tishchenko, J. R. Gour, W. Li, J. J. Lutz, P. Piecuch, and D. G. Truhlar, J. Phys. Chem. A 113, 5786 (2009). [243] Y. Zhao, N. E. Schultz, and D. G. Truhlar, J. Chem. Theory Comput. 2, 364 (2006). [244] S. G. V. Ornum, R. M. Champeau, and R. Pariza, Chem. Rev. 106, 2990 (2006). [245] K. B. Wiberg and J. M. Lavanish, J. Am. Chem. Soc. 88, 5272 (1966). [246] G. L. Closs and P. E. Pfeffer, J. Am. Chem. Soc. 90, 2452 (1968). [247] R. B. Woodward and R. Hoffman, The Conservation of Orbital Symmetry (Academic, New York, 1970). 204 [248] K. A. Nguyen and M. S. Gordon, J. Am. Chem. Soc. 117, 3835 (1995). [249] P. B. Shevlin and M. L. Mckee, J. Am. Chem. Soc. 110, 1666 (1988). [250] M. J. S. Dewar and S. Kirschner, J. Am. Chem. Soc. 97, 2931 (1975). [251] R. Srinivasan, A. Levi, and I. Haller, J. Phys. Chem. 58, 1775 (1965). [252] K. B. Wiberg and R. A. Fenoglio, J. Am. Chem. Soc. 90, 3395 (1968). [253] D. A. Mazziotti, J. Phys. Chem. A 112, 13684 (2008). [254] R. Berner and A. Luchow, J. Phys. Chem. A 114, 13222 (2010). [255] J. Shen and P. Piecuch, Chem. Phys., submitted (2011). [256] W.L. Glab, J. Chem. Phys. 107, 5979 (1997). [257] P. Farmanara, O. Steinkellner, M.T. Wick, M. Wittmann, G. Korn, V. Stert, and W. Radloff, J. Chem. Phys. 111, 6264 (1999). [258] M.S. Child and W.L. Glab, J. Chem. Phys. 112, 3754 (2000). [259] J.H. Fillion, R. van Harrevelt, J. Ruiz, M. Castillejo, A.H. Zanganeh, J.-L. Lemaire, M.C. van Hemert, and F. Rostas, J. Phys. Chem. A 105, 11414 (2001). [260] W.L. Glab, J. Chem. Phys. 117, 9316 (2002). [261] F. Edery and A. Kanaev, Eur. Phys. J. D 23, 257 (2003). [262] B.-M. Cheng, C.-Y. Chung, M. Bahou, Y.-P. Lee, L.C. Lee, R. van Harrevelt, and M.C. van Hemert, J. Chem. Phys. 120, 224 (2004). [263] J.H. Fillion, J. Ruiz, X.-F. Yang, M. Castillejo, F. Rostas, and J.-L. Lemaire, J. Chem. Phys. 120, 6531 (2004). [264] O. Steinkellner, F. Noack, H.-H. Ritze, W. Radloff, and V. Hertel, J. Chem. Phys. 121, 1765 (2004). 205 [265] D.M. Hirst and M.S. Child, Molec. Phys. 77, 463 (1992). [266] M. von Dirke, B. Heumann, K. Kuhl, T. Schroder, and R. Schinke, J. Chem. Phys. 101, 2051 (1994). [267] F. Schneider, F. Di Giacomo, and F.A. Gianturco, J. Chem. Phys. 104, 5153 (1996). [268] R. van Harrelvelt and M.C. van Hemert, J. Chem. Phys. 112, 5777 (2000). [269] R. van Harrelvelt and M.C. van Hemert, J. Chem. Phys. 114, 9453 (2001). [270] V. Engel, V. Staemmler, R.L. Vander Wal et al., J. Phys. Chem. 96, 3201 (1992). [271] X. Li and J. Paldus, J. Chem. Phys. 133, 024102 (2010). [272] C. Tanner, C. Manca, and S. Leutwyler, Science 302, 1736 (2003). [273] D. Bruhwiler, G. Calzaferri, T. Torres, J. H. Ramm, N. Gartmann, L. Q. Dieu, I. Lopez-Duarte, and M. V. Martinez-Diaz, J. Mater. Chem. 19, 8040 (2009). [274] F. E. Hernandez, S. Yu, M. Garcia, and A. D. Campiglia, J. Phys. Chem. B 109 9499 (2005). [275] J. M. Goldberg, S. Batjargal, and E. J. Petersson, J. Am. Chem. Soc. 132, 14719 (2010). [276] M. Thut, C. Tanner, A. Steinlin, and S Leutwyler, J. Phys. Chem. A 112, 5566 (2008). [277] T. A. Wesolowski and A. Warshel, J. Phys. Chem. 97, 8050 (1993). [278] T. A. Wesolowski, In Computational Chemistry: Reviews of Current Trends; edited by J. Leszczy´ski, Vol. 10, pages 1–82, (World Scientific, Singapore, 2006). n [279] T. A. Wesolowski, Phys. Rev. A 77, 012504 (2008). [280] K. Pernal and T. A. Wesolowski, Int. J. Quantum Chem. 109, 2520 (2009). [281] T. Wesolowski, J. Am. Chem. Soc. 126, 11444 (2004). 206 [282] W. J. Hehre, R. Ditchfield, and J. A. Pople, 56, 2257 (1972). [283] P. C. Hariharan and J. A. Pople, Theor. Chim. Acta 28, 213 (1973). [284] T. Clark, J. Chandrasekhar, G. W. Spitznagel, and P. v. R. Schleyer, J. Comput. Chem. 4, 294 (1983). [285] R. Krishnan, J. S. Binkley, R. Seeger, and J. A. Pople, J. Chem. Phys. 72, 650 (1980). [286] A. J. Sadlej Coll. Czech. Chem. Commun. 53, 1995 (1988). [287] S. Coussan, Y. Ferro, A. Trivella, P. Roubin, R. Wieczorek, C. Manca, P. Piecuch, K. Kowalski, M Wloch, S. A. Kucharski, and M. Musial, J. Phys. Chem. A 110, 3920 (2006). [288] K. Kowalski, S. Krishnamoorthy, O. Villa, J. R. Hammond, N. Govind, J. Chem. Phys. 132, 154103 (2010). [289] M. Ehara, J. R. Gour, and P. Piecuch, Mol. Phys. 107, 871 (2009). [290] J. A. Hansen, P. Piecuch, J. J. Lutz, and J. R. Gour, Phys. Scr. 84, 028110 (2011). [291] T. H. Dunning Jr., J. Chem. Phys. 53, 2823 (1970). [292] T. Nakajima and H. Nakatsuji, Chem. Phys. Lett. 79, 280 (1997). [293] M. Ishida, K. Toyota, M. Ehara, and H. Nakatsuji, 347 493 (2001). [294] M. Ishida, K. Toyota, M. Ehara, H. Nakatsuji, and M. J. Frisch, J. Chem. Phys., 120 2593 (2004). [295] G. Herzberg, Molecular Spectra and Molecular Structure III: Electronic Spectra and Electronic Structure of Polyatomic Molecules (Van Nostrand, London, 1967). [296] A. J. Merer and D. N. Travis, Can. J. Phys. 44, 353 (1966). [297] A. J. Merer and D. N. Travis, Can. J. Phys. 43, 1795 (1966). [298] Y. Ohtsuka, P. Piecuch, J. R. Gour, M. Ehara, H. Nakatsuji, J. Chem. Phys. 126, 164111 (2007). 207 [299] P. D. Rajendra and P. Chandra, 114, 1589 (2001). [300] J. J. Lutz and P. Piecuch, In Nuclei and Mesoscopic Physics, Workshop on Nuclei and Mesoscopic Physics, WNMP 2007, edited by P. Danielewicz, P. Piecuch, and V. Zelevinsky AIP Conferenc Proceedings, Vol. 995, pages 62-71, (American Institute of Physics, Melville, NY, 2008). [301] A. J. C. Varandas, Chem. Phys. Lett. 443, 398 (2007). [302] A. J. C. Varandas, J. Chem. Phys. 127, 114316 (2007). [303] K. Hirao and H. Nakatsuji, J. Comput. Phys. 45, 246 (1982). [304] E. R. Davidson, J. Comput. Phys. 17, 87 (1975). 208