From:	<BETTENS \\at// MPS.OHIO-STATE.EDU>
Date:	Fri, 16 Dec 1994 14:10:00 -0500 (EST)
Subject:	Spin contamination & AM1 "ROHF" versus UHF

Dear Netters,

I posted a number of questions to the CCL about spin contamination
on the 14th November.  Below is the original questions and a summary
of the responses.  I also posted on the 7th of December questions
relating to the 5th question below.  The responses to this posting
are also summarized here.

======================================================================
                         14th November Posting
======================================================================

I have a number of questions regarding the effects of spin
contamination on total electronic energies and structures.  My
understanding of spin contamination is that unrestricted Hartree-Fock
(UHF) wave functions are not eigenfunctions of the total spin
operator, so the electronic wave function of interest can be
contaminated by functions corresponding to states of higher spin
multiplicity.  This brings me to several of my questions:

1.  Given that, (a) we have performed an ab initio study on an open
shell molecule using a single-determinant wave function (i.e.,
variational), (b) we have found its theoretical equilibrium geometry
for the ground electronic state, and (c) have projected out ALL
contaminating higher spin multiplicity states along the way to the
minimum energy geometric configuration.  Will the total calculated
electronic energy be the lowest possible calculated energy for the
given basis set and Hamiltonian?

2.  What can be said about 1, above, regarding the total calculated
electronic energy when a perturbation treatment to the configuration
interaction is introduced, e.g., MP4(SDTQ)?

3.  What can be said regarding the total calculated electronic energy,
if in the case of 2, above, condition (1c) is not fully met, i.e.,
some spin contamination remains while going to the optimum geometry.

4.  Given that different spin states correspond to different
equilibrium geometries, is an optimized structure, which had
significant spin contamination in its electronic wave function all the
way down to its "minimum" energy geometrical configuration, some kind
of mixture of structures involving the state of interest and the
different contaminating higher spin multiplicity states?

5.  Regarding semiempirical calculations, in the paper of Novoa et
al., Inorganica Chemica Acta, Vol. 198-200 (1992) 133, the heats of
formation of some very large carbon clusters were calculated.  In
their paper the authors state:  "For an odd-membered linear C_n with
13 <= n < 20, the AM1 calculations predict the triplet state to be
more stable than the singlet, due mainly to the spin contamination of
the UHF calculations."  The calculated stabilities of these larger
linear carbon clusters are not expected based on what is known for the
smaller linear clusters where, for odd n (n > 1), the singlet states
are more stable than the triplet states.  My question is this, is it
possible that the authors are not correct regarding their reason for
the for greater stabilities of the triplet states?  (It is not my
intention to attack the above authors work.  I merely wish to evaluate
the quoted heats of formation because I require some kind of estimate
for the heats of formation for these species.)

Regards,

Ryan Bettens,
OSU Physics Department,
BETTENS -AatT- MPS.OHIO-STATE.edu

                       ________________________
----------------------| Responses to Questions |----------------------
                       ------------------------

From:  Christopher J. Cramer,
       University of Minnesota,
       Department of Chemistry

   Re your queries on spin contamination. I can not speak to all of
your questions, but one point I think worth making is that the UHF SCF
procedure is variational in that it produces the unique orthogonal set
of MO's that give rise to the lowest energy wavefunction for whatever
basis set has been chosen. However, those MO's derive from the UHF
calculation itself, which includes all spin contamination. Any
projection operator formalism applied thereafter can provide you some
(significant) improvement in energy because of the annihilation of
higher spin states, BUT that procedure does NOT reoptimize the MO's.
So, even if you optimize a geometry by hand using PUHF (or PMPn)
energies (since gradients for these methods are not, to my knowledge,
available) you will still suffer from potentially poor wavefunctions.

   Alternatives include ROHF (which I personally don't like because we
tend to be interested in spin-polarization, which you can't get at
this level), spin-adapted MCSCF treatments, and, something you might
consider, spin-polarized DFT. The latter is in essence a UHF-like
calculation using an approximate scheme for exchange/correlation (as
all DFT does), but it has been shown to be far less prone to spin
contamination.


             ----------------------------------------

From:  Dave Ewing,
       John Carroll University,

See warren J. Hehre, et. al., Ab initio Molecular Orbital Theory (John
wiley & Sons, 1986), pp. 203-4.

My own experience has been that structures are OK as long as the spin
contamination isn't too gross, e.g. no more than =1.0 for a
doublet.


             ----------------------------------------

From:  I. Shavitt,
       Ohio State University,
       Department of Chemistry

1.  The wave function obtained by projecting out all spin
    contamination AFTER having calculated a UHF solution
    is not the lowest-energy solution for the spin-state
    in question.  This is so because the orbitals have been
    optimized for the spin-contaminated function, before
    projection, not for the projected function.

