UFDC Home  myUFDC Home  Help 



Full Text  
ON ELECTRONIC REPRESENTATIONS IN MOLECULAR REACTION DYNAMICS By BENJAMIN J. KILLIAN A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY UNIVERSITY OF FLORIDA 2005 Copyright 2005 by Benjamin J. Killian This work is dedicated to Breana, Antonia, and Blake Jones. May the love and wonder of science always be near. ACKNOWLEDGMENTS I am very grateful for the environment of the Quantum Theory Project and its members. The QTP is a wonderful place of education where continual discussion and strong science foster intellect and curiosity. I wish to express tremendous gratitude to Dr. Remigio CabreraTrujillo for his hours of discussion, correction, and teaching throughout my graduate education, as well as the help and friendship offered by the other END Research Group members with whom I collaborated, Dr. Mauricio CoutinhoNeto, Dr. Anatol Blass, Dr. David Masiello, Dr. Denis Jacquemin, Dr. Svetlana Malinovskya, and Mr. Virg Fermo. I offer special acknowledgment to Prof. Yngve Ohrn and Dr. Erik Deumens for their patient and masterful daytoday teaching and mentoring in the field of quantum dynamics. My grateful thanks are given to my mother and father, Jane and James Killian, for always encouraging me to strive in every endeavor and for support throughout the long education process. Without a parent's motivation and loving support, very little can be truly accomplished in a person's life. Finally, and most importantly, I wish to express my utmost love and gratitude to my wife, Donna, for her iii.1'i i.,iiiir. compassion, and emotional support throughout our life together. Without her devoted love and encouragement, (melPh.D.) 0. TABLE OF CONTENTS page ACKNOWLEDGMENTS ................... ...... iv LIST OF TABLES ...................... ......... vii LIST OF FIGURES ................... ......... viii ABSTRACT ....................... ........... xi CHAPTERS 1 INTRODUCTION ........................... 1 1.1 Pricis of Quantum Dynamics ...... ......... 3 1.2 Solving the Schridinger Equation ...... ...... ... 4 1.3 TimeDependent HartreeFock (TDHF) ............ 6 1.4 ElectronNuclear Dynamics (END) ....... ........ 8 1.4.1 The END Equations of Motion ............. 9 1.4.2 Minimal ElectronNuclear Dynamics .......... 10 2 THEORY OF COLLISIONS .......... ............. 12 2.1 Scattering Theory ........... ..... ....... 12 2.1.1 Deflection Functions and Scattering Angles . 13 2.1.2 Cross Sections ................. . .. 15 2.1.3 Identical Particles .............. .. .. 17 2.1.4 Reference Frame Transformations . .... 18 2.2 Quantum Mechanical Treatment of Scattering Phenomena 22 2.2.1 The Integral Equation and its Relation to the Scat tering Amplitude .................. .. 22 2.2.2 The Born Series ............ .. .. .. 24 2.3 SemiClassical Treatment of Scattering Phenomena: The Schiff Approximation .. . . ... . 26 2.3.1 Schiff Scattering Amplitude for Large Angles . 27 2.3.2 Schiff Scattering Amplitude for Small Angles . 34 3 BASIS SETS FOR DYNAMICAL HARTREEFOCK CALCULA TIONS .................. ....... ..... 41 3.1 The HartreeFock Approximation . . ..... 42 3.1.1 Partitioning of the Molecular Wave Function . 42 3.1.2 The HartreeFock Wave Function . . ... 45 3.1.3 Solving the HF Equations: Basis Set Expansions 49 3.2 General Forms and Properties of Basis Sets . ... 52 3.2.1 SlaterType and GaussianType Orbitals . ... 52 3.2.2 The Structure of Basis Sets . . .... 56 3.3 Method for Constructing Basis Sets Consistent with Dynam ical Calculations . . ..... ..... 63 3.3.1 Basis Set Properties for Dynamical Calculations 64 3.3.2 Physical Justification for the Basis Set Construction Method .......... .... .. ....... 66 3.3.3 Construction of the Basis Set . . 71 3.4 Comparative Results ................... .. 82 3.4.1 Atomic Energetics .............. .. .. 83 3.4.2 Charge Transfer Results . . . 87 3.4.3 Properties of Diatomic and Triatomic Molecules 98 4 VECTOR HARTREEFOCK: A MULTICONFIGURATIONAL REPRESENTATION IN ELECTRONNUCLEAR DYNAMICS 104 4.1 Introduction of the Lagrangian and Verification of the Equa tions of Motion . . .......... ... 104 4.2 Parameterization of the State Vector . . 106 4.3 The VHF Equations of Motion . . 114 4.3.1 Derivation of the Equations of Motion . ... 115 4.3.2 Evaluating the Equations of Motion . ... 120 4.4 Implementation of the Vector HartreeFock Method . 140 5 CONCLUSION.... ....... ........ .......... .. 147 APPENDIX BASIS SET LIBRARY ................. .. .. 149 REFERENCES ................... ............. 168 BIOGRAPHICAL SKETCH .................. ......... .. 173 LIST OF TABLES Table page 3.1 Notation employed in this review. ............. . 42 3.2 Slater Exponents and Coefficients for He, Li, and Be. . ... 79 3.3 Slater Exponents and Coefficients for B, C, and N. . .... 80 3.4 Slater Exponents and Coefficients for O, F, and Ne. . ... 81 3.5 Atomic energies and electronic excitations in Helium . ... 84 3.6 Atomic energies and electronic excitations in Lithium . ... 85 3.7 Molecular properties of the nitrogen molecule. .. . ..... 99 3.8 Molecular properties of the carbon monoxide molecule. . ... 100 3.9 Molecular properties of the water molecule. . . 102 LIST OF FIGURES Figure page 2.1 Diagram of a classical collision system. .... . ... 13 2.2 Three possible deflection angles. A demonstrates a positive deflection angle. B demonstrates a negative deflection angle. C demonstrates a deflection angle with magnitude greater than 27. . ... 14 2.3 Collision of distinguishable particles. ................ .. 17 2.4 Collision of identical particles. ................ ..... 17 2.5 Relation of collision pair and centerofmass to origin. . ... 19 2.6 Relation of the scattering angles between reference frames. . 19 3.1 Comparison of STO and GTO representations of the radial part of the hydrogen Is orbital. Top: A single GTO function fit to a sin gle STO function. Bottom: A linear combination of six GTO func tions fit to a single STO function. ................. .. 55 3.2 Top: Plot of the radial distribution function for the Is (), 2s (  ), 3s (. .), and 4s ( ) orbitals of the hydrogen atom. Bottom: Plot of the radial distribution function for the 2p (), 3p ( ), and 4p (. .) orbitals of the hydrogen atom. . . 69 3.3 Top: Plot of the radial distribution function for the Is (), 2s (  ), and 3s (. .) orbitals of the argon atom. The Isorbital is scaled by a factor of twothirds. Bottom: Plot of the radial distribution function for the 2p () and 3p ( ) orbitals of the argon atom. .70 3.4 Plot of the Is orbital exponent for the atoms through Kr as a func tion of atomic number. The data are from Clementi and Raimondi (+). ...... .. .. .. .. ... .. .. .... ... ..... 73 3.5 Plot of the 2s and 2p orbital exponents for the atoms through Kr as a function of atomic number. Top: The 2s orbital exponents. Bot tom: The 2p orbital exponents. The data are from Clementi and Raimondi (+) and from the present work (o). . . ... 74 3.6 Plot of the 3s and 3p orbital exponents for the atoms through Kr as a function of atomic number. Top: The 3s orbital exponents. Bot tom: The 3p orbital exponents. The data are from Clementi and Raimondi (+) and from the present work (o). . . ... 75 3.7 Plot of the 4s and 4p orbital exponents for the atoms through Kr as a function of atomic number. Top: The 4s orbital exponents. Bot tom: The 4p orbital exponents. The data are from Clementi and Raimondi (+) and from the present work (o). . . ... 76 3.8 Comparison of the probability for nearresonant charge transfer be tween H+ and Li at 10 keV collision energy. The following Li basis sets are used: 631G (red), 631B (blue), BJK01 (purple), BJK02 (black). The BJK01 basis set for H was used for each run. . 89 3.9 Comparison of the probability for nearresonant charge transfer be tween H+ and Li at 1 keV collision energy. The following Li basis sets are used: 631G (red), 631B (blue), BJK01 (purple), BJK02 (black). The BJKO1 basis set for H was used for each run. . 90 3.10 Comparison of the differential cross section for nearresonant charge transfer between H+ and Li at 10 keV collision energy. The follow ing Li basis sets are used: 631G (red), 631B (blue), BJKO1 (pur ple), BJK02 (black). The BJKO1 basis set for H was used for each run .............. ................ .. 91 3.11 Comparison of the total cross section for nearresonant charge trans fer between H+ and Li as a function of energy. The experimental data are: Varghese et al. (+) and Aumayr et al. (x). The theoret ical data are: Allan et al () and Fritsch and Lin ( ). For the END data the following Li basis sets are used: 631G (0), 631B (o), BJK01 (A), BJK02 (0). The BJK01 basis set for H was used for each run. . . . . .. ... . 92 3.12 Comparison of the probability for resonant charge transfer between Li+ and Li at 1 keV collision energy. The following Li basis sets are used: 631G (red), BJK02 (blue). ............... .. 93 3.13 Comparison of the probability for resonant charge transfer between Li+ and Li at 10 keV collision energy. The following Li basis sets are used: 631G (red), BJK02 (blue). ............... .. 94 3.14 Comparison of the differential cross section for resonant charge trans fer between Li+ and Li at 1 keV collision energy. The following Li basis sets are used: 631G (red), BJK02 (blue). . . ... 95 3.15 Comparison of the total cross section for resonant charge transfer be tween Li+ and Li at as a function of collision energy. The experi mental data are form Lorents et al. (x). The theoretical data are from Sakabe and Izawa (). For the END data the following Li basis sets are used: 631G (0), 631B (o), BJK02 (A). . 96 3.16 Comparison of the probability for resonant charge transfer between He+ and He at 5 keV collision energy. The following He basis sets are used: 631G** (red), 631B** (blue), and BJKO1 (purple). 97 3.17 Comparison of the differential cross section for resonant charge trans fer between He+ and He at 5 keV collision energy. The experimen tal data are from Gao et al. (*). For the END data the following He basis sets are used: 631G** (red), 631B** (blue), and BJKO1 (purple). . . . . . .. . 98 4.1 Structural flow chart for minimal END electronic calculations. . 141 4.2 Structural flow chart for Vector HartreeFock END electronic calcula tions. ........ .. .. .. .. ... .. .. .. .. ... 142 Abstract of Dissertation Presented to the Graduate School of the University of Florida in Partial Fulfillment of the Requirements for the Degree of Doctor of Philosophy ON ELECTRONIC REPRESENTATIONS IN MOLECULAR REACTION DYNAMICS By Benjamin J. Killian August 2005 Chair: Nils Y. Ohrn M.,. ir Department: Chemistry For many decades, the field of chemical reaction dynamics has utilized compu tational methods that rely on potential energy surfaces that are constructed using stationarystate calculations. These methods are typically devoid of dynamical cou plings between the electronic and nuclear degrees of freedom, a fact that can result in incorrect descriptions of dynamical processes. Often, nonadiabatic coupling expressions are included in these methodologies. The ElectronNuclear Dynamics (END) formalism, in contrast, circumvents these deficiencies by calculating all intermolecular forces directly at each time step in the dynamics and by explicitly maintaining all electronicnuclear couplings. The purpose of this work is to offer two new frameworks for implementing electronic representations in dynamical calculations. Firstly, a new schema is proposed for developing atomic basis sets that are consistent with dynamical calculations. Traditionally, basis sets have been designed for use in stationarystate calculations of the structures and properties of molecules in their ground states. As a consequence of common construction techniques that utilize energy optimization methods, the unoccupied orbitals bear little resemblance to 1.r, i .,l virtual atomic orbitals. We develop and implement a method for basis set construction that relies upon 1ir, i, .,l properties of atomic orbitals and that results in meaningful virtual orbitals. These basis sets are shown to provide a significant improvement in the accuracy of calculated dynamical properties such as charge transfer probabilities. Secondly, the theoretical framework of END is expanded to incorporate a multiconfigurational representation for electrons. This formalism, named Vector HartreeFock, is based in the theory of vector coherent states and utilizes a complete active space electronic representation. The Vector HartreeFock method is fully disclosed, with derivation of the equations of motion. The expressions for the equation of motion are derived in full and a plan for implementing the Vector HartreeFock formalism within the current ENDyne computer code is given. CHAPTER 1 INTRODUCTION Chemistry is the study of how atoms and molecules change with time, and fundamentally requires a dynamical description to correctly describe these changes. Experimental ]li, i .,1 chemistry is concerned with measurement of chemical changes and the enumeration of the properties associated with these changes. Quantum chemistry, on the other hand, seeks to describe the mechanism through which chemical changes occur and to calculate properties associated with these changes as a function of fundamental properties of the molecule itself. Virtually all of chemistry is described in the vocabulary of dynamics. Chemists speak of electron transfers, transition states, equilibria. The most fundamental measurements made in a 1ir, i .,1 chemistry laboratory involve transference of heat over time, the number of vibrations a molecule undergoes in a unit of time, the inversion of electron populations as a function of time. Chemistry is dynamic. Despite all of this, nondynamical approaches are still most often employed when performing chemical calculations. These approximations are often sufficient for obtaining average properties, but fail to offer a deeper understanding of chemical processes that comes through dynamical methods. Traditionally, chemists have been interested in answering two general classes of questions regarding chemical reactions. The first, "To what extent will a reaction occur?" is purely thermodynamic in nature. Provided the difference in the free energy between the product species and the reactant species is negative, the reaction has the potential to occur spontaneously to some extent. Calculations of this type are wellsuited to timeindependent quantum mechanics; stationary state energy values can be calculated for various ground and excited states, which can yield (through the application of statistical thermodynamical principles) the free energies of the products and the reactants. The second class of question that a chemist might ask is, "How fast will the reaction occur?" This type of question is kinetic in nature and relies upon some (usually detailed) knowledge of the dynamics of the reaction. This type of question is, in generally, only successfully answered using timedependent quantum mechanical treatments. To correctly describe the kinetics of a reaction, one must first describe the interactions through and responses to inter and intramolecular forces as a function of the timescale of the reaction. The focus of this dissertation is to investigate methods of improving represen tations of electrons for use in timedependent dynamical calculations. Specifically, this work will investigate two sides of the very same coin. Firstly, a discussion will be made toward basis set expansion for use in construction of the electronic wave function. A new method will be proposed for the construction of dynamically con sistent basis sets. Two particular features of the proposed construction method are the pli, .1l basis that underlies the method and its ease of application. Secondly, discussion will be made toward the construction of a multiconfigurational wave function for use with the ElectronNuclear Dynamics formalism. The full set of equations of motion are derived for the Vector HartreeFock implementation of ElectronNuclear Dynamics in terms of a general expansion of atomic basis func tions. Additionally, discussion is made as to a possible scheme for implementation. This chapter provides a basic review of quantum molecular dynamics. Subse quently, a brief discussion is made of the general structure of the TimeDependent HartreeFock method. Finally the ElectronNuclear Dynamics formalism in its simplest form is introduced and discussed and compared to TimeDependent HartreeFock. 1.1 Pr6cis of Quantum Dynamics The lr, li, .,1 description of any chemical object is predicated on the concepts of a quantum mechanical state and a quantum mechanical configuration [1]. A quantum mechanical configuration is the set of all descriptive variables in the dynamical phase space (e.g., momentum and position) and/or in the electronic Hilbert space (e.g., orbital and spin angular momentum) which describes (within the completeness of the set) a given chemical object at any given instant in time, t = r [1]. The configuration is usually represented as a ket in Dirac notation, written as 1, t = r), where 9 is a representative of the required set of descrip tive variables. If r is specifically defined, then the configuration may be more conveniently expressed as 1). A quantum mechanical state, or wavefunction, can also be expressed in Dirac notation, 1, t), where t indicates the general dependence on time. The state can be defined as any subset of the complete set of configurations that contains a given configuration, 90, t = r), as well as all configurations that result from the dynamical evolution of 1, t = r) and all configurations that evolve into 91, t = r) with the passage of time (assuming no perturbation of the system by outside influences) [1]. Some descriptions must include the time dependence explicitly in the wave function, while others allow for the timedependence to be treated as separable factor, depending on the appearance of time in the Schridinger Equation [2]. For a state to be completely described, one must know the relation between all variables that determine the state at a specific time, usually denoted as time t = 0. Investigation of chemical processes requires the knowledge of how a given state evolves with time how the state changes dynamically. The principle tool in the descriptions of quantum chemical dynamics is the Schridinger Wave Equation, a timedependent differential equation of second order spatially and first order temporally. The Schridinger Equation defines the equations of motion for a given quantum mechanical state; once the initial state is completely described, then the time evolution of the state as defined by the statespecific Schrodinger Equation allows one to correctly predict the future properties of the state based on its initial conditions. The statespecific Schrodinger Equation can be expressed in function notation as N V 21tt. + t) 1 (1.1) i= 1 On the left hand side of Equation (1.1), the first term is the kinetic energy of the state; the summation is over the N particles comprising the state, M, is the mass of the ith particle, and V2 is the Laplacian of the ith particle. The second term on the left is the potential energy of the state. The potential energy of the state is generally a complicated function of several or all of the dynamical variables and is of particular importance to the description of quantum dynamics. One can combine these two terms and define the action of these terms by an operator, 4r, called the Hamiltonian operator. Thus, Equation (1.1) can also be written in operator form, as \0 t} = h t). (1.2) Equation (1.2)is the set of equations of motion required to calculate the dynamics of a quantum mechanical state. 1.2 Solving the Schr6dinger Equation The bulk of quantum dynamical research is involved with devising and refining methods to solve Equation (1.2). A molecular Hamiltonian operator will take the form (in atomic units, where h = c m 1) K N K K K N N N 4 IV'772 7v72 Y y Za yyZa z z' ,(1.3) a=1 i=1 a=1 3 % +1 = a,13 I 1 i ,a j1 i=j+1 i' where the indices a and 3 refer to elements of the K nuclei and i and j label elements of the N electrons. Furthermore, Za is the charge of nucleus a and r,/p, rF,i, and rij are the distance between nuclei a and 3, the distance between nucleus a and electron i, and the distance between electrons i and j, respectively. The complexity of the molecular Hamiltonian operator renders the Schr6dinger equation analytically intractable for all but the simplest systems. For this reason, certain approximations must be made to simplify the system of equations. Many methods treat the nuclear and electronic degree of freedom separately, in what is usually referred to as the BornOppenheimer approximation. Born and Oppen heimer [4] first introduced the concept of an adiabatic separation, in which the electronic and nuclear degrees of freedom are decoupled. In a simplistic description, the mass of the electrons, m, is so small in comparison to the nuclei, M, that the electronic motion is significantly faster than the motion of the nuclei, allowing one to perform a separate electronic calculation at each nuclear configuration. The re sult of this approximation is that the electronic potential energy can be constructed at a large number of nuclear configurations, resulting in a potential energy surface that is a function of the nuclear coordinates. At this point, the nuclear Hamilto nian is cinl''v,d in which the potential energy term is now a sum of the Coulomb repulsion of the nuclei and the specific electronic potential energy that correspond to the given nuclear configuration. This approximation gives rise to such methods as molecular dynamics (or quantum molecular dynamics if the nuclei are quantum). By epinll,ving excited state surfaces, multisurface methods arise, such as surface hopping methods. One immediate deficiency of these BornOppenheimertype methods is the lack of dynamical coupling between the electrons and the nuclei. As a result, it is desirable to discuss methods in which no potential energy surfaces are utilized. This requires direct calculation of the potential energy between nuclei and elec trons at each timestep. Such methods will be the focus of the remainder of this dissertation. Two specific methods will be discussed. First is the TimeDependent Hartree Fock method. Second is the ElectronNuclear Dynamics formalism. Both methods are timedependent methods that do not require the calculation of potential energy surfaces. 1.3 TimeDependent HartreeFock (TDHF) The TDHF method was first proposed by Dirac in 1930 [5]. The principal assumption of the TDHF method is that the selfconsistent central field approx imation holds. In this approximation (which is equivalent to the HartreeFock approximation) the total electronic wave function is composed of a product of wave functions for each electron. Furthermore, each electron moves in an average potential. The reference state for TDHF methods is the single determinant HartreeFock ground state, represented as a Slater determinant (using second quantization), ID) =a...a>vac). (1.4) The TDHF state vector has the form I) = eiF,," A(t a lD) (1.5) which has a corresponding Dirac density operator of the form P e x (rs t)a ', 'i I, . ty r,, (t)a~ (1.6) The density operator has a corresponding matrix form that is idempotent (F2 F) and the trace is equal to the number of electrons in the system [6]. The quantum mechanical Lagrangian takes the form = [(ih W ) + (lih l) ] ('WF P '), (1.7) 2 Ot Ot where # is the Fock operator, which has matrix elements of the form rs = hrs Y (rs pq)F,p (1.8) pq where hr, is an element of the oneelectron Hamiltonian matrix. The equations of motion for the TDHF method are derived from the fact that the Lagrangian should be stationary with respect to small variations in the TDHF state vector, such that "') > t 6A,,at aI), (1.9) r'S as well as a corresponding variation in the density operator. By requiring the Lagrangian to be stationary under such an arbitrary variation, one obtains the TDHF equations of motion i = [F, r], (1.10) where the dot indicates the time derivative and the brackets denote the commuta tion relation [a, b] = ab ba. At this point, it is assumed that the timeevolution will act as a small pertur bation to the Fock operator. The oneelectron Hamiltonian will be perturbed such that h,s h, + 6hs, (1.11) where ho, is the oneelectron Hamiltonian that corresponds to the unperturbed HF stationary ground state. The second term contains the timedependence. Secondly, the twoelectron component is perturbed through the density, such that p F + 6F, (1.12) where, again, the subscript indicates the HF stationary ground state solution and the timedependence is carried by the perturbative term. The above perturbations are then substituted into the TDHF equations of motion and, traditionally, terms of order two or greater in the perturbations are truncated. This implementation results in a method that is linearized and limits the electronnuclear couplings within TDHF methods [7]. Furthermore, as both components of the Fock operator depend explicitly upon the electron density, the equations of motion must be solved iteratively at each timestep. Again, this is traditionally not done. In most TDHF methods, the reference state is determined as either a configuration interaction (CI) expansion to the level of single excitations (resulting in the TammDankoff Approximation) [8] or as a perturbation expansion to the level of single and some double excitations (resulting in the Random Phase Approximation) [9]. By allowing the expansion coefficients to be timedependent, the timeevolution of the individual orbitals is now no longer required and the single reference is used throughout the dynam ics [7]. This is also a limitation to the TDHF method, as the state is not permitted dynamically outside of the limit of the single configuration included in the HF reference state. 1.4 ElectronNuclear Dynamics (END) The END theory is a nonadiabatic formulation allowing a complete dynamical treatment of electrons and nuclei that eliminates the need for potential energy surfaces that is equivalent to a generalized TDHF approximation. Electron Nuclear Dynamics is derived from application of the timedependent variational principle (TDVP) on a family of approximate state vectors parameterized in terms of Thouless coefficients. The END formalism and application of END to varied pllvi, .,l problems has been explained in detail in previous works [10, 11, 12, 13]. A brief account of the derivation of the END equations of motion will be provided in this section. 1.4.1 The END Equations of Motion The derivation begins with the quantum mechanical action functional A dt, (1.13) where Y is the quantum mechanical Lagrangian, defined as (O (1.14) Here ( defines a complete set of timedependent parameters that describe the particular choice of state vector and A is the quantum mechanical Hamiltonian operator for the system. By requiring that the action remain stationary under variations of the parameters (, one obtains the EulerLagrange equations d OY OY (1.15) dtA aQ1' for the set of dynamical variables {(Q} [14]. Application of the principle of least action results in a coupled set of firstorder differential equations of motion o iC* ( 0 ( ) (1.16) iC 0 * where E((*, () (1.17) is the energy of the system and acts as the generator of infinitesimal time transla tions and where C is an invertible metric matrix with elements defined as C 2 In S ( C 3 0(* (1.18) The metric matrices define the couplings between the dynamical variables. The argument of the logarithm in Equation (1.18) is the overlap, defined as S = ((I). The equations of motion given in Equation (1.16) are exact; any approximations to the END formulation arise due to the choice of dynamical parameters ( and the completeness of the electronic basis set. 1.4.2 Minimal ElectronNuclear Dynamics The minimal END theory is implemented in the ENDyne computer code. In the current version, nuclei are treated classically and electrons are represented using single Thouless determinants that are complex and single valued [15]. All electronnuclear couplings are retained. The dynamical parameters are chosen as (={Rk, Pk, Zph, Zph}, where Rk and Pk are the average position and momentum of the kth nucleus and Zph and z*h (spanned by a set of atomic orbitals) are the Thouless coefficient corresponding to the pth atomic orbital of the hth spin orbital and its complex conjugate, respectively. Choosing the state to be represented as a product coherent state vector I) = z; R IR, P) } \z)\ (1.19) where the nuclear wave function is expressed as a product of traveling Gaussians )J exp ( R )2 k(Xk Rk)j (1.20) k taken in the narrow wave packet limit (,, 0, for all ,, ), and where the electronic state is represented with a single determinant electronic wave function 4z = det{xh(Xp)}, (1.21) allows for a consistent description of the electronnuclear dynamics. The molecular spin orbitals in Equation (1.21), Xh, are spanned by a set of basis functions, {uh(x)}, that are Gaussiantype atomic orbitals (GTOs) centered on the average nuclear positions. By iC 0 iCp 0 iC* iC iC  where the dynamic m making these assumptions, Equation (1.16) takes iCR iCp z E/Oz iCh Cp z* oE/Oz* CRR I +CRpP E/OR I+ CpR Cpp P \ E/OP etric elements of Equation (1.22) are a2 In S (CxY)ij;kl 2 Im 0 OXikOiYj R'=R the form (1.22) (1.23) (1.24) 82 In S (Cxk)ph n S Ozd, Xik 1 R' =R and Chg 02 nS (1.25) OZ iOZqg The coupling is explicit in the metric elements given in Equations (1.23)  (1.25), nuclearnuclear coupling, nonadiabatic nuclearelectronic coupling, and pure electronicelectronic ('liii, respectively. Minimal ElectronNuclear Dynamics is a generalization of the TDHF method [10]. The electronic and nuclear interactions are still regulated through the Fock matrix, but there are some differences in the dynamical evolution of the state. Firstly, the variation of the onedensity in the END formalism is achieved through a general variation of the Thouless coefficients and is not truncated at any particular order. As a result, END is equivalent to fully nonlinearized TDHF. Secondly, the reference state is permitted to change as warranted by the dynamics of the system. This provides additional flexibility to the dynamical evolution, as the reference state is not limited to a stationary state as calculated previous to the dynamical evolution. CHAPTER 2 THEORY OF COLLISIONS Many chemical reaction principles and results can be elucidated and obtained by virtue of socalled " II iiii, exp. i iiil I.". The experiments most generally involve a beam of projectile species incident upon a reaction cell containing target species or upon a second beam of target species. The reactions and energy transferences occur in the volume of the reaction cell or the union volume of the crossed beams. Theoretical descriptions of such scattering phenomena are contained in the classical, quantal, and semiclassical realms of collision theory. 2.1 Scattering Theory While experimental results are obtained for bulk phase reactions in general, excellent theoretical descriptions of these scattering processes can be made by considering the interaction of a single projectile particle with a single target particle. At it's simplest, classical scattering theory involves the reduction of a twobody problem into a reduced onebody problem. This transformation results in a description of the motion in the centerofmass reference frame. This type of analysis is similar to the analysis of unbound Kepler motion, however, the nature of the central force is considerably different (though the functional form may be very similar) and all information about the orbits is lost (only the incident momentum, the final momentum, and the angle between the two is observed) [16]. In this section, we will discuss the fundamentals of classical scattering theory and methods to extend these fundamentals into the language of atomic and molecular collisions, quantum scattering theory. 2.1.1 Deflection Functions and Scattering Angles One begins by considering a classical collision system consisting of two point masses in the laboratory frame with some undefined potential acting between them that is a function only of the separation between the particles. It will be assumed that, initially, one particle is in motion (the projectile) and one is a rest (the target). Figure 2.1 demonstrates such a collision system. The scattering axis is defined as the axis parallel to the incident motion of the projectile and containing the location of the target, as demonstrated by A in Figure 2.1. The impact parameter, labeled as b, is defined as the distance at which the projectile is initially located in a direction perpendicular to the scattering axis. As one can assume a spherically symmetric potential, V(r), this can be in any direction perpendicular to the scattering axis. The other principle feature of Figure 2.1 is the deflection angle, labeled as T. The deflection angle is the angle through which the projectile is deflected by the potential, as measured in a counterclockwise manner from the positive scattering axis. The deflection angle can be positive due to a repulsive potential (as demon strated in A in Figure 2.2), or it can be negative due to an attractive potential (see B in Figure 2.2). Furthermore, under certain conditions, the projectile can orbit Projectile b /A Target Figure 2.1: Diagram of a classical collision system. ....."^ A. > "" B. Figure 2.2: Three possible deflection angles. A demonstrates a positive deflection angle. B demonstrates a negative deflection angle. C demonstrates a deflection angle with magnitude greater than 2w. the target for one or more periods, resulting in a deflection with magnitude greater than 2r radians (see C in Figure 2.2). It should be noted that, experimentally, it is impossible to distinguish between the three deflections shown in Figure 2.2. Experiment can only determine the angle at at which the projectile is scattered, relative to the scattering axis, and cannot elucidate whether a particle is scattered through a positive deflection, a negative deflection, or through an orbiting deflection with magnitude greater than 2r radians. Thus, we must define a parameter that corresponds to what is l1ieri .'ily measured; a parameter called the scattering angle. The scattering angle, denoted 0, is defined as 0 = I mod 2I1. (2.1) For any scattering system there exists a mapping that relates the impact parameter to a resulting scattering. This mapping is called the deflection function, O(b). It is obvious from the previous paragraph that the deflection function is not an injective Im.,i iiinJ.. as more than one impact parameter may result in the same scattering angle. Furthermore, the deflection function is not a surjective iii..''iir:. as there is no guarantee that the projectile will be scattered into every angle between 0 and r. 2.1.2 Cross Sections While reaction rates are usually the quantities desired from , .ii I ig exp. iinii. ii. the fundamental observable of such experiments is the differential cross section [17]. The differential cross section (da/dQ) is defined as da scattered current per unit solid angle ( (2.2) df incident current per unit area At this point it becomes necessary to assume some limited quantum nature of the scattering system. Without further justification at this point in the discussion, we will introduce the general form of the scattering wave function. It is customary to choose the ansatz wave function for the target beam by considering the ..i 'mp totic regions where the effect of the scattering potential produced by the target (assumed to be finite) is negligible. In the postscattering region, the scattered wave function will be a linear combination of incident plane wave and scattered spherical waves [18], eikfr e e""ik i) f(, ikr (2.3) where ki and kf are the magnitudes of the initial and final moment of the pro jectile, respectively. Both the numerator and the denominator in Equation (2.2) require evaluation of the probability current density, j, for the wave function, defined as 1 j (y*V V *) (2.4) 2mi where m is the mass of the projectile. It should be noted that vectors are denoted using boldfaced letters, while respective magnitudes are represented using non boldfaced letters. The numerator of Equation (2.2) is found from the current density of the outgoing spherical component by the expression jot dA, where dA is the unit differential area normal to the solid angle subtended by dQ. For this case, dA = r2 dr. It should be noted that, due to the area, only the r component of jot is required, greatly simplifying the calculation. The denominator of Equation (2.2) is likewise found considering the current density of the incoming plane wave, ji,. As all incident particles are considered, the denominator is just the incoming current density. Thus, the the differential cross section can now be expressed as da jot dA dQ Jii (2.5) kIf(0, 0) 2 ki Thus, the differential cross section depends only on the square of the amplitude of the scattered spherical wave. It should be furthermore noted that a statetostate differential cross section may be obtained if a probability amplitude for transition is included in the scattering amplitude The second quantity of import in ', ., I iiin exp. I iii, iii is the total cross section, a, defined as a= sin Ododo. (2.6) The total cross section depends only upon the relative kinetic energy of the colliding particles and is an effective area of scattering or reaction. If the projectile strikes within this effective area it will be scattered or reacted, however, because A B e e is not equivalent to B B B A Figure 2.3: Collision of distinguishable particles. of the avn.,.ii,:. the specific information contained within the differential cross section is lost. 2.1.3 Identical Particles Consider the collision of two particles in the center of mass frame (Figure 2.3). To simplify the explanation, it will be assumed that axial symmetry exists, eliminating a dependence on 9. Because each particle is distinguishable throughout the entire collision process, each can be determined to scatter through either and angle of 0 or an angle of 7 0. However, a interesting problem arises when identical particles are considered. Even if the particles are distinguishable at some time of separation, once the particles enter the interaction region they are no longer distinguishable. The particle scattered through angle 0 is indistinguishable from the particle scattered through angle r 0 (Figure 2.4). A A A A equivalent to A A Figure 2.4: Collision of identical particles. When considering identical particles, the ansatz wave function for the projec tile now takes the form [19] eikfr Se &ir + [f(0) f( 0)] eik (2.7) where the sign is determined by whether FermiDirac or BoseEinstein statistics are employed. As the particles are identical, momentum is conserved in the collision. From the same argument provided in Section 2.1.2, the differential cross section for identical particles in the centerofmass frame is = If(0)f( ) 2. (2.8) 2.1.4 Reference Frame Transformations One benefit of the END formalism is that calculations need not be performed in the centerofmass frame of reference. This is beneficial as experimental results are generally reported in the lab frame. However, it was demonstrated in the previous section that the calculation of the differential cross section for the reaction of identical particles is best handled in the centerofmass frame. Thus, a schema for transformation between lab and centerofmass frames must be developed. Consider the collision of two particles, as depicted in Figure 2.5. The first particle (projectile) has a mass, mi, the second particle (target) has mass, m2. In the lab frame, the particles have positions, ri and r2, and moment, kl = mn l and k2 = r22 (where the dot indicates differentiation with respect to time). In the frame relative to the centerofmass, the particles are defined by positions, sl and S2, and moment, pi = ms1 and p2 = m22. The centerofmass for the collision system is defined as having mass, M = mi + m2, position, R, and momentum, P = MR, relative to the origin. It is evident from the definition of the system that the position of a particle in the lab frame differs from the position in the centerofmass frame, at any time S r 2/ Origin Figure 2.5: Relation of collision pair and centerofmass to origin. during the dynamics, by the position of the centerofmass, ri s1 + R. (2.9) Differentiation of the above equation and multiplication by mi yields the relation amongst moment, k P mi kl = pl + P. M (2.10) If one considers the moment for the projectile in the two reference frames after the collision, one may obtain the relationship between the scattering angles in the two reference frames. The projectile possesses a momentum k{ in the lab frame P ___i^_p /, ^^ Figure 2.6: Relation of the scattering angles between reference frames. and pl in the centerofmass frame. The superscript f indicates a final condition. The projectile is scattered through an angle of a in the lab frame and an angle of 0 in the centerofmass frame (Figure 2.6). It is evident that the transverse components of the moment are equal, p{ sin(0) kf sin(a), (2.11) and that the longitudinal components differ by the momentum of the centerof mass, p{ cos(0) f cos(a) P. (2.12) M Transformations between the scattering angles are obtained by taking the appropriate ratios of Equations (2.11) and (2.12). The centerofmass scattering angle is obtained from the expression sin(a) tan(0) si (a) (2.13) cos(a) J This generalized equation can be much simplified under two assumptions. If the particles are identical then r = If the target is initially stationary in the lab frame, then the magnitude of the centerofmass momentum is equal to the magnitude of the initial momentum of the projectile (k ) at any time during the collision. Thus, the simplified expression takes the form sin(a) tan(0) os(a) (2.14) cos(a) 7 where 7 Likewise, the reverse relation can also be obtained by inverting the ratio. The lab frame scattering angle can be obtained from the general expression sin(0) tan(a) sin (2.15) cos(0) + 7 p17 Again, the expression can be simplified. If the particles are assumed to be identical, then (as before) 4 = 1 Additionally, the collision is elastic in the centerofmass frame, requiring that pl = P at all times in the course of the dynamics. Thus, sin(0) tan(a) cos() (2.16) which indicates that a = 10. Once the scattering angle is converted to the centerofmass frame, the deflection function and the differential cross section can be calculated. The center ofmass frame differential cross section must then be transformed back to the lab frame. This is accomplished by realizing that fact that the number of particles scattered into a given solid angle must be conserved between the two frames [20], such that d sin(a) da =( sin(0) do (2.17) " lab dQ CM or (da dc \ d[cos(a)] )dlab dJOM d[cos(a)] ( 21) The multiplicative factor in Equation (2.18) is evaluated by first applying the Law of Cosines to Figure 2.6 and substitution into Equation (2.12), yielding cos(O) + cos(a) cos( + (2.19) V1 + e2 +2 cos(0)' where = ~. Now the derivative can be taken with respect to cos(0), d[cos(a)] cos() + 1 (2.20) d[cos(O)] (l+ 2 2 cos(0))3 2' As the numerator above never goes to zero in the domain 0 = [7, r], Equation (2.20) can be inverted to yield the required factor, d[cos(O)] (1+ 2 + 2 cos(0))3/2 d[cos(a cos() (2.21) d[cos(a)] ( cos(0) + 1 The expression in Equation (2.21) is completely general to a collision pair. As before, simplifications can be made by considering the particles to be identical. As seen earlier, this assumption requires that = 1. Under this assumption, the transformation from the centerofmass frame to the lab frame takes the form d 4 cos (2.22) lab CM 2.2 Quantum Mechanical Treatment of Scattering Phenomena To this point the description of scattering phenomena has been tacitly quantal. The scattering wave function ansatz was introduced and a discussion of the identical particle problem was made, but no reference to the Schrodinger Equation has yet been made. In this section, the differential form of the scattering Schr6dinger Equation will be introduced, a derivation of the integral form of the Schr6dinger Equation will be made, and this integral form will be related to the scattering amplitude. The Born Series and the various Born Approximations, a set of selfconsistent solutions to the integral form of the Schr6dinger Equation is then introduced. Finally, the Schiff Approximations for large and small scattering angles, which are derived from the Born Series, will be fully derived and discussed. 2.2.1 The Integral Equation and its Relation to the Scattering Ampli tude The timeindependent Schrodinger Equation hV2(r)+ V(r)2(r) E (r) (2.23) 2m can be rewritten to take the form [V2 + k2] () (r)(r), (2.24) where k = /(2mE/h2) is the magnitude of the momentum vector for the scattered particle and where U(r) = 2mV(r)/h2 is the scattering potential. Equation (2.24) is an inhomogeneous differential equation with a solution [17, 18] y(r) o o(r) + G(r, ro)U(ro)(ro)dro, (2.25) where yo(r) and G(r, ro) are, respectively, a general solution and the Green's Function that corresponds to the homogeneous counterpart of Equation (2.24), which is of the form [V2 + k2] o(r) 0. (2.26) The above homogeneous differential equation is nothing more than the free particle Schrodinger Equation, and the solution is yo(r) = eik. Additionally, the Equation (2.26) is in the form of the Helmholtz Equation and therefore the corresponding Green's Function that satisfies Equation (2.25) takes the form [14] ikrrol G(r, ro) > ro (2.27) 47r ro Thus, one finds that Equation (2.25) now becomes, upon substitution, 1 eikrrol y(r) ek I U(ro)(ro)dro, (2.28) which is the integral form of the Schrodinger Equation with scattering potential located at position ro. It is now desirable to find an expression that relates the scattering amplitude to the integral form of the Schrodinger Equation. To accomplish this, one returns to the postscattering conditions in which the ansatz of Equation (2.3) is assumed to hold, namely that r > ro. Under this requirement, the angle between r and ro approaches zero and Ir rol r ro r. This approximation allows one to express the integral form of the Schrodinger Equation as 1 eikr. (2.29) (r) eikr e ikro U(ro) (ro)dro. (2.29) 47 r One can then immediately compare the modified integral form of the Schr6dinger Equation as given above with Equation (2.3), which demonstrates clearly that f(0) = f(k) eikroU(ro)(ro)dro. (2.30) 47 2.2.2 The Born Series Though one has now seen the solution to the scattering wave function in Equation (2.28), the solution is analytically intractable, as the integrand itself depends on y. Solution of the scattering wave function requires numerical integra tion to selfconsistency. This selfconsistent method yields the Born Series, with truncations yielding the Born Approximations [19, 21] of various orders. The First Born Approximation assumes that the scattering potential has a negligible effect on the incoming plane wave, thus the scattered wave function takes the form = eki'r, (2.31) The subscript on r is for indexing purposes only. This term can be substituted into Equation (2.30), and the scattering amplitude in the First Born Approximation then becomes f(k, kf) = / e ikfrU(rl)eiki.rldri. (2.32) 47 In some cases conditions are such that the First Born Approximation is sufficient for analysis, most notably when the scattering potential is not very strong and its effective range is quite small [19]. However, the First Born Approximation rarely provides adequate accuracy. The Second Born Approximation is begun by iterating on the wave function. The scattered wave function for the second order of approximation now becomes Seiki'1 + f G(ri r2)U(r2)2(r2)dr2. (2.33) Substitution of Equation (2.31) into the above equation yields the wave function in the Second Born Approximation Seikir1 + G(ri r2)U(r2)eikr2dr2 (2.34) which generates the second order scattering amplitude f(k, kf) = [ eikr l(rl)eikildr (2.35) + ek*rlU )G(r) r U r2(r2, k r2dr, rl This pattern can then be continued, and one can make successive iteration to an arbitrary degree and obtain the nth order wave function Seikfr G(ri r2)U(r2)eikir2 dr2 + f I G(r r2)U(r2)G(r2 r3)U(r3)e ik 3drdr +...+ (2.36) .+ G(r, r2)U(r2)G(r2 r3)U(r3) X ... x G(r,_1 rn)U(rn)eiki.rndrn ... dr2. The infinite Born Series is obtained by letting n oo (with careful reordering of the indices). The infinite Born Series leads to a scattering amplitude with the general form f (k, kf) = eikf.U(rn)G(rn r,_l)U(r,_l) rn=1 (2.37) x .. x G(r2 rl)U(rl)eiki' dr... drl. Now, while Equation (2.37) is exact, it is very unwieldy to solve and still depends upon the scattering potential function, a property not utilized within nor obtained from END. Furthermore, the Born Series is often afflicted with slow or no convergence [2]. So it is desirable for a number of reasons to express the Born Series in a simplified form, in particular a form that does not depend explicitly upon the scattering potential. 2.3 SemiClassical Treatment of Scattering Phenomena: The Schiff Approximation The work in this section was first introduced by L. I. Schiff in 1956 [22]. This discussion is based on that paper. The scattering potential and the projectile particle can be characterized by several pir, i. .,1 properties. The scattering potential can be described by V and R, which are rough indications of the strength of the potential and its effective range, respectively. Likewise, the projectile is characterized by its kinetic energy, T; its wave number (magnitude of the momentum), k; its speed, v; and its scattering angle, 0. These characteristics provide a qualitative means by which to discuss the scattering process. For example, the First Born Approximation is valid under conditions when the scattered wave is insignificantly perturbed by the scattering potential (c.f. Equation (2.31)). This can be expressed in a qualitative manner by the magnitude of V being very small compared to the collision energy (weak scattering potential) or by the magnitude of R being very small (short range scattering potential). Schiff provides a general condition under which the First Born Approximation is valid as given by the expression (IVIR)/(hv) << 1 [22]. The Schiff Approximation provides an approximate scattering amplitude for large collision energy collisions. Specifically, the assumptions imposed are that IV/T << 1, that 0 be very large or very small in comparison to 1/(R), and that the scattering potential be slowing varying when compared to the incoming wavelength. In contrast with the First Born Approximation, the Schiff Approximation is valid for any magnitude of (IVIR)/(hv) [22]. 2.3.1 Schiff Scattering Amplitude for Large Angles The Schiff Approximation consists of approximately representing each of the terms in the Born series using the stationary state approximation [22, 23, 24]. One begins the derivation with the infinite Born series, Equation (2.37), rewritten slightly to the form f (k, kf) ... ei(kf r~krl)U(r)G(r r,_i)U(r,_i) n=1 (2.38) x x x G(r2 rl)U(rl)dr, dri. A change of variables is made, such that p, = r r,_l, so that r, p,, _+r,_1. The argument of the exponential term can be expanded iteratively as kf r, + ki ri kf (p,_ + r,_) + k (r2 pi) kf p,_1 kf (p~2 + r_2) + ki (r3 P2) ki pi Skf p,_ kf p_2 kf p,_3 ... kf *(p, + rm) (2.39) + ki (r, Pm) ... ki p2 k *i pi = kf P_ kf Pn2 kf p_ ... kf p, + q r, ki pi ... ki p2 ki pi, where q = ki kf is the momentum transfer for the collision process. To complete the transformation of variables in Equation (2.38), one must calculate the Jacobian of the transformation, specifically pi '_ '' 'P 'P 2 L Opl 'P '' 'P I 'P 2 L S_ Or2 Or2 Or2 Or2 r2 Or2 (2.40) Opi 'p_ 'P L 'P L 'P 2 'P L Dr.n Drn Drn Orn Orrn rrn, noting that the r, row and column are omitted. The derivatives take the general form 1 forf = g apf 1 for f = g (2.41) 0 otherwise From this, it can be justified that the Jacobian matrix has the values of 1 along the diagonal, the values of 1 for each element immediately below the diagonal, and the value 0 elsewhere. From this is can be seen that the determinant takes the value J = 1. Therefore, dp,_l...dp = dr,...drm+ldrm _...dri, (2.42) and the limits of integration do not change. The change of variables may now be accomplished, transforming Equation (2.38) into the following form f (ki,k) =  / / e' ikp 1...k*Pm+1+qrmk'Pm 1...kiP .n lrn 1 x U(rm Pn+ p + Pn2 + . + Pm)G(Pni) x U(rm + Pn2 + P3 P+m)G(pn2) x ...X x U(rm + pm)G(pm)U(rm)G(pm )U(rm Pm1)G(pm2) (2.43) X ...X x U(rm Pmi Pm2 .. P3)G(P2) x U(rm Pmi Pm2 .. P2)G(Pl) X U(rm Pmi Pm2 Pi) x drmdpn1 dp, ...dp2dPl. where the fact that r im+pf_l+pf_i+...+P, for f >m rf rm 1Pf2 f > (2.44) rm Pm1 Pm2 P1 for f < 7 was used to indicate the p dependence of the potential terms. The second summa tion (over m) arises from the the stationary phase approximation. The majority of the integral will come from the regions where the phase is stationary, specifically where the derivative of the phase with respect to the wave vector is zero [24]. Thus, one must sum over all of these stationary phase points to completely enumer ate the integral. Now, recognizing that the transformed Green's function for the scattering amplitude equation now takes the form eikfP G(r, r,_l)  G(p,_) e p (2.45) 47p the above equation can be partitioned into the form Sn nl x p ei(kfPn1kf'Pn k g(p n)d 1 J (2.46) X Pn_2ei(kP 2k p 2)g(pn_2)dp_2 X ... SP21ei(kfP2ki'P2)g(P2)dP2 x f Pei(kpki pl)g(Pi)dl In Equation (2.46), the terms g(pj) represent a product of all U terms that depend explicitly upon pj. Thus, for an elastic scattering process, the scattering amplitude has been reduced to a product of integrals of the form I pg(p)ei(kk)dp (2.47) *2j7 To evaluate this integral, one needs to transform to spherical polar coordinates (p, 0, q) where 0 measures the angle of k with respect to p, I= j d sin Od pg(p, 0, )e(kkcos)p, (2.48) where the scalar product is geometrically evaluated. One can now make the familiar change of variable p = cos 0, where dp = sin OdO. Substitution (with appropriate change of limits) yields I= d0 dp pg(p, 0)eikp(p)dp. (2.49) One can now integrate by parts with respect to pt, resulting in the expression I= f r pd i p ( )eikp() ikp(1p)d} (2.50) Integration of the second term in the braces will bring down another factor of 1/(kp). By applying the limits on the first term, one can write I = du pdp i [g(p, 1,) g(p, 1, )eikp] + (k2). (2.51) o J kp As the collision energy is assumed to be large (k >> 1), one can reasonably expect the magnitude of kp to be somewhat greater than 1. Thusly, the exponential will be highly oscillatory and subsequent integration over p will result in a negligible contribution to I. As a result, one can omit this term, as well as the terms of order k2 and higher. Furthermore, if one assumes a spherically symmetric scattering potential, the value of 9 becomes arbitrary and integrates out to a factor of 27. Thus, I i j g(kp)dp. (2.52) k Jo The kp arises due to the facts that p = 1 (0 = 0) and that 9 is arbitrary, and so the potential term now depends only on the magnitude of p in the direction of the scattered momentum. So, the integration over the angles serves to transform the dependence on the vector p in the scattering potential into a dependence upon the magnitude of p in the direction of the momentum. Now, returning to Equation (2.46) one finds that integration of the angular dependence results in the equation i oo n ( / n1 coo f(k, kf) 4= k eiq'r drm/ g(kfpn1)dpn1 n=lI m=/ k J (2.53) x g(kfpn2)dp.2 x ... x g(kiP2)dP2 g(kip1)dpi}. The individual g terms can be factored back to their respective U terms, with the p dependence replaced by kp dependence. Specifically, one finds that f.k<,f +) {( i) n1j drn dP... dner oo n n1n2 x U (r + kf(pn_1 P+ n2 ** Pmn)) X (2.54) x U rm + kf Pm U (rm) U (rm kipmi) x ... xU (rn ki(Pml + Pm2 +... + P) At this point, a second change of variables can be made, such that Pj + Pj+ + ... + pPm2 Pm,1 for j < m sj P= m + Pm+i +  + Pj2 + ji1 for j > m (2.55) Pm for j = m There are two recursion relations that arise, specifically sj = p + sj+l for j < m and sj = sj1 + pj for j > m. It is clear from following the previous analysis that the Jacobian is unity, but the lower limits of integration must be transformed. When integrating over pj for j < m, one finds that the lower limit of pj = 0 transforms to sj+l. Likewise, when one integrates over pj for j > m, one finds that the lower limit transforms to sj1. Therefore, the change of variables may be imposed, transforming Equation (2.54) into the form i n ( n1 o f(k kf) = { ( i /drmei r ds1U (rm m i1s) 47 1 mnl I 1 2 SdS2U (rm iS) X ... X j dSm r m iSm2) X j dSm1U (rm kism,l) dsmU (rm kfSm) 0 0 0 X dS+iU (rm+ kfSm+l) x ... x dSn2 (rm+ kfSn2) JS "Sn3 X ds,lU rm+ fSn1) s } Sn2 (2.56) At this point, the above equation can be partitioned into two parts, the part explicitly depending upon the scattered direction (kf) and the part depending explicitly upon the incident direction (k ). If one considers the product depending only on the incident direction, then yet another substitution can be made, such that W1j U (rrn kiSj 1_ S (2.57) and Wo U (rm kismi) dsm1. (2.58) Thus, calling the product in question K, one finds that K j dsiU (ri sis j dS2U (m iS2) ... K S2 J S300 X dSm2U rm kim2) dSmlU (mr kiSm1) (2.59) SWml_ 1 0 d m2 .. W/ dH0W _2 X ... X W dW2 / di. 0 0 1 1 The expression for K can now be integrated, with the following results: rWo rWm fW3 rW2 Jo Jo Jo Jo rWo rWm 1W? K dW, dW, x ... x dWW2 2WO d,_1 id. X.. 2W J O Jo Jo 2 x 3 d1WMl dWm_2 x ... x I dW4W4 2 x 3./n .n .O (2.60) S2x3x (  3) dW m 2x .." m )J Jo dWm2W,7rn3 1 fWO 2x3x... x (m 3) x (m 2) O l l 1 2 x 3 x ...x (m 3) x (m 2) x (n 1) r/IW m1  [(m l)!] [ U (rm kis ds A similar substitution can be made for the product involving the direction of scattered momentum, allowing for the evaluation of these n m integrals, j dSmU (rm kfSm) j dSm+U (rmn + kfSm+l) .. X dSn2U (rm +kfSn2) X X dSnU (rm + fS n1 Sn3 n2 [(n m)!]1 [ U (r + kfs) ds] The scattering amplitude can then be expressed using these terms, . f(ki, kf)  drme 1 e U (rm) n 1 rn 1 x [(m 1)!] [ U (rm is) ds] S[(n m)!]1 U (rm+ k) ds] nm} (2.61) (2.62) At this point, the above equation can be further simplified, but cannot be made devoid of its dependence upon the scattering potential, which was the intention. However, a second approximation can be made which will allow a general expression that does not depend explicitly upon the scattering potential. 2.3.2 Schiff Scattering Amplitude for Small Angles The above expression for the scattering amplitude is valid under conditions that k is large, so that the second term in Equation (2.51) is negligible and that the scattering angle is large, requiring that the n different stationary phase points are distinct. For this discussion, one begins with the assumption that the scattering angle is very small, such that 0 << (kR)1/2. This results in two specific requirements. First, the n distinct stationary phase points now converge to a single point, thus eliminating the requirement for summation over m. Secondly, the scalar product in the exponential of Equation (2.62) can be written in its Cartesian components, eiq'rm ei(q9m+qyym+qzzm) (2.63) Since the scattering angle is assumed to be so small, the angle between ki and kf is nearly zero, ensuring that q, is negligible. This can be explicitly demonstrated if one requires that ki coincide with the zaxis. Thus, q q ki S(k kf) k, k kikf ' = ki kf cos (2.64) = kf (1 cos 6) kf 02 2 In the above demonstration, the fifth line was obtained by virtue of the fact that ki = kf for elastic scattering. Furthermore, the sixth line was obtained by nipl,ving the Maclaurin series representation of cos 0 with truncation at the second term. It becomes clear that, because 0 is so small, the longitudinal component of the momentum transfer is negligible, as would be expected ]li, i, .illy. This is not true of the transverse components of the momentum transfer. Again, one must consider the momentum transfer, qAq q = + q, (2.65) where the longitudinal component has been omitted. Using the definition of the momentum transfer, one finds that 2 2 )2 q1 + qY (ki kf)2 = +k k+ 2kikf cos0 = 2k 22k cos0 (2.66) 2k (1 cos0) 2k2 (02) kf 2 Thus, it is clear that the transverse components of the momentum transfer are both of the order of ~ kO and cannot be neglected. Converting to Cartesian components, the scattering amplitude now becomes f(k kf) { f /" x do 1 ,I I n=l x o i(qx+q~y) foo dzm x dzmU (xm, [m, zm) [(m 1)!]1 U (xrm, ym, z) dz (2.67) JOO [J d x [( )!]1 U (r, yz) dz ,zm  To obtain the above form, it must be recognized that the first integral over s in Equation (2.62) represents the prescattering trajectory while the second integration over s represents the postscattering trajectory. As the scattering angle is so small, these to wavevectors are effectively the same. Therefore, instead of integrating over two different trajectories one can integrate a single trajectory for both regions (oo to Zm and then Zm to oo). The integration over Zm is of principle interest, so one can define it as M, such that M f dZmU (xm,, Ym, Zm) [(m 1)!]1 U (xm, Y, z) dz m (2.68) x [(n m)J!]}~ U (xm", m, z) dz One can then introduce the notation that w U (xm, Ym, z) dz (2.69) 00 and r00 a U (xm, m,, z) dz. (2.70) oO In order to be able to integrate over w, one must evaluate the derivative of w with respect to m,, dw d fzm dz dz / U (Xm, ym, z) dm dZm oo (2.71) SU (Xm, yn, Zm) . Therefore, it is clear that dw = dzmU (xm, y,, Zm). Finally, the new limits of integration must be evaluated, leading one to find that j U (Xm, ym z,) = 0 for zm = 0 w U = (2.72) f J U (,m' Y z) ) a for zm = oc Thus, the integral over Zm now can be rewritten as M [(m 1)! (n m)!] 11 a ) mdw. (2.73) This equation must now be integrated by parts through m 1 iterations. For the first integration, one can let u = w"1 and dv = (a m)""dw. Thus, M j W (a w)ndw S__m__ wnm _ (m )! (n m)! n (m )(2. 1 1 )Wjm(a ( w) w )) dw" M I 12 (m 2)!(n (m 1))! o For the second integration, one can let u = w"2 and let dv = (a w)"(l)1dw, such that M = ja W 2(a w)n(m 1)dw (m 2)!(n (m 1))! o 1 1 a m2(a W)n(m2) (m 2)!(n (m 1))! n (m2) 2 (2.75) n (2 (%) ( W)n(m))d n (m 2) 1 3 (m3)!(n (m 2))!(a 2 Following this pattern, one finds that after m 1 integration, M takes the form M (a w)ldw. (2.76) (n 1)! f The final integration can now be performed, yielding the following expression, M = (a w)ndw (n 1)! S (a w) (2.77) (n 1)! n n! The scattering Thus, amplitude is no longer dependent upon m, as should be expected. f (k, k) = J ax 1!) a 47 t 2k k n! 24rf J x 2k) n! k m "I 2 17 00 2k n! 2 1 2 J dx dyei(q +qy) (e i) ik ( r"  Ix .hpi,. +qyy) 1 exp U (x, U, z) dz . 27x 1 2k ) (2.78) At this point, it should be noted that the integration over x and y is taken over the plane perpendicular to the scattering axis. In the limit of an axially sym metric scattering potential, x and y can be represented by the impact parameter, b = /x2 2, and the azimuthal angle, 0, such that x = b cos 0 and y b sin 0. Thus, the transformation to cylindrical coordinates renders the scattering ampli tude into the form f(ki, kf) = d(1 e(qb cos+qb sin)bdb I exp U (b, z) dz] 27 Jo 2koo (2.79) Now, as the scattering potential is spherically symmetric, the momentum transfer will always be parallel with b. As a result, the expression qb cos 0 + qyb sin 0 will have the same value regardless of the value of 0. Therefore, one can choose any value of 5, and then sum all of the contributions from 0 to 27. Specifically, if one allows = 0' = 0, then q, = q and qy 0. Consequently, qb cos 0 + qb sin = qb cos 0'. (2.80) The integral over 0 now takes the form of the zerothorder Bessel function [14], 2 eqbcosdo 27Jo(qb). (2.81) The small angle Schiff Approximation to the scattering amplitude can now be expressed as f(k, kf ) ik Jo(qb) exp U (b, z) dz bdb. (2.82) Unfortunately, the scattering amplitude is still an explicit function of the scattering potential, a term that is not obtained during an END calculation. This expression can, however, be transformed into one that does not depend on the scattering potential [23]. If one assumes that r is the position of the projectile at any point during the trajectory, then the fact that the scattering angle is very small requires that r r (z + b) (z + b), (2.83) where z is the position of the projectile along the scattering axis. As the impact parameter is perpendicular to the scattering axis, it is clear that r2 = 2 + b2. From this, one finds that dr dz = (2.84) [1 (b)]1/2 Therefore, the argument of the exponential now becomes I 1 f U(r) i U (b, z) f U(r) ( )d 2k o 1 2 (b)]1/2 o [(2.85) 1 p 1U(r) k [1 (b)]1/2 where the fact that the limits of integration remain the same when the transforma tion from z to r is made. Furthermore, the symmetry of the scattering potential was utilized in the second step. The above expression is the MasseyMohr ap proximation for the semiclassical phase shift, 6(b), in the small scattering angle limit [21, 25, 26] which is the semiclassical analog of the Kennard approximation for classical scattering processes [23, 27]. Specifically, the M..  Ml 1ir approxima tion is given the form [23] 6(b) k 1 2 (2 12 dr. (2.86) In the .i'. mptotic limit, r >> b, therefore, U (b, z) 26(b). (2.87) 2k Finally, this means that the Schiff Approximation scattering amplitude for small angle scattering takes the form [23, 28] f (ki, kf) ik Jo(qb) (1 exp [2i(b)]) bdb. (2.88) While Equation (2.88) no longer depends explicitly upon the potential, it has introduced the semiclassical phase shift, which itself is not a value pro vided by END calculations. This is easily remedied, though, by the fact that the semiclassical phase shift is directly related to the deflection function by the expression [19, 21, 23] (b) (b) (2.89) k db The Schiff Approximation provides an effective method for calculating the scattering amplitude without any prior knowledge of the scattering potential. The principle power of the Schiff Approximation is the fact that it explicitly contains all of the terms in the infinite Born series. Furthermore, the Schiff Approximation provides good differential cross sections in the region of small scattering angle, which the region most often reported by experimental studies. CHAPTER 3 BASIS SETS FOR DYNAMICAL HARTREEFOCK CALCULATIONS Quantum mechanical calculations on molecular systems are traditionally divided into two general categories, the localized valence bondtype (VB) calcula tions and the delocalized molecular orbitaltype (\10) calculations. In many ways, these two descriptions serve as complements to each other. Where VB methods, in general, provide more visual 1l., 1i .,1 descriptions, such as molecular geome tries, they generally fail with important properties, such as paramagnetism and bond strengths. The MO method tends to provide a more rigorous quantum me chanical description of molecules, such as electron delocalization over the entirety of a molecule and variations in bond strengths, but generally offers little in the way of 1ir, i' .illy intuitive descriptions of the molecular system. Both classes of methods are formulated from rather coarse approximations, and, as a consequence, each possesses deficiencies, in particular with relation to electron correlation. The VB methods tend to overestimate electron correlation, while MO methods tend to underestimate correlation effects [29]. While each has its strengths and weaknesses, and while each converges to a common molecular wave function rep resentation in the absence of their respective generalized approximations, a large community of quantum chemists has chosen to employ MO methods for quantum chemical calculations. The MO method, and its application in the HartreeFock approximation, will be the focus of this chapter. Some address must be made toward notation. The literature of ab initio quantum chemistry, as with most other fields, is not completely consistent. While it is not crucial to describe a set of canonical variables, clarification of symbolic Term Table 3.1: Notation employed in this review. I Meaning x a,o3 dV dr Uppercase Roman indices Lowercase Roman indices Lowercase Greek indices Exact electronic wave function Generally approximated electronic wave function Electronic spinorbital Spatial electronic factor Spin electronic factor General electronic basis function Element of spatial volume (no spin included) Element of volume (spin included) Nuclear indices Molecular orbital expansion coefficients Atomic orbital expansion coefficients notation up front is the most lucid manner of handling the issue. The notation of Szabo [30] will be most closely followed, with minor adjustments as needed. Table 3.1 provides an exhaustive list of notation used in this chapter. 3.1 The HartreeFock Approximation The MO method in quantum chemical calculations of molecular systems has become synonymous with the HartreeFock (HF) approximation. In this section, the HF approximation and the HF wave functions will be introduced. Additionally, the concept of electron correlation and the correlation effects will be addressed. Finally, an in depth discussion of HF basis sets will be made. 3.1.1 Partitioning of the Molecular Wave Function In Section 1.2, the nonrelativistic molecular Hamiltonian operator was introduced, with the explicit form K 1 72 K NK K K N N Y KV72 K ZAZB >1>1 ZA K ,(3.1) A=l i= A=1B i= A+1a j i=j+l i, where the indices A and B refer to elements of the K nuclei and i and j label elements of the N electrons. Furthermore, ZA is the charge of nucleus A and TA,B, TA,i, and rij are the distance between nuclei A and B, the distance between nucleus A and electron i, and the distance between electrons i and j, respectively. The complicated coupling between the nuclear and the electronic degrees of freedom cause the Schrodinger equation to be practically insoluble by direct methods. To simplify the problem, the BornOppenheimer (BO) approximation was introduced, in which the nuclear degrees of freedom are assumed to vary on a much larger time scale than the electronic degrees of freedom. In essence, the BO approximation allows for a "clamped nucleus" model, in which the nuclei are fixed in space and the electronic degrees of freedom are allowed to vary on this fixed nuclear framework. As the nuclei are motionless, the nuclear kinetic energy terms (the first summation in Equation (3.1), above) become identically zero. This assumption decouples to nuclear degrees of freedom from the nuclear degrees of freedom. The resulting form of the Schr6dinger equation is now called the electronic Schr6dinger equation and takes the form N K K K NN N ZC za ZAZ Z f ZA YEY (3.2) i=1 A=l B=A+1 rAB Al i=1 r j= i=j+1 i It is the solutions to this equation with which the majority of quantum chemistry is focused. Following the example of L6wdin [31], one can partition the electronic Hamil tonian operator into three parts, N N N ea = ho + h' j ih. (3.3) i=1 i=1 j=i+1 In the above equation, the term ho is the nuclear repulsion (or zeroelectron Hamiltonian) term, which depends only upon the fixed nuclear degrees of freedom and is defined as K K h A ZAZB (3.4) A=1 B=A+1 The second term in Equation 3.3, h}, is the oneelectron Hamiltonian. The oneelectron Hamiltonian includes the additive kinetic energies of the individual electrons as well as instantaneous attractive potentials between the N electrons and the K nuclei and takes the form h' v. Z (3.5) A=1 'ia The final term, hi, is the twoelectron Hamiltonian, which generates the in stantaneous repulsive forces between each individual electron and the remaining molecular electrons. The twoelectron Hamiltonian takes the form h = (3.6) The terms denoted by ho result in a constant for a given nuclear configuration, and therefore will not have any bearing on our choice of molecular wave function. Rather, the total energy will only be increased by the fixed potential energy values associated with h. The above partitioning of the molecular Hamiltonian operator in Equation (3.3) allows one to write the molecular wave function as a product wave function of the form ed ({fi}; {(A}) I ({fi; {rA})I)2({ri}; { A}), (3.7) where the factor wave functions are eigenfunctions of the operators in Equation 3.3, specifically, 4K1b = 1i, (3.8) ,42 eA2, (3.9) where e1, and e2 are the energy eigenvalues. The total energy for the given nuclear configuration, 8, becomes the sum of the energy eigenvalues and the nuclear potential energy defined by h, S= Si + e2 + Vn. (3.10) The major focus of quantum chemistry involves the solution of Equations (3.3), (3.7), and (3.10), and we will use these equations as the starting point for the discussion of the HartreeFock approximation. 3.1.2 The HartreeFock Wave Function The HF approximation begins with the assumption that the total electronic wave function can be approximated by a product of oneelectron wave functions. Furthermore, one must assume that the potential experienced by a given electron is an average of the potentials produced by the remaining electrons [29]. In this approximation, the term ho and the N different h' terms are maintained, but the (N 1) twoelectron terms are replaced by N additional one electron terms. By replacing the twoelectron terms with oneelectron terms, the HF approximation does not explicitly treat the instantaneous interaction of individual electrons. Rather, each electron is treated as if it were influenced by an average field produced by the other electrons in the molecule. This approximation allows for relatively accurate quantum chemical calculations despite the gross approximations imposed, though several important 1.r, i, .,l descriptions are omitted. These correlation effects will not be discussed in this chapter. The new electronic Hamiltonian can be rewritten in the form N i [ =v[h HF] (3.11) i= 1 where h1 is the oneelectron Hamiltonian (as previously defined) and VHF is the new HF oneelectron potential energy term (to be defined later in this section). Because the zeroelectron Hamiltonian only adds a constant factor to the electronic energy, it is omitted. The most striking benefit of this reformulation of the electronic Hamiltonian is the fact that, because the Hamiltonian operator is now a sum of oneelectron terms, one can approximate a single Nelectron wave function as a product of N oneelectron wave functions (called a Hartree product), N HP ({xi; {A}) j( ; FA}), (3.12) j=1 which is an eigenfunction of the oneelectron Hamiltonian, j7lHP= EHP. (3.13) As before, the eigenvalue is the energy of the system. The wave function defined in Equation (3.12) is a product of a set of N spin orbitals, Xy, which are themselves eigenfunctions of the respective oneelectron Hamiltonians, [h + VHF] X = X, (3.14) where ce is the energy eigenvalue of the jth orbital. The total energy, E, is a sum of the orbital energies, N E = e, (3.15) j=1 The electronic Hamiltonian in Equation (3.11) does not have a spin dependence, therefore a transformation from spatial electronic coordinates, {ri}, to spatialspin coordinates, {xi} by employing a product form of the spinorbitals does not change the energy eigenvalues of Equation (3.13). Thus, in the absence of relativistic effects, the N spinorbitals can be represented by a product of a spatialdependent factor, y)({ri}), and a spindependent factor, either "spin up", a(wc), or "spin down", P(cw), such that Xji = ( }) ) (3.16) Y j({rt})X(1c) The variable ac represents the general spin variable. It should be noted that the product wave functions indicated in the above equation have the same spatial factors with varying spin factors. This is consistent with restricted HF method. Allowing each spinorbital to have a unique spin AND spatial factor results in the unrestricted HF method. [30]. The partitioning of the spinorbitals given by Equation (3.16) satisfies the Pauli exclusion principle by allowing two spinorbitals to have the exact same spatial factor, but opposite spins (that is, one molecular orbital with a paired set of electrons). However, a further requirement of the molecular electronic wave function is that it must obey the FermiDirac statistics, in particular that the electronic wave function must be antisymmetric with respect to exchange of electron indices. The Hartree product given by Equation (3.12) does not satisfy the FermiDirac statistics, being that the sign of IJHP remains unchanged if two indices are exchanged. Adherence to the FermiDirac statistics traditionally required detailed group algebra. This was greatly simplified by Slater, who circumvented group theoretical descriptions by introduction of the spinorbital function directly into a determinant that would later bear his name [32]. Slater exploited the property of matrix algebra that, given a matrix, inter changing any two columns of the matrix will change the sign of the determinant of the matrix. Thus, Slater constructed a matrix in which the spinorbitals are placed as the columns and the occupying electrons are placed as the rows of the matrix. The determinant of this matrix is the most general antisymmetrized product wave function. Mathematically, one finds that, Xi(xi) X2(X) .. XN (X) sater (N!)1 (N!) det {Xj(xi)}. (3.17) X1(xN) X2(XN) ... XN(XN) At this point some discussion of the approximated oneelectron potential in Equation (3.11) can be made. It is beyond the scope of this review to explicitly derive the form of the term HF, so it must suffice to say that the potential form can be obtained by employing Lagrange's method of undetermined multipliers to the energy eigenvalue equation [30] esl3tSlater = E TSlater. (3.18) The resulting HF Hamiltonian (variously referred to as the Fockian and denoted F) is a sum of N Fock operators, which satisfy the eigenvalue equations of the form N N hiJ+ A Xi = iXi (3.19) ji ji The term in brackets in the above equation is the Fock operator for orbital i. It can be seen that the approximated oneelectron potential has been split into two components. The first is the coulomb operator, which defines the interaction of electrons with an average potential. The action of the coulomb operator is to provide an average repulsive potential felt by an electron at the position xl that arises from an electron in a second orbital. The coulomb operator has the inverse r form of a coulomb interaction, weighted by the probability density of the orbital to be averaged, specifically, J(xi)x(xi) [/dx2,X(x2) 2rxi, j(X1). (3.20) The coulomb operator arises as a consequence of the assumption that the electronic Hamiltonian is a sum of oneelectron operators only. The second operator results from the antisymmetrization of the wave function through the use of a Slater determinant. This operator, the exchange operator, results in the exchange of two electrons and produces a oneelectron potential that is dependent upon the value of the orbital in question throughout all space [30]. The form of the exchange operator is J(xi)xj(xi) dx2Xj(x2)r 1x(x2) Xj(X). (3.21) The effect of the coulomb operator is localized and has a classically intuitive interpretation, however, the action of the exchange operator is nonlocal and depends upon the location of these two orbitals in the spinorbital space. 3.1.3 Solving the HF Equations: Basis Set Expansions In Section 3.1.2, we introduced the HartreeFock integraldifferential equations (now rewritten using Dirac notation), j QSlater) E IISlater (3.22) by using Equations (3.18) (3.21). Traditionally, one makes consistent use of Slater determinants for the molecular wave function, and superscripts on T will consequently be dropped for the remainder of the chapter. The solution of these equations is still a nontrivial task. The first methods of solving these equations, specifically for small atomic systems, was through numerical integration [29]. A considerable breakthrough was introduced by Roothaan [33] in 1951. The computation routine of Roothaan involved expanding the molecular wave function in a basis of atomic spinorbitals with the general form of a linear combination of spatial atomic orbital basis functions multiplied by the appropriate spin function. This expansion allowed the HartreeFock differential equations to be written as a set of algebraic matrix equations, which could be solved using available linear algebraic techniques. The general form of the basis set expansion of the ith spatial component of the wave function takes the form K 1 = (3.23) P11 where the basis set is composed of K functions, {q,L = 1, 2,..., K}. The term c1i represents the expansion coefficient for the ptth basis function to form the ith spatial component of the wave function. The size of the expansion (K) is in general not limited to a specific number. In fact, an infinite expansion would be desirable, as this would correspond to the full HartreeFock wave function. However, this is practically impossible due to computational limitations. The size of K should ideally be large enough to offer the best descriptions of the molecular orbitals without becoming computationally inefficient. This, along with other specific criteria that must be considered will be addressed in the next section. At this point it will suffice to assume that some general function form and expansion size has been decided upon. The HartreeFock equation for a given spatial component of the wave function can now be written as .) = e .), (3.24) where ce is the orbital energy of the ith spatial orbital. The basis set expansion for '. can now be introduced. The result becomes K K Y,#E ^ K (3.25) p=1 p=1 At this point, one can then multiply through on the left of the equation by an arbitrary basis component ({w1 and integrate. The ith orbital energy is a number, and thus can be extracted from the integration, as can the expansion coefficients. Thus, the equation now becomes K K > (9wci?qV)cIO )Cpw. (3.26) p~ 1 p~ 1 At this point, one must introduce two matrices with elements that are related to the terms in the above equation. The first is the Fock matrix, F. The Fock matrix has elements defined as F = (   ). (3.27) The second matrix is the overlap matrix, S. The elements of the overlap matrix are given as S1 = (< ,1) (3.28) and arises from the fact that the basis functions are not necessarily orthogonal. The overlap matrix provides a measure of the linear dependence of the set of basis functions. Due to assumed normalization, the diagonal elements of the overlap matrix all have a magnitude of one. The offdiagonal elements will range in magnitude between zero and one. Elements approaching one will demonstrate a strong linear dependence between two basis functions, while a value approaching zero indicates a strong linear independence. The HF equation can thusly be expressed using the newly defined matrix elements, K K F F.,pci = Scici. (3.29) /p=1 p1=1 It is clear at this point that the above equation is an element of one single matrix equation of the form FC = SCe, (3.30) where F and S are as defined above, C is the K x K matrix of expansion coeffi cients, and c is a K x K diagonal matrix with the orbital energies as the diagonal elements. Equation (3.30) is commonly referred to as the RoothaanHall equa tion [34]. A companion equation exists for unrestricted determinants, called the PopleNesbet equation [30, 35]. The PopleNesbet equations have individual matrix equations for each set of spinorbitals, as each pair of spinorbitals has a different spatial component in the unrestricted formalism. As the overlap matrix and the Fock matrix are both Hermitian (and in many cases real and symmetric), relatively simple solution techniques are available for solving the RoothaanHall equations [36, 37]. The most common method involves diagonalization of the matrices by use of a unitary transformation. This method presents the eigenvalues (elements of c) and eigenfunctions (elements of C) for the matrix equation. The rudiments of solving the RoothaanHall or PopleNesbet equations will not be discussed any further in this dissertation. Rather, in the next two sections specific interest will be placed on the general form and construction of basis sets for use with these matrix equations. 3.2 General Forms and Properties of Basis Sets To this point, no mention has been made as to the functional form that the basis set should take. In general, any functional form is possible, but certain properties are desirable. Specifically, the wave function for the system must be singlevalued, finite, continuous, and squareintegrable [38]. It is thus desirable that one choose basis functions that possess these characteristics. A set of atomic orbitals (AOs) is an immediate choice for a given basis set, as AOs satisfy the above criteria as well as offer a chemical intuitiveness lacking in other choices. 3.2.1 SlaterType and GaussianType Orbitals A first choice of basis sets was a set of spatial orbitals with the same functional form as the hydrogenic orbitals, (2( ))2n+'1 1/2 )STO [ (2) I e2 +l ). (3.31) Equation (3.31) is called a Slatertype orbital (STO) [29, 34]. The equation involves the terms (, which is the orbital exponent; n which is the principle quantum number; 1, which is the azimuthal quantum number; and m, which is the magnetic quantum number. The term YQ(0, 9) is the spherical harmonic. As the functional form of hydrogenic orbitals are derived from specific linear combinations of STO basis functions, one may reasonably expect that they will be good approximations to orbitals in multielectron atoms. In fact, STO functions have correct functional forms at small r values and also at very large r values. As can be seen from Equation (3.31), STO functions do not have any nodal structure, and therefore linear combinations must be constructed to properly mimic the nodal structure of atomic orbitals. However, basis sets built from relatively small linear combinations of STO basis functions have been quite successfully employed in quantum chemical calculations [29]. Despite the desirable functional form of the STOs, the principle drawback is numerical [34]. The large portion of computational time in the HF method is the calculation of the manycenter twoelectron integrals. Slatertype orbital basis functions do not admit simple analytical expression for such twocenter integrals, and must therefore be numerically integrated [39]. This is a time consuming process and, as a result, for any system with more than a few atoms the accuracy obtained by using STO bases is outweighed by the severe decrease in computational efficiency. This limitation was circumvented by Boys, who ii.'., ed using Gaussian type functions instead of exponential functions [40]. The functional form of a Gaussiantype orbital (GTO) is [41] [G1(2/)1/2 () 1/2 (GTO ) 2 n 1 e2C (,2(0) YrO [ (2n 1)!! r (3.32) where n, 1, and m are the same quantum numbers as given in Equation (3.31). The parameter a is a different orbital exponent specific to the GTO basis. Note that angularly, STO and GTO basis functions have the same functional form, the only difference lies in the radial factors. The primary deficiency of the GTO basis is immediately clear. A GTO function does not provide the correct functional form as r 0 or at r oo. In particular, an STO function will have dysTo/dr > 0 at r = 0, whereas a GTO function will always have dGTO/dr = 0 at r = 0. Furthermore, the GTO basis function has a much faster dropoff in the tail of the orbital as r oc than does the STO function. This is a fairly severe limitation to the accuracy of GTO bases in computations, as the functional form of the STO basis sets leads to superior accuracy over GTO bases. The top panel in Figure 3.1 demonstrates the considerable difference between the radial part of the hydrogen Is orbital represented as a single STO function and a single GTO function. This red curve is the hydrogenic Is STO function and the blue line is a single GTO function in which the exponent has been optimized to provide the best leastsquares fit to the STO [41]. One can improve the structure by building a linear combination of GTO basis functions (with optimized exponents) for each STO function cupllvd. Both the set of exponents for each primitive Gaussian orbital and the set of contraction coefficients must be optimized to fit the Slater in question. The bottom panel in Figure 3.1 shows that a linear combination of six GTO functions (blue) provides a much better fit for the hydrogen Is orbital as represented by a single STO function (red), but that it is still not particularly good at the cusp of the orbital. In general, a very large number of GTO functions with large exponents would be needed to correctly mimic the cusp of an STO function, but no matter how many terms where included in the linear combination, the derivative of the GTO function would still be zero at r = 0. Despite this cusp deficiency, quite good accuracy can be obtained from a good sized linear combination of GTO basis functions per Slater. Yet, one must ask why a linear combination of six or more GTO functions is favorable over a single STO 18 16 14 12 12  08  06 04 02 0 ...........  0 1 2 3 4 5 Radial Distance (a u) 2 18 16 14 12 E S08  06 \ 04 02  0  0 1 2 3 4 5 Radial Distance (a u ) Figure 3.1: Comparison of STO and GTO representations of the radial part of the hydrogen Is orbital. Top: A single GTO function fit to a single STO function. Bottom: A linear combination of six GTO functions fit to a single STO function. function, particularly when the accuracy of the STO basis is superior to the larger GTO basis. The answer lies in the fact that when a product of two GTOs is taken, the result is a third GTO [30]. This reduces a multicenter twoelectron integral to a considerably simpler analytic form [39]. As a consequence, a calculation utilizing a larger GTO basis is much more efficient that a calculation using a smaller STO basis. This substantially greater computational speed has lead to the fact that most calculations of polyatomic systems have traditionally cmiql..d GTO basis sets [42]. Other basis sets have been i'pl,vid in various studies, such plane wave basis sets. However, these are still not as widely used as STO bases and, particularly, GTO bases. The remainder of this chapter will specifically focus on Slater and Gaussian basis sets. 3.2.2 The Structure of Basis Sets Now that the general form of an orbital basis function has been chosen, either STO or GTO, the full basis set for an atom must be constructed. While there are (theoretically) no limitations on the construction of basis sets, practicality has restricted quantum chemists to certain accepted forms of contraction. Minimal basis sets Any atom must have at least enough orbitals available to completely contain the required number of electrons. Any basis set designed to completely represent only the ground state orbital structure of an atom is referred to a minimal basis set [30]. For example, a minimal basis set for a Mg atom would contain no less than 6 orbitals, the Is, 2s, 2p,, 2py, 2p,, and the 3s orbital. It would therefore be sufficient to build a minimal basis of 6 STO basis functions, three with n = 0 and three different values of ( as well as three others with n = 1 and a fourth value of Furthermore, it would be possible to build a minimal basis set of six single, uncontracted GTO basis functions in the same manner. However, as the previous discussion indicated, the accuracy would be considerably less than the for the comparable STO basis. So, to remedy this, a specific number of Gaussian primitive functions (single GTO basis functions) are generally fit to each STO function in the minimal set. This describes the formalism to construct a minimal STONG basis set [43]. The title of this basis set indicates that N Gaussian primitives are optimized to fit a single STO function. While any number of primitives may be fit to a given STO function, in practice the STO6G is the largest Gaussian minimal basis that is eiplv1d. Without exception, the STO3G is the most widely used minimal basis set. Minimal basis sets are notoriously inaccurate basis sets. In general STONG basis sets offer qualitatively correct descriptions of fundamental chemical proper ties such as bonding and can be "ipl1,vd for initial guesses and for calculations involving very large molecules where more complete basis sets would be computa tionally inefficient. However, the small size of the STONG basis set is prohibitive to the use of minimal bases for calculations in which even moderate accuracy is required [30, 44]. As a final point concerning STONG minimal bases it should be noted that for certain atoms these bases are not truly minimal. Many times chemical bonding is not correctly mimicked using a truly minimal basis set. For this reason, the Group 1A and Group 2A metals, as well as the first two rows of transition metals, will include the lowlying p orbitals even though they are unoccupied in the unbound atom [34]. The s and p exponents for a given principal energy level are identical. This is computationally more efficient than allowing for separate s and p exponents as each set of s and p orbitals will then have the same radial behavior and can, consequently, be integrated together [30]. Doublezeta and splitvalence basis sets One of the main limitations of the minimal basis set is the fact that there is no flexibility for the generated orbitals to change size under the influence of intramolecular surroundings. Each orbital has a single set of exponents that control the size and shape of the orbital, and while the amplitude of the orbital can be adjusted through the HF coefficients, the spatial size cannot be changed. To remedy this, one would desire to include more than a single linear combination of basis functions for a given orbital. This idea gives rise to the next level of basis sets, the doublezeta and the split valence basis sets. Much of this discussion will be made in terms of STO bases. It should be remembered that these ideas can be translated directly to GTO bases by requiring that each STO function be constructed as a linear combination of GTO functions. Specifically, the doublezeta functions allow for each orbital to be a composed of a linear combination of two STO basis functions, each with a different exponent. The larger exponent (a tighter function) is, in general, slightly larger than the optimal exponent for the single zeta function, while the smaller exponent (a more diffuse function) is slightly smaller [30]. Furthermore, as the doublezeta basis is more complete than a comparable singlezeta basis, a correctly optimized linear combination of doublezeta functions will be a better representation of the ]li,, i. .1 orbital than a singlezeta function. This results in an improvement in the ground state energy, as demonstrated by Clementi and Roetti in their seminal paper on doublezeta functions for atoms [45]. A simplification of the doublezeta basis can be made by realizing that, during a chemical process, the size and shape of core atomic orbitals will not change significantly. Therefore, under most any conditions, a single welloptimized STO function will provide a sufficient representation of a core orbital. One can then allow for the valence orbitals to be represented using a combination of two STO basis functions, as in the doublezeta basis. This is the formula for constructing splitvalence basis sets [34]. Splitvalence basis sets provide ground state energies that show improvements over minimal basis sets, but that are not as good as doublezeta bases. Again, this is due to the fact that a more complete description of core orbitals is obtained with the doublezeta functions. However, this energy difference is small when compared to the computational efficiency gained through using splitvalence bases [30]. The most common splitvalence basis sets that are employed in computational chemistry are the l:mnG basis sets. These are combinations of GTO functions such that the core orbitals for the atom are represented as a contraction of I GTO basis functions. A given valence orbital is represented as a linear combination of two basis functions, one composed of m GTO primitives and the second composed of n GTO primitives [34]. The only exceptions occur with H and He, in which there are no core shells and only the valence shell structure is used (and as a consequence the 4:31G bases for H and He are identical to the 6:31G bases). The most common examples of splitvalence bases are the 3:21G [46], 4:31G [47], and 6:31G [48] bases. This pattern can be extended to larger basis sets such as the 6:311G basis [44]. Doublezeta and splitvalence basis sets improve the electronic representation in several ways. First, as mentioned above, the doublezeta bases allow for a better description of atomic orbitals by virtue of the increased completeness of the basis set. Also, the valence orbitals are now flexible enough to change size during a chemical process. In particular, this allows for better descriptions of bonding and anisotropic chemical processes, such as the anisotropy of the bonding porbitals when forming a and 7 bonds in systems with bond orders greater than unity [34]. One last feature of splitvalance bases is that, with proper optimization, the more diffuse valence functions can behave as virtual orbitals, a property not available in minimal basis sets. This feature will become increasingly important when one desires to investigate dynamical processes. Polarization basis sets Many ]li, i. .,l processes require not only a change in the size of atomic orbitals over time, but also a change in the shape of the orbital. Examples include the behavior of the orbitals for an atom subject to an external electric field or the orbitals of an atom which has some nonzero momentum. Both of these processes result in a polarization of the atomic orbitals. This polarization causes a net increase of electron density in one area offcenter from the nuclear center and a corresponding net decrease in electron density in the region of space immediately opposite the nuclear center. While splitvalence basis sets allow for the size of the orbitals to fluctuate, they do not permit the shape of the orbital change. This cannot be accomplished by merely changing the size of the orbital exponents. In order to change the shape of atomic orbitals, the basis must flexible enough to allow combinations of basis functions that represent occupied atomic orbitals with higher angular momentum basis functions[30, 34, 44]. The most common methods of polarization involve the addition of basis functions that mimic a dorbital to the elements from Li to Ar. This level of polarization is denoted using a single asterisk. The most common example is the 6:31G* basis [49], which is the basic splitvalence 6:31G basis described above with the addition of a single dsymmetry basis function [34] (or fsymmetry basis function to transition metals) [44]. The second form includes the d(f)orbitals for heavy atoms as well as an basis function with psymmetry to H and He. This level is denoted with two asterisks (such as 6:31G**)[49]. Again, this pattern can be "'aplvd using larger splitvalence bases, resulting in such combinations as 6:311G**, a basis that is commonly employed for correlated calculations [50]. The effect of including polarization functions has traditionally been observed in structural properties, particularly in constrained systems where the electron density is shifted away from the nuclear centers [34], and in systems subjected to external electric fields [30]. In particularly, atoms which can be multiply bonded will, in particular, require a greater degree of polarization. This is evident from the valencebond description of chemical bonding, in which the formation of "hybrid" atomic orbitals (which are nothing more than polarized atomic orbitals) is the underlying principle of chemical bond formation [51]. For this reason, it has long been accepted that it is more important to include polarization basis sets on main group elements rather than singly bonding species, such as Group 1 elements, as indicated by the strong statement of Szabo and Ostlund [30] that it Ii.i' been empirically determined that adding polarization functions to the heavy atoms is more important than adding polarization functions to hydrogen." While this is generally true for structure calculations, it cannot be accepted when dynamical calculations are being performed. In particular, any atom that possesses nonzero momentum will experience a polarization of its electronic orbitals due to the motion of the atom. This effect will be present in all atoms, including H and He. For this reason, it is of crucial importance to include polarizing pfunctions on H atoms for dynamical calculations. Diffuse basis sets The previously mentioned basis set structures do a good job of describing various chemical processes, however, all of them locate the electronic density relatively close to the nuclear centers. The splitvalence structure allows for increasing the size of orbitals, but the exponents are always close in magnitude to the exponent in a comparable singlevalence basis function. This limits the ability of the orbital constructed from a splitvalence basis to expand beyond small fluctuations around the size of the orbital constructed from a singlevalence basis. Additionally, the polarization functions allow for shifting of the electron density away from the nuclear center, say to a chemical bond. Again, this shift in the density is not large. As a consequence, systems with large electron densities that are located a significant distance from the nuclear center (such as anions and systems involving Rydberg states) are not properly modeled using minimal, splitvalence, or polarization basis sets. To properly describe such systems, diffuse basis sets must be employed [34, 44]. Diffuse basis sets are structured in a manner very similar to splitvalence basis sets. A minimal (or splitvalence or polarization) basis set in constructed and additional basis functions are included to provide for the diffuse atomic orbitals. However, the exponents of these diffuse basis functions are much smaller than for the valence basis functions, resulting in an electron density that is located much further away from the nuclear center. In general, the diffuse functions are of the same angular momentum as the valence basis functions. This means that a carbon atom would incorporate one additional s and one additional pfunction. This diffuse structure is denoted using a + symbol [34]. If a single diffuse sfunction is added to a hydrogen or helium atom, then this is denoted by two plus signs. Thus, one can now begin to evpl']' a virtual alphabet soup of such combinations as 3:21+G*, 6:31++G*, or 6:311+G**. Eventempered and universal eventempered basis sets Further advances in the building of basis sets were made when it was realized that, as a basis set got larger, the orbital exponents within a given angular momen tum converged to a geometric sequence. This geometric sequence takes the form (i = CoP, (3.33) where Q( is the ith exponent in the sequence, (o is the largest exponent, and 3 is a constant that is specific to the angular momentum. A basis set that is constructed using this type of method is called an eventempered basis set [39, 42]. The general feature of an eventempered basis set is that it limits the number of parameters that must be optimized. Furthermore, eventempering ensures that a GTO expansion of an STO is wellspanned, with no regions in which the representation is particularly poor. There is in general a small energy price that must be paid, but this usually is on the order of several hundredths of a Hartree [42]. Furthermore, it has been postulated that, if enough eventempered exponents are include in a contraction, then this set of exponents would eventually become identical over an entire row of the periodic table. This leads to a universal even tempered basis set [39, 42]. Other basis set structures In this section, a brief discussion has been made of the general structure of basis sets, with details given about those basis sets most commonly employed in quantum chemical calculations. This is only a very small sampling of the basis sets available for computation, however, most of these basis sets include the basic principles listed above. It is beyond the scope of this work to provide an in depth discussion of the different types of basis sets used in calculations. The comparisons made in the forthcoming sections will, in general, be related to the types of basis sets reviewed in this section. 3.3 Method for Constructing Basis Sets Consistent with Dynamical Calculations The structure of a basis set is heavily dependent upon the types of 1pr, i cal properties that one desires to calculate. In some cases diffuse functions are required, in others polarization functions are required. In most cases, some bal anced combination of all of the properties is needed. Because of this, many have viewed basis set construction as an art (or black magic in some cases). Yet, no matter what form the basis set takes, most have one trait in common: with very few exceptions, basis sets must be optimized. Usually basis sets are optimized with respect to the energy of the ground state of the system by means of the vari ational principle. However, the HF equations are nonlinear, and therefore any optimization process becomes computationally costly. An additional feature that is common to most basis sets is that they are built for use in stationary state calculations of the ground state of a given system. While most of these basis sets do provide representations of unoccupied orbitals either due to constructing a splitvalence or diffuse basis set, these unoccupied orbitals bear little or no resemblance to the virtual orbitals in the system. In this section a new method for the construction of basis sets will be intro duced. This method has a simple ].li,i .,1 underlying justification. The method does not require expensive and complex energy optimizations (and does not, in fact, require any optimizations at all). Finally, this method allows for the construc tion of 1.1r, i, .,lly meaningful virtual orbitals, a necessity for the computation of a wide variety of dynamical properties. 3.3.1 Basis Set Properties for Dynamical Calculations In stationary state calculations minimal basis sets are rarely, if ever, sufficient for the description of the chemical species in question. As outlined in the previous section, a variety of extra basis functions must be included to improve the descrip tion. In this case any set of functions added to the minimal basis set generally demonstrates no 1ir, i .,1 resemblance to atomic orbitals in the system. Rather, they just act to provide a more complete spanning of the electronic Hilbert space, increasing the flexibility of the basis set. These extra basis functions serve to lower the ground state energy (due to an increase in the accuracy of the representation of the occupied atomic orbitals and in some cases by partially accounting for the correlation energy in the system) but they do little else. For dynamical calculations, particularly charge transfer processes, the basis set must be flexible enough to allow for correct description of electronic transi tions between atomic orbitals, either within a single atom or molecule or between the collision pair. One specific aspect of this requirement is the fact that virtual orbitals (atomic or molecular) must be available for occupation throughout the dy namical processes. To this point, little effort has been made by the computational community to construct atomic basis sets that properly describe virtual atomic orbitals, mainly due the fact that most basis sets are optimized with respect to the total ground state energy of the system rather than optimized with respect to the individual atomic orbital energies. For this reason, and with very few exceptions, stock basis sets that are most commonly "mplvid in computational chemistry are inadequate for use in dynamical calculations. Methods for improving stock basis sets traditionally follow along the lines of increasing the size of the basis sets by the inclusion of more and more uncontracted diffuse primitives. This method will, in the limit of infinite expansion sizes, lead to correct virtual energy levels by virtue of the HylleraasUnsold separation the orem [52]. However, this brute force method is extremely inefficient for basis set construction. As the number of basis function increases, so does computation time. While this is extremely limiting in the area of structure theory, it is virtually im possible in dynamical methods such and END, where a large number of calculations must be made per trajectory (with many trajectories required for a single collision energy). For this reason, a new method must be devised for building basis sets that include correct representations of virtual orbitals. The construction of dynamically meaningful virtual orbitals is dependent most strongly upon two properties of the atomic orbitals, the energetic of the orbitals and the shapes of the orbitals. The energetic of the orbitals are the most obvious concern. If the orbital energies are not correct, then the energy required for electronic transitions within the basis set will not properly model the energy differences in nature. Too small of an energy gap will result in increased transfer probability while too large of a energy difference will have the opposite effect. In a single determinantal treatment of the electrons, the orbital energies should mimic the energetic of the system. Therefore, the correct energetic of the atomic orbitals should be a measure of the ability to correctly simulate dynamical transitions. The shape of the orbital wave function is also important for dynamical calculations. The nodal structure of the wave function determines regions in which electron density is zero or where it is nonzero. Again, a basis set must correctly model the electron density in an atom. More discussion about these properties will be made in the next sections. However, at this point it will suffice to say that both of these properties can be addressed quite effectively through the use of STO basis sets. Specifically, the energetic of an orbital is largely dependent upon the structure of the tail of the orbital wave function. As was mentioned in the previous section, STO bases correctly describe the tail of hydrogenic orbitals and, likewise, reproduce the orbital tails in manyelectron atoms quite well. Additionally, while single STO functions do not contain any nodal structure, linear combinations of STO functions can if carefully built. For this reason, it seems most reasonable to construct basis sets using STO functions, at least initially. 3.3.2 Physical Justification for the Basis Set Construction Method As a part of the pli, i. .,1 justification of the proposed method for basis set construction, one must first return to the radial factor of a Slatertype orbital basis function, given by the form S(2)2n+rl1 1/2 RsO r lec r. (3.34) (2n)! j When one compares this form to the hydrogenic orbital functions, one finds that the orbital exponent is related directly to the nuclear charge of the hydrogenic atom in question, specifically z, (3.35) where Z is the nuclear charge and n is the principal quantum number. In the case of a hydrogenic system, only one electron is associated with the system and therefore the electron will always feel the full nuclear charge (that is, there is no nuclear shielding due to the presence of other electrons). This idea can be extended to the construction of a wave function for an orbital in any chosen atom. It should be noted, however, that for a given electron in an arbitrary orbital in a manyelectron atom the nuclear charge felt by that electron will not be the full nuclear charge. Rather, the full nuclear charge will be shielded by electron density located between the orbital in question and the nucleus. This gives rise to the concept of an effective nuclear charge, Zeff. The orbital exponent now takes the form S Zeff (3.36) The effective nuclear charge will vary as a function of the principal quantum number, in general an orbital with a smaller principal quantum number will have a larger effective charge. Zener provided values for these effective charges based on variational calculations [53]. Slater [54] determined an empirical method for calculating the effective nuclear charge for any arbitrary orbital and presented the equation Zff = Z s, (3.37) where Z in the full nuclear charge and s is a screening constant that is a function of the orbital and the number of electrons. For s and porbitals, the screening constant was defined to be s = 0.35N, + 0.85N,_1 + 1.OON,,2, (3.38) where N, is the number of additional electrons in the same principal level, N,_1 is the number of electrons in the principal level immediately lower, and N,_2 is the number of all remaining electrons in lower principal levels [54]. Having defined the effective charge, Slater then II..i r construction of atomic orbitals as single STO functions of the same form as Equation (3.31), with the exception that the orbital exponent takes the form of Equation (3.36) and the principal quantum number n is replaced by an effective quantum number, n*. The effective quantum number deviates from the true quantum number only for n > 3 [54]. This, the form of the wave function consists of a single STO for each set of n and 1 quantum numbers, with the orbital exponent equal to the effective nuclear charge felt by the corresponding atomic orbital divided by the effective principal quantum number associated with the orbital. Having discussed a specific method for the construction of basis sets for many electron atoms, it is now time to consider some 1ir, i, .1l aspects of atomic orbitals in a bit more detail. One property that a basis set should properly model is the radial distribution of the electronic orbital, defined as [55] D,(r) = r2 [Rn(r)2 (3.39) In the above equation, r is the radial distance from the nucleus, n is the principal quantum number, I is the azimuthal quantum number, and Rni is the radial factor of the wave function for a given n and 1. The top panel of Figure 3.2 demonstrates the radial distribution functions for the first four sorbitals in the hydrogen atom. The humps in the radial density function for a given orbital indicate regions in which the probability for the electron to exist is greatest. This feature demonstrates that the electron density corresponding to an orbital takes the form of concentric shells of electron density [29, 55]. From this it becomes clear that an atomic orbital with principal quantum number n will posses n regions of electron density, each becoming increasingly closer to the nucleus but with smaller probability. This fact is mirrored in the nuclear screening expression advanced by Slater and given in Equation (3.38). The shielding due to electrons in the next lowest principal level do not completely shield the nucleus, rather they have an only an e.i' effective shielding due to the penetration of the higher principal levels. This ]1li,i .l1 feature of atomic orbitals is not limited to sorbitals, nor is it limited to the description of the hydrogen atomic orbitals. In Figure 3.2, the 06 05 04 03 02 01 02 0 5 10 15 20 25 30 Radial Distance (a u) 02 0 18 0 16 014 0 12 01 0 08  006 004 002 0 "   . 0 5 10 15 20 25 30 Radial Distance (a u ) Figure 3.2: Top: Plot of the radial distribution function for the Is (), 2s ( ), 3s (. .), and 4s ( ) orbitals of the hydrogen atom. Bottom: Plot of the radial distribution function for the 2p (), 3p ( ), and 4p (. ) orbitals of the hydro gen atom. bottom panel demonstrates the radial density functions for the 2p, 3p, and 4porbitals of the hydrogen atom, where a similar shell structure is observed. Furthermore, Figure 3.3 shows the radial distribution functions for the occupied s and porbitals in the Ar atom, in the top and bottom panel, respectively. In both plots, the orbitals wave functions that are plotted are the doublezeta wave functions of Clementi and Roetti [45]. Again, the shell structure is clearly evident. At this point, the most significant pil,vi, .1l attribute of these shells of electron density becomes apparent. In the case of s and porbitals in the atoms of the first 35 3 25 O 15 05 0 05 1 15 2 25 3 Radial Distance (a u ) 25 2 15  05 / "" 05 0 05 1 15 2 25 3 Radial Distance (a u ) Figure 3.3: Top: Plot of the radial distribution function for the Is (), 2s ( ), and 3s (. .) orbitals of the argon atom. The Isorbital is scaled by a factor of twothirds. Bottom: Plot of the radial distribution function for the 2p () and 3p ( ) orbitals of the argon atom. few rows of a periodic table, the radial location of the shells is largely independent of the principal quantum number associated with an orbital. In other words, all sorbitals have a shell of electron density that has roughly the same radial location as the shell of electron density due to the Is orbital. Likewise, all sorbitals with principle quantum number n have n1 shells that have roughly the same radial location as the n1 lower energy sorbitals. This is also true for porbitals. The consequence is that any orbital of s symmetry possesses a partial character of all of the lower energy sorbitals. Further, each of these characteristic shells exists within the same region of effective nuclear charge that is specific to that shell of electron density. Now, relating the regions of effective nuclear charge back to the idea of the Slater orbital exponent, this now means that a given orbital can be constructed as a linear combination of all of the previous orbital basis functions, each with a specific orbital exponent that relates to the effective charge region experienced by the corresponding shell of electron density. Specifically, the orbital wave function can be written as Sn, = cNN,,e (3.40) i=1 where N,,I is the appropriate normalization factor for the STO in question, the ci's are the expansion coefficients, and the & is the orbital exponent (effective nuclear charge) for the ith shell. On the surface, this is nothing new. Relating the orbital exponent to an effective charge was proposed by Slater and the linear combination is nothing more than a restatement of the superposition principle [1] through which the HF method determines the HF eigenstates as a linear combination of the basis vectors [29]. However, this plr, i, .,l insight does serve as an important undertone for the basis set construction proposed in this work. 3.3.3 Construction of the Basis Set One begins construction of the basis set by defining a linear combination of STO functions, each with orbital exponents derived from the effective charge experienced by the orbital in question. The form of the effective charge may be determined in any number of ways, the simplest of which is to employ Slater's formulation for the screening constant [54]. However, Slater's screening constants are empirically modeled and may not be as accurate as those determined by other methods. Instead, one may consider the work of Clementi and Raimondi [56] as an extension of the earlier work of Slater and Zener. Clementi and Raimondi made a study of the elements through Kr, representing each atomic orbital as a single STO function (that is a minimal basis set) with variable orbital exponents. The exponents where then optimized (with respect to the ground state energy) using an SCF procedure. At this point, the proposed method has not deviated from older methods of basis set construction. By eiply,ving this method, one can construct a minimal basis set for the ground state of the atom in question, however, no recourse is available for construction of virtual orbitals. The optimization process could be extended to determine the orbital coefficients for the virtual states, but this would require calculations beyond the HF level to do so, as the transition energies would be dependent upon electron correlation. The minimum level of theory that could be eipl1v d would be configuration interaction. This would increase the computational effort required to optimize the virtual orbital exponents. To remedy this, a new approach to determining these virtual orbital exponents is proposed in this work; a method that is extremely simplistic in its application, yet has proven to be very powerful. The method begins through the investigation of the behavior of Clementi's shielded orbital exponents as a function of the atomic number. Figure 3.4 demonstrates the functional behavior for the Is orbital exponents through Kr. As can be seen from the figure, the orbital exponents exhibit a very linear dependence on the atomic number. The top panel of Figure 3.5 shows the dependence of the 2s orbital exponents on the atomic number. There are two data sets in this figure, the data of Clementi and Raimondi [56], which are denoted using the plus symbols, and the data from the present work, which are denoted using the open circles. Clementi's data still demonstrates a linear dependence of the exponent on the atomic number, however, it becomes clear at this point that more than one linear region is observed. As an example, the orbital exponents for the elements from Li through Ne have a 4 0    i  i  i  i  i  40 35 35 + + 5+0 + + + 0 5 10 15 20 25 30 35 40 + ++ 10 + 15 Atomic Number Figure 3.4: Plot of the is orbital exponent for the atoms through Kr as a function of atomic number. The data are from Clementi and Raimondi (+). slightly different slope and intercept than for the remaining elements. Each block of elements will have a slightly different slope. The data points current to this work will be discussed at length later. The data in the bottom panel of Figure 3.5 show the same relation for the 2p orbital exponents. Likewise, the plots in Figure 3.6 demonstrate the dependence for the 3s (top) and 3p (bottom) orbital exponents and those in Figure 3.7 present the data for the 4s (top) and 4p (bottom) orbitals. The most striking feature when comparing all of the previous plots is that, while the individual linear regions become more distinct from one another as the principal and azimuthal quantum numbers increase, the (local) linear dependence of the orbital exponent on the atomic number is still quite strong. It is this feature that defines the proposed method for virtual orbital construction. Discussion must now be made with regards to the remaining sets of data points, those denoted with the open circles. These data points are new to this work and are derived from the data of Clementi and Raimondi. Specifically, these points are orbital exponents corresponding to virtual orbitals for the atoms in question. are orbital exponents corresponding to virtual orbitals for the atoms in question. 14 12 10 O+ 8 o + IQ + 0 0 5 10 15 20 25 30 35 40 Atomic Number 18 14 12 + + 10 + 0 5 10 15 20 25 30 35 40 Atomic Number 16 +  + + Figure 3.5: Plot of the 2s and 2p orbital exponents for the atoms through Kr as a function of atomic number. Top: The 2s orbital exponents. Bottom: The 2p + orbitalexponents. he data are rom lementi and aimondi () and rom the present work (o). These virtual orbital exponents are determined by first considering the hydrogen atom. As the hydrogen atom has only a single electron, there will never be any nuclear shielding for that atom. This the electron will always experience the same effective nuclear charge (of unit magnitude) no matter in which orbital the electron has probability for existing. This means that Zeff = 1 always, and the orbital exponent for any orbital in the H atom is just equal to the reciprocal of the principal quantum number. At this point, one makes reference to the linear behavior of the occupied orbital exponents. The virtual orbital exponents are then 12 + 0 I 7  65 4 + + ++ 0 5 10 15 20 25 30 35 40 7+ Atomrc Number 6 4 + 2  + o o o + ao 0 5 10 15 20 25 30 35 40 Atomic Number 7 5 + + 4 + I+ S+ + + + Figure 3.6: Plot of the 3s and 3p orbital exponents for the atoms through Kr as a function of atomic number. Top: The 3s orbital exponents. Bottom: The 3p orbital exponents. The data are from Clementi and Raimondi (+) and from the present work (o). determined by making a linear interpolation between the H atom virtual orbital exponent and the exponent corresponding to the first available occupied orbital in that symmetry. For example, in the case of the 3p orbitals, the interpolation is made between the exponent corresponding to the H atom 3p orbital and the exponent that corresponds to the 3p orbital of Al (atomic number = 13). This method for determining the exponents that correspond to virtual atomic orbitals relies on welldocumented trends exhibited by a parameter that is related directly to a 1]li, i. .1l property, namely the trend in the regions of effective charge + 25 + 25 + 15 0 1  0 oo 0 5 10 15 20 25 30 35 40 Atomic Number 25 + + + + 15 ++ 0o a function of atomic number. Top: The 4s orbital exponents. Bottom: The 4p orbital exponents. The data are from Clementi and Raimondi (+) and from the present work (o). as demonstrated earlier in this section. However, it is important to note that a number of severe assumptions have been made. Perhaps the two strongest assumptions are that the virtual orbital exponents demonstrate the same linear behavior as do the occupied orbital exponents and the assumption that there is not a large change in the magnitude of the orbital exponent as one transitions from the occupied orbital exponents to the unoccupied orbital exponents. While neither of these assumptions can be tested without the construction of energyoptimized virtual orbital wave functions, the severity of the assumptions is tolerated in lieu of the ease of application. And, in spite of these assumptions, the calculations using basis sets constructed from this starting point have yielded surprisingly good results for both stationary state and dynamical calculations, as will be demonstrated in the next section. Once the interpolations have been accomplished, the next step is to construct the atomic orbital wave functions. This is begun with the wave function for the Is orbital, which is represented by a single Slatertype orbital, is = Nls((ls)eClsr. (3.41) The expansion coefficient in just the normalization coefficient for an STO basis function with n = 1 and with orbital exponent (Ps. From this, the wave function of the 2s orbital can then be constructed as a linear combination of the Is wave function (providing for the cusp) and a single STO basis function with n = 2 that is used to represent the tail portion of the orbital. The orbital wave function takes the form 02s = Cls 1 1(C,)er c2 2N2~(,)re2 (3.42) Again, the terms N1, and N2s are the normalization coefficients for the specified STOs. In Equation (3.42), two expansion coefficients must be determined. This requires the simultaneous solution of a set of two equations. In this case the two equations are the normalization condition for the 2s orbital ((2s i2s) = 1) and the orthogonii.lil' of the Is orbital with the 2s orbital ((1is 2s) = 0). It becomes clear at this point that the above process can be iterated over as many sorbitals as are required for an atomic basis set, be they occupied or virtual. Specifically, the ns orbital wave function is defined as '.. csNs((is)rl e('sr. (3.43) i= 1 As before, the wave function has n undetermined coefficients and therefore requires the solution of a set of n simultaneous equations. These equations take the form the normalization of the ns wave function ({(Q, Qns,) = 1) and the orthogon. lil v of the ns wave function with the other wave functions ((91.. ..) = ( '.. '.) =. ({ ._)s'..'. = 0). This forms a set of n simultaneous equations that can be used to determine a set of expansion coefficients. The construction of the porbital wave functions follows the same schema. This proposed construction method is employed through a template designed for any commercial computational package, such as Maple, into which the orbital exponents are input. The program then calculates the expansion coefficients by solving the orthonormality conditions for each set of orbitals. The power of this method is that it does not require any expensive nonlinear energy optimizations. Furthermore, it allows for a general construction of any atomic orbital, either occu pied or virtual, provided that the orbital exponent is known or can be interpolated using the above mentioned method. Wave functions have been built for the atoms from He through Ne. Table 3.2 presents the STO exponents and coefficients for He, Li, and Be. These exponents and coefficients were determined using the previously described method. Table 3.3 lists the wave function parameters for B, C, and N. Lastly, Table 3.4 contains the exponents and coefficients for 0, F, and Ne. As a final step, the STO wave functions must then be expanded in a basis GTO functions to allow for computations to be performed. There are two main methods by which this is accomplished. The first method is through use of a linear leastsquares fitting program which will determine the best set of GTO exponents and expansion coefficients. The program that has been ei'plvid for some of the results reported in this work used a variation of the Amoeba program from Numerical Recipes [37]. In this method the GTO orbital exponents were Table 3.2: Slater Exponents and Coefficients for He, Li, and Be. Atom: Helium Configuration: He : 1s22s02p03s03po4so4po sFunctions pFunctions Exponents Exponents is 2s 3s 4s 2p 3p 4p 1.6875 0.5698 0.3836 0.2847 0.6777 0.4185 0.2935 Coefficients Coefficients is 1.00000 0.00000 0.00000 0.00000 2s 0.29929 1.04383 0.00000 0.00000 2p 1.00000 0.00000 0.00000 3s 0.22486 0.97897 1.37212 0.00000 3p 0.75715 1.25430 0.00000 4s 0.17763 0.8.''.O 1.84082 1.72886 4p 0.58192 1.40551 1.52534 Atom: Lithium Configuration: Li : 1s22s12p03so3po4so4p0 sFunctions pFunctions Exponents Exponents Is 2s 3s 4s 2p 3p 4p 2.6906 0.6396 0.4338 0.3193 0.S..4 0.5036 0.3370 Coefficients Coefficients is 1.00000 0.00000 0.00000 0.00000 2s 0.16487 1.01350 0.00000 0.00000 2p 1.00000 0.00000 0.00000 3s 0.13471 0.90398 1.36057 0.00000 3p 0.69166 1.21589 0.00000 4s 0.04890 0.38438 1.11657 1.56430 4p 0.48801 1.20388 1.42094 Atom: Beryllium Configuration: Be : 1s22s22p03so3po4so4p0 sFunctions pFunctions Exponents Exponents is 2s 3s 4s 2p 3p 4p 3.6848 0.9560 0.4841 0.3540 1.0330 0.5888 0.3805 Coefficients Coefficients is 1.00000 0.00000 0.00000 0.00000 2s 0.18228 1.01832 0.00000 0.00000 2p 1.00000 0.00000 0.00000 3s 0.09188 0.54041 1.13214 0.00000 3p 0.65076 1.19310 0.00000 4s 0.03755 0.23713 0.87282 1.46603 4p 0.43111 1.08063 1.35846 Table 3.3: Slater Exponents and Coefficients for B, C, and N. Atom: Boron Configuration: B : 1s22s22p13so3po4so4p sFunctions pFunctions Exponents Exponents is 2s 3s 4s 2p 3p 4p 4.6795 1.2881 0.5343 0.3886 1.2107 0.6740 0.4241 Coefficients Coefficients is 1.00000 0.00000 0.00000 0.00000 2s 0.21294 1.02242 0.00000 0.00000 2p 1.00000 0.00000 0.00000 3s 0.07096 0.37100 1.06382 0.00000 3p 0.62253 1.17794 0.00000 4s 0.03208 0.17493 0.7.1 1.42357 4p 0.39304 0.99801 1.31754 Atom: Carbon Configuration: C : 1s22s22p23so3po4so4p sFunctions pFunctions Exponents Exponents Is 2s 3s 4s 2p 3p 4p 5.6727 1.6083 0.5846 0.4233 1.5679 0.7591 0.4676 Coefficients Coefficients is 1.00000 0.00000 0.00000 0.00000 2s 0.22393 1.02477 0.00000 0.00000 2p 1.00000 0.00000 0.00000 3s 0.05798 0.28600 1.i 's ,, 0.00000 3p 0.47',, 1.10861 0.00000 4s 0.02788 0.14169 0.75459 1.40332 4p 0.29066 O.,.',) 1.27113 Atom: Nitrogen Configuration: N : 1s22s22p33so3po4so4p sFunctions pFunctions Exponents Exponents is 2s 3s 4s 2p 3p 4p 6.6651 1.9237 0.6348 0.4580 1.9170 0.8443 0.5111 Coefficients Coefficients is 1.00000 0.00000 0.00000 0.00000 2s 0.23080 1.02629 0.00000 0.00000 2p 1.00000 0.00000 0.00000 3s 0.04925 0.23459 1.02580 0.00000 3p 0.39933 1.07678 0.00000 4s O.'2.1 0.14023 0.90396 1.52899 4p 0.23578 0.79242 1.24296 Table 3.4: Slater Exponents and Coefficients for 0, F, and Ne. Atom: Configuration: Oxygen 0 : s22s22p43s03p04s04po sFunctions pFunctions Exponents Exponents is 2s 3s 4s 2p 3p 4p 7.6579 2.2458 0.6851 0.4926 2.2266 0.9294 0.5546 Coefficients Coefficients is 1.00000 0.00000 0.00000 0.00000 2s 0.23710 1.02772 0.00000 0.00000 2p 1.00000 0.00000 0.00000 3s 0.04298 0.19877 1.01854 0.00000 3p 0.35976 1.06274 0.00000 4s 0.01173 0.05606 0.46221 1.27638 4p 0.20750 0.75025 1.22490 Atom: Fluorine Configuration: F : 1s22s22p53s03p04s4p0 sFunctions pFunctions Exponents Exponents Is 2s 3s 4s 2p 3p 4p 8.6501 2.5638 0.7353 0.5273 2.5500 1.0146 0.5981 Coefficients Coefficients is 1.00000 0.00000 0.00000 0.00000 2s 0.24136 1.02871 0.00000 0.00000 2p 1.00000 0.00000 0.00000 3s 0.03835 0.17388 1.01419 0.00000 3p 0.32735 1.05222 0.00000 4s 0.02018 0.09298 0.71660 1.37637 4p 0.18516 0.71747 1.21094 Atom: Neon Configuration: Ne : ls22s22p63s03po4s04p0 sFunctions pFunctions Exponents Exponents is 2s 3s 4s 2p 3p 4p 9.6421 2.8792 0.77 0.5619 2.8792 1.0997 0.6416 Coefficients Coefficients is 1.00000 0.00000 0.00000 0.00000 2s 0.24439 1.02943 0.00000 0.00000 2p 1.00000 0.00000 0.00000 3s 0.03483 0.15570 1.01138 0.00000 3p 0.30169 1.04452 0.00000 4s 0.01858 0.08429 0.71010 1.37115 4p 0.16788 0.69204 1.20000 constrained to be eventempered. Depending on the use of the basis set, any number of GTO functions (up to ten) were used for each STO fit. While the first fitting method employs a fairly simple linear optimization to determine exponents and coefficients for GTO expansions of STO orbitals, the second method does not require any such optimizations. This method uses the results of Stewart [41]. Stewart performed least squares calculations to determine the best set of parameters that allows for the construction of an expansion of up to six GTO functions for a specified single STO function. As previously mentioned, this method does not require tedious optimizations, as they have already been performed. This allows for a simple template to be constructed in a numerical spreadsheet program. The orbital exponents are input into this spreadsheet template, along with the required parameters from Stewart's paper, and the GTO orbital exponents and expansion coefficients are the resulting output. While the method for basis set construction outlined in this section is very simplistic in its application, the ]li',i .,l underpinning of the construction is quite strong and maintained throughout the method. Furthermore, this formalism is easy to inrplv, with no costly or timeconsuming optimizations required if engineered correctly. In the following section, results that have been obtained with basis sets constructed using this method are presented for comparison with several common stock basis sets that were built using energy optimization methods. 3.4 Comparative Results This section will provide comparisons between basis sets constructed using the method proposed in the previous section and with some of the more commonly used stock basis sets, such as the 321G, the 631G, and the 631G* basis sets. The comparisons will be made such that a computation performed with, for example, a 631G basis set will be compared with a newly constructed basis set with the same size parameters (same number of expansions per orbital). Comparisons will be made for several types of ]li, i, .,1 properties, specifically for the excitation energies within the specific atom, for charge transfer probabilities and cross sections, and for the vibrational and electronic properties of several diatomic and triatomic molecules. These first two comparisons are for properties that are important to dynamical processes. The third comparison relates properties that are traditionally calculated based on energy (such as diagonalization of the Hessian to obtain vibrational frequencies). This will offer comparison between energy optimized basis sets with the basis sets constructed from the current method (with no energy optimizations). 3.4.1 Atomic Energetics In this section, the relative accuracy of electronic excitations within atoms will be compared. Six types of basis sets will be used in this comparison, including three stock basis sets (the 321G, 631G, and 631G** basis sets) and three comparable basis sets built using the previously proposed method (here called the 321B, 631B, and 631B** basis sets). The new basis sets have been constructed such that the same number of primitive Gaussians are used for each orbital contraction. This will help to ensure that meaningful comparisons can be made between the energy optimized stock basis sets and the newly constructed basis sets. The structure of the new, dynamically consistent basis sets can be found in the basis set library located in the Appendix. The data in this section are comprised of excitation energies as calculated using the Gaussian 98 computational suite [57]. The excitation energies are calculated as the absolute energy difference between the two states in question. The absolute energies are calculated using the multiconfigurational capabilities of Gaussian 98, as found within the CASSCF routine [44]. The complete active space is defined to be the set of valence shell electrons and the valence shell and all Table 3.5: Atomic energies and electronic excitations in Helium Stock Basis Sets Atomic HF Energy (Experimental: 2.9031840 a.u.) 321G 631G 631G** (a.u.) (a.u.) (a.u.) 2.8505767 2.8701621 2.8873650 Excitation Energies Expt. 321G % 631G % 631G** % Transition Excitation Excitation Error Excitation Error Excitation Error (cm1) (cm1) (cm1) (cm1) Is2(1S)+ s2s(3S) 159843.3 442932.5 177 322814.3 102 326589.9 104 1s2(1S) ls2s(1S) 166264.7 557807.8 235 421708.2 154 425452.4 156 1s2(1S) ls2p(3P) 169074.1 N/A N/A 483411.3 183 12(1S) ls2p(1P) 171122.2 N/A N/A 558949.8 227 Dynamically Consistent Basis Sets Atomic HF Energy (Experimental: 2.9031840 a.u.) 321B 631B 631B** (a.u.) (a.u.) (a.u.) 2.7123567 2.8118598 2.8119307 Excitation Energies Expt. 321B % 631B % 631B** % Transition Excitation Excitation Error Excitation Error Excitation Error (cm1) (cm1) (cm1) (cm1) Is2(1S) ls2s(3S) 159843.3 146507.5 8.34 155196.7 2.91 155212.3 2.90 1s2(1S) ls2s(1S) 166264.7 155066.0 6.74 162821.9 2.07 162778.0 2.10 1s2(1S) ls2p(3P) 169074.1 N/A N/A 169004.0 0.04 1s2(1S) ls2p(1P) 171122.2 N/A N/A 171592.9 0.28 Table 3.6: Atomic energies and electronic excitations in Lithium Stock Basis Sets Atomic HF Energy (Experimental: 7.8848995 a.u.) 321G 631G 631G** (a.u.) (a.u.) (a.u.) 7.3815132 7.4312358 7.4313723 Excitation Energies Expt. 321G % 631G % 631G** % Transition Excitation Excitation Error Excitation Error Excitation Error (cm1) (cm1) (cm1) (cm1) 2s'(2S)+ 2pl(2P) 14903.8 14620.3 1.90 15626.9 4.85 15647.1 4.99 2s'(2S)+ 3sl (2S) 27205.8 43415.8 59.6 46379.9 70.5 48098.6 76.8 2s1(2S) 3p (2P) 30925.9 43662.6 41.2 49549.6 60.2 46408.9 50.1 2s1(2S) 3d1(2D) 31283.2 N/A N/A 112671.0 260 Dynamically Consistent Basis Sets Atomic HF Energy (Experimental: 7.8848995 a.u.) 321B 631B 631B** (a.u.) (a.u.) (a.u.) 7.3252176 7.4161167 7.4161167 Excitation Energies Expt. 321B % 631B % 631B** % Transition Excitation Excitation Error Excitation Error Excitation Error (cm1) (cm1) (cm1) (cm1) 2s' (2S) 2pl(2P) 14903.8 14810.7 0.62 15678.5 5.20 15678.5 5.20 2s (2S)+ 3 s (2S) 27205.8 26096.7 4.08 26795.2 1.51 26795.2 1.51 2s1(2S) 3p (2P) 30925.9 43823.9 41.7 44431.2 43.7 44431.2 43.7 2s1(2S)+ 3d (2D) 31283.2 N/A N/A 35503.5 13.5 virtual orbitals in the given basis set. By using the CASSCF methodology, essential correlation can be accounted for in the basis set comparisons. The data in Table 3.5 demonstrates the excitations for the helium atom. The top part of the table relates the experimental ground state atomic energy and the excitation energies using the stock basis sets. The bottom half of the table shows the same values calculated in the newly constructed and dynamically consistent basis sets. It becomes clear immediately upon inspection of the data that the stock basis sets provide a much more accurate ground state atomic energy. However, this is to be expected, as the stock basis sets were constructed in such a way as to minimize the HartreeFock atomic energy and the dynamically consistent basis sets were not. It should be noted that, even in the worst case, the dynamically consistent basis sets are no more than 6.57 percent in error with the experimental energy, taken from Moore [58]. The most striking feature of the table is the excitation energy comparisons. The experimental excitation energies are taken from Bacher and Goudsmit [59]. The stock basis sets do not represent the orbitaltoorbital excitation energies well. One's attention can first be directed to the percent errors of the calculated excitation energies (using stock basis sets) with respect to the experimental values. The smallest error is about 100 percent. Even more important is the fact that all of the excitations correspond to virtual orbitals that are bound states in the He atom. The excitation energies are all greater in magnitude than the ionization energy of the atom (about 197000 cm1) [59]. This fact is denoted by the red coloration of the excitation energy values. By contrast, the dynamically consistent basis sets constructed in this work show a remarkable improvement in the excitation energetic. The percent errors are reduced from a minimum of 102 percent in the stock basis sets to a maximum of 8.34 percent in the newly constructed basis sets. Moreover, all of the newly constructed orbitals represent bound states in the atom. Table 3.6 relates the same data, though for the lithium atom. The experimen tal atomic energy value is calculated in this case by considering the correlation energy obtained by Eggarter and Eggarter using secondorder perturbation meth ods [60], and the excitation energies are from Bacher and Goudsmit [59]. The atomic energies demonstrate the same pattern as the helium atom: the stock basis sets offer good representations of the atomic energy, while the newly constructed basis sets are not as good (though now the largest percent error is only 7.1 per cent). The excitation energy data offers a considerably different comparison than in the case of helium. Particularly one finds that in this case the stock basis sets do allow for a good representation of the 2p virtual orbital excitation which are, in fact, better than the orbitals arising from the new basis sets. However, the 3s, 3p, and 3d orbitals are all unbound. The dynamically consistent basis sets all provide good 2p, 3s, and 3d virtual orbitals, yet perform very poorly in the description of the 3p orbitals. This can be attributed to either of the principal assumptions made in the interpolation method (c.f. Section 3.3.3). 3.4.2 Charge Transfer Results From the preceding data, it becomes clear that a basis set can be constructed that has a relatively small size and 1.li,i .,lly meaningful excitations into virtual orbitals. While this is a ].li,i .,1 characteristic that is very important in many dynamical properties, one must specifically investigate the behavior of the basis sets when calculating such dynamical properties. In this section, charge transfer probabilities are calculated using the 631G stock basis set and several basis sets constructed using the newly proposed methodology. The contraction schemes for each of the basis sets are given in the Appendix. Three charge transfer systems will be investigated. The first is the near resonant charge transfer in the H+/Li collision system. In this process, the n = 2 orbitals are in nearresonance with the 2p orbitals in Li (with orbital energies of 0.125 a.u. and 0.130 a.u., respectively) [61]. Likewise, the n = 3 orbitals of H are in nearresonance with the 3p orbitals in Li (with orbital energies of 0.055 a.u. and 0.057 a.u., respectively) [61]. These nearresonances should result in a few regions of large probability for transfer, as the electron can be excited into one of the virtual orbitals in Li and then transfer over to the H atom. The second v1 inII is the resonant charge transfer in the Li+/Li collision system. In this case, the orbitals in both collision species have the same energies, promoting strong resonances in the charge transfer process. Finally, the resonant transfer between He+ and He is investigated. The minimal END formalism has been applied successfully to investigations of resonant charge transfer processes, particularly the collision of H+/H [62]. The same methods are employed in this investigation. Specifically, the probability for electron transfer is calculated to be the difference in Mulliken population between the incident projectile and the fastest particle after collision, the differential cross section is calculated using the distinguishable or identical particle scattering amplitudes (as the case warrants), the scattering amplitude is calculated from the small angle Schiff Approximation, and the total cross section for resonant transfer is calculated using the semiclassical formula oRT = bP(b)db, (3.44) where b is the impact parameter and P(b) is the transfer probability. All of the calculations in this section were performed using ENDyne, version 5 [63]. The first comparison to be made is the transfer probability. Figure 3.8 demonstrates the nearresonant charge transfer probability in the H+/Li collision 