Semiempirical Molecular Orbital Methods Applied to Polymeric Systems

by

Paul M. Lahti

Lederle Graduate Research Tower
Department of Chemistry
University of Massachusetts
Box 34510, Amherst, MA 01003-4510

Introduction


Scope of this article: Both industrially and academically, the field of polymer chemistry has become almost intimidatingly large. The ability to "fine-tune" desirable properties of known polymers, or to select structures appropriate to give new and novel properties, still rests in great part upon the efforts of teams of scientists working to synthesize large arrays of new molecules for evaluation purposes. Any rational strategy which limits the number of choices in an extensive structure-property relationship study would greatly reduce the cost of designing a material, by comparison to a brute force strategy based on instinct, sheer hard work, and a certain degree of luck. Small companies and small research groups cannot realistically compete at the brute force level with larger combines, and even large research teams find it desirable to reduce the effort in designing new materials.
Computational chemistry has reached a level of predictability and availability that makes it a valuable tool in the design of new molecular and polymeric materials. Far from being merely effective for production of attractive graphical presentations, modern computational chemistry allows a knowledgeable scientist to explore the effect upon properties caused by various substituent and structural changes in known materials, and in some cases even allows evaluation of completely novel materials. Like a multistep synthesis in a "wet" laboratory, a well-done computational materials study requires a specific definition of a problem to solve, a well-designed and precedented methodology to attack the problem, and an analysis of the limits of capability of the method to distinguish variations in the properties evaluated. It is thus not inapt to speak of carrying out "computational experiments".
This chapter is aimed at helping the laboratory chemist to become familiar with semiempirical molecular orbital computational-aided design (CAD) of polymers and materials. A tremendous volume of work has been carried out by many groups to develop and calibrate the various semiempirical methods discussed herein. It is not our intent to review this work in entirety, hence we refer the reader to the various citations for the background of methodological development and testing. Rather, it is our intent to describe the qualitative theoretical underpinning of semiempirical MO methods, to outline some of their strengths and weaknesses, and to demonstrate their application to some specific problems. We do not aim to give a complete review of the literature of semiempirical MO methods used to explore polymer-related problems, but to give a practical primer for a laboratory chemist who wishes to evaluate such methodology for application in real-world situations.


Brief Background for Molecular Orbital Calculations : Molecular orbital (MO) methods are based upon the Schroedinger hamiltonian expression for a multi-electron molecule (equation 1). This expression eludes exact solution, hence a variety of schemes have been made to obtain approximate solutions. For the hamiltonian H , a set of wavefunctions exists that gives discrete energy solutions E for the molecular system. This is a classic eigenvector-eigenvalue problem, where the MO wavefunction eigenvectors correspond to the MO energy eigenvalues.


The form of the wavefunction varies with the level of approximation used. It is very common for the linear combination of atomic orbitals (LCAO) approximation to be used, where all are made by combinations of AOs from the constituent atoms of the molecule. The set of AOs used to make up the MOs is called the basis set. Linear combinations of the AOs give a number of MOs equal to the number of basis set orbitals, where the MO eigenvectors form an orthonormal set according to the equations (2)-(4).



Basis Sets and Their Effect on Predictive Computational Capability : A wide variety of AO basis sets have been used. For ab initio MO computations, the minimal level of basis set (termed single-zeta) uses both core and valence AOs. For hydrogen, the single-zeta basis set is the 1s orbital, for first row elements it is the 1s/2s/2px,y,z set of orbitals. For ease of computational integration, almost all modern ab initio computations approximate AOs as summations of gaussian type orbital (GTO) functions (equation 5). For higher level work, complex basis sets have been devised, using two or more shells composed of summations of gaussian functions in order to simulate each occupied shell of an atom (and often even the higher-lying empty shells). The reader is referred elsewhere for further descriptions of ab initio basis sets.[1]



The rationale behind using complex basis sets for ab initio computations is that any approximate set of MO eigenvectors i will yield a molecular energy that lies above the "true" energy. This is due to the variational theorem, which states that Eapprox > Eactual in eqn 1 for nonapproximate hamiltonian expressions H (in this case the nonrelativistic, time-independent hamiltonian is appropriate). The greater the flexibility of the basis set, the greater the flexibility in the approximate MOs i, and the closer Eapprox will come to Eactual. The cost for this greater level of accuracy is an increase in the time required to run a computation, and an increased complexity in interpreting the final result. These time constraints can be very substantial for either medium to large molecules, or for large basis set computations. Therefore, ab initio theory is practically usable only for certain types of problems in materials chemistry, even with the present state of the art of fast programming algorithms and ever-faster computers to run them. A previous chapter in this book describes the application of ab initio theory to polymer and materials chemistry problems.

For semiempirical MO calculations, truncation of the AO basis set is typical. Usually a valence shell basis set is used: 1s for the H atom, 2s/2px,y,z for the second row elements. In the Pariser-Parr-Pople method[2-4] to be described later, a -only basis set is used, and all saturated sites are ignored as minor perturbations. In many cases, Slater-type orbitals (STOs) are used in the basis set, with the general form given by equation 6.[5] STOs are not conveniently integrated in analytical form, but are more similar in behavior to true AOs than are the more-easily integrated gaussian functions. Since semiempirical methods use other approximations that speed computation, the use of STOs does not appreciable slow a semiempirical calculation.




Any of the basis set approximations described up to this point may be used with an ab initio quality hamiltonian in the Schroedinger equation. A major portion of the computational time spent in ab initio work is the evaluation of the many in the hamiltonian. We will now outline the various hamiltonian terms, in order to understand the approximations that are made in semiempirical MO methodology.


Choosing Approximations to Compute Properties for the Multielectron Molecule : Although the Schroedinger equation cannot be solved analytically for a multielectron molecule, the interactions in H can be separated as shown in equation 7.