2.  After the addition of a correlation treatment, like
    MP4(SDTQ), the effects of spin contamination would
    usually be diminished, though not eliminated.  The
    effect on the energy cannot be predicted with confidence,
    because there are diverse factors determining how the
    choice of orbitals will affect the correlation energy.
    In fact, orbitals which give the lowest Hartree-Fock
    energy are not guaranteed to produce the lowest
    correlated energy, though they usually do.  Coupled
    cluster methods, such as CCSDT, are less sensitive to
    the choice of orbitals, and therefore less likely to
    suffer from the use of nonoptimal orbitals.

3.  See 2.

4.  The effect of spin contamination on optimized geometries
    depends both on the degree of spin contamination and on
    the relative energies and characters of the higher-spin
    states.  I don't think I can give a general answer to this
    question.

5.  Personally, I am very skeptical about the ability of
    methods like AM1 to give definitive answers to questions
    concerning the relative energy of different electronic
    states.  But the same is true, even more strongly, for
    UHF calculations.  For some molecules it is very difficult
    to determine which structure or multiplicity is lower in
    energy.  An example is C_4, for which a linear triplet state
    and a cyclic singlet are almost of the same energy, and
    different levels of even high-quality calculations give
    conflicting answers.  My impression was that odd-membered
    carbon chains have a singlet ground state, and it would take
    much more than AM1 calculations to convince me otherwise.


             ----------------------------------------

From:  Ab Initio Molecular Orbital Calculations for Chemists, 2nd ed.,
       W.G. Richards and D.L. Cooper,
       Clarendon Press, Oxford, 1983.

I found this comment in the above book after posting the above
questions.  The comment essentially answers my question 4.  I quote
this book here (from page 59) as it may be of convenience to people.

     The average interaction interaction for alpha-spin electrons
     and beta-spin electrons may be different in open-shell
     systems so that it is not unreasonable that orbitals
     differing only in their spin quantum number should have
     different spatial functions.  This is the origin of the
     unrestricted Hartree-Fock (UHF) method implemented in some
     programs.  Unfortunately, UHF wavefunctions are not always
     eigenfunctions of S^2.  If for example we carried out
     calculations on a molecule with a doublet state and a
     quartet state that were close in energy, then the
     equilibrium geometry of the doublet state would be partly
     characteristic of that for the quartet state.
     
     Although there are methods for projecting out spin
     eigenfunctions at the end of the UHF calculation, many
     theoreticians prefer the restricted Hartree-Fock (RHF)
     method and this is implemented in many open-shell SCF
     programs.  In the RHF method, orbitals differing only in
     spin quantum number . . . have identical spatial parts.

======================================================================
                         7th December Posting
======================================================================

Thanks to those of you who responded to my last posting on the 14th of
November regarding spin contamination.  Before summarizing these
responses I would like to ask one further question which is related to
my last question in the previous posting.  I will then compile all
responses and post a summary.  My question concerns ROHF versus UHF
calculations using the AM1 semiempirical model.  Is anybody aware of
any published work which recommends the use of either in calculating
heats of formation, specifically for hydrocarbon open shell species?

I have examined the original paper of Dewar, M.J.S.; Zoebisch, E.G.;
Healy, E.F. and Stewart, J.J.P., J. Am. Chem. Soc., 107 (1985) 3902
where AM1 was introduced.  In this work the heats of formation of 6
hydrocarbon radicals (all doublets) where compared with the
experimental values.  I could find no indication as to whether the
reported calculations used the UHF or ROHF approaches.  This is
important because the calculated heats for each approach produces
significantly different results.  I consequently repeated the
calculations on the 6 hydrocarbon radicals given in the above paper
and found the following heats of formation (kcal mol^{-1}).  I also
give the values of S^2 and the reported heats as well as the current
experimental values (from J. Phys. Chem. Ref. Data, 17 (1988) Sup. 1)
for comparison below.

           AM1           AM1    Reported
           UHF    S^2    ROHF      AM1     Exp.     Exp. - Calc.(ROHF)
CH3       29.95  0.7613  31.25    31.25   34.8(3)       3.6
C2H3      60.46  0.8589  64.78    64.78   63.4(10)     -1.4
C2H5      15.49  0.7619  18.14    15.48   28           10
CH2CHCH2  30.20  0.9300  38.58    38.58   39            0
(CH3)2CH   3.61  0.7622   6.80    10.07   22.3(6)      15.5
(CH3)3C   -6.14  0.7623  -2.66    -2.66   11.0(6)      13.7

