edited by
JAMES P. BIRK Arizona State Univenily Tempe. AZ 85281
computer series, 139 The Gaussian Programs as a Teaching Tool A Case Study on Molecular Hydrogen Calculations Brian J. Duke School of Chemistry and Earth Sciences Northern Territory University Box 40146 Casuarina, Northern Territory,Australia
This study demonstrates how modern ab initio calculations using Gaussian basis sets and the GAUSSIAN program ( I )can be used to make the historical calculations on molecular hydrogen come alive. (These calculations are discussed in many textbooks.) At the same time the students increase their understanding of these basis sets and the use of this ab initio program. We suggest that students be exposed as early as possible to ab initio quantum chemistry programs, such as GAUSSIAN, that are now readily available. A valuable part of many courses in theoretical chemistry has been the case study of molecular hydrogen in which students compare the various predictions for the total energy and the equilibthat were obtained using rium internuclear distance (Re) various methods. This popular case study is a suitable vehicle for gaining experience with these computer programs. We anticipate using this case study with senior undergraduate students in a course that emphasizes the general form of the hydrogen wave functions rather than the details of evaluating energies. I t has been used with success in an earlier form with final-year undergraduates students in Australia and first-year graduate students in America. Suitable reference texts are "Atoms and Molecules" by Karplus and Porter (2) and "Coulson's Valence" (3).Both references give a table of energies and bond distances for
Table 1. Results of Some Simple Molecular Orbital Calculations for the Hydrogen Molecule
Minimum basis MO (no scaling: z = 1.0)(Coulson, 1937)
Total Energy (H) -1.0990
The Hydrogen Molecule I t has now been over fiRy years since the first ab initio molecular orbital calculations were camed out independentlv bv Coulson (1937) and Hellman (1937) (4). In 1927. loch (5)proposed a n "iknite orbitaln b ex&&* the delo: calization of electrons res~onsiblefor conduction in metals. In 1928, in the hands of ~ u l l i k e n(6)and LennardJones (71,this concept gave rise to the fundamental idea of the linear combination of atomic orbitals (LCAO) on which the whole structure of molecular orbital theory is based. Coulson believed t h a t this approximation offered a mathematically more convenient framework than that provided by the valence-bond approach, which was then universally used. He supported this view by carrying out the series of molecnlar orbital calculations on the hydrogen molecule that are shown in Tables 1and 2. Choosine the hvdrorten molecule as our svstem of stndv allowed u s to develop a hands-on approaci to theoretic& chemistrv. Bv ~erformine a calculation that is very similar to the histohial calcnla%on, the students gain fresh insight into that earlier caIculation. At the same time, students develop an understanding ofbasis sets. They also get an opportunity to see the effect that improving the basis set can have on such calculations. ~
and 5. A hlank iine must follow the control line. Then a title line and another blank line follow. The title can he anything you wish. .Line 6. The line that contains only '0 1" specifies two things: The molecule is neutral, that is, has a charge of zero; and the state is a singlet. Lines 7and 8. The next two lines are a very simple example of a Z matrix. Atom 1is H. Atom 2 is H, which is bonded to atom 1at distance R. 'Line 9. The Z matrix is followed by a blank line. .Line 10. This line specifies that the initial value of R, to start the geometry optimization,is 1.74A. .Line 11.There is a final hlank line. To alter the method and the basis set, only the control line must be changed. (The title line may also be changed if desired.) Specifying the Basis Set A calculation comparable to Coulson's free-atom hydrogen-orbital calculation can be achieved in several ways. The easiest choice is also the one that gives the students
530
Journal of Chemical Education
Results for these basis sets are shown in Table 1.Both basis sets, but particularly the STO-6G(z = 1.0) basis set, give reasonable agreement with Coulson's result. The ST03G(z = 1.0) is simply a rather poor fit to the free-atom hydrogen orbital. Scaling Orbitals by a Direct Method These two basis sets can be scaled most easily by a simple direct method. Expecting that the optimumz value will be close to the Coulson result near 1.19, we repeat the above calculations with the scale factor set to 1.18, 1.19, and 1.20. After confirming that these results are consistent with a minimum energy in this range, we can fit the three enereies to a auadratic eauation and find the value of z that &es the minimum or the quadratic. The calculation can then be r e ~ e a t e dat this value ofz to obtain the enerm and bond lengih.
Table 3. Simple Configuration-Interaction Results
Method
Table 4. Appmachlng the Best ConfigurationInteraction Results for Hydrogen
Total Energy (M Re /A
Minimum basis Ci (no scaling: z = 1 .O) (Weinbaum,1933)
-1 .I187
0.884
GAUSSIAN Basis
H basis
Total Energy (H)
Re (A)
-
3-21 G
21
-1.1479
0.750
STO-3G CID
631G
31
-1.1517
0.746
STMG CID
6-31 1G
311
-1.1535
0.747
STWG(z=1.0)CID
631W
31(P)
-1.1652
0.739
STMG(z=1.0)CID
6-311G"
~WP)
-1.1683
0.743
6311G(d,2p)
3ll(2P)
-1 .I 708
0.742
6311 G(d3p)
311(3P)
-1.1710
0.741
6-311G(d,3pd)
311(3pd)
-1.125
0.742
-1.1 731
(at exp Re)
Minimum basis CI (scaled:z= 1.193) (Weinbaum,1933)
-1.1480
0.757
STC-3G(zOpl) CID (z= 1.194)
-1.1390
0.759
STMG(zopt)CID (z = 1.194)
-1.1476
0.757
orbital and three virtual orbitals: on,, ag,and on. This means that a split-valence basis set supports configuration-interaction calculations with double excitations (CID) between four configurations. In addition, because one of the virtual molecular orbitals has the same symmetry as the occupied orbital, calculations can also be carried out using the configuration-interaction method with single and double excitations (CISD). Thus, it is not possible to directly compare the CID or CISD results with those obtained by Weinbaum because he used a smaller CI expansion. This size difference results from the method of scaling. Weinbaum uses a basis of two scaled orbitals, that is, the GAUSSIAN programs scale with an extra contracted orbital added to the basis set on each center. The resulting extra configurations give a better result than Weinbaum's two-configuration expansion. They are considered in the next section. Appmach to the Best Configuration-lnteractlonResult Table 4 shows the result of CISD calculations using better than minimum basis sets compared against the result obtained by Hagstrom (1963). Hagstrom's calculatiou is close to the best configuration-interactioncalculation that can be obtained. To do better than this, you must include configurations for nonbouud states or--even better-include the electron repulsion term rlz directly in the wave function. Here students can see the effect of increasing the size of the CI expansion by addition of extra sets of p and d functions. Thus, adding a set of p functions to the 31 and 311 basis a t CISD lowers the energy by 0.0135 and 0.0143 H, respectively. Using the best huilbin basis set gives a result of -1.17% H, which is not that far from Hagstom's value of -1.1731 H. The addition of p and d functions at this level lowers the CISD energy by 0.0190 H. Students can compare this with the MO case in which the addition of p and d functions lowers the energy by a mere 0.0051 H. Including polarization functions causes a more pronounced improvement in the CI results than in the MO results. This is to be expected because MO calculations can not take full advantage of polarization functions. Consider, for example, the effect of placing three p functions on each hydrogen atom. In an MO calculation, two of these functions (those perpendicular to the molecular axis) will lie in the virtual space. Thus, the basis here has limited flexibility because it is supplemented by one p function, not three. This p function points along the molecular axis. In a CI calculation, however, the basis has full flexibility because all three p functions are brought into the
532
Journal of Chemical Education
Best CI without r12 terms (Hagstrom,1963)
1. The note belowTable 2 aDDlies here eauallv well. 2 All ca ma! ons are LS ng'l'reClSDformalrsm. 3. For tne larger bass sels. C ca c.lat~ons are rather expensve A cneaper n l t a r n a t , ~w o r d oa to one of me oenJroalan methods for n1rw.cn.a a .- Lnc - - ~ correlation such as MP2. MP3, or MP4. in panicular, MP2 is much cheapertiian Cl. However MP results are not variation& and these methods may be more confusing to Some groups of siudents. Tables, comparable to Table 4, for MP2, MP3, and MP4 results are available directly from the authors. ~
~
~
calculation. A similar argument can be made for d and f functions. Running the Program Before the students use the GAUSSIAN program, they should become familiar with their particular local implementation. This affects the control lines for calling up the program. However, the main data is standard. For a run to optimize the bond length with the MO method and the STO-3G basis set, the data was given earlier. To run a CISD calculation with a 6311G** basis set (311(p)for hydrogen), the folowing line is the control line. # RHF
CISD 6-311** OPT
Before asking students to do particular calculations, the instructor should give careful consideration to the time required for particular calculations and to the resources available. Although most calculations take only a minute or so on a machine such as a MicroVAX, some larger basis sets, particularly a t CISD level, take considerably longer. Conclusion The approach outlined above is an economical way to illustrate the historical calculations on the hydrogen molecule and to introduce students to the GAUSSIAN programs. At the MO level it is interesting to note that, using the GAUSSIAN programs' best built-in basis, it takes a mere 17 minutes on a MimVAX to obtain a result that is near the best HartreeFock limit calculationthat has been published. An improved description of the atomic cusp by the addition of further s-type Gaussian functions to the basis would give further improvement. Additional p and d functions have the meatest effect at the CI level where the energy is loweredby anamount some fourtimesgreaterthan at the MO level. To reach the CI limit would require the student to include many different configurations. Although clearlv our CISD calculation bas not reached the h i t , we fairly close with a relatively small basis. can After completing this case study, students and instructors should find more-advanced projects feasible. The book by Hehre et al. (9)is a fruitful source of ideas for such projects.
Literature - ~ ~ - Cited -- ~~ ~
.
~
-
~
.
I T,.ch. >I J Head Gordon. M lhrlu. O W . Fwcman. J 5: Srhlcgel, H 8 : ~ n ~ h. h b~b M ~A B~N ~ W ~C Y J hs .(ionrn~ex.c ~ ~ .~ ~ k. r ~~ eJ rF. ~~ . D J .\Vh#tpadr. H A . % ? f i r . R . M c l ~ u s . C . F : S s k ~ r ..JM a n ! n . H I.. K.4hn.l.R . 31cworl J J I,:nplo~ S . P+.. .I A c ~ . ~ ~ ~ I . A Nand oI, scralon3 Gaussian:Pittsburgh, 19m. 2. Kmplus, M.;Pater, R. MAfom ond Mokulos;W . A Benjamin: NY, 1970. 3. McWeeny, R. Coulsonk Vale-; W o r d University, 1979. 4. A full bibliography of ealeulations on the hydmpn molecule to 1980 is given in: M c k a n , A. D.; Weisa, A,: Yoshimhe, M. hu.Mod. Phys. 1960,32,211. Reference to individualeahlations (by author and year)cited in t h e k t are homthir bibliograph~ 5. Sloeh, F 2.Physik 1928,52,555. 6. Mulliken, R. S.Phya. h v . 1928,32,186,388,761:1029,33,607,730. 7. LennardJones, J. E. Pam. Fomd. Sae 1029.25.668. 8. Kauzmann, W.Quonfum Chomistry;Academie: New York, 1957:Chapter 12. 9. Hehre. W J.; Radom, L.;Schleyer, P. v. R.; Pople, J. A.Ab inirio mokulorarbitd I h r y : John Wlley and Sans:New Yo* 1986.
Teaching Molecular Modeling An Introductory Course for Chemists, Implemented at the Universite de Montreal Hermann Dugas Un vers re de Montreal Monrrbal. Oueoec. Canada H3C 3J7
General Approach in Molecular Modeling With the outcome of more sophisticated theoretical methods of calculation and the accessibilitv of com~utational resources, computational chemistry has now been accepted as a useful and powerful tool in research for industry as well as for the academic wmmunity. Molecular modeline of a molecule bv com~utereenerallv implies a graphics'presentation of t6e geometry or the c&fi&ation of the atoms of the molecule followed, if no data from crystallography are available, by the application of a theoretical method, like ab initio quantum mechanics, semiempirical quant;m mechanics or molecular mechanics, to evaluate the physical properties of the molecule under study. Computational chemistry is thus a discipline that is based on concepts of theoretical chemistry and computer science and ax such became in recent years a usefulcomplement to classical experimentation. This approach does hot intend to replace experiments but is aimed a t analyzing results more critically, at proposing a n interpretation of physiwchemical phenomena, and a t suggesting new experiments. Any simulation of a process is called modeling and the manifestation of such simulation is called a model. As t simulate some such. molecular modeline is an a t t e m ~to molecular properties in order to nndektand them better. In other terms. it correswnds to the numerical determination of structural properties, energetics and/or dynamics of a molecular svstem through a model based on the ensemble of the geometric confi&ations ohtained from a digital calculation. One of the objectives of molecular modeling is to reproduce and predict &alitatively, and possibly q&ntitatively, the physicochemical behavior of molecules and macromolecules. Molecular simulation should, for example, allow the evaluation of the difference in energy between conformers, the effect of solvent on conformations or on a reaction, the difference in energy between molecular complexes, ete. Not o d y theoretical chemists but also crystallographers, spectroscopists, inorganic chemists, organic i s t ~ chemists, biochemists and b i ~ o r ~ a n i c ~ h e mpolymer and pharmacologists have reasons to be interested in molecular modeling. The following enumeration gives some of the objectives of molecular modeling in chemishy.
Some Objectives of Molecular Modeling in Chemistry Computer graphics viauslization and desi~mof rnulecules from structural data using rhe facility of n computer Use eumouterircd srruc~urnldata banks to rdentlh, molecular systems with common points. Use empirical force fields to determine molecular properties as well as atomic distances. Correlate molecular properties with a given electronic structure. .Gain information on dynamic molecular movements and their energies. Apply Wle design of molecules by computer to molecular recognition in organic, bioorganic, and medicinal chemistry as well as material sciences. For this, different techniques of visualization, manipulation. analysis and calculation of molecular structures are known. ~ i h e m a t i c a l lone ~ , distinguishes the techniques of molecular graphics allowing a representation on a screen of the 2D or 3D structure of a molecule, its manipulation (rotations, translations, superpositions, changes in conformation, etc.) in an interactive way and the possibility to analyze its molecular coordinates (evaluation of geometric such as distances, angles, accessible surfaces, etc.). These structures can be obtained by constmction from elements (atoms. emuns. nucleic acid or peptidie residues) taken frbm a ~r&r~m'library. In some cases, they are extracted from crvstalloera~hicdata banks: CamS D rookh ) haven Protein Data bridge Structural D ~ ~ ~ C & Bank (PDB). In other cases, they are constructed from experimental results, that is from X-rays or NMR data. To pass from an initial model to a fmal one, many steps of optimization must be overcome with the use of different techniques of calculation like molecular mechanics and/or molecuiar dynamics. These methods constitute a semnd aspect of molecular modeling. They are based on the use of different software to evaluate energies, to do simulations, and to simulate molecular motions Purpose and Goal of the Course I t is with these objectives, just outlined above, that a course in molecular modeline was develo~edrecentlv a t the Universitb de Montrbal. it is an introductory co&e where the ma1 is to allow our eraduate students in chemistry, biochimistry, and pharmacology to familiarize themselves with computational chemistry. The aim is to lay down in a simple way the concepts of ionic and molecular interactions through methods of molecular simulation applied to organic and bioorganic molecules. Particular attention is given to molecular mechanics calculations that provide useful information concerning variation in potential energy as a function of variation in molecular stmcture. The focus here is mainly with organic and bioorganic molecules. I t is hoped also that through such exercise, the students will benefit from a better understanding of stereochemistry and conformational analysis. Materials for the Course and Evaluation of Students' Performance This one-semester course is given for a period of three hours ner week where the students can mani~ulatedifferent sokware on the screen. They have access to either an Evans and Sntherland PS-390, wnnected to a micro-VAX 3500, (for MacroModel) or four Silicon Graphics Personal IRIS (for MacroModel and Chiron). Four IBM-PC (for PCModel and REACCS) and two Macintosh MAC11 (for REACCS and Chemprotekt) are also available. The wurse is given mainly to graduate students in chemistry but students from other departments also are welcome. In addition, it is offered as an optional course to fmal-year undergraduate students in chemistry with a particularly good background in organic chemistry. Relevant textbooks on Volume 69 Number 7 July 1992
533