Under the Born-Oppenheimer approximation, we can separate electronic and nuclear terms to a very high degree of accuracy, due to the small mass of electrons. Enuclear-nuclear is readily obtained as a coulombic repulsion force among the stationary nuclei of known nuclear charges. One may then add the attraction of Eelectron-nuclear to obtain the zeroth order hamiltonian used by such approximate methods as Hückel theory. Up to this point, no account is taken of repulsion or correlated motion between the moving electrons. While a variety of qualitative molecular properties may be obtained at this level of theory (such as the nodal properties of MOs), the crude approximation involved renders quantitative prediction very unlikely. Until the advent of computers, this level of theory was extremely attractive, since Enuclear-nuclear and Eelectron-nuclear may be evaluated fairly easily.

The evaluation of interelectronic repulsion, Eelectron-electron, is a substantially harder problem, for reasons that are well described in texts such as the one by Karplus and Porter.[6] By appropriate manipulation of the hamiltonian H , one may obtain an expression for the eigenvalue energies E in terms of the zeroth order terms plus a correction.[7]



The final result is given in equation 8, where (Hii)core are the zeroth-order hamiltonian terms described earlier above, Jij are interelectronic coulombic interactions, and Kij is the so-called exchange energy of the system. While J is intuitively understandable, K is a term derived from the quantum theoretical requirement that electrons be interchangable. Gross molecular properties are largely determined by the zeroth-order energy, but subtle variations properties from one molecule to another are determined by the variation in J and K, hence these energy terms must be considered to get any reasonable quantitative approximation of property variation as a function of structure. We may define semiempirical MO methods as being MO computational procedures that use various means to approximate the terms in eqn 8.

Various good references exist describing how the simple treatment described above may be extended to larger polyatomic, multielectron systems.7 In most MO procedures, the Hartree-Fock (HF) approximation F to the hamiltonian expression H is used, such that determinants of one-electron functions may be set up to represent the LCAO-derived MOs. For closed-shell molecules, MOs must contain two or zero electrons (the restricted Hartree-Fock, or RHF method); for molecules with unpaired electrons, - and -spin electrons may be treated separately in one computation by the unrestricted (UHF[7]) method. We will not concern ourselves herein with the derivation or computational details of usage of the Fock operator F, save to note that the HF method is standard usage today.


Caveat About the Use of Semiempirical MO Theory : Semiempirical methods generally neglect a variety of coulombic and exchange interactions, and parameterize the remainder of interactions in a manner to produce good predicted values of experimentally known properties. Once parameterized functions are found that fit a known set of molecules, one assumes that new, unknown molecules will be well approximated by the same parameterization. Much work has gone into finding functions and parameters that allow computational modeling of a variety of different properties. Some approximations give excellent predictions of one molecular property, but poor predictions of other properties. It is incumbent upon the researcher to become aware of the strengths and weaknesses of any semiempirical method used for predictive purposes.

In addition, the variational principle does not apply to the approximate hamiltonians used by semiempirical methods. As a result, one cannot be assured of coming closer to reality by doing more complex semiempirical computations, as often can be done in ab initio work. Each semiempirical method follows a strict procedure with known approximations -- the procedure and parameters appropriate to each method must be followed rigidly, or predictive capability becomes impossible to guess. This is (and should be) highly discomfiting to anyone doing semiempirical MO computations.
Modern programs automatically incorporate the appropriate approximations and parameters, with easy-to-use menus of choices for the nonexpert researcher. It is fairly easy to use these programs, hence even easier to misuse them and obtain meaningless results. Fortunately, it is straightforward in many cases to evaluate the usefulness of a semiempirical procedure by using common sense procedures. In the following section, we outline levels of approximation, strengths, and weaknesses of some semi-empirical methods that may be applied to polymer and materials problems.


The Computational Instrument 

Preamble : The computational methods in this section are described in various publications, and have been tested versus experimental results. Several are implemented in computer programs that are available for a nominal cost through the Quantum Chemistry Program Exchange (QCPE) at Indiana University, or through the original publishing authors. A number have also been built into graphically-interfaced chemistry CAD modeling programs that are available for various computer hardware platforms. Appropriate references are given to allow the interested researcher to check the applicability of these methods to various problems, and to obtain the programs if desired.

The Pariser-Parr-Pople (PPP) Method : This method[3,4] was one of the earliest semiempirical schemes, representing a fairly modest step beyond the zeroth order approximation of methods such as Hückel theory. It is a -electron only theory, but allows variation in properties with bond lengths and heteroatom substitution. Like Hückel theory, the PPP method is best suited to qualitative comparison of MO nodal properties, but it has proven well-suited to prediction of electronic transitions. This justifies its inclusion in our discussions, despite the severe limitations of its use for quantitative property computations.

In PPP computations, all sigma and lone pair electrons are ignored. All - atomic centers are assigned parameterized atomic energies (one-center coulomb energies), depending both on the atomic number and hybridization of the atom involved; for instance, a pyrrole nitrogen has two electrons in its z-orbital, and so must be differently parameterized to a pyridine nitrogen with one -electron. Resonance integral interactions ij are parameterized for each diatomic connection. All exchange interactions are ignored, but some coulombic terms are parameterized according to various schemes. Beveridge and Hinze[2] parameterized coulombic interactions ij by using a simple empirical equation[8,9] (not a theoretically justified one). This method gives good predictions for UV-VIS spectral transitions of -conjugated molecules and radicals[10,11]. We shall exemplify its use in the case study section. Key features to note are: (1) conjugated -electrons alone are considered, (2) variation in properties with geometry is allowed, (3) parameterization is based upon specific experimental data, and approximation is fairly extreme, especially in the elimination of exchange terms.


Neglect of Differential Overlap -- CNDO and INDO : Pople and Beveridge[7] have described the Complete Neglect of Differential Overlap and Intermediate Neglect of Differential Overlap methods, CNDO and INDO, respectively. These are among the earliest full valence shell semiempirical schemes that have been used extensively. While they have been mostly supplanted by the next generation of NDO methods (MNDO, AM1, PM3), they are worth a brief description as examples of how semiempirical methods are developed. An excellent description of computational molecular orbital theory background, as well as the basis and foundations for the CNDO and INDO methods, is given in the book by Pople and Beveridge.[7]