For C2H5 the authors seem to have reported the UHF value.  I was not
able to reproduce the reported result for (CH3)2CH with either the
ROHF or the UHF result.  I was able to reproduce the calculated heat
of formation for this species given by Dewar, M.J.S. and Thiel, W., J.
Am. Chem. Soc., 99 (1977) 4907, using the MNDO Hamiltonian.  In this
earlier work the authors explicitly stated that they used ROHF heats
of formation calculated with MNDO.  It is noted that for the remaining
species the reported AM1 heats of formation also appear to be the ROHF
values.  While agreement with the experimental values are not terrific
and appear particularly bad for the larger saturated hydrocarbon
species, they are adequate for my purposes given that these errors can
be taken as typical for calculated heats of formation of analogous but
larger open shell hydrocarbon radicals.

Does anybody know of any further studies, using the AM1 Hamiltonian,
on radicals where the calculated heats of formation have been compared
with experiment?  It seems to me, from the above table, that the UHF
results are not useful in comparison with the ROHF results, and that
the extent of spin contamination tends to increase with molecular size
and unsaturation, being close to 0.75 for the most saturated species.
For instance, I repeated the calculation on triplet linear C19 given
by Novoa et al., Inorganica Chemica Acta 198-200 (1992) 133, using UHF
and ROHF.  I obtained the same heat of formation as these authors with
no assumed symmetry in the linear chain and using the UHF AM1
Hamiltonian.  However, I obtained an S^2 of 4.356 (should be 2) and a
ROHF heat of formation of 678.72 kcal mol^{-1}, which is greater than
the UHF value by 24.17 kcal mol^{-1}!  I would very much appreciate
any helpful comments or literature references which deals with any of
the above issues I've raised here about AM1 and its applicability of
open shell systems.

Regards,

Ryan Bettens,
OSU Physics Department,
BETTENS : at : MPS.OHIO-STATE.edu

                       ________________________
----------------------| Responses to Questions |----------------------
                       ------------------------

From:  Andrew Holder,
       University of Missouri,
       Department of Chemistry

   You have indeed found some the reporting errors in the original AM1
paper.  The problem with your larger systems is that the description
of an open shell system of this size will usually require more than
one determinant.  (That is what the bad spin contamination values
indicate.)  A better way to handle this is to do everything at the
same level of configuration interaction (CI).  Be sure that you
compare everything in a reaction profile or an analogous series
computed at the same level of theory.


             ----------------------------------------

From:  John M. McKelvey,
       Eastman Kodak Company, Rochester NY
       Computational Science Laboratory

There is no proper ROHF in QCPE MOPAC or AMPAC codes!  It is a half-
electron method supplemented by a 3x3 CI.  I THINK there is proper
ROHF for AM1 and PM3 and MNDO in GAMESS.


             ----------------------------------------
             
Comment to J.M. McKelvey's response.

It seems I had too loosely used the ROHF terminology in my posting.
When I stated that Dewar and Thiel, ". . . explicitly stated that they
used ROHF heats of formation calculated with MNDO" I ignorantly
mistook the statement actually made by the authors about the method
they used, i.e., the half-electron method, for ROHF (which from what I
gather from various modern books on ab initio theory is Roothaan's
open shell method).  Therefore, in the above posting wherever ROHF
appears, read: half-electron method.

John McKelvey pointed me toward the publication:  Dewar, M.J.S. and
S. Olivella, J. Chem. Soc. Faraday Trans. II 75 (1979) 829, which
compared the calculated heats of formation and geometries for the
half-electron (HE) method with what is now meant by ROHF using
MINDO/3.  Their results showed that the ROHF heats of formation were
"systematically more negative than the HE ones".  They also found that
for the molecules studied the r.m.s. differences in the heats of
formation for open shell doublets, triplets and singlets were 2.2, 2.5
and 3.7 kcal mol^{-1} respectively.

For those interested, in the ROHF method the wavefunction is (a)
variationally optimized and (b) is an eigenfunction of the S^2
operator, unlike UHF wavefunctions.  In the HE method the wavefunction
is not variationally optimized and is the same for analogous open
shell states of different multiplicity.  A non-mathematical
description of the HE method can be found in the introduction of the
Dewar and Olivella paper, above.



Similar Messages
08/01/1995:  Spin contamination, effect on energy and structure.
12/09/1994:  Open shell hydrocarbons, AM1 ROHF versus UHF
11/14/1994:  Spin contamination, effect on energy and structure.
04/16/1997:  Re: CCL:Chem Topic: Spin Contamination
04/08/1997:  Chem Topic: Spin Contamination
08/01/1996:  Re: CCL:M:Heat of formation calculation using MOPAC.
06/15/1993:  Responses to restricted vs. unrestricted semi-empirical question
03/07/1997:  Chem Topic: Excited States
11/01/1996:  G:SUMMARY: Wavefunction Instability
12/11/1995:  Summary: Koopmans' Theorem and Neglect of Bond in CAChe MOPAC Input


Raw Message Text