The CNDO level of approximation starts by parameterizing zeroth order energy terms (U ii) and two-center resonance integrals ij. All exchange terms in the Fock operator are ignored, but two-center coulombic repulsions ij were parameterized by Pople as a function of overlap between atomic 2s functions. Other schemes have been used elsewhere for evaluation of ij (especially in spectral parameterizations). Core electrons are ignored, but all valence orbitals are treated. Variation of interactions as a function of interatomic distance and orbital overlap is incorporated into all the parameterized functions. Due to the lack of exchange terms, open-shell molecules with the same orbital occupancies but different spin multiplicity cannot be properly distinguished. Hence, CNDO is not the best method for calculating excited state properties.

The INDO level of approximation uses the basic terms and parameter-izations employed in CNDO, and adds a parameterized set of exchange terms between orbitals on the same atom. This level of exchange is the minimum necessary to differentiate molecular states of different spin multiplicity, a task for which INDO has proved well suited. The amount of computational time increases only by a negligible amount over CNDO. As a result, INDO and its variants have enjoyed considerable popularity since its inception.

Both the CNDO and INDO methods are available in the program CNINDO, using both closed-shell RHF and open-shell UHF methods. Several variants of CNINDO[12] are available through the QCPE for different computers. The key features to note in CNINDO are:

Bond lengths and angles are predicted quite well in both closed and open-shell molecules, but ionization potentials are substantially overestimated. The most effective use of the INDO method is for prediction of spin density distributions in radicals or other open-shell molecules. Also, CNDO and INDO have been reparameterized for use with configuration interaction to predict UV-VIS spectral transitions. We shall describe this usage in a section below.


Extension of NDO Methods -- MNDO, AM1, and PM3 : The work of Dewar, Stewart with their coworkers has resulted in a considerable improvement of the NDO semiempirical methodology. Originally, a modification of the INDO method was carried out (the so-called MINDO, or Modified INDO method[13]), but this was further evolved into the presently much-used MNDO[14] and AM1[15] methods. A new parameterization of MNDO has been introduced by Stewart, the PM3 method,[16] which has not yet been as extensively tested as the other two, but which shows much promise. We will give a superficial overview of these methods, along with some of their strengths and weaknesses.

Unlike CNDO/INDO/MINDO, MNDO/AM1/PM3 use only monatomic parameters, so far fewer individual parameters must reside in the program than when all possible diatomic terms must be parameterized. STOs, rather than GTOs, are used. One-center atomic integrals and resonance integrals are separately parameterized for -type and -type AO interactions. The main differences between MNDO, AM1, and PM3 reside in the differing methodology used to parameterize the various coulombic and exchange terms that are retained. These differences have been well summarized by Stewart.[17] MNDO and AM1 have been used sufficiently to show the superiority of AM1 for a variety of problems, such as prediction of enthalpies of formation and of transition state geometries. Like MNDO, AM1 does not treat hydrogen bonding very well in all cases.[18,19] The recent PM3 reformulation of MNDO uses a new fitting procedure that incorporates experimental data for an unusually large number of molecules into its parameterization. PM3 shows a considerable improvement in treatment of hypervalent atoms, but the method is so new that further evaluation is needed to compare it properly to AM1.

MNDO, AM1, and PM3 are presently the state of the art of semiempirical MO methods, and are capable of reproducing a wide variety of geometric and electronic molecular properties for a wide variety of molecules, to give good agreement with experiment in the majority of cases. Their implementation in the programs MOPAC[20] and AMPAC[21] allows them to be easily used for a variety of property computations, simply by specifying appropriate choices of keywords in the programs. Geometric energy minimization (optimization with or without constraint of some structural parameters), reaction coordinate following, transition state location, vibration mode prediction, population analysis, optimization and/or energy comparison of various CI corrected states, dipole and hyperpolarizability predictions: all of these and more may be carried out by using MOPAC and AMPAC, which may be obtained from the QCPE.[22]

There are failures of these methodologies for some types of problems,[17-19] but on the whole they are the methods of choice for optimization and computational exploration of large systems. The size of systems that may be treated depends upon editing a defaults file upon installation of the program, and upon the hardware used. At a fraction of the cost and CPU time expenditure of ab initio MO computations, semiempirical MO methods can give agreement with experiment that is as almost as good, given that ab initio methods must usually make assumptions of molecular conformation, lack of solvent effects, and limitation of basis set. More importantly, semiempirical methods allow the exploration of problems of such size that ab initio methods would be impossible or at least impractical. This is achieved at a the cost of a reduction in the theoretical rigor that would be expected in an ab initio computation, hence all semiempirical results should be initially viewed with caution, and in the context of related results both computational and experimental. In the case study section given below, we will describe ways to minimize the likelihood of obtaining unreasonably inaccurate results from semiempirical MO methods (or at least how to recognize when such a result has been obtained).


Configuration Interaction Methods in Semiempirical Computations : Among the most important properties predictable from MO theory are UV-VIS spectroscopic transitions. The qualitative concept of electron excitation from occupied to unoccupied MOs is common knowledge to chemists. Quantitative prediction of UV-VIS transitions is rather more difficult, because of the failure of simple HF-MO methods to give a minimally correct description of many excited states. Unlike molecular ground states, excited states are seldom well described by a single electronic configuration. This weakness can be offset by the use of configuration interaction (CI).
CI is in some ways a computational analogue to the organic chemist's use of resonance theory. In both cases, the best electronic model for a molecule lies somewhere between two or more well-defined but simplistic models. Benzene, for instance, corresponds to neither of its Kekulé cyclohexatriene resonance structures, but is electronically and structurally a mixture of these.

Likewise, most open-shell states are often best described as mixtures of well-defined electronic configurations, represented by placement of zero, one, or two electrons in the various MOs computed for the ground state of a molecule. CI is a computational algorithm to mix these basis configurations to achieve a final MO-CI description of the molecular energy of a state.
Consider the set of MOs shown in Figure 1. The ground state is well-described by the configuration GS with all electrons paired (closed-shell). One singly-excited state is well-described by wavefunction ES1 with unpaired electrons in the HOMO and LUMO (open-shell). If the HOMO-LUMO gap is small (for instance, in a highly conjugated molecule or polymer), doubly-excited configurations with two electrons placed into virtual orbitals also become energetically plausible. In the limit where the HOMO-LUMO gap approaches zero, one cannot easily describe a closed-shell singlet state, because two different configurations, GS and ES, will have similar energies. A CI computation will mix these two configurations to give two new states


1 = (GS - ES)
and
2 = (GS + ES)


positive and negative linear combinations of the configurations, where and are normalization coefficients. 1 will be lowered in energy relative to GS, and 2 will increase in energy relative to ES. The energy difference between 1 and 2 will be a reasonable approximation of the UV-VIS vertical transition energy, whereas the HOMO-LUMO gap would not.





Figure 1 : Single and double excitations from a closed-shell ground state configuration wavefunction.

A fuller treatment of how CI algorithms function is beyond the scope of this chapter. There are important considerations of point group symmetry and MO overlap that affect which configurations mix most strongly to produce a molecular state. Modern computer programs incorporate such algorithms, making automatic selection of configurations for CI. Since the total number of possible excited state configurations grows extremely rapidly as one allows excitation of more electrons into more unoccupied virtual orbitals, realistic systems must be treated by limited (or truncated) CI methods, with excitations of a limited number of electrons at any one time, within a limited selection of occupied and unoccupied MOs. Figure 1 shows how one can "freeze" some occupied (core) MOs and some unoccupied virtual orbitals such that electrons are not promoted from or into them. Unfortunately, it is possible to get unrealistic results from limited CI treatments, since configurations important to a given state may be missed.

For semiempirical methods, CI is often limited to -orbitals, or to a subset of MOs within a fairly narrow energy band of the frontier orbital region. Typically only single (CIS) or single+double (CISD) excitations are used, although the MOPAC and AMPAC programs generate all configurations within a subset of orbitals chosen by the programmer, and then (in the default setup) use the 100 configurations that are estimated to contribute most strongly to a set of CI wavefunctions. Any such limited CI computation will necessarily be approximate (though still likely to be better than a simple HF-SCF calculation). When quoting results for a CI computation, a researcher should describe the algorithm by which configurations were generated, to allow others to reproduce the results if desired . It is often useful to compare a small scale CI computation to experimental results, and then to carry out CI computations of increasing complexity, to see if UV-VIS transition energies or other properties of interest cease to change appreciably. If this occurs, one can be reasonably confident that the CI strategy is as stable and reproducible as it will get without substantial changes in methodology (basis set, inclusion of higher levels of excitation, etc .)

Semiempirical methods that use CI can be specially parameterized to reproduce UV-VIS spectral results. Beveridge and Hinze have described a PPP-CI method[2] that uses single excitations for -conjugated systems. Good agreement has been achieved between observed and experimental results for quite large conjugated neutral systems by PPP-CI, as well as for large conjugated radical[10,11] cations and anions. We shall demonstrate the use of this method in a case study given later in this chapter.

Del Bene and Jaffé[23,24
If a CI scheme is used that is different from that originally published, it is critical to double-check computed vs. experimentally known results, since the programmed parameter set is specifically meant for certain types of CI generation schemes. CNDO/S and INDO/S are examples of methods that are calibrated to predict specific properties (in this case, UV-VIS transitions and ionizatin potentials), and which are not suited to prediction of other tasks. The user should be aware of the limits of any computational tools employed .

MOPAC and AMPAC have CI algorithms that differ from those used by ZINDO, but which are easy to use with the MNDO, AM1 and PM3 methods. Both programs allow optimization of a given state as part of a CI computation. All three methods were parameterized without the inclusion of CI, hence their use with CI guarantees some double inclusion of correlation effects, which are implicitly included during parameterization. Therefore, the use of CI with MNDO, AM1, and PM3 is normally only desirable in dealing with problems where electron correlation effects are critically important, such as finding transition state energies for reaction coordinates, looking at relative state energies of open-shell molecules, or evaluating the electronic structure and geometry of excited states.
SCF-MO plus CI methods are greatly important in the study of electronic properties of molecules, despite the fact that they consume much more CPU time than simple SCF-MO methods. The researcher who must understand excited state electronic structure or predict UV-VIS spectral properties for large p-conjugated systems, will be well rewarded by investing the time to understand how to use semiempirical MO plus CI methodology for the programs that are presently available to carry out such computations.


Band Structure Methods : Band structure methods have been developed to treat periodic systems such as perfectly regular polymers, rather than finite systems. We shall not treat the mathematics of how these methods work, save to note that they give results somewhat different from those from molecular MO methods. Instead of discrete MO eigenvalues, these methods produce bands associated with particular orbital nodal properties. The bands have discernable energy widths, and may even overlap with other bands to produce regions of high electron density in the electronic structure of polymer.

Figure 2 illustrates how one may conceptualize the evolution of a set of MOs for monomeric system into a polymer band structure diagram. Ethylene has a well-defined HOMO and LUMO in its -MO set. As one forms larger polyenes, the HOMO-LUMO gap decreases, but does not tend to zero as N -> . A plot of the HOMO and LUMO energies of these systems as a function of increasing monomer number would be roughly equivalent to plotting the HOMO and LUMO bands of polyacetylene. The range of values for the HOMO would constitute the band width of the HOMO band, and likewise for the LUMO band. The energy spacing from the top of the HOMO band to the bottom of the LUMO band is the band gap Eg for the system, and is nominally equivalent to the first UV-VIS transition in a polymer. In reality, the band gap in a semiempirical computation is seldom equal to the first UV-VIS transition energy, although the two values may be correlated within set of similar molecules, as we shall see in the Case Studies section. There is no theoretical reason to expect that the LUMO energy of a ground state computation is directly correlated with the energy of an excited state, hence it is best to consider band structure Eg values as semiquantitative, until direct experimental comparisons are available as benchmark for some members of a set of computations.




Figure 2 : The progression of the HOMO-LUMO gap in polyenes into the frontier band structure in polyacetylene.



For an infinitely extended polymeric chain, solid state physics can be applied to describe the electronic band structure. A substantial amount of such work has been accomplished using fairly simple approaches, such as Hückel type approximations. The popular Su-Schrieffer-Heeger model is based on such a -electron approach.[28-30] Whangbo[31] has also described some of the ideas and work in this area in a brief review article.

The MOPAC program can generate one-dimensional oligomer input decks for a polymer unit through the MERS=n keyword, and compute band structures for polymeric molecules. MOPAC uses a single point approach to band structure computations, rather than the solid state Brillouin zone approach used by some other programs. Because of this, it is important to have a unit cell length of at least 10 Å, so the MERS=n keyword can be used to generate an oligomeric monomer for the calculation. Qualitative band structure diagrams and density-of-states (DOS) diagrams can be obtained by using the K=(m,n) keyword described later below. These diagrams are useful in identifying materials that are effective electronic conductors (or insulators), in predicting the nature of electronic transitions and of different mechanisms of conduction, and in simulation of X-ray photoelectron spectra (XPS). Uncorrected polymer band gap or ionization potential predictions by AM1 and PM3 may not agree particularly well with experimental values, but trends among structurally similar molecules are often useful.

The valence-effective hamiltonian (VEH) method has been described by Brédas, Nicolas and Durand,[32,33] and elaborated by Brédas and coworkers[34-41] into a method described by the latter[37] as being "structured to reproduce double zeta quality ab initio calculations". Brédas' formulation replaces the usual one-electron Fock operator by simulating a set of atomic potentials within polymer unit cells in a band structure computational format. The critical atomic potentials are obtained through computational model studies by ab initio methods. The result is a methodology which is not semiempirical since experimental data is not explicitly included in the atomic potentials, but which is approximate rather than ab initio, and so is included in this chapter. Input geometries for VEH computation must be obtained by other means.[41] Ionization potentials obtained by VEH require correction to give reasonable agreement with experiment,[35] hence a de facto semiempirical correction factor is used here. Remarkably, ->* type band gaps predicted by VEH are in quite good agreement with experiment for -conjugated systems, although * bands are often predicted to be anomalously low in energy. The band gap capability predictive is surprising, since no explicit account is taken of excited state natures in the computation of the VEH atomic potentials.
While less widely available and applicable than some other methods, VEH has proven remarkably effective in predicting band gaps for a variety of -conjugated polymers that are of interest as precursors to doped conducting polymers, and so is of considerable interest to researchers working in the area of highly conjugated materials. The reader is referred to the original references for further description of this method.



Case Studies

Testing a Computational Procedure : While by no means all-inclusive, the array of semiempirical methods described above may already seem confusing to the computational novice. It is a good policy to trust (but verify) when undertaking a semiempirical MO study. While the approximate nature of all these methods may seem to make trust impossible, one may gain substantial confidence in a particular computational procedure by a standard sequence of tests, a sort of "computational checklist" that should always be followed.


Estimating the Band Gap of a Conjugated Polymer : In many applications it is desirable to estimate the UV-VIS spectrum of a conjugated system as a function of structure. For example, this capability is useful in evaluating potential doped conducting polymers, light-blocking substances, and very large dyes. SCF-MO-CI methods are specifically intended to solve this sort of problem. We will demonstrate this by investigating polyacetylene.

First, a model geometry must be obtained. This could be done by using experimental model structures, by molecular mechanics, or by MO computation. While the polyacetylene structure is well-established, for instructive purposes we demonstrate how to obtain a model geometry using MOPAC. One may generate an oligomeric polyene structure (N=8-10) using any of a variety of graphically based programs that are commercially available to generate a MOPAC input deck automatically. This method is quick and easy, but does not assure perfect repeat unit symmetry. It is adequate for most quick survey purposes, but should be checked for strong asymmetry in the oligomeric units, which should not occur save possibly at the end group monomer units.
It is preferable to use the Z-matrix format to create a periodic polymer geometry for input into the MOPAC cluster computation algorithm.[42] In Z-matrix format, atomic positions are referenced to other atomic positions in terms of bond lengths, angles, and teratomic torsions. Initially, a selection must be made of three reference atoms to define an origin, an axis, and a plane. From this point, all atoms are referenced to the positions of previously defined atoms. Clark's book gives an excellent description of this format.[43] In a polymer cluster computation, we must input a unit cell repeat distance, along with the atomic positions of the monomer. In MOPAC, the repeat unit is a translation vector (symbol Tv) across which one would translate the original monomer to attach the next unit. Often use is made of "dummy atoms" (symbol XX), as imaginary reference points to define the




Figure 3 : Z-matrix description for polyacetylene.



positions of atoms. Use of a translation vector assures repeat unit symmetry in oligomers generated using the keyword MERS=n. Figure 3 shows a Z-matrix definition of polyacetylene.

A MOPAC Z-matrix input deck POLYAC.DAT for polyacetylene is given in Appendix A. We refer the reader to the MOPAC manual for details of the format. The polymer is assumed to have a planar, all trans-transoid conformation, and is optimized within those limits. Running the deck will produce both an output file POLYAC.LOG and an archive file POLYAC.ARC. The output file contains details of the computation, including interatomic distances, the computed heats of formation of intermediate structures tested along the way to the final structure, and any other special information that is requested by use of appropriate keywords (MO eigenvectors, electronic charge and spin distributions, etc.). The archive file for the MERS=6 case is shown in Appendix B. It contains a summary of the computation, along with the final geometry, heat of formation, etc. This file may be edited and used as input for another

Table I : AM1 optimized bonding parameters for trans-transoid polyacetylene.

Parameter

r C-C

r C=C r C-H C-C-C H-C-C

Translation Vector (Tv)

Hf /C2H2

AM1
value

1.346 Å

1.441 Å 1.104 Å 123.0 120.5

2.46 Å

13.6
kcal/mol

Experimental values

1.38a
1.36b Å

1.43a
1.44b Å

122a

2.43 Åc

aRef 44. bRef 45. cRef 17. dCarried out with a hexamer or large unit cell. See ref 17.


MOPAC computation. Table I shows some AM1 results for trans-transoid polyacetylene by comparison to experimental results. Hf and Tv are given per monomeric C2H2 unit, and the Hf value is computed using an extended size "monomer" unit cell that contains several C2H2 units. The reason for the extended monomer unit is described elsewhere,[42] and is the same reason that an extended polyacetylene monomer unit is used in the MOPAC band structure computation that is described below.

Once the POLYAC.ARC file is obtained, we use it with the MERS keyword to generate oligomers N=3-10 allowed by the program limits. It is convenient to convert the Z-matrix decks to Cartesian coordinates for use with other programs. This can be easily done by running each input deck through the MOPAC program with the 1SCF keyword, and editing the output file to obtain the coordinates for the MO computations that yield the UV-VIS predictions.

We demonstrate band gap prediction using the PPP MO-CI oligomer extrapolation method.[46] A more detailed theoretic treatment along related lines is given by Soos and coworkers.[47-49] In the simplest approach, we predict the UV-VIS spectral bands for oligomers N=3-10, take the longest wavelength allowed transition as being the band gap Eg, and plot Eg as a function of 1/N. The y-intercept at (1/N)=0 will correspond to the band gap of the N = polymer. A similar extrapolation method may be used to predict the ionization potential (IP) of a polymer, by plotting the INDO/S or AM1/PM3 HOMO energies for a set of oligomers vs. (1/N). For IP and Eg predictions, we must recall that solid state packing interactions may cause the extrapolation assumption to fail, since all computations are for isolated gas phase systems. In such a case, an empirical correction may yield good correlation with experiment for a particular property and methodology, such as the -1.9 eV correction applied to VEH-computed IPs.[35]

Many variants of the PPP-CI program exist, for which specific input decks must be generated. A portion of the output file from a PPP-CI program[50] using the Beveridge-Hinze[2] parameterization is shown in Appendix C, for the N=4 oligomer of polyacetylene. Up to 16 frontier MOs were allowed for single-excitations configuration generation (8 occupied, 8 virtual), and up to 256 configurations. The output includes predicted band positions and x,y,z-transition intensities. This data can be used for simulating UV-VIS spectral band positions and relative intensities in a conjugated system, using the relationship log(e) = log (OS) + 4.5, where OS = oscillator strength10. We need concentrate here only on the longest wavelength allowed transition (3.87 eV). The graph in Figure 4 shows the results for polyene oligomers N=3-10. The extrapolated polyacetylene Eg(N = ) = 2.68 eV. The experimental band gap is 1.5-2.0 eV for trans-transoid polyacetylene.[35] The more flexible PPP-CI option within the program ZINDO was also employed,[27] using Mataga integrals and the same MO subspace for configuration generation as was used in the first



Figure 4 : Prediction of Eg for polyacetylene at AM1 optimized geometry by Hinze-Beveridge PPP-CI (open squares, ref 2) and ZINDO-PPP-CI (colored squares, ref 27).

calculation: up to 65 configurations were kept in this calculation after an energy selection criterion was applied. Using the ZINDO PPP-CI algorithm, Eg(N = ) = 2.19 eV using the AM1 geometries. The difference in results reflects different parameterization methods used in the two programs.
By using experimentally modeled geometries, Eg predictions can sometimes be improved further. Table II reproduces experimental vs. Beveridge-Hinze type PPP-CI band gaps for a number of highly conjugated polymers.[46] The agreement of theory and experiment is typically within about 0.5 eV by this method, although this grows less satisfactory as the experimental band gap becomes small.


Table II : Comparison of PPP-CI and experimental band gaps (Eg) for some conjugated polymers. See ref 46.

Polymer

Eg PPP / eV

Eg expt / eV

Polymer

Eg PPP / eV

Eg expt / eV

2.2

1.9

3.0

3.0

3.6

3.4

3.1

2.7

4.6

4.9

2.5

2.2

3.3

3.0

1.8

~1.5


The procedure of extrapolating oligomer properties to N = is quite simple to use. Small variations in the input geometry ( 0.01 Å in bond length, 2-3 in bond angles) do not greatly change computed band gaps or ionization potentials by comparison to the inherent precision of semiempirical methods. Perfect repeat unit symmetry of the oligomers is not required by the oligomer extrapolation method to get an estimate of Eg, though it is much preferred for procedural rigor. Experimental geometries can often be substituted for those computed by semiempirical or molecular mechanics methods. For many conjugated systems, semiempirical MO-CI band gaps and UV-VIS predictions are in good accord with experiment. We refer the reader elsewhere for additional examples of applying this approach.[38,48,49]

A Sample Band Structure Calculation : The oligomer decks generated above can be used to produce electronic and phonon band structure diagrams, using the MOPAC cluster calculation algorithm. As mentioned earlier, one must have a sizable oligomer repeat unit length to obtain proper results. We can use the hexamer archive deck shown earlier, with Tv= 14.7 Å, modified to produce band structure output by adding the keyword K=(0.01,4). The value 0.01 is a mesh size for the band calculation, while the value 4 refers to a four-atom



Figure 5: AM1 band structure diagram and density of states for polyacetylene.

electronic monomer unit, C2H2. Figure 5 shows graphs[51] representing the qualitative band structure and density of states for polyacetylene. The figure identifies the bands as -HOMO, -LUMO, and types. Most of the bands for the virtual orbitals are not shown, although they are included in the actual computation output. The band gap predicted by this calculation is too large at 5.55 eV, based on the minimum HOMO-LUMO gap. The ionization potential from the top of the HOMO is overestimated at 7.75 eV compared to the experimental value of 4.7 eV, a tendency shared by the uncorrected VEH method.[35]

Band structure methods are fairly qualitative in nature, and are of most use to gain insight into the electronic structure for highly conjugated polymers, and for prediction of XPS spectral shapes. Compounds with large conducting bandwidths and modest sized band gaps (by comparison to other structures) are likely to be of interest as potential dopable conducting polymers and as electro-optic materials. Quantitative use of band structure methodology is less satisfactory at present, and limited to specialized applications to specific properties, as in the VEH method.


A Model Structure for a Neutral Soliton : Our final test case will be to model the neutral soliton, a spin-bearing defect that has been of interest as a potential electron conduction carrier in polyacetylene. We chose a (=C-H)15 oligomeric unit, since the soliton has been previously predicted to be of approximately this size,[52] and since an odd number of conjugated carbons is needed to induce a neutral soliton in polyacetylene. Since the soliton is a radical, we optimze it with AM1 in MOPAC using either the UHF DOUBLET keywords, or the MECI MS=0.5 C.I.=(5,1) OPEN(1,1) keyword sequence. Appendix D gives the archive file SOLITON.ARC for the the UHF and MECI calculations, and Figure 6 shows the geometric results for the MECI computation. No symmetry contraints were placed on these calculations, except both were set to be planar, all trans-transoid starting points. The PRECISE keyword was also used, which requires MOPAC to




Figure 6 : MECI(UHF) optimized bonding lengths in ångstroms for (=C=H)15. Lower portion shows a contour diagram of the HOMO of the system.



use a more selective criterion than is the default for deciding when the geometry optimization is sufficiently converged to stop. The use of PRECISE is highly recommended[17] for systems where interpreting small variations in geometric parameters would be of importance.

The lower part of Figure 6 shows a contour diagram of the HOMO of (C-H)15, prepared using the MOPAC companion programs DENSITY[53] and DRAW.[54] The contour is a two-dimensional slice obtained at 0.005 Å above the chain. Since the HOMO contains the unpaired electron, this diagram is a useful approximation of the unpaired -spin distribution in the soliton. More sophisticated displays of charge and spin density, as well as MO contours, may be obtained from a variety of commercially available programs that incorporate or interface with the MOPAC program or the AM1/PM3 MO methods.

Both the UHF and MECI computations yield structures in which central portion of the chain shows C=C delocalized bond lengths. Delocalization decreases in a symmetrical fashion as we move toward either end of the chain from the middle. This behavior is consistent with that predicted in other computational studies,[55] which show that the neutral soliton will have a finite chain length in a polyacetylene strand. ENDOR has suggested a soliton delocalization length rather larger[30] than computation suggests, hence the comparison between experiment and theory is still under study.
The use of AM1-CI and other open-shell and CI methods to study spin-bearing sites in polymers and model oligomers is an area of substantial research activity due to recent interest in conducting and ferromagnetic materials. As a result, such computations are of much value for qualitative and semiquantitative predictions, such as descriptions of bond delocalization given above, or the energetic spacing of states with different spin multiplicity, or the effect of geometric distortions upon interactions between unpaired electrons.[56]


Conclusion

Polymers present special problems in computational modeling. They are extended molecular systems with properties that may arise from secondary and tertiary structure effects, as well as from the primary microstructure. Their size and structural complexity render them difficult to treat, even without worrying about interchain packing effects. Molecular mechanics is well suited to study geometric and structural effects that are due to steric interactions, but give no insight into electronic structures. Ab initio molecular orbital computations allow considerable rigor and a high degree of trustworthiness, but only for the limited size systems that may be so treated. Semiempirical molecular orbital computations are an excellent compromise between these approaches; although yielding some theoretical rigor due to approximation, they allow the treatment of realistically large systems by electronic structure methods. J. J. P. Stewart has aptly commented that "semiempirical methods are not particularly good at any one thing (but they) have become extremely versatile". Because of their availability, their widespread use, and their demonstrated ability to reproduce experimental properties to reasonable levels of accuracy, semiempirical molecular orbital methods are and should remain for some time a powerful tool for the investigation of polymer electronic structural problems. Appendix E lists some of the programs that this author has found of use in such investigations, as well as the source of these programs.


References

(1) See for instance Hehre, W. J.; Radom, L.; Schleyer, P. v. R.; Pople, J. A. Ab Initio Molecular Orbital Theory , Wiley, New York, NY, 1986.

(2) Beveridge, D. L.; Hinze, J. J. Am. Chem. Soc. 1971 , 93 , 3107.

(3) Pariser, R.; Parr, R. G. J. Chem. Phys. 1953 , 21 , 466, 767.

(4) Pople, J. A. Proc. Phys. Soc., Ser. A. 1955 , 68 , 81.

(5) Slater, J. C. Phys. Rev. 1930 , 36 , 57.

(6) Karplus, M.; Porter, R. N. Atoms and Molecules: An Introduction for Student of Physical Chemistry. ; W. A. Benjamin: Menlo Park, CA, 1970, pp p. 282 ff.

(7) Pople, J. A.; Beveridge, D. A. Approximate Molecular Orbital Theory ; McGraw-Hill: New York, NY, 1970.

(8) Nishimoto, K.; Mataga, N. Z. Phys. Chem. 1957 , 12 , 335.

(9) Nishimoto, K.; Mataga, N. Z. Phys. Chem. 1957 , 13 , 140.

(10) Cársky, P.; Zaharadník, R. In Progress in Phys. Org. Chem. ; Streitweiser, A., Jr.; Taft, R. A.; Eds.; 1973, Vol. 10, 327.

(11) Cársky, P.; Zaharadník, R. Fortschr. Chem. Forsch. 1973 , 43 , 2.

(12) Dobosh, P. A. QCPE Program 141.

(13) Bingham, R. C.; Dewar, M. J. S.; Lo, D. H. J. Am. Chem. Soc. 1975 , 97 , 1294.

(14) Dewar, M. J. S.; Thiel, W. J. Am. Chem. Soc. 1977 , 99 , 4899.

(15) Dewar, M. J. S.; Zoebisch, E. F.; Healy, E. F.; Stewart, J. J. P. J. Am. Chem. Soc. 1985 , 107 , 3902.

(16) Stewart, J. J. P. J. Comput. Chem. 1989 , 10 , 209, 221.

(17) Stewart, J. J. P. J. Computed-Aided Design 1990 , 4 , 1.

(18) A review describing the use of AM1 to model hydrogen bonding is given by Dannenberg, J. J.; Evleth, E. M. Int. J. Quant. Chem ., 1992 , 44 , 869.

(19) For examples of cases where AM1 does not appear to treat hydrogen-bonding well, see the following: Dado, G. P.; Gellman, S. H. J. Am. Chem. Soc. 1992 , 114 , 3138. Baldridge, K. K.; Siegel, J. S. J. Am. Chem. Soc. 1993 , 115 , 10782.

(20) Stewart, J. J. P. QCPE Program 455.

(21) Dewar, M. J. S.; Stewart, J. J. P. QCPE Program 506.

(22) The Quantum Chemistry Program Exchange, Department of Chemistry, Indiana University, Bloomington, IN 47405.

(23) Del Bene, J.; Jaffé, H. H. J. Chem. Phys. 1968 , 48 , 1807.

(24) Del Bene, J. J. Chem. Phys. 1969 , 50 , 1126.

(25) Zerner, M. C.; Ridley, J. E. Theor. Chim. Acta 1973 , 32 , 11.

(26) Zerner, M. C.; Ridley, J. E. Theor. Chim. Acta 1976 , 42 , 233.

(27) Zerner, M. C.; et al . Program ZINDO, Florida Quantum Theory Project, University of Florida, Gainesville, FL.

(28) Su, W. P.; Schrieffer, J. R.; Heeger, A. J. Phys. Rev. Lett. 1979 , 42 , 1698.

(29) Su, W. P.; Schrieffer, J. R.; Heeger, A. J. Phys. Rev. B 1980 , 22 , 2099.

(30) Campbell, D. K.; Bishop, A. R.; Rice, M. J. In Handbook of Conducting Polymers ; Skotheim, T. A., Ed.; Marcel Dekker: New York, NY, 1986; Vol. 2; pp 937.

(31) Whangbo, M.-H. Acc. Chem. Res. 1983 , 16 , 95.

(32) André, J. M.; Burke, J. A.; Delhalle, J.; Nicolas, G.; Durand, P. Int. J. Quant. Chem. 1979 , 13, 283.

(33) Nicolas, G.; Durand, P. J. Chem. Phys. 1980 , 72 , 453.

(34) Brédas, J.-L.; Chance, R. R.; Silbey, R.; Nicolas, G.; Durand, P. J. Chem. Phys. 1981 , 75 , 255.

(35) Brédas, J.-L.; Chance, R. R.; Baughman, R. L.; Silbey, R. J. Chem. Phys. 1982 , 76 , 3673.

(36) Brédas, J.-L.; Thémans, B.; André, J. M. J. Chem. Phys. 1983 , 78 , 6137.

(37) Brédas, J.-L.; Elsenbaumer, R. L.; Chance, R. R.; Silbey, R. J. Chem. Phys. 1983 , 78 , 5656.

(38) Brédas, J.-L.; Silbey, R.; Boudreaux, D. S.; Chance, R. R. J. Am. Chem. Soc. 1983 , 105 , 6555.

(39) Brédas, J.-L.; Baughman, R. H. J. Chem. Phys. 1985 , 83 , 1316.

(40) Brédas, J.-L. In Handbook of Conducting Polymers ; Skotheim, T. A., Ed.; Marcel Dekker: New York, NY, 1986; Vol. 2; pp 860.

(41) Brédas, J.-L.; Quattrocchi, C.; Libert, J.; MacDiarmid, A. G.; Ginder, J. M.; Epstein, A. J. Phys. Rev. B 1991 , 44 , 6002.

(42) Stewart, J. J. P. New Polym. Mater. 1987 , 1 , 53.

(43) Clark, T. A Handbook of Computational Chemistry ; Wiley: New York, NY, 1985.

(44) Chien, J. C. W. Polyacetylene: Chemistry, Physics, and Material Science ; Academic Press: Orlando, FL, 1984; pp 107, 115-117.

(45) Clarke, T. C.; Scott, J. C. In Handbook of Conducting Polymers ; Skotheim, T. A., Ed.; Marcel Dekker: New York, NY, 1986; Vol. 2; pp 1128.

(46) Lahti, P. M.; Obrzut, J.; Karasz, F. E. Macromolecules 1987 , 20 , 2023.

(47) Ramasesha, S.; Soos, Z. G. J. Chem. Phys. 1983 , 80 , 3278.

(48) Soos, Z. G.; Ramasesha, S. Phys. Rev. B 1984 , 29 , 5410.

(49) Soos, Z. G.; Ramasesha, S. Mol. Cryst. Liq. Cryst. 1985 , 118 , 31.

(50) Based on a program written by F. van Cataledge and adapted for use on Silicon Graphics Indigo R4000. Similar programs are Bloor, J. E.; Gilson, B. R. QCPE 71 (SCFCIO), Janiszewski, T. QCPE 76 (POPLE PI), and Bloor, J. E.; Gilson, B. R. QCPE 77(SCFOPEN).

(51) These graphs were output with the aid of programs BANDRD.FOR and MOPDRAW.BAS written by P. M. Lahti. BANDRD.FOR scans output from MOPAC v6.0 and produces a file containing band structure and density of states data. These are then plotted to a Hewlett Packard 7475 plotter by MOPDRAW.BAS, written in Turbo Basic.

(52) Ref 44, pp 196-198 and references therein.

(53) Storch, D. M. QCPE Program 492.

(54) Storch, D. M. QCPE Program 493 (version 2.0).

(55) See the discussion in Chance, R. R.; Boudreaux, D. S.; Brédas, J.-L.; Silbey, R. In Handbook of Conducting Polymers ; Skotheim, T. A., Ed.; Marcel Dekker: New York, NY, 1986; Vol. 2; pp 825.

(56) See Lahti, P. M.; Ichimura, A. S. J. Org. Chem . 1991 , 56 , 3030 and references therein.