Coupled-cluster based methods for excited state energies and gradients

MISSING IMAGE

Material Information

Title:
Coupled-cluster based methods for excited state energies and gradients
Physical Description:
vi, 140 leaves : ill. ; 29 cm.
Language:
English
Creator:
Gwaltney, Steven Ray, 1970-
Publication Date:

Subjects

Subjects / Keywords:
Cluster theory (Nuclear physics)   ( lcsh )
Quantum chemistry   ( lcsh )
Excited state chemistry   ( lcsh )
Chemistry thesis, Ph.D   ( lcsh )
Dissertations, Academic -- Chemistry -- UF   ( lcsh )
Genre:
bibliography   ( marcgt )
non-fiction   ( marcgt )

Notes

Thesis:
Thesis (Ph.D.)--University of Florida, 1997.
Bibliography:
Includes bibliographical references (leaves 129-138).
Statement of Responsibility:
by Steven Ray Gwaltney.
General Note:
Typescript.
General Note:
Vita.

Record Information

Source Institution:
University of Florida
Rights Management:
All applicable rights reserved by the source institution and holding location.
Resource Identifier:
aleph - 028628193
oclc - 38863556
System ID:
AA00017654:00001


This item is only available as the following downloads:


Full Text










COUPLED-CLUSTER BASED METHODS FOR
EXCITED STATE ENERGIES AND GRADIENTS










By

STEVEN RAY GWALTNEY


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


1997











ACKNOWLEDGMENTS


Soli Deo Gloria, to God alone be the glory.

I thank God for creating a universe so full of wonder and mystery and for creating

me and allowing me to explore it. I also thank my wife, Charity, for her constant support

and encouragement, even when it meant that she had to do more than her fair share. I

love you and will always remember how much you sacrificed for me and for this degree.

I want to thank my parents for never letting me settle for being less than I could be,

and I thank them for pushing me when I was growing up and then letting me go far

away when it was time.

I would like to thank my research advisor, Prof. Rodney Bartlett, from whom I have

learned so much. I would also like to thank Prof. Ernest Davidson. Who knows if I would

be doing quantum chemistry if he had not allowed me to do undergraduate research under

him, even though I did not know what I was doing. Thanks go out to all the members of

the Quantum Theory Project, to the members of Dr. Bartlett's group, and especially to

my office mates through the years. The environment of cooperation here is wonderful.

Thanks to Prof. Michael Zerner for help understanding the free base porphin problem.

Finally, special thanks go to Prof. Marcel Nooijen. I learned much in our discussions.

The work presented here has been funded by the National Science Foundation through

a Graduate Research Fellowship and by the United States Air Force Office of Scientific

Research through grant number F49620-95-1-0130 and through AASERT grant number

F49620-95-I-0421.











Computer time for the porphin calculations was provided by the CEWES Depart-

ment of Defense Major Shared Resource Center. The CEWES computer center is grate-

fully acknowledged for their help in getting these calculations done. The EOM-CCSD

porphin calculations were performed in conjunction with Dr. Howard Prichard at Cray

Research/SGI.











TABLE OF CONTENTS

ACKNOWLEDGMENTS ........... ...................... ii

ABSTRACT .... ..................................... vi

CHAPTERS

1. AN OVERVIEW OF CURRENT METHODS FOR EXCITED STATES .... 1
Coupled-Cluster M ethods ....................... ........ 1
Other Current Excited State Methods ............... ....... 3
Gradients for Excited State Methods ......................... 8

2. PARTITIONED EQUATION-OF-MOTION COUPLED-CLUSTER METHODS
FOR EXCITED STATE ENERGIES ........................... 9
Introduction .................... .................. 9
Theory .................... .. ....................13
Partitioned EOM-CC for Excited States ..... .............. 13
MBPT[2] Ground State ............................. 15
Calculations ................... ................... 17
Be Atom ................... .................. 17
Example M olecules ............................... 18
Conclusions ................... ................... 21

3. GRADIENTS FOR THE PARTITIONED EQUATION-OF-MOTION MBPT(2)
METHOD FOR EXCITED STATES .................. ........ 28
Introduction ................... ................... 28
Theory ..................... .................... 30
EOM-CC Gradients ............................... 30
P-EOM-MBPT(2) Gradients ...........................41
Application to Diatomic Molecules ................ .. .... ..46
Vertical and Adiabatic Excitation Energies .................. 48
Bond Lengths and Vibrational Frequencies .... .......... 50
Application to Polyatomic Molecules ............. .... ...... 51
The S1 State of Ammonia ............... ............ 51
Trans-bent and Vinylidenic Isomers on the S1 Surface of C2H2 53
Simple Carbonyls ................................ 54
Conclusions ................... ................... 55





iv












4. THE SIMILARITY TRANSFORMED EQUATION-OF-MOTION
COUPLED-CLUSTER METHOD .......................... 69
Introduction .... .......... ................. ....... 69


The STEOM-CC Method ...........
Discussion About STEOM-CC .

5. THE SPECTRUM OF FREE BASE PORPHIN
Introduction ...................
Computational Details ...............
STEOM-CC ................
Basis Set and Geometry .......
Results . . .
Ionized and Electron Attached States ..
Excited States ...............
Discussion ....................
Conclusions ...................


S73
. 78

. 83
. 83
. 85
. 85
. 86
S88
. 88
. 89
.91
.94


6. GRADIENTS FOR THE SIMILARITY TRANSFORMED
EQUATION-OF-MOTION COUPLED-CLUSTER METHOD .......... 110
Theory ..................... .................. 110
A applications .. ... . .. .. .. .. 120

7. A FINAL WORD ................... ............... 127

BIBLIOGRAPHY ......... ............................. 129

BIOGRAPHICAL SKETCH ....... ....................... 139














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

COUPLED-CLUSTER BASED METHODS FOR
EXCITED STATE ENERGIES AND GRADIENTS

By

Steven Ray Gwaltney

December 1997

Chairman: Rodney J. Bartlett
Major Department: Chemistry

Coupled-cluster based methods have become widely used for calculating energies and
properties of molecular ground states. In this dissertation, I present a new method, the
partitioned equation-of-motion coupled-cluster method based on a second order many-
body perturbation theory ground state (P-EOM-MBPT(2)) for calculating energies of
molecular electronic excited states. I also describe how to calculate gradients for
two excited state methods, P-EOM-MBPT(2) and the similarity transformed equation-
of-motion coupled-cluster singles and doubles (STEOM-CCSD) method. These new
methods are then applied to several example molecules.














CHAPTER 1
AN OVERVIEW OF CURRENT METHODS FOR EXCITED STATES


Coupled-Cluster Methods


Coupled-cluster theory'-14 provides a routinely applicable method for calculating

energies and properties of molecules in the lowest state of a given spin and symmetry.

By adding triple,15, 16 quadruple,17, 18 and higher excitations, it is possible to converge on

the exact answer within a given basis set. The perturbative inclusion of triple excitation

effects19-21 is now routine, meaning that energies and properties can approach near

chemical accuracy for many molecules. Unfortunately, these methods do not, in general,

work for electronic excited states.

The primary problem with describing excited states with normal coupled-cluster (CC)

theory is its inherent single reference nature. The wavefunction for the state can be

expressed as a linear combination of determinants, with the reference determinant having

a coefficient of one. Although in principle coupled-cluster theory can be made exact, the

cost of such a calculation would increase in cost at the same rate as full configuration

interaction (CI), i.e. as the size of the system factorial.22 Therefore, the set of excitations

used in the calculation are truncated. In such a truncated CC, a second determinant with a

coefficient of more than 0.2 to 0.3 can cause the accuracy of the computed results to drop

significantly. However, for an open-shell singlet excited state, the second determinant

must have a coefficient of one.

Coupled-cluster theory has been extended in order to be able to describe such

multireference states. Multireference coupled-cluster theories can be divided into two









2

main categories. They are the Hilbert space approaches23-31 and the Fock spaces

approaches.32-45 A very simplified view of Hilbert space multireference coupled-cluster

theory is that each determinant in the active space gets its own set of T amplitudes,

and a small effective Hamiltonian is diagonalized over the space of these wavefunctions.

In the Fock space coupled-cluster method, a single exponential operator is used for all

determinants in the active space, but this exponential contains operators that change the

number of electrons in the system. Hilbert space coupled-cluster theory is more suited

to studying one or a few states, while Fock space coupled-cluster theory is more suited

to studying energy differences between many states and to calculating spectra.

An alternative method for calculating excitation energies within a coupled-cluster

framework is to follow the response of the system to an external time dependent

perturbation. This was introduced by Monkhorst,46 and has been developed into the

coupled-cluster linear response (CCLR) method.47-54 The same equations can also be

derived from an equations-of-motion55 framework based on a coupled-cluster ground

state. This gives the equation-of-motion coupled-cluster (EOM-CC) method.56-59 In

EOM-CC (or CCLR) the ground state coupled-cluster operator is used to perform a

similarity transformation of the Hamiltonian, which is then diagonalized over a space

of excited determinants. EOM-CC and CCLR are identical for excited state energies,

but the approximations made by the two methods differ for oscillator strengths and

polarizabilities. Typically, the ground state operator and the space of excited determinants

is truncated to include only single and double excitations. This gives the EOM-CCSD58 59

(or CCSDLR)53 54 method. Models have also been developed which partially include

the effects of triple excitations.60-63









3

Recently, Nooijen and Bartlett have developed a new method for calculating excited

state energies and properties. It is the similarity transformed equation-of-motion coupled-

cluster (STEOM-CC) method.64-66 For singly excited states, STEOM-CC is closely

related to the Fock space coupled-cluster method,66 but conceptually they are very

different. In STEOM-CC, ionization potential EOM-CCSD67 and electron attachment

EOM-CCSD68 calculations are done. Information about the states with one fewer and

one more electron is then used to perform a second similarity transformation on the

Hamiltonian69 such that the terms which directly couple single and double excitations

are set to zero. A detailed description of STEOM-CC will be given in Chapter 4.


Other Current Excited State Methods

An attempt to catalog all ab initio methods currently used for calculating excited

states would be futile and will not be attempted. Instead, several of the methods which

are used most often will be reviewed. But first it will be useful to take a more general

view of excited states and break down how the various methods describe excited states.

In this view, the excited state energy is divided into several components. They are the

zeroth order description of the wavefunction, the orbital relaxation, and the dynamic

correlation in the excited state.

Properly describing the zeroth order description simply means being able to include

all of the determinants which play a dominant role in the excited state. For example,

this means being able to give the two determinants in an open-shell singlet the same

weight. The orbital relaxation refers to how much the orbitals are allowed to rearrange

themselves in response to the excitation.









4

The dynamic correlation for the excited state is further broken down into the dynamic

correlation for the ground state and the differential dynamic correlation between the

ground and excited states. The reason for this division of the dynamic correlation is that,

especially for the single reference theories, first a ground state is calculated and then the

excited states are built upon it. The justification for this division is that most electrons

do not take part in the excitation and therefore still interact with each other in the same

way in the excited state as in the ground state.

These definitions are not unique, especially when considering orbital relaxation versus

differential dynamic correlation. The coupled-cluster viewpoint will be taken here. The

key to the coupled-cluster viewpoint is Thouless' theorem,70 which states that any two

non-orthogonal determinants can be written as


|)i) =eT J|), (1-1)

where Ti is a single excitation operator. Therefore, any set of single excitations or any

product of single excitations will be considered orbital relaxation, and dynamic correlation

will be reserved for "connected" double and higher excitations.13

Table 1-1 provides a convenient representation for how the selected methods break

down into the four terms. A "+" means that the method treats this term adequately.

Two pluses mean that the method treats the term very well. A "-" means that the term

is included, but poorly. Finally, a "0" means that the term is completely neglected.

The question mark means that how well the term in question is treated varies based on

specifics of the calculation.

The first method to be considered is probably the simplest of all ab initio excited state

methods. It is the Tamm-Dancoff approximation,71, 72 also known as the configuration









5
Table 1-1: An overview of excited state methods.

h o r o ground state differential
zeroth order orbital da
dynamic dynamic
description relaxation correlation
correlation correlation
CIS + 0 0
RPA + + 0 0
MCSCF ++ (?) ++
CASPT2 ++ (?) ++ + +
CIS(D) + + +
EOM-CCSD ++ + ++ +
STEOM-
SD + + ++ ++ (?)
CCSD
MRCI ++ (?) ++ (?) ++ (?) ++ (?)

interaction singles (CIS)73-75 method. The names will be used interchangeably. In

CIS the excited state is determined from diagonalizing a Hamiltonian matrix made up

of single excitations with respect to the Hartree-Fock reference determinant. Since, by

Brillouin's theorem,22 single excitations cannot couple with the Hartree-Fock determinant,

this diagonalization leaves the ground state unchanged and gives excited states. Because

all possible single excitations are included, any excited state which corresponds to a single

electron being promoted can be described. However, the orbital relaxation is extremely

limited and the method contains no dynamic correlation.

The next method is the random phase approximation (RPA).76, 77 Although RPA can

be derived several ways, the one most useful here is to consider RPA as the response

of a wavefunction written as12, 13 1 ) = (1 + Ti + T2) TI2o). Clearly, then, RPA also

contains no dynamic correlation, but it does contain more orbital relaxation than CIS.

The next method to be considered is multiconfigurational SCF (MCSCF).78, 79 In

MCSCF, a set of configurations is chosen. Then the coefficient for each configuration









6

is variationally optimized while simultaneously optimizing the orbitals from which the

configurations are built. Since the procedure involves variationally choosing the optimal

orbitals, orbital relaxation can be treated quite well. Because of the flexibility to choose

any configurations wanted, any excited state, no matter its excitation level, can be

described. Thus the zeroth order description can also be very good if the configurations

are chosen properly. Because of the multideterminantal nature of the final wavefunction,

some dynamic correlation is recovered, but normally the amount is small. A special way

of choosing the configurations is to choose a set of orbitals and to create all possible

configurations with a given number of electrons distributed among that set of orbitals.

This gives rise to the complete active space SCF (CASSCF)80 variant of MCSCF.

The CASSCF wavefunction can be used as the zeroth order wavefunction in a

multireference perturbation theory expansion. The most popular of the methods derived

this way is the CASPT2 method.81, 82 Since it is based on a CASSCF wavefunction, the

arguments above about the zeroth order description and the orbital relaxation still apply.

The perturbation theory then allows a treatment of the dynamic correlation. Overall,

this method works well, but because of issues involving properly choosing the active

space and a potential intruder state problem, it is sometimes difficult to get consistently

good results.

Starting with the CIS wavefunction being used as the zeroth order wavefunction in a

perturbation expansion, it is possible to derive the CIS(D) method.83 CIS(D) shares CIS's

weakness in term of orbital relaxation, but it does allow for a second order perturbation

description of both the ground state dynamic correlation and the differential dynamic

correlation.









7

The next method to be considered is EOM-CCSD,58 59 discussed above. Since the

excited state wavefunction consists of both single and double excitations, it can have a

proper zeroth order description of not only singly excited states, but also doubly excited

states and states which are mixtures of singles and doubles. On the other hand, since the

reference state determines the orbitals, they are fixed, and the excited state wavefunction

has only limited flexibility to relax them through the other singles and through the doubles

acting as a single times a single. This linear excited state also limits to what extent

differential dynamic correlation is included. More is probably included than for any other

method considered so far, but in EOM-CCSD, unlike CASPT2, the differential dynamic

correlation needs to compensate for the orbital relaxation. The ground state dynamic

correlation is very well described through the ground state coupled-cluster calculation,

and therefore for states without large orbital relaxation and differential correlation, EOM-

CCSD works well.

In STEOM-CCSD,64, 66 the use of a singles only excited state operator limits the

zeroth order description to those kinds of states describable in CIS. The orbital relaxation

is also somewhat limited, but some orbital relaxation comes indirectly from the con-

sideration of the ionized and electron-attached states. The ground state coupled-cluster

calculation again describes the ground state dynamic correlation very well. The differen-

tial dynamic correlation is handled through the second similarity transformation, which is

built from an active space of ionized and electron-attached states. How well this works

in practice for a large variety of molecules is not clear yet, but the amount of differential

correlation included will depend heavily on the completeness of the active space.66









8

Finally, there is multireference CI (MRCI).84 Because of the flexibility offered by

a MRCI calculation, it is possible to achieve extremely good accuracy by properly

accounting for all four parts of the excited state energy. But MRCI has two very

significant drawbacks. It often takes an experienced practitioner to properly build the

trial wavefunction, and the cost of the method is such that it can only be applied to very

small systems. It is not size extensive, either, and its lack of size extensivity will cause

the error in the calculation to grow as the molecule grows.


Gradients for Excited State Methods

One of the most important advances in quantum chemistry was Pulay's derivation

of gradients for SCF wavefunctions.85 This allowed for the convenient computation

of the derivative of the energy with respect to a general perturbation. By taking the

perturbation to be moving an atom, his discovery made practical the exploration of

ground state potential energy surfaces.86 People could readily find minima and transition

states and calculate vibrational frequencies. Pople et al.87 extended gradient technology

to correlated methods with the development of gradients for MBPT(2). For a review of

current analytical gradient techniques, see Ref. 88.

In terms of excited states, gradients were first introduced for the MCSCF method.86

Although not in general use, gradients are also available for MRCI wavefunctions.89

Gradients for CIS were first developed in 1992.75 Then came gradients for EOM-CCSD

in 1993,90-92 followed by CIS(D) gradients.93 For a general overview of gradients for

coupled-cluster based methods for excited states, see Ref. 94.













CHAPTER 2
PARTITIONED EQUATION-OF-MOTION COUPLED-CLUSTER
METHODS FOR EXCITED STATE ENERGIES

Introduction

The equation-of-motion coupled-cluster (EOM-CC) method57-59 is a conceptually

single reference, generally applicable, unambiguous approach for the description of

excited,57-59 electron-attached68 or ionized states.67 All follow from simple consideration

of the Schrodinger equation for two states, a reference state ITO) (not necessarily the

ground state), and an excited (electron-attached or ionized) state I',,). Considering H to

be in second quantization, where the number of particles is irrelevant, the Schr6dinger

equations for the two states are


HI o) = Eol o), (2-1)

and

HIT,) = EI4,). (2-2)

By choosing to represent the excited state eigenfunction as


|T.) = 7Z[lo), (2-3)

it is readily obtained that


[H, Rn]|o}) = R, Io) (2-4)

for w, = E, Eo, from subtracting Eq. (2-1) from Eq. (2-2) after left multiplication

by 7K.









10

The choice of 7Z defines the particular EOM with i, j, k, ... indicating occupied

orbital indices and operators, while a, b, c, ... are unoccupied orbitals and operators.

Also, p, q, r, refer to orbitals and operators of either occupation. For electronic

excited states,


EE = r() + E r ({ai + + E r (s){ a ibtj+ (2-5)
i,a i

for electron attachment


IZEA a E a(){a} E rKb(){atjb} --' (2-6)
a a,b,j

and for ionization


RIP = f) {ri){} + 2J E .(\)fi j + f ., (2-7)
i ij,a

Coupled-cluster theory is introduced by choosing


|Wo) = eTIO) (2-8)

where T is the usual excitation operator, and 10) represents some independent particle

model reference. Since [7,, T] = 0 for any of the above choices, we can commute the

operators to give


[H, J]|O) = (HZ)710) = wj,,Z10), (2-9)


where

H = e-THeT, (2-10)

and (HR 1,) indicates the open, connected terms that remain in the commutator. Two

operators being connected means that they share at least one index. In this way, all of









11

the ground state CC information is contained in H, now generalized to have three- and

four-body terms.58 Note, H is also non-Hermitian, necessitating that both its left C, and

right R, eigenvectors, which form a biorthogonal set, LiR.j = 6ij, be considered for a

treatment of properties. In matrix form,

(HR,)c = Ruw.
(2-11)
LH = Lw,.

Consequently, EOM-CC reduces to a CI-like equation that provides the relevant excitation

energies directly.

Besides the obvious approximations, such as T = Ti + T2, and R, being limited to

single and double excitations, which defines EOM-CCSD, and various triple excitation

extensions, EOM-CCSD(T),61 EOM-CCSDT-1,61, EOM-CCSD(T),63 and EOM-CCSDT-

3,63 one can conceive of many other approximations to the basic EOM-CC structure.

Instead of Eq. (2-8), for the reference state a perturbation approximation might be

chosen, such as


m))= (1 + T(1)+ T(1)T(1) + (2) +. ..) (2-12)
S2 (2-12)

where one truncates at a particular order, m. (See Ref. 95 for the coupled-cluster, non-

Hartree-Fock definition of the various orders.) Alternatively, rather than retaining the

full perturbation approximation to Io(m)), we could truncate H itself to some order.96, 97

It is possible to conceive of perturbative approximations to R,, too. The latter are,

perhaps, most easily viewed from the partitioning approach to perturbation theory.98 That

is, Eq. (2-12) can be partitioned into the spaces P and Q, where P represents the principal









12
configuration space (of dimension p) and Q (of dimension q) represents its orthogonal

complement. Then it is well known that we can consider an effective Hamiltonian,


Hpp(w ) = Hpp + HpQ(wK HQQ)-Hp, (2-13)


whose eigenvectors are solely defined in the P space


Hppcp = cpwu (2-14)


for the first several eigenvalues, w,. Expanding the inverse in Eq. (2-13) provides a

series of perturbative approximations to Hpp or to the eigenvectors cp in Eq. (2-14).98

With P chosen to be the space of single excitations, and Q that of double excitations, such

partitioned EOM-CCSD results, first presented in 198957 have been shown to retain most

of the accuracy of the full EOM-CCSD method. Other approximations can be made. To

introduce triple and quadruple excitations, selective double excitations could be retained

in P and at least diagonal approximations could be made for the triple and quadruple

excitation blocks of HQQ.

In this chapter such partitioned and perturbation-based approximations to the full

EOM-CCSD method will be reconsidered. It will be demonstrated that very good

accuracy may be obtained within a much less expensive computational structure. In

particular, whereas the full EOM-CCSD is proportional to noNAirt oc n6, for n basis

functions, a n cNc rt cx n procedure will be developed that shows promise for large

molecules.









13
Theory

Partitioned EOM-CC for Excited States

In a typical EOM-CCSD calculation, the excitation energy and excited state properties

are calculated via diagonalization of the non-symmetric matrix,58

() [ ss HSD (2-15)
()= [HDS HDD (2-15)

where Hss stands for the singles-singles block of the matrix, etc. The major step in an

iterative diagonalization99 of the matrix is multiplying the matrix by a trial vector C. In

the EE case, the equations for this multiplication are as follows:58' 60


[IssC]a = FacCF FmiCi + WamieC, (2-16)
e m cm

[HSDC]? = FmC, + i WamifC WmniCa, (2-17)
em me f rnne

[HDSC/]j =P(ab) E WmaijCb, + P(ij) E WabejCf
m C

+ P(ab) Wbn fe ) (2-18)
f \me /

P(Cj) (nW Iie
n me
abriaa 1


[HDDC]j =P(ab) Fb P(j) FmC + b WafC'
e m ef

+ WmnijCn + P(ab)P(ij) E WbmjeC
mn emn

2 rP(abE) EWnc f C ea (2-19)
S\mnec

+ P(ij)E EWnfC, )III
n \w-ef









14
The permutation operator P(qr) is defined as


P(qr) (...q ... r...) = ...q.. .r...) (...r...q...) (2-20)

Fpq and Wpqrs are presented elsewhere.58

In the partitioned scheme, HDD is replaced by Ho1n, where H0 is the usual M0ller-

Plesset unperturbed Hamiltonian from many-body perturbation theory.13 In other words,

the unfolded matrix analogous to Eq. (2-13) to be diagonalized is approximated by

IT = /SS HSD \
(H= F HoR ) (2-21)
\HDS HODD/

and Eq. (2-19) becomes


[HODDC]b = P(ab) E fE ,Cje P(ij) E ,,m a, (2-22)
C m

Here, fpq is the Fock operator. In the case when the reference function is composed

of canonical or semi-canonical Hartree-Fock orbitals, HooD becomes diagonal with

differences of orbital energies on the diagonal. From examining the equations, it is

clear that the method is formally an iterative n5 method, and the results are invariant to

rotations among occupied or unoccupied orbitals. This partitioning scheme differs from

that presented in Ref. 57 in two ways. Geertsen et al.57 included the full H elements

on the diagonal instead of just including Ho elements in the doubles-doubles block, thus

destroying the orbital invariance. Also, Geertsen et al.57 did not include the three body

terms (the last two terms in Eq. (2-18)) in their work.

After a matrix in the form of Eq. (2-13) has been diagonalized, the eigenvectors

consisting of single excitations and double excitations are built. Therefore, the same

techniques used to calculate the properties of a full EOM-CCSD wavefunction58 can









15

be used to calculate properties of the partitioned EOM-CCSD (i.e. P-EOM-CCSD)

wavefunction.

Since H is not Hermitian, its left and right hand eigenvalues are the same, but the

eigenvectors differ.58 The equations for the left hand side are similar to those presented

above.

MBPT[2] Ground State

The first logical perturbative approximation to

H =e-THeT = (HeT)c
(2-23)
= Fpq {p q} + 1 Wpqrs{ptqtsr} + higher order terms (2-2
pq pqrs
is obtained by keeping only terms through second order. Such an approximation defines

an EOM-MBPT[2] method for excited states (i.e. an EOM-CC calculation based on a

MBPT[2] ground state instead of a coupled-cluster ground state).97 If we insist upon a

method correct through second-order for the non-Hartree-Fock case, the explicit equations

for H(2 are


Fai = 0, (2-24)

Fab = fab + tafmb t (mallfb) ta (mnllbe), (2-25)
m frm, emn

Fj = fij + tfie + t,(imije) Z te(im ef), (2-26)
e em cfm

Fia = fia + te (iml ae), (2-27)
em
Wabij = 0, (2-28)

Wijkl = (ijkl) + P(kl) ft(ike) + ~{(ij Ief), (2-29)
e oef









16

Wacd = (ab lcd) P(ab) t (amllcd) + t+n (mnllcd), (2-30)
m mn
Waibc = (aillbc) ta(millbc), (2-31)

Wijka (jkllia) + E t |(jkllea), (2-32)
C'
ijab = (ijlab), (2-33)

Waijb = (ailjb) + t t(ailleb) t (milljb) tjr(milleb), (2-34)
e m em





1
Wabci =(abllci) P(ab) t ( er mra (2-35)
+ tImfe + E t (abllce) P(ab) E trn(mbllci),
m e m

Wijk =(ia P(ij) T k(im\je) + yt(ia||ef)
en ef (2-36)
+ Etf + (imjk) P(ij) E t(ia I ek).
e m e
In all of these equations tq and tO refer to their values through first order. That is,

tl1] = 1 and t[lab (abllij)
i ~ ,- ij ~ ,+j-6a-h'

Since approximating the excited state with a partitioned EOM-CC calculation and ap-

proximating the ground state with a MBPT(2) calculation are independent approximations,

we can obviously combine them into a partitioned EOM-MBPT(2) (P-EOM-MBPT(2))

method. When combined, the H elements Wijkl, Wijab, and the most numerous, Wabcd,

are not needed. Therefore, in the Hartree-Fock case, the (abllcd) integrals never need to

be calculated. The (abllcd) integrals would contribute to the Wabci H element multiplied

by a T11, but in the Hartree-Fock case all TI ]'s are zero. A few n6 terms still remain

in the calculation of the H elements, but these terms only need to be calculated once.

In most calculations the cost of the iterative n5 step in the excited state calculation will

dominate over the n6 step involved in calculating the H elements.









17

Calculations

Be Atom

Table 2-1 lists energies calculated using the four methods described above (EOM-

CCSD, P-EOM-CCSD, EOM-MBPT(2), and P-EOM-MBPT(2)) for the first few singlet

excited states in beryllium. The basis set from Ref. 100 was used in the calculations.

From the mean absolute errors, it would appear that only the EOM-CCSD calculation is

able to reproduce the full CI excitation energies (which are exact within the basis set used).

However, if the 11D state is excluded from the determination of the mean absolute error,

the errors become 0.010 eV for EOM-CCSD, 0.151 eV for P-EOM-CCSD, 0.431 eV for

EOM-MBPT(2), and 0.300 eV for P-EOM-MBPT(2). With an approximate excitation

level (AEL) value of 1.60, the I'D state is dominated by double excitation character.

The AEL is a measure of the number of electrons excited in an excitation.58 A value of

1.00 is a pure single excitation, while 2.00 corresponds to a pure double excitation. Since

the doubles-doubles block of the EOM wavefunction is approximated in this partitioned

scheme, it is reasonable to expect any state with appreciable double excitation character

will be poorly described by this partitioned EOM calculation.

Even with the I'D excluded, the calculations based on a MBPT(2) ground state are

still very poor. This poor behavior can be ascribed to the inadequacy of the MBPT(2)

wavefunction in describing the ground state of Be. In a CCSD calculation the largest

T2 amplitude is 0.065, while the largest TT1 amplitude in the MBPT(2) calculation is

only 0.034. Also, the MBPT(2) calculation only obtains 73% of the correlation energy

recovered by the CCSD calculation. Clearly, MBPT(2) is not adequate to describe the









18
ground state, and any excited state calculation based on that MBPT(2) ground state will

suffer accordingly.

In Table 2-2 the first few triplet excited states of beryllium starting from the singlet

ground state are presented. Since all of the calculated states are single excitations, the

EOM-CCSD method does exceedingly well, while the P-EOM-CCSD also does a good

job of describing the states.

Example Molecules

While comparisons with full CI calculations are useful, it is also informative to look

at the performance of the methods for more common molecules. These calculations

will also provide an opportunity to compare the current methods with other single

reference methods used today. Tables 2-3 to 2-6 present calculations on four molecules:

formaldehyde, acetaldehyde, ethylene, and butadiene. The calculations are performed at

the MP2/6-31G* geometries given in Refs. 101-103. A 6-311(2+,2+)G** basis setl01

is used for formaldehyde and ethylene, with a 6-311(2+)G* basis set102 being used for

acetaldehyde and butadiene. All electrons are correlated except in butadiene, where the

first four core orbitals are dropped.

All of the methods presented, except for the CIS-MP2,75 can be viewed as approxi-

mations to the full EOM-CCSD method. TDA (or CIS)104 is the crudest of the methods in

that the excited state is given as a linear combination of single excitations out of a single

reference ground state. This method includes no dynamic correlation. CIS(D)83 provides

a non-iterative n5 perturbation correction to the CIS energy. CIS-MP2 also provides a

correction to the CIS energy, but the method is not size-consistent and scales as n6.83

The partitioned methods provide an iterative n5 excited state method based on either a









19

n5 ground state (P-EOM-MBPT(2)) or an iterative n6 ground state (P-EOM-CCSD). The

full EOM-CCSD method is an iterative n6 method.

Except for two differences, our assignment of states agrees with the assignments

of Wiberg et al. and with the corrected assignments of Head-Gordon et al.,105. The

states have been ordered based on their experimental excitation energies, or, where not

available, their EOM-CCSD excitation energies, instead of their CIS excitation energies.

Also, Head-Gordon etal.83 had incorrectly assigned the 'A1 EOM-CCSD state at 9.27 eV

in formaldehyde to the valence 41A, state. Based upon the state's properties, including

an value of 54.4 au2 (compared to the ground state value of 20.5 au2), the

'Al EOM-CCSD state has been reassigned to the Rydberg 3'A, state. The 4'A, state

(the 7r*+--r state), with an value of 32.6 au2, has an EOM-CCSD excitation energy

of 10.00 eV.

Considering the basis sets and the relatively inexpensive methods used here, these

calculations are not meant to be definitive. For previous results on these molecules

see the references in Refs. 101-103. Instead, the mean absolute errors of the methods

compared to experiment and compared to EOM-CCSD will be used to access the quality

of the methods. The mean absolute error for each of the methods with respect to the

experimental values given is 0.66 eV for CIS, 0.41 eV for CIS-MP2, 0.32 eV for CIS(D),

0.17 eV for P-EOM-MBPT(2), 0.20 eV for P-EOM-CCSD, and 0.14 eV for the full EOM-

CCSD method. Compared to the more complete EOM-CCSD method, the mean absolute

errors of the various methods are 0.70 eV for CIS, 0.39 eV for CIS-MP2, 0.35 for CIS(D),

0.12 eV for P-EOM-MBPT(2), and 0.16 eV for P-EOM-CCSD.









20

As should be expected from its simple nature, TDA performs poorly. The average

errors for CIS-MP2 and CIS(D) are similar, although, as previously noted,83 the CIS-MP2

energies are more erratic. The P-EOM-MBPT(2) and P-EOM-CCSD energies are also

similar, seldom differing by 0.1 eV. This suggests that for these cases where MBPT(2) is

able to well represent the ground state, P-EOM-MBPT(2) should be able to well describe

singly excited states. A balance argument between the ground and excited state would

also tend to favor the P-EOM-MBPT(2) method, since the partitioning, as discussed

above, can be viewed as a second-order in H perturbation expansion for the excited state.

Looking only at the valence states, as designated by Wiberg et al.,101-103 the mean

absolute errors from the EOM-CCSD energies are 0.45 eV for CIS, 0.68 eV for CIS-MP2,

0.33 eV for CIS(D), 0.23 for P-EOM-MBPT(2), and 0.38 eV for P-EOM-CCSD. This is

compared to the Rydberg states, where the mean absolute errors from the EOM-CCSD

energies are 0.76 eV for CIS, 0.32 eV for CIS-MP2, 0.35 eV for CIS(D), 0.10 for P-

EOM-MBPT(2), and 0.11 eV for P-EOM-CCSD. It appears that the partitioned methods

on average do not describe valence states as well as they describe Rydberg states. Since

orbital relaxation is often greater for valence states, it is possible the reduced relaxation

available by approximating the doubles vector might cause larger errors for valence states.

These calculations also give a good chance to measure the time savings of the

partitioned method. For example, for butadiene the full EOM-CCSD calculation took

25905 seconds on an IBM RS/6000 model 590 to calculate eighteen excited states (both

right and left hand sides). To calculate the same states with the P-EOM-CCSD method

took 3998 seconds. These times only include the time for the excited states. P-EOM-









21

MBPT(2) also saves time in the calculation of the ground state and in the formation of

H, as well as saves disk space.


Conclusions

In this chapter a formalism is presented and results are given for partitioned methods

based on the equation-of-motion coupled-cluster theory, where the ground state can be

described by either a CCSD or an MBPT(2) wavefunction. The partitioned methods

provide an iterative n5 method (plus an n6 step for forming H elements) for excited states.

When the ground state of the system is well described by an MBPT(2) wavefunction,

the P-EOM-MBPT(2) method provides an inexpensive way to accurately calculate the

energies and properties of singly excited states. For systems less well described by

a MBPT(2) wavefunction, the P-EOM-CCSD method is a generally accurate, but more

economical, alternative to a full EOM-CCSD calculation. For vertical excitation energies,

P-EOM-MBPT(2) is a superior n5 method to CIS(D).











Table 2-1: Be singlet excitation energies (in eV)


EOM-
CCSD
5.323
6.773

7.139
7.468
8.055

8.084
8.309
8.548
8.583d
8.700


P-EOM-
CCSDb

5.666
6.985
10.579
7.653
7.892
8.203
8.432
8.548
8.694
8.792


0.485
(0.151)e


a)Ref. 60. The AEL is for the EOM-CCSD wavefunction.


b)Because of a different implementation of the partitioning, these numbers are slightly

different than those in Ref. 57.

C)Ref. 100.

d)Ref. 60.


e)Average error without 11D double excited state.


state
liP
2'S
I'D
21P
2'D
31S
31P

3'D
4'S
4'p


AELa
1.07
1.06

1.60
1.06
1.21
1.04
1.05
1.09
1.04
1.04


EOM-
MBPT(2)
4.869
6.329
6.682
7.023
7.618
7.648
7.873
8.115


8.267


0.428
(0.431)e


P-EOM-
MBPT(2)
5.228
6.576
10.170
7.224
7.477
7.788
8.011

8.132
8.278
8.374


full CIc
5.318
6.765
7.089
7.462
8.034
8.076
8.302

8.536
8.600
8.693


mean
abs. error


0.014
(0.010)e


0.578
(0.300)e











Table 2-2: Be


EOM-
CCSD

2.729

6.447

7.301

7.748

7.991

8.278

8.456

8.58d

8.70d



0.008


23
triplet excitation energies (in eV)


P-EOM-
CCSDb

2.819

6.583

7.424

7.866
8.089

8.366

8.539

8.647

8.766



0.104


EOM-
MBPT(2)

2.271

5.997
6.868

7.316

7.557


8.025







0.436


P-EOM-
MBPT(2)

2.366
6.140

7.006

7.448

7.668

7.949
8.121

8.228

8.348


a)The AEL is for the EOM-CCSD wavefunction.

b)See footnote b, Table 2-1.

c)Ref. 100.


d)Ref. 59.


AELa

1.01

1.04

1.04

1.04

1.04

1.03

1.03


state
13P

23S

23P
13D

33S

43P
23D

43S

53P


full CI

2.733

6.444

7.295

7.741

7.985

8.272

8.449

8.560

8.686


mean
abs. error


0.321


--









24
Table 2-3: Formaldehyde (energies in eV)

CIS- P-EOM- P-EOM- EOM-
CISa MP2a CIS(D)b MBPT(2)c CCSDC CCSD Exp.a
11A2(V) 4.48 4.58 3.98 4.31 4.41 3.95b 4.07
11B2(R) 8.63 6.85 6.44 7.10 7.10 7.06b 7.11
21B2(R) 9.36 7.66 7.26 7.97 7.95 7.89b 7.97
2'A1(R) 9.66 8.47 8.12 8.02 8.01 8.00b 8.14
21A2(R) 9.78 7.83 7.50 8.25 8.23 8.23b 8.37
31B2(R) 10.61 8.46 8.21 8.99 8.97 9.07b 8.88
I'BI(V) 9.66 9.97 9.37 9.61 9.68 9.26b
31Ai(R) 10.88 8.75 8.52 9.24 9.21 9.27c
41B2(R) 10.86 8.94 8.63 9.39 9.38 9.40b
4'1A(V) 9.45 9.19 8.80 10.08 10.24 10.00c


a)Ref. 102.

b)Ref. 83.


C)Present work.









25
Table 2-4: Acetaldehyde (energies in eV)

CIS- P-EOM- P-EOM- EOM-
CISa MP2a CIS(D)b MBPT(2)C CCSDc CCSD Exp.a
I'A"(V) 4.89 5.27 4.28 4.65 4.71 4.26b 4.28
21A'(R) 8.51 6.71 6.13 6.88 6.84 6.78b 6.82
3'A'(R) 9.22 7.57 7.04 7.57 7.52 7.49b 7.46
2'A"(R) 9.37 7.37 6.90 7.70 7.64 7.64b
41A'(R) 9.30 8.00 7.42 7.77 7.70 7.68b 7.75
51A'(R) 10.19 8.09 7.70 8.42 8.35 8.39b 8.43
61A'(R) 10.26 8.08 7.70 8.53 8.47 8.51b 8.69
3'A"(R) 10.31 8.10 7.74 8.58 8.51 8.57b
4'A"(V) 9.78 10.34 9.34 9.61 9.65 9.23c
71A'(V) 9.73 9.07 8.50 9.58 10.07 9.44C


a)Ref. 102.

b)Ref. 83.


C)Present work.











Table 2-5: Ethylene (energies in eV)
CIS- P-EOM- P-EOM- EOM-
CISa MP2a CIS(D)b MBPT(2)c CCSDc CCSD Exp.a
llB3u(R) 7.13 7.52 7.21 7.45 7.51 7.31b 7.11
11BluV) 7.74 8.39 8.04 8.20 8.36 8.14b 7.60
l1Big(R) 7.71 8.14 7.84 8.08 8.14 7.96b 7.80
11B2g(R) 7.86 8.12 7.86 8.12 8.18 7.99b 8.01
21Ag(R) 8.09 8.42 8.18 8.43 8.48 8.34b 8.29
2'B3u(R) 8.63 8.92 8.69 8.95 9.00 8.86b 8.62
I'Au(R) 8.77 9.00 8.80 9.07 9.12 9.01b
3'B3u(R) 8.93 9.14 8.96 9.23 9.28 9.18b
2'B1u(R) 9.09 9.38 9.18 9.42 9.51 9.39C 9.33
21Big(R) 9.09 9.31 9.12 9.42 9.48 9.38c 9.34


a)Ref. 101.

b)Ref. 83.


c)Present work.











Table 2-6: Butadiene (energies in eV)
CIS- P-EOM- P-EOM- EOM-
CISa MP2a CIS(D)b MBPT(2)c CCSDc CCSD Exp.a
l'Bu(V) 6.21 7.00 6.29 6.52 6.63 6.42b 5.91
1'Bg(R) 6.11 6.73 6.11 6.40 6.39 6.20b 6.22
11Au(R) 6.45 7.03 6.44 6.72 6.73 6.53b
21Au(R) 6.61 7.11 6.55 6.85 6.84 6.67b 6.66
21Bu(R) 6.99 7.58 7.03 7.29 7.34 7.17b 7.07
21Bg(R) 7.22 7.66 7.17 7.47 7.47 7.31b 7.36
2'Ag(R) 7.19 7.74 7.19 7.43 7.44 7.10b 7.4
31Bg(R) 7.25 7.74 7.24 7.54 7.53 7.39b 7.62
41Bg(R) 7.39 7.87 7.40 7.70 7.69 7.55b 7.72
3'Ag(R) 7.45 7.88 7.44 7.74 7.73 7.61b
31Bu(R) 8.05 8.40 8.01 8.32 8.31 8.21c 8.00
3'Au(R) 7.78 6.75 7.73 8.07 8.07 7.92c 8.18
41Au(R) 7.92 7.66 7.86 8.19 8.18 8.06c 8.21


a)Ref. 103.

b)Ref. 83.


C)Present work.














CHAPTER 3
GRADIENTS FOR THE PARTITIONED EQUATION-OF-
MOTION MBPT(2) METHOD FOR EXCITED STATES


Introduction

The equation-of-motion coupled-cluster method with single and double excitations

for excitation energies (EOM-CCSD)58 59 provides an accurate method for calculating

the energy and properties of many excited states of molecules. With the development of

gradients for EOM-CCSD90-92 it has become possible to study the excited state potential

energy surfaces, just as the development of coupled-cluster gradients for ground state

methods106-109 made routine the study of ground state potential energy surfaces. More

recently, Stanton and Gauss have developed second derivatives for EOM-CCSD110 and

gradients for ionization potential EOM-CCSD (IP-EOM-CCSD).67 The key to all coupled-

cluster gradients has been the recognition of how the interchange theorem of Dalgarno

and Stewart,"1 initially used by Adamowicz et al.106 and Bartlett107 for coupled-cluster

gradients and sometimes known in quantum chemistry as the Z-vector method,112 could

be used to avoid computing the response of the ground state T amplitudes to each

perturbation.

One of the primary problems with the EOM-CCSD method is its cost. An EOM-

CCSD gradient calculation involves four iterative steps which scale as nccNvirt, where

nocc stands for the number of orbitals occupied in the reference determinant, and Nvirt

stands for the number of orbitals unoccupied in the reference determinant These

iterative no~N4irt steps are calculating the ground state T amplitudes, calculating the









29
right hand eigenvector 7, calculating the left hand eigenvector C, and inverting H as

part of calculating the Z vector.

In the philosophy of replacing ground state cluster amplitudes with second-order

perturbation approaches,113, 96 Stanton and Gauss97 developed an EOM-CC method based

on a H that was expanded through second order. For calculations based on a Hartree-

Fock reference, expanding H through second order and replacing the CCSD amplitudes

with the MBPT(2) amplitudes in the EOM-CC equations are functionally equivalent.

This truncation of H reduced the cost of calculating the T amplitudes to an non-iterative

nocc(noc + Nvirt)4 step and reduced the cost of inverting H to n,,.ocir.97 However,

the iterative norcNirt steps for calculating R and L still remained. Since this did not

significantly reduce the cost, the method was not very useful in practice.

In order to further reduce the cost and to create a more practical method geared

toward providing a robust second-order-like treatment of excited states, in analogy with

the robust MBPT[2] for ground states,114 we proposed partitioning away the doubles-

doubles block of H. It is in the doubles-doubles block of H that all of the iterative

N.ccNirt steps in calculating R and arise. In practice this means replacing H in

the doubles-doubles block with Ho, which is diagonal in the Hartree-Fock case. Thus,

the iterative nccNirt steps involved in calculating R and are replaced with iterative

ccNirt steps. The overall effect is that the cost is significantly reduced, with little loss
in accuracy in vertical excitation energies.114

Now, however, the non-iterative noccNvirt step involved in calculating H elements

can also become significant. The point at which the cost of calculating the one n' NaN

step in H becomes dominant over the cost of calculating the three iterative n Cirt steps









30

in the excited state can be estimated. The two will have approximately equal cost when

the number of virtual orbitals is three times the number of iterations. Since it typically

takes around fifteen iterations to solve for each R and each L, this would imply that

the cost of calculating H would become significant when the number of virtual orbitals

is more than about ninety per excited state being calculated. In practice, some of the

terms not included in this estimate also contribute significantly to the cost, increasing

the break-even point substantially. Also, when calculating vertical excitation energies,

multiple excited states can be calculated while forming H only once.

Theory

EOM-CC Gradients

Before we beginning with an overview of EOM-CC theory and EOM-CC gradients,

it will be helpful to define some notation. Given a set of spin orbitals and a reference

determinant 10), the labels i,j, k,. represent spin orbitals occupied in 10), the labels a,

b, c, .represent spin orbitals unoccupied in 10), and the labels p, q, r, represent

spin orbitals whose occupation is not specified. The reference determinant is usually,

but not necessarily, the SCF determinant. The space of all possible determinants with

n electrons formed from the spin orbitals is represented by jh). The space Ih) is then

divided into Ih) = Ip) D Iq). Here, 1p) represents the space spanned by the operators

C and R, where it is assumed that T, L, and R all have the same rank. Finally, Ip) is

divided into lp) = 10) Ig). For EOM-CCSD |g) consists of all determinants singly and

doubly excited with respect to 10). This is the same notation used by Stanton.90

The ground state coupled-cluster equations can be written as


(01H 0) = Ecc,


(3-1)











and

(gIHrO) = 0. (3-2)

In pure CC methods (e.g. CCSD or CCSDT), H can be written as


S= e-THe = (He (3-3)

where the subscript c indicates that the terms are connected; they share at least one index

in common.


T = tai + {a tibtj} +... (3-4)
i,a '.,
a,b
is an excitation operator which accounts for the ground state correlation. In the iterative

triples methods such as the CCSDT-n methods,15 or in many-body perturbation theory,

H does not have such a concise form. Therefore, Eqs. (3-1) and (3-2) will be taken

as the defining equations of H.

It is now possible to write the EOM-CC equations58 as


(0oJHlp) = E(O|llp), (3-5)

and

(plHR I0) = E(pIIIO). (3-6)

Here, E is the total energy for the excited state,


Si= lo + {i'ajtb} + ... (3-7)
a,i a,b
ij
is a left eigenvector of Ht, and


S= ro + Z ai} + I EI afibtj} +... (3-8)
ia ,,j
.,b









32
is a right eigenvector of H. It is also possible to include the equation-of-motion in the

7Z equation. Eq. (3-6) then becomes59


(p| [H, Z] 10) = w(plR|0), (3-9)

where

w = E ECC (3-10)

is the excitation energy.

It can be seen that in EOM-CC the excited states are the eigenvectors of H, with

the corresponding eigenvalues being the energies. Since H is formed from a similarity

transformation, it is not Hermitian, but its eigenvectors do form a biorthogonal set.

Choosing the norm of the excited state vectors to be one gives that


(01l7Zk0) = 6jk (3-11)

for all states j and k. Therefore, from either Eq. (3-5) or (3-6) we can get an expression

for the excited state energy


E = (0|1 T H |0). (3-12)


When Stanton first derived the equations for EOM-CC gradients,90 he followed the

same development as that of Refs. 106, 107, 109 by taking the derivative of Eq. (3-12)

with respect to a general external perturbation. In the process, he had to introduce the

perturbation independent quantity zeta (Z) to account for the effect of the excited state

perturbation on the ground state T amplitudes. Since L and R are eigenvectors for the

excited state, they are variationally optimum and therefore stationary with respect to the









33
perturbation. On the other hand, T is never stationary (even for the ground state), and

another operator must be included to account for its response.

Szalay94 then used the approach of a general functional analogous to the A functional

in ground state coupled-cluster theory14, 116 by introducing Z to force T to be stationary

with respect to the external perturbation. Following this approach, here is introduced

an EOM-CC functional which has the property that all quantities are stationary, and

therefore the functional will satisfy the generalized Hellmann-Feynman theorem.117, 118

Such a functional is


F = (OICH1Z0) + (01ZH 0) + E(1 (01L 10)). (3-13)


By construction


Z = ('ita + iaC aj-b + (3-14)
a,i a,b

is a pure deexcitation operator of rank equal to T. Since the last two terms of Eq. (3-13)

do not contribute to the value of the functional, the value of the functional is just the

energy. By taking the derivative of each of the quantities on the right hand side of

Eq. (3-13) and setting them equal to zero, the EOM-CC gradient equations will be

recovered. For example,

OF
DE = 0 = 1 (01OR 0) (3-15)

is a restatement of Eq. (3-11), the normalization condition for L and R. Taking the

derivative of the functional with respect to R (C) gives Eq. (3-5) (Eq. (3-6)), which is

the equation for C (R). Taking the derivative of the functional with respect to Z gives

Eq. (3-2), which is the ground state T equation.









34
However, taking the derivative of the functional with respect to T is much more

complicated. Since the form of H, and therefore the equations for T, vary between

different CC methods, it should not be surprising that the equations for Z also vary

between different EOM-CC methods. Here will be assumed a pure EOM-CC method

(e.g. EOM-CCSD, EOM-CCSDT, ...). When the derivative with respect to one of the

operators is taken, L, R, T, or Z, what is really being done is separately taking the

derivative with respect to every coefficient in the operator. For example, -, really

means 9 for all t amplitudes in T. This still leaves the excitation or deexcitation part of

the operator. So when jr is taken as a part of getting 4, what results is
8H
= -Te-THeT + e-THeTQT
OT (3-16)
= [fH,QT].
Here

T ti + {atibtj} ... (3-17)
i'a '>j
a,b
is the pure excitation operator without any coefficients. Therefore,

OF
0T = 0 = (ol0[H, QT]n IO) + (OIZ[H, QT] 0). (3-18)

After significant simplification, Eq. (3-18) reduces to

(0|ZIg) = (01Lr|q)(qjlZg) ((gl(Eccl H) g))-, (3-19)

where 1 stands for a unit matrix of rank Ig). This is exactly equivalent to the equation

for Z in Ref. 90.

Finally, since the functional satisfies the generalized Hellmann-Feynman theorem,

for any perturbation X

OF OE
ox= = (OlHxZIO) + (OlZHrXlO). (3-20)









35
Here, Hf represents the derivative of the bare Hamiltonian elements with respect to the

perturbation X within the similarity transformed Hamiltonian H. For the pure EOM-CC

methods this means


Hx = e-T He. T(3-21)

Eq. (3-20) is the same as Eq. (36) of Ref. 90. Although the form of the functional was

simply postulated, the fact that it produces equations (and therefore density matrices)

identical to those produced through straightforward derivation of the energy expression90

proves that the functional is valid. For a discussion of the properties of the density

matrices resulting from these equations, see Ref. 119.

It is also possible to derive an energy expression from Eq. (3-9). Projecting on the

left by (0Lp) gives the expression


E = (OI[H, R] I0) + (I0H O), (3-22)

where the second term, Ecc, is added in to give a total energy instead of just the

excitation energy. The functional which corresponds to this energy expression is


F = (0[H, ]] 10) + (0oHI0) + (OIZrI0) + w(1 (0IrL10)). (3-23)


Superficially, Eqs. (3-13) and (3-23) appear similar. In fact, Eq. (3-22) can be straight-

forwardly derived from Eq. (3-12), and both functionals have the same value, the total

energy E. However, these energy functionals contain several important differences. The

first is that w, the excitation energy, appears on the right instead of E. More importantly,

though, the equations for and value of Z have changed. This will seen later, but first it

is important to see that the EOM-CC equations are still contained in the new functional.









36
As before, the derivative of each of the quantities on the right hand side of the

functional will be taken and set to zero. First,

OF
= 0 = 1 (070) (3-24)

again gives the normalization condition for C and R.

OF
S= 0 = (glHI0) (3-25)

again gives the coupled-cluster equations.

Next, taking the derivative of the functional with respect to L gives

OF
OL = 0 = (pl [H, 7] 10) w(pIR.0), (3-26)

which is equivalent to Eq. (3-9). Expanding the commutator gives


w(PlIZ0) = (plrIZ0) (plRHI0). (3-27)

Inserting a resolution of the identity Ih)(hl into the last term gives

w(pnlZIO) = (pIHRIO) (pjlZq)(qlH0O)

(pIR7g)(glHj0) (pIRo0)(0IHIO) (3-28)

= (plHRIO) (pR0IO)Ecc,
or

(Ecc + w)(plZIO) = (pIHR0O). (3-29)

The third term in the first equality of Eq. (3-28) disappears because of the solution of the

coupled-cluster equations. The second term in the first equality of Eq. (3-28) disappears

because R, being an excitation operator, can never decrease the excitation level. Thus









37
it cannot have any non-vanishing terms between (p| and 1q). Note that Eq. (3-29) is

equivalent to Eq. (3-6).

Taking the derivative of the functional with respect to R gives
OF
S= = 0 = I(O[f a] 10) w(oILOc R1). (3-30)

In analogy to QT above, R is the excitation part of the operator R without any

coefficients. The difference is that fT acting on 10) gives 1g), while Q7 acting on |0)

gives |p), since Qf- contains a constant part and QT does not. Rewriting Eq. (3-30) gives

w(0oCp) = (0lCHft o|0) (01Loc rH0). (3-31)

Inserting a resolution of the identity into the last term gives
wo(01lP) = (0oIHIp) (OILnzO0)(OIHI0) (01LS lp)(plHo0)

(OIl0 CQq)(q|lR0) (3-32)

= (0oLCHp) (OIClp)Ecc,
or

(w + Ecc)(0O Ip) = (0|CHflp). (3-33)

This is simply Eq. (3-5). The third term in the top equality of Eq. (3-32) disappears

because of the solution of the coupled-cluster equations. The last term in top equality

of Eq. (3-32) disappears because Qn acting on |q) gives back 1q), and cannot give

non-vanishing terms between (01 and |q).

For the Zeta equations a pure EOM-CC method will again be assumed. Taking the

derivative of the functional with respect to T gives
OF
T = =(0[[r, r],7Z] i0) + (01 [H, QT] 0) + (OIZ[H ,T] 10)

= (OI[H, RT]) 10) (01CZ[H, QT] I0) + (0|[, QT] 10) (3-34)

+ (01Z [r, T] 10).








38
Notice that the first and last terms are the same as Eq. (3-18). Expanding the commutators
and inserting resolutions of the identity gives
0 =(0o CHq)(qlTRZlO) + (0|IHlp)(plTRZI0) (OI|Tjq)(qlHRI0)

(0L7TIp)(pIHR10) (I0|7Zq)(q[H, rT] 0)

(0|jlZg)(gl [f, T] |0) (I7ZI0)(01 [f, T] 0) + (01 [H, T] 10) (3-35)

+ (0|Zlq)(ql rT0j) + (0|Z|g)(glHnTI0) + (0IZl0)(0|HfTI0)

(0 q)(q- (IZTlqqlI (OIZnTrg)(gljH0) (0IZQrTO)(0Ho10).
The third, fifth, and twelfth term disappear for the same reason that the (OILQ q)

disappeared above. Since the norm is chosen to be one, the seventh and eighth terms
cancel. The ninth and eleventh terms disappear since Z cannot connect (01 to either 10)

or Iq). The thirteenth term disappears because of the solution of the coupled-cluster
equations. Eqs. (3-5) and (3-6) can be used to simplify the second and fourth terms.

Doing this gives
0 =(0oJH q)(qjlZjg) + E(OILIp)(p)ORTZI) E(0)ITIp)(p RI0)

(01LRg)(g Hf2T 0) + (Oj0C(g)(g sT q)(qjHI0)
(3-36)
+ (0[7Zlg)(glrTIg)(glHI0) + (OjLR7g)(glTIO0)(OjHI0)

+ (OIZIg)(glHlg) (0IZIg)Ecc.
The fact that R and QT, both being excitation operators, commute has been used. Since
R acting on 10) only connects with (pl and since acting on (01 only connects with Ip),

the second and third terms cancel. The fifth term disappears since QT cannot connect

(gl to Iq). The sixth term disappears because of the solution of the coupled-cluster
equations. Finally, the Zeta equation becomes
(OjZjg)(gj(Eccl H) )g) =(0LfIq)(qj7Zg)
(3-37)
+ (0 7Z g)(g (Ecc1 R) Ig),











or

(0OIZg) = (ol0rq)(q|Rl g)((gl(Ecc1 f) g))-1 + (0olRg). (3-38)

The first term in Eq. (3-38) is the same as Eq. (3-19). However, the second term is new.

The expression for the energy derivative corresponding to this functional is

OE
O = (OIL[Hx, R] 10) + (olI0XI) + (OIZHXI0). (3-39)

In fact, this equation and Eq. (3-20) produce the same density matrices. In the first

term of Eq. (3-20) Hx and RZ are not required to be connected, as they are in Eq.

(3-39). Those disconnected term can be divided up into two categories. In the first, C

is connected to only R. When all such terms are added up, the norm can be factored

out, leaving the second term of Eq. (3-39). The other category is when C is connected

to Hx. These terms lead to the extra term in Eq. (3-38).

It is important to note that the two different functionals have different equations for

Z, but since they give the same density matrices, either may be used. Therefore, the

operator Z is not unique, but only has meaning within the context of its functional.

To better understand the differences between the functionals, let us examine how the

ground state fits into each. In EOM-CC theory, the ground state is just an eigenstate of

H with the special property that R = 1, and C = 1 + A, where A, a pure deexcitation

operator, is the usual A from coupled-cluster theory.107 Substituting these values into

the first functional gives


F = (01(1 + A)HIO),


(3-40)









40
the usual ground state functional. In this functional Z for the ground state is zero. This

can be seen from Eq. (3-19), where the constant operator 1 cannot connect (ql and 1g).

These give for the energy derivative, Eq. (3-20),

9E
= (01(1 + A)HX |0), (3-41)


again the usual coupled-cluster expression.

For the second functional the picture is somewhat more complicated. Substituting the

appropriate values for R and C into Eq. (3-38) gives that (0|IZg) = (0(1 + A)l1g) =

(0A|g). In other words, Z for the ground state in the second functional is A. But in

this case, R for the ground state, a constant operator, gives zero for the commutators in

the first terms of Eqs. (3-23) and (3-39). So the usual coupled-cluster expressions, Eqs.

(3-40) and (3-41), are still obtained.

The above equations derived from the second functional have a small computational

advantage over those derived from the first functional. Specifically, the second Zeta

equation suggests a simpler way to calculate some of the terms in the density matrix

equations. Stanton and Gauss92 suggest combining Z and roL into a composite operator,

and using this new operator in a ground state coupled-cluster gradient code to calculate

many of the excited state density matrix contributions. The approach presented here

suggests including all of (0|C7Ig) into the composite operator, reducing even further the

number of terms which would need to be calculated separately. In practice, this should

make very little difference in terms of how long a calculation takes. However, it should

make programming somewhat easier.








41
It is possible to derive a third Z equation. Starting with the first line of Eq. (3-34)
and not breaking the commutators gives

(0|Zig) = ((0[[H., ], T] 10)+ (0 [H,s,] 10))((gl(Eccl Hr)g))-l. (3-42)

In the first term H must be connected to both R and QrT, while in the second term H
is connected to just Tr, but Z and L cannot appear. In diagrammatic language saying
that H is connected to OT is equivalent to saying that at least one of the lines from H
must be dangling from the bottom.

Since the second and third Z equations were derived starting from the same place,
they must be equivalent. For EOM-CCSD the third Z equation has two more diagrams
than the first Z equation. None of the diagrams where the two differ would change the
overall computational scaling of the method.

P-EOM-MBPT(2) Gradients

The defining equations for P-EOM-MBPT(2) are identical to Eqs. (3-4) and (3-5), but
with a modified H. Detailed equations for H, Z, and have been given previously.114
For (non-partitioned) EOM-MBPT(2) the functional is97

F=(0 c [(Ho] + H1+ H1 + ]T[11R) 0)

+ 0 H0o + H1 + (H[']11 Tl) 0)
(3-43)
+ (0 Z[(H[] Eo)T[11 + HXI] o)

+ E(1-(0 R o)).

The c means that the operators are connected; they share at least one index. Here, the
usual many-body perturbation theory partitioning of the Hamiltonian is assumed. H[n]








42
is the nth order term of the Hamiltionian, and TIn] is the nth order T amplitude from
many-body perturbation theory. For P-EOM-MBPT(2) the functional becomes

F =( (1 + 2)[(H[o+] + H[11 + H[l]T[1I) I] 0)

+ (0 [(H[] + H[11 + HI1]T[l])72] 0)

+ (0 C2(H[oZ2) 0)
(3-44)
+ (0 H[ol + H11l + (Hl]T[1I) 0) O

+ (0 Z[(H[] E[])T[1] + H['I] 0)

+ E(1- (0 I 0)).

This form of the Zeta operator comes from the defining equation for the T[11 amplitudes,
which is

(g (H[] E[])T1l + H 1 0) = 0. (3-45)

For the rest of the paper it will be assumed that the reference determinant is the
Hartree-Fock solution. In this case,

H[ol = EZ {iti+ EZ{aa}, (3-46)
i a
and
H1 = H H[o
S (pq rs) ptq sr (3-47)
pqrs
Also,

T[I= _1[ a ibtj } (3-48)
nfh
where
ab1l] (ab ij)
tia= (3-49)
bl z + tj a -- b,








43
Since T111 is now a double excitation operator, (0 1 [HIITL1' 2] 0) must vanish,
and the functional is slightly simpler. Another consequence of T['] being a pure double
excitation is that Z consists only of double deexcitation operators.

Taking derivatives of the (simplified) Eq. (3-44) gives the equations for the gradients.
It can be shown that taking derivatives with respect to the various R and C operators
recovers the energy expression. Again, by construction, taking the derivative with respect
to Z recovers the T[1] equation. Finally, taking the derivative with respect to T[11 gives

OF 0 =(o0 I (H[llR,,) 0)+(0 HI[l] d)

+ (0 Z (H Elo]) d).

Here, Qd = 1 C {afibtj} is the pure doubles excitation operator. The fact that R
i,j,a,b
commutes with Qd has been used. Solving for Z gives

(o0 Z d) ={(o I(H[ ']Z d)c 0) (0 H[1 d)}

x((d (Eil H[l) d))(3-51)

Previously, solving for Z involved inverting the full H, but now only the inverse
of (d (E[o] H[]) d is needed, which, in the Hartree-Fock case, is just the
denominator from MBPT(2). It is now convenient to introduce a new operator, E, which
is the inhomogeneous part of the Z equation. For P-EOM-MBPT(2)


ai (3-52)
=(0 1L(H['I zl d) 0o) (0 H111 d).
Computationally, to solve for Z, E must first be calculated and then the denominator
applied.









44
Now the derivative can be rewritten as

OE OF OHi (l 2[O ] OH[1] ____
9E aF i90 all['] a98HI
-O 0 ( + + +-- + rT[1] 0
9Oxx QX ax

S [(H[o] H9Il 9T9H11 ) \
C
O X -,- -X + x T[1 R2] 0/

+(0 o 2 (HOIz2) 0) (3-53)


\ 9 + 5 x 10
+(o H[o] 9 H[' + TH[1 T+I1)

0 Z T x ) 0.

This expression only contains quantities which are independent of the perturbation plus

derivatives of the Hamiltonian elements. Thus, the derivative can be recast as one- and

two-particle effective density matrices times derivatives of integrals;

E Ofq q + pq( rs) (3-54)
Ox = x rs ax

Explicit spin orbital equations for the effective density matrix elements in terms of R,

L, T[1], and Z are

P, = I ;i je 1 / imtef
S= 2 1E e j. + ij, (3-55)
e m,e,f

P= lbm + 2 e "ae (3-56)
b m,n,a
SImn ftea lmnracf. V7 m aca
p= -- l ef i mn 2 lrfn cmi + c rim, (3-57)
m,n,e,f m,n,e,f m,e
Pa = 0, (3-58)

1ij
Pka -(ijr)P (kl)ik, (3-59)

p = --4 1 r (3-60)
ff












a k ea m r (8 + e E f tjk
S= + ) 42lr a + 8 1a
e nm,e,f m,e,f

8 k-() Yme, f
m,e,f


Pb 1 P_(ij)P_(ab)lra,

Pb 4E 1 e64


+ (ab) m l,"i m n
m,n,e


mn b .ea ba
"c ro tinn ar f tan eTba
m,,ec mn
m,n,e.


pi = IE ,m a
Pbc mia m
m
P2 -i _(ij) E\ lmretab P_ ) r "'a Teb
ab -' e i mj -4 2 e m rnj
m,e m,e
1 1
+ P-(ij)P-(ab) I, m raet +
4e ri I -4 4 ,
m,e
pa = 0
cd 0,


pb = 0.
Pzj


The explicit spin orbital equation for 5 is


iJb = P(ab) E lmr(ij Ieb) P_(ij)
m ,e
+ P(ij)P-(ab) lr (rnjIeb)
m,c
P- (ij)P-(ab) e e(m )
m,n,e
+ P(ij)P-(ab) E 'r' (fj leb)
m,e, r
+ P-(ab) 1 i'r e(fmllbe)
m,e,f


- P(ij) r (jnjjne) 1 r'(fmllba)
m,n,e m,ef)

- (ien) + (ijjllab).
m,n,e


(3-61)


(3-62)


(3-63)


(3-64)




(3-65)


(3-66)


(3-67)


m1rm (mj lab)
me


(3-68)









46
In these equations P.(pq) stands for [1 P(pq)], where P(pq) means permute p and q.

The equations for p and for ( are subsets of the corresponding equations for EOM-

MBPT(2) gradients,97 and the P-EOM-MBPT(2) gradients have been implemented in the

ACES II program system* starting from Stanton and Gauss's EOM-MBPT(2) gradient

code.97

One final note should be made about the cost of the P-EOM-MBPT(2) gradients.

Previously it was argued that the iterative n Ncctirt steps should dominate the time over

the non-iterative noccNirt step for normal vertical excitation energy calculations. This

was partially because the cost of the no cANirt step can be amortized over multiple

excited states. However, the gradient for only one excited state can be calculated at a

time. Therefore, this amortization cannot occur. More importantly, there are now many

more noccNirt steps, which, when added together will also contribute significantly to

the cost of the calculation. Therefore, P-EOM-MBPT(2) gradients is most accurately

considered a non-iterative n2 N 4irt method.


Application to Diatomic Molecules


In order to assess the performance of the CIS(D)83 method for excited state potential

energy surfaces, Stanton et al.93 studied singlet excited states for six diatomic molecules.




* The ACES H program is a product of the Quantum Theory Project, University of
Florida. Authors: J. F. Stanton, J. Gauss, J. D. Watts, M. Nooijen, N. Oliphant, S. A.
Perera, P. G. Szalay, W. J. Lauderdale, S. R. Gwaltney, S. Beck, A. Balkovi, D. E.
Bernholdt, K.-K. Baeck, H. Sekino, P. Rozyczko, C. Huber, and R. J. Bartlett. Integrals
packages included are VMOL (J. Alml6f and P. R. Taylor), VPROPS (P. R. Taylor), and
a modified version of the ABACUS integral derivative package (T. U. Helgaker, H. J.
Aa. Jensen, J. Olsen, P. Jorgensen, and P. R. Taylor).









47
In order to compare this new work with the methods presented there, those systems have

also been studied here. The six diatomics are CO, N2, C2, H2, BH, and BF. Except for

N2, the studied states are the lowest singlet excited state of the molecule. For N2 the

11'll state was studied.

The basis sets used were 6-31G*,120 aug-cc-pVDZ,121 and aug-cc-pVTZ.121t From

attempting to match the reported numbers, it would appear that for the 6-31G* and

aug-cc-pVDZ basis sets all six Cartesian d functions were used previously, and for the

aug-cc-pVTZ basis set only the spherical d and f functions were used. The same was

done here.

These states have been carefully studied with very accurate calculations. However,

the point of this work is to compare the current method with other single reference

excited state methods. Therefore, an attempt to review the volumous literature on these

molecules will not be made. Instead the three inexpensive methods (CIS, CIS(D), and P-

EOM-MBPT(2)) will be compared to the more expensive and more complete EOM-CCSD

and to experiment.122 All of the other methods can be viewed as various approximations

to EOM-CCSD. For valence excited states there are frequently important contributions

from triples to consider.60-63





t Basis sets were obtained from the Extensible Computational Chemistry Environment
Basis Set Database, Version 1.0, as developed and distributed by the Molecular Sciences
Laboratory which is part of the Pacific Northwest Laboratory, P. O. Box 999, Richland,
Washington 99352, USA, and funded by the U.S. Department of Energy. The Pacific
Northwest Laboratory is a multi-program laboratory operated by Battelle Memorial
Institute for the U.S. Department of Energy under contract DE-AC06-76RLO 1830.
Contact David Feller, Karen Schuchardt, or Don Jones for further information.









48

Results for the diatomic molecules are presented in Tables 3-1 to 3-5. There are

some small differences between the numbers reported here and the values given in Ref.

93. Several values given in the previous paper were in error. Those errors have been

corrected here.

Vertical and Adiabatic Excitation Energies

The vertical excitation energies are given for the six molecules, for the four methods,

and for the three basis sets in Table 3-1. The adiabatic excitation energies are given in

Table 3-2. The ground state minimum geometries used are listed in Table 3-3, with the

excited state minimum geometries in Table 3-4. The appropriate ground state minimum

for CIS is the SCF minimum. The appropriate ground state minimum for CIS(D) and P-

EOM-MBPT(2) is the MBPT(2) minimum. The appropriate minimum for EOM-CCSD

is the CCSD minimum.

As has been noted before,93 the CIS answers are too poor to be reliable CIS has

several errors of greater than 1 eV, and for C2 CIS even predicts the wrong sign for

the excitation energy. The mean absolute deviation with respect to experiment for the

adiabatic excitation energies with the aug-cc-pVTZ basis is 0.724 eV. These failings are

directly attributable to the method's complete lack of dynamic correlation.

CIS(D) improves on the CIS in all cases, except those four where the CIS is

accidentally already in good agreement. In all cases, the CIS(D) energies are within

0.2 eV of the EOM-CCSD energies, and in all of the cases with the best basis set, the

error with respect to experiment is less than 0.3 eV. The mean absolute deviation with

respect to experiment for the adiabatic excitation energies with the aug-cc-pVTZ basis is









49
0.112 eV. This is partially due to the choice of excited states considered and is actually

better accuracy than was initially reported for the method.83

The P-EOM-MBPT(2) method, on the other hand, actually performs worse than

originally reported."14 The vertical excitation energies are slightly worse than the CIS(D)

energies. The adiabatic energies are slightly better to somewhat worse than the CIS(D)

energies, with their mean absolute deviation with respect to experiment for the aug-cc-

pVTZ basis being 0.279 eV.

The fact that the P-EOM-MBPT(2) energies are worse than the CIS(D) energies seems

to be caused by the choice of states and by a curious feature of CIS(D). In the previous

study,114 it was noted that P-EOM-MBPT(2) performed better on Rydberg states than

on valence states. This can be understood, since Rydberg excitations essentially involve

pulling an electron out of the valence region and putting it in a very diffuse orbital, while

valence excitations put the electron back into the valence space in a different arrange-

ment. Thus valence excitations should involve more orbital relaxation and differential

correlation compared to Rydberg excitations. Therefore, it is reasonable that a method

as simple as P-EOM-MBPT(2) could have more trouble accurately describing valence

excitations than Rydberg excitations. EOM-CCSD is also more accurate for Rydberg

states than valence states, where it takes triples to partially correct the problem.123

CIS(D), on the other hand, performs better for valence states than for Rydberg

states.114 This is curious, since the method does not allow the states to relax in the

presence of electron correlation. All of the states studied here are valence states, which

is probably why the CIS(D) energies are better than the P-EOM-MBPT(2) energies.

Also, these states are well separated from other states of the same symmetry, meaning









50
that mixing of states, which the zeroth order wavefunction, the CIS, would have trouble

handling, does not occur in these problems.

The EOM-CCSD aug-cc-pVTZ results agree quite well with experiment, which is not

unexpected for such relatively simple states. The mean absolute deviation with respect

to experiment for the adiabatic excitation energies is 0.106 eV. Actually, there are not

significant differences between the basis sets, except for H2, where 6-31G* does not

contain any polarization functions.

Bond Lengths and Vibrational Frequencies

Excited state equilibrium bond lengths are presented in Table 3-4 and excited state

vibrational frequencies are presented in Table 3-5. Certain trends hold true. CIS

bond lengths are always too short and vibrational frequencies are all too high. This

is completely analogous to the situation for Hartree-Fock for the ground state. However,

like Hartree-Fock, the bond lengths and vibrational frequencies are reasonable. The

CIS(D) and P-EOM-MBPT(2) bond lengths and vibrational frequencies are quite similar.

And, as should be expected, the EOM-CCSD is superior to the other methods.

For CIS, the mean absolute deviation in the bond length is 0.028 A, and the mean

absolute deviation in the frequency is 188 cm-1. For comparison, the mean absolute

deviation in the ground state bond lengths for the SCF was 0.015 A. The mean absolute

deviation in the bond length drops to 0.015 A for CIS(D), with the mean absolute deviation

in the vibrational frequencies being 108 cm-1. This compares to a mean absolute deviation

of 0.011 A for the P-EOM-MBPT(2) bond lengths, 106 cm-1 for the P-EOM-MBPT(2)

vibrational frequencies, and 0.009 A for the MBPT(2) ground state bond lengths. The

EOM-CCSD bond lengths have a mean absolute error of 0.008A, and the vibrational









51

frequencies have a mean absolute error of 67 cm-1. The ground state CCSD has a

mean absolute deviation of 0.004 A for the bond length. All comparisons were made

between the aug-cc-pVTZ results and the experimental results. In general, the ground

state methods were slightly better than the excited state methods, and EOM-CCSD was

better than P-EOM-MBPT(2), which was slightly better than CIS(D).


Application to Polyatomic Molecules


The P-EOM-MBPT(2) method has also been applied to some prototypical polyatomic

molecules. Because of the much more complicated chemistry possible, polyatomics can

provide a much more demanding test of a method.

The S1 State of Ammonia

The Si state of ammonia has been studied extensively.124 125 As a methods test,

one of its attractive features is its predissociative nature. The character of the state

changes from Rydberg like near its D3h minimum to valence like near the C2v transition

state.126 Therefore for a method to properly describe the barrier height, it must equally

treat the valence and Rydberg parts of the potential energy surface. Studying the Si state

of ammonia will also provide an opportunity to compare the P-EOM-MBPT(2) with the

(non-partitioned) EOM-MBPT(2) method97 in order to assess the effect of the partitioning

apart from the effect of replacing the CCSD ground state with the MBPT(2) ground state.

The basis set used for this study was the "A" basis set of Ref. 127. This basis set

contains 65 contracted functions and includes multiple diffuse functions at each atom. It

is flexible enough to give a reasonable description of the entire area of interest of the

potential energy surface. Table 3-6 gives the geometry and vibrational spectrum for the









52

D3h minimum. It is believed that the N-H bond distance is 1.055 0.008 A.128 The

EOM-CCSD bond distance falls within these error bars. However, a larger basis set would

tend to shorten it.129 The CIS bond length is once again too short and the frequencies

are too large. The CIS(D) P-EOM-MBPT(2), and EOM-MBPT(2) bond lengths are

all similar and slightly shorter than the EOM-CCSD bond length. For the vibrational

frequencies, the CIS(D) results are slightly more erratic compared to the EOM-CCSD

than the are P-EOM-MBPT(2) and EOM-MBPT(2), which are similar.

The infrared intensities show some interesting patterns. In every case CIS(D)

substantially overestimates the correlation correction to the intensity. The P-EOM-

MBPT(2) is only partially able to correct that behavior. Most of the error in the P-EOM-

MBPT(2) can be traced to the very approximate treatment of the doubles-doubles block,

since the EOM-MBPT(2) intensities are much closer to the EOM-CCSD intensities. The

intensities provide a very sensitive measure of the quality of the wavefunction, and these

large errors for the P-EOM-MBPT(2) and especially the CIS(D) methods suggest that

their descriptions of the wavefunction are not as good as the bond lengths and vibrational

frequencies suggest.

The structure and vibrational spectrum of the transition state are presented in Table

3-7. The CIS(D) bond length for the hydrogen being extracted is much too long. This

bond length for the P-EOM-MBPT(2) method is also too long, although shorter than

for CIS(D). CIS actually manages to get this bond length right compared to EOM-

CCSD. These elongated bonds for CIS(D) and P-EOM-MBPT(2) cause w6 to be too

small93 with the CIS(D) actually predicting that the transition state is not C2v. Stanton

et al.93 suggested that the long bond length was caused by the tendency of CIS(D) to










53

underestimate the energy of Rydberg states,83 thus causing a "late" transition state. For

P-EOM-MBPT(2) the problem is the opposite. It tends to overestimate the valence states.

However, that will also lead to a "late" transition state and a longer C-H bond.

Once again, CIS(D) and P-EOM-MBPT(2) have problems with the intensities, with 12,

13, and 15 having significant errors with respect to EOM-CCSD and EOM-MBPT(2). On

the other hand, the CIS(D) and P-EOM-MBPT(2) dipole moments are pretty reasonable.

The barrier heights, though, are much too large, with the CIS(D) barrier height being

much worse than the P-EOM-MBPT(2) barrier height. Even the EOM-CCSD barrier

height is well above the experimental barrier height, estimated to be about 2100 cm-'.128

Trans-bent and Vinylidenic Isomers on the Si Surface of C2H2

The S1 surface of acetylene provides a nice test for excited state methods. The

C2h trans-bent minimum for the S1 state was described in the 1950's.130-132 However,

in a recent theoretical paper, it was suggested that the global minimum was a C2v vinyli-

denic structure with both hydrogens bonded to one of the carbons.133 The vinylidenic

structure was calculated to be 13 kcal/mol lower than the trans-bent acetylenic structure.

This prediction that the vinylidenic state was lower was confirmed with multireference

CI, approximate coupled-pair functional,134 and averaged quadratic coupled-cluster'35

calculations. On the other hand, both CIS and CIS(D) predict the wrong ordering of

the isomers.93

The structures, vibrational frequencies, and energies of the two isomers are presented

in Tables 3-8 and 3-9. The calculations were performed with a TZ2P basis set.136 Not

including zero point energy corrections, CIS predicts the acetylenic structure to be more

stable by 2.1 kcal/mol. The CIS(D) predicts the acetylenic structure to be more stable









54

by 0.9 kcal/mol. The P-EOM-MBPT(2) also has the wrong order, with the acetylenic

structure being 2.9 kcal/mol lower in energy than the vinylidenic structure.

Simple Carbonyls

As the simplest of the carbonyls, formaldehyde's spectrum has been studied

extensively.137' 138 The first excited state, n--+7r*, has two distinctive geometrical fea-

tures. The C-O bond lengthens and the molecule becomes pyramidal. CIS severely

underestimates the bond lengthening,102 while CIS(D) severely overestimates the bond

lengthening.93 Also, CIS(D) predicts the molecule to be almost flat.93

Table 3-10 presents the P-EOM-MBPT(2) results. Other than predicting a C-O bond

that is 0.02 A too long, the P-EOM-MBPT(2) geometry agrees well with experiment,

especially considering the relatively poor basis set used (6-31G*,120 the same that was

used in Ref. 93). The frequencies are also reasonable, except for w4. The dipole is

somewhat too large, but is still in better agreement with experiment than the CIS(D)

dipole.

Table 3-11 presents a comparison of the formaldehyde, acetaldehyde, and acetone

geometries, all with the 6-31G* basis set.120 Hydrogen "a" is a hydrogen attached to

the carbonyl carbon. Hydrogen "b" is the "in plane" hydrogen. Hydrogen "c" is the

hydrogen pointing into the pyramid. Finally, hydrogen "d" is the hydrogen pointing

away from the pyramid. Since acetaldehyde and acetone are pyramidal in their first

excited state like formaldehyde is, there is no true in plane hydrogen. However, in both

cases there is a hydrogen with a dihedral angle only about five degrees out of planarity.

The difference with the ground state, however, is that the in plane hydrogenn points away

from the oxygen instead of being hydrogen bound to it.139 This can be attributed to the









55
removal of one electron from the in plane lone pair on oxygen, letting the steric effects

dominate over the weakened hydrogen bonding.

Zuckermann et al.140 carefully measured the spectrum of acetone around the origin

of the n--+7r* band. Their analysis led to the positive analysis of three vibrational

bands, with uncertainty about two others. Their assignments, along with the calculated

vibrational frequencies, are presented in Table 3-12. Zuckermann et al. could not

definitely distinguish between the first possibility and the second possibility; though

they preferred the first. Determining which is correct depends upon the frequency of

v19. The computed frequency of 368.4 cm-1 is well below the 465.4 cm-1 for the second

possibility. On the other hand the P-EOM-MBPT(2) is low on all of the other frequencies

except for v12, considering that the calculated frequencies are harmonic frequencies as

opposed to the measured fundamental frequencies. Also, in acetaldehyde the P-EOM-

MBPT(2) gets very good agreement between the calculated vibrational frequency (368.4

cm-1) and the measured frequency (370 cm-')139 for the corresponding mode. The next

calculated mode in acetone is at 811.1 cm-1. Hence, the calculated spectrum supports

the second possibility over the first possibility. It would take better calculations, for

example an EOM-CCSD calculation with a better basis set and possibly with triples, to

confirm which one it is.


Conclusions

An alternative and simpler derivation for EOM-CCSD gradients based upon an

excited state EOM-CC functional is presented. Gradients for the P-EOM-MBPT(2)

method have been derived and their cost discussed. Then, by studying a series of









56
diatomics and prototypical polyatomics the performance of P-EOM-MBPT(2) for low-

lying valence potential energy surfaces was studied versus other single reference methods

for excited states for which gradients have been derived.

The performance of CIS, CIS(D), and EOM-CCSD has been discussed previously93,

and so they will only briefly be discussed here. The CIS energies are too poor to be

even qualitatively reliable. However, the CIS geometries are normally quite reasonable.

CIS(D) has the opposite behavior. Its energies for valence states often has errors of 0.3

eV or less, but its geometries sometime fail dramatically. Even worse, there seems to

be no indication of when the CIS(D) geometries will fail. EOM-CCSD performed quite

well in all of the tests presented here, both for the energies and for the geometries.

Overall, the performance of the P-EOM-MBPT(2) method was mixed. The energies

were reasonable but were not as good as had been reported previously.114 The geometries

were also reasonable but were noticeably worse than the EOM-CCSD geometries. Over-

all, in these tests P-EOM-MBPT(2) performed very similarly to CIS(D) when CIS(D) did

not fail and was qualitatively correct for those cases for which CIS(D) did fail.

The discrepancy between the P-EOM-MBPT(2) and the EOM-CCSD geometries

is attributable primarily to the partitioning of the doubles-doubles block in P-EOM-

MBPT(2), since the (non-partitioned) EOM-MBPT(2) performed similarly to the EOM-

CCSD for ammonia. However, it is this partitioning that makes the method attractive,

since it is the partitioning that reduces the cost from an iterative n2c rt to an iterative

nocNirt plus a non-iterative n Nirt










57
Table 3-1: Vertical excitation energies (in


eV) for the excited singlet states


CISa
CIS(D)a
P-EOM-MBPT(2)
EOM-CCSDa
CISa
CIS(D)a
P-EOM-MBPT(2)
EOM-CCSDa
CISa
CIS(D)a
P-EOM-MBPT(2)
EOM-CCSDa
CISa
CIS(D)a
P-EOM-MBPT(2)
EOM-CCSDa
CISa
CIS(D)a
P-EOM-MBPT(2)
EOM-CCSDa
CISa
CIS(D)a
P-EOM-MBPT(2)
EOM-CCSDa


a)Ref. 93.


6-31G*
15.354
15.341
15.349
15.257
3.029
2.991
2.991
3.117
9.385
8.911
9.196
8.833
10.378
9.444
9.883
9.425
6.888
6.865
6.935
6.898
-1.180
1.310
1.423
1.526


aug-cc-pVDZ
12.639
12.574
12.574
12.484
2.849
2.864
2.875
2.972
9.272
8.768
9.047
8.638
10.360
9.389
9.805
9.339
6.559
6.481
6.560
6.482
-1.258
1.151
1.261
1.283


aug-cc-pVTZ
12.714
12.832
12.835
12.717
2.852
2.810
2.824
2.914
9.330
8.767
9.076
8.666
10.550
9.532
9.980
9.514
6.601
6.455
6.545
6.454
-1.229
1.222
1.328
1.308


--


~I


~I











Table 3-2: Adiabatic excitation energies (in eV) for the excited singlet states

6-31* aug-cc- aug-cc- Expt.a
pVDZ pVTZ
H2 CISb 12.765 11.269 11.352

CIS(D)b 12.980 11.290 11.429
P-EOM-MBPT(2) 12.943 11.274 11.408
EOM-CCSDb 13.110 11.228 11.353 11.3694
BH CISb 3.024 2.845 2.849

CIS(D)b 2.989 2.861 2.807
P-EOM-MBPT(2) 2.989 2.872 2.821
EOM-CCSDb 3.117 2.972 2.913 2.8685
CO CISb 8.784 8.739 8.802

CIS(D)b 8.269 8.219 8.246
P-EOM-MBPT(2) 8.566 8.527 8.580
EOM-CCSDb 8.330 8.229 8.256 8.0684
N2 CISb 9.413 9.432 9.582

CIS(D)b 8.802 8.775 8.864
P-EOM-MBPT(2) 9.264 9.231 9.346
EOM-CCSDb 8.760 8.714 8.839 8.5900
BF CISb 6.806 6.519 6.563

CIS(D)b 6.760 6.428 6.404
P-EOM-MBPT(2) 6.822 6.507 6.495
EOM-CCSDb 6.811 6.444 6.417 6.3427
C2 CISb -1.275 -1.350 -1.319
CIS(D)b 1.143 0.997 1.075
P-EOM-MBPT(2) 1.293 1.138 1.208
EOM-CCSDb 1.301 1.066 1.103 1.0404


a)Ref. 122.

b)Ref. 93.











Table 3-3: Equilibrium distances (in A) for the ground states

aug-cc- aug-cc-
pVDZ pVTZ
H2 SCFb 0.7300 0.7481 0.7345
MBPT(2)b 0.7375 0.7549 0.7374
CCSDb 0.7462 0.7617 0.7431 0.7414
BH SCFb 1.225 1.233 1.221

MBPT(2)b 1.233 1.241 1.216
CCSDb 1.244 1.249 1.220 1.2324
CO SCFb 1.114 1.110 1.104
MBPT(2)b 1.150 1.147 1.134
CCSDb 1.141 1.138 1.124 1.1283
N2 SCFb 1.078 1.078 1.067

MBPT(2)b 1.130 1.131 1.110
CCSDb 1.113 1.113 1.093 1.0977
BF SCFb 1.260 1.270 1.249

MBPT(2)b 1.279 1.294 1.264
CCSDb 1.281 1.296 1.263 1.2626
C2 SCFb 1.245 1.253 1.241
MBPT(2)b 1.264 1.276 1.254
CCSDb 1.252 1.265 1.241 1.2425


a)Ref. 122.

b)Ref. 93.










60
Table 3-4: Equilibrium distances (in


A) for the excited singlet states


6-31G* aug-cc- aug-cc- Expt.a
pVDZ pVTZ
H2 CISb 1.544 1.239 1.239

CIS(D)b 1.599 1.256 1.273
P-EOM-MBPT(2) 1.626 1.273 1.295
EOM-CCSDb 1.616 1.267 1.283 1.2928
BH CISb 1.204 1.214 1.204

CIS(D)b 1.219 1.223 1.199
P-EOM-MBPT(2) 1.216 1.223 1.199
EOM-CCSDb 1.241 1.242 1.211 1.2186
CO CISb 1.228 1.220 1.213

CIS(D)b 1.295 1.286 1.263
P-EOM-MBPT(2) 1.293 1.282 1.259
EOM-CCSDb 1.252 1.242 1.224 1.2353
N2 CISb 1.200 1.200 1.192

CIS(D)b 1.250 1.249 1.231
P-EOM-MBPT(2) 1.245 1.244 1.227
EOM-CCSDb 1.221 1.220 1.202 1.2203
BF CISb 1.316 1.312 1.287

CIS(D)b 1.350 1.349 1.312
P-EOM-MBPT(2) 1.353 1.349 1.312
EOM-CCSDb 1.345 1.342 1.304 1.3038
C2 CISb 1.293 1.301 1.289
CIS(D)b 1.333 1.346 1.320
P-EOM-MBPT(2) 1.325 1.337 1.312
EOM-CCSDb 1.334 1.347 1.318 1.3184


a)Ref. 122.

b)Ref. 93.











Table 3-5: Harmonic vibrational frequencies (in cm-') for the excited singlet states

6-31* aug-cc- aug-cc-pt.a
pVDZ pVTZ
H2 CISb 1495 1610 1589

CIS(D)b 1438 1485 1422
P-EOM-MBPT(2) 1407 1410 1337
EOM-CCSDb 1428 1439 1368 1358.09
BH CISb 2576 2536 2544

CIS(D)b 2426 2441 2517
P-EOM-MBPT(2) 2445 2443 2515
EOM-CCSDb 2180 2243 2372 2251.0
CO CISb 1646 1615 1633

CIS(D)b 1282 1238 1323
P-EOM-MBPT(2) 1281 1244 1330
EOM-CCSDb 1559 1517 1592 1518.2
N2 CISb 1939 1909 1897

CIS(D)b 1612 1583 1615
P-EOM-MBPT(2) 1638 1612 1635
EOM-CCSDb 1856 1825 1854 1694.21
BF CISb 1339 1273 1367

CIS(D)b 1173 1097 1239
P-EOM-MBPT(2) 1160 1095 1241
EOM-CCSDb 1199 1130 1279 1264.9
C2 CISb 1830 1797 1794
CIS(D)b 1618 1577 1631
P-EOM-MBPT(2) 1674 1635 1687
EOM-CCSDb 1605 1563 1630 1608.35


a)Ref. 122.

b)Ref. 93.











Table 3-6: Geometries (in A), harmonic vibrational frequencies (in cm-1), infrared intensities
(in km/mol), and energies (in Hartrees) for the D3h equilibrium geometry of the S1 state of NH3


CIS(D)a


1.0441
2814.0
736.5
3277.2
1378.6
57.6
1538.5
981.9
35 -56.240 3


P-EOM-
MBPT(2)
1.0420
2971.4
767.2
3185.5
1370.6
48.3
3620.2
464.4


EOM-
MBPT(2)b
1.0487
3180.8
769.6
3021.0
1331.1
14.7
4622.1
293.0


376 -56.225 152 -56.234 687


EOM-
CCSDb
1.0512
2993.1
741.2
2997.5
1335.2
9.7
4447.1
376.0
-56.246 539


CISa


rNH
wI(al')
w1(ai')
w2(a2")
w3(e')
w4(e')
12
13
I4
Energy


1.0213
3180.8
842.4
3356.7
1517.9
0.1
5803.2
19.4
-55.968 6


a)Ref. 93.

b)Ref. 97.











Table 3-7: Geometries (in A and degrees), harmonic vibrational frequencies (in
cm-1), infrared intensities (in km/mol), dipole moments (in Debyes), energies (in
Hartrees), and barrier heights (in cm-1) for the C2v predissociative transition state of
the Si state of NH3. The symmetry unique hydrogen is denoted by an asterisk.


CISa


1.3497
1.0106
123.43
3538.5
1568.5
1446.6i
1025.3
3737.2
498.8
293.7
5.4
463.7
50.2
2.0
73.0
3.750
-55.952 723

3492


CIS(D)a


1.4371
1.0421
126.78
3133.7
1539.1
1510.3i
1276.0
3313.9
104.8i
640.4
0.3
5.5
88.4
6.7
66.5
2.951
-56.214 426


rNH*
rNH
0(H*NH)
wl(al)
w;2(al)
w3(al)
w4(bl)
w5(b2)
W6(b2)
I1
12
13
14
15
16
/t'
Energy
Barrier
Height


P-EOM-
MBPT(2)
1.4108
1.0370
125.82
3203.6
1525.9
1543.2i
1191.1
3384.0
129.9
590.5
0.0
112.1
68.0
3.7
103.4
2.985
-56.205 848

4237


EOM-
MBPT(2)b
1.3219
1.0423
124.51
3069.0
1478.5
1977.3i
1042.0
3274.9
396.5
1081.2
48.9
2301.1
43.9
58.1
158.0
2.631
-56.223 512

2452


EOM-
CCSDb
1.3421
1.0441
124.27
3051.6
1456.0
1897.3i
1023.1
3255.2
419.6
1010.6
41.0
1921.3
69.9
56.5
167.5
2.727
-56.234 405

2663


5695


a)Ref. 93.

b)Ref. 97.











Table 3-8: Geometries (in A and degrees), harmonic vibrational
frequencies (in cm-'), infrared intensities (in km/mol), and energies
(in Hartrees) for the trans-bent isomer of the S1 state of acetylene

P-EOM-
CISa CIS(D)a EOM-CCSDb
MBPT(2)
rCH 1.0788 1.0889 1.0890 1.0907
rcc 1.3521 1.3724 1.3717 1.3575
0(HCC) 124.55 122.06 121.86 123.64
w1(ag) 3293.5 3144.7 3141.7 3107.7
w2(ag) 1545.1 1423.6 1436.4 1471.4
w3(ag) 1162.4 1108.2 1114.1 1106.6
w4(au) 1230.8 687.2 907.8 614.6
w5(bu) 3280.3 3131.6 3126.2 3083.7
w6(bu) 805.7 882.4 887.5 745.8
14 2.8 111.1 40.7 100.1
15 2.5 3.6 3.9 6.9
16 342.4 281.6 295.6 405.9
Energy -76.681 20 -76.969 01 -76.966 77 -76.981 76


a)Ref. 93.

b)Ref. 133











Table 3-9: Geometries (in A and degrees), harmonic vibrational
frequencies (in cm-1), infrared intensities (in km/mol), and energies
(in Hartrees) for the C2v vinylidenic isomer of the S1 state of acetylene

P-EOM-
CISa CIS(D)a ) EOM-CCSDb
MBPT(2)
rcH 1.0847 1.0881 1.0879 1.0890
rcc 1.3877 1.4225 1.4174 1.4289
O(HCC) 122.94 122.15 122.38 122.28
wl(al) 3173.5 3086.1 3089.2 3060.1
w2(al) 1618.8 1464.2 1476.7 1486.5
w3(al) 1381.5 1225.7 1246.5 1213.1
w4(b2) 3235.5 3172.0 3171.0 3138.9
ws(b2) 1015.5 894.6 903.9 927.1
w6(a2) 977.1 702.3 708.6 763.2
I1 17.1 0.1 0.2 3.2
12 1.3 29.1 24.5 14.4
13 0.3 57.0 41.7 31.3
14 9.7 2.2 2.1 7.0
15 4.7 3.9 3.5 5.3
16 1.6 1.1 0.4 2.1

Dipole 3.817 2.926 3.029 2.875
moment
Energy -76.677 85 -76.967 63 -76.962 08 -77.003 07


a)Ref. 93.

b)Ref. 133











Table 3-10: Geometries (in A and degrees), harmonic vibrational frequencies (in cm- ), infrared
intensities (in km/mol), dipole moments (in Debyes), and energies (in Hartrees) for the equilibrium
geometry of the A 'A" state of formaldehyde. Experimental frequencies are fundamentals.

P-EOM- EOM-
CISa CIS(D)b P-EOM- EOM Expt.a
cIS CS(D)b MBPT(2) CCSDb Expt
rco 1.258 1.384 1.340 1.324 1.321
rcH 1.085 1.088 1.093 1.096 1.097
0(OCH) 117.3 117.8 117.0 115.8 118.0
7(HOCH) 148.4 171.5 153.0 145.9 148.3
wi (a') 3209.5 3154.1 3102.3 3072.9 2847
w2(a') 1437.9 1420.5 1412.7 1407.9 1290
w3(a') 1660.8 895.7 1187.2 1248.5 1173
w4(a') 561.6 172.6 564.9 719.8 683
w5(a") 3299.6 3295.8 3223.6 3189.5 2968
w6(a") 993.5 914.0 932.9 934.3 898
11 29.6 2.6 0.7 2.8
12 52.0 6.3 0.8 3.5
13 2.1 109.6 66.8 93.7
14 108.2 38.3 36.0 30.8
15 0.0 1.6 2.9 2.8
16 3.4 3.6 3.4 4.0
Dipole 1.507 2.374 2.094 1.776 1.56(7)c)
moment
Energy -113.669 47 -114.043 84 -114.028 60 -114.052 75


a)Ref. 102.

b)Ref. 93.

c)Ref. 141.











Table 3-11: A comparison of the geometries of the n--47r* states in formaldehyde,
acetaldehyde, and acetone. Geometries are in Angstroms and degrees.


Formaldehyde
rcn 1.340


rCHa
rCHb
rcHc
rCHd
0(oCC)
0(OCHa)
0(CCHb)
0(CCHc)
O(CCHd)
T(HaCOH)
T(HaCOC)
r(CCOC)
T(HbCCO)
r(HcCCO)
T(HdCCO)


1.093







117.0


Acetaldehyde
1.359
1.494
1.095
1.091
1.099
1.094
116.0
114.3
110.2
110.6
110.2


153.0


149.9


-177.5
61.6
-57.1


Acetone
1.376
1.494


1.092
1.101
1.094
113.4


110.3
110.2
110.6


146.9
-174.5
65.0
-53.7










68
Table 3-12: A comparison of two possible assignments of the low frequency vibrations in
the SI state of acetone to the calculated frequencies. All frequencies are in cm-1.


v12 torsion (antigearing)
V24 torsion (gearing)
V8 C-C-C bend
v23 C=O out-of-plane wagging
v19 C=O in-plane wagging


first
possibility
155.5
172.5
373


177.5


a)Ref. 140.


second
possibility
155.5
172.5
373
333
465.4


P-EOM-
MBPT(2)
206.8
190.1
353.8
329.7
368.4














CHAPTER 4
THE SIMILARITY TRANSFORMED EQUATION-
OF-MOTION COUPLED-CLUSTER METHOD


Introduction

The equation-of-motion coupled-cluster (EOM-CC) method,56 57 especially when

based on a coupled-cluster singles and doubles (CCSD)9 ground state, offers an attractive

and unified formalism to extend single reference coupled-cluster methods to excitation

energies (EE-EOM-CCSD, sometimes just EOM-CCSD),58 59 ionization potentials (IP-

EOM-CCSD),67 and electron affinities (EA-EOM-CCSD).68 The method is conceptually

single reference, in that the entire calculation is built upon one reference determinant.

Therefore, the only choices necessary in the method are the choice of basis set, the choice

of reference determinant, usually the Hartree-Fock solution, and the choice of excitation

level to include in the ground and excited state operators.

Several applications58, 142-145, 123 of EE-EOM-CCSD and its twin, coupled-cluster

singles and doubles linear response,53, 54 have been made. The method has also been

extended to include various effects of triple excitations.6-63 But EOM-CCSD, and

especially EE-EOM-CCSD, has some drawbacks. The first is the cost of the method.

Every EE-EOM-CCSD excited state is approximately as expensive as a ground state

CCSD calculation. Approximations can be made to reduce the cost,114 with varying

degrees of success. The second drawback is the accuracy of the method. For ethylene,

butadiene, and cyclopentadiene, the average error for the Rydberg states was 0.17 eV,

but the average error for the valence states was 0.26 eV.123









70

The key to EOM-CC theory is the similarity transformation of the Hamiltonian. After

the ground state coupled-cluster equations have been solved, the cluster amplitudes are

used to transform the Hamiltonian. A similarity transformation changes the eigenvectors

of an operator without changing its eigenvalues. Hence, the energies of the excited states,

electron attached states, and ionized states, which are all eigenvalues of the Hamiltonian

matrix in terms of the appropriate bases, have not changed.

So far nothing has been accomplished. To calculate the energies of these states, it is

still necessary to diagonalize a matrix that is, in principal, infinite. After a finite basis set

has been introduced, the matrix is now finite, but its size grows as the size of the basis set

factorial.22 Therefore, approximations must be introduced. The typical approximations

involve truncating the Hamiltonian matrix to include only a subset of possible excitations.

In EE-EOM-CCSD the matrix is limited to only single and double excitations, i.e.

determinants where one or two virtual orbitals replace one or two occupied orbitals. For

EA-EOM-CCSD the space is limited to determinants where an electron has been added

to a single virtual orbital and determinants where along with the electron being added

to a virtual orbital, an occupied orbital is replaced with a virtual orbital. These are the

so called lp and 2plh states. For IP-EOM-CCSD the space consists of the lh and 2hlp

determinants, which are those determinants with one electron removed from an occupied

orbital and those determinants with one electron removed and one electron excited.

After the truncation of the matrix is when the similarity transformation becomes

important. The truncation will change the eigenvalue spectrum. Many eigenvalues

will no longer appear, and those that remain will be shifted by an amount which will

depend on how much of the true eigenvector lies in the omitted space. What gives the










71
similarity transformation its power is that by choosing the transformation so that the

important eigenvectors are contained, to the greatest extent possible, in the space kept

after the truncation, the effect of the truncation on the energies of the states of interest

can be minimized. Note that an similarity transformation is needed, since a unitary

transformation leaves the eigenvectors unchanged. An equivalent picture is to choose

the transformation such that the parts of the matrix that couple the kept space with the

excluded space are minimized, since these parts of the matrix determine the extent to

which eigenvectors span both spaces.

For singly excited states, the most important terms deleted in the truncation to singles

and doubles are the triple excitations. In the untransformed Hamiltonian matrix, the terms

that couple single excitations to triple excitations are the two electron integrals (abllij),

where a, b, c, and d are unoccupied in the reference determinant and i, j, k, and I are

occupied in the reference determinant. These are precisely the terms which are set to

zero by the similarity transformation. The transformation does add new coupling terms,

but they are of the form of a t2 amplitude times a different two electron integral. Since

for a normal single reference problem, the largest t amplitudes are of the order of 0.1, the

new coupling terms are at least an order of magnitude smaller than the former coupling

terms. Therefore the "true" eigenvector for a singly excited state is contained more

fully within the space of single and double excitations, and the eigenvalue for the state

with the truncated similarity transformed Hamiltonian is much closer to the eigenvalue

of the untruncated Hamiltonian than is the eigenvalue for that state with the truncated

untransformed Hamiltonian.









72

In the similarity transformed equation-of-motion coupled-cluster (STEOM-CC)

theory64-66 the blocking of the matrix by similarity transformations is carried one

step farther. Starting with the EOM-CC transformed Hamiltonian, a second similar-

ity transformation69 is performed in order to minimize the coupling between the singles

space and the doubles space. The dominant terms that couple these two blocks in the

EOM-CC Hamiltonian are the (ablIci) and (ialljk) two electron integrals. The most

important of these coupling terms are set to zero in the second similarity transformation,

leaving as the main coupling term between the singles and doubles blocks again a t2

amplitude times a two electron integral. Once again, the coupling is reduced by at least

an order of magnitude.

Since the second similarity transformation is done in such a way as to preserve

the reduction of the singles-triples coupling from the first similarity transformation, the

eigenvectors for states dominated by single excitations are almost exclusively composed

of only single excitations, meaning that the truncation needs to leave only the space of

singles. Now the final diagonalization to get the eigenvalues and eigenvectors only scales

as the fourth power of the basis set, as opposed to the sixth power scaling for EOM-

CCSD. The trade-off is having to perform the second similarity transformation, but the

second transformation only scales as the fifth power of the size of the basis set.

If the final diagonalization is over the space where an electron has been removed

from an occupied orbital and placed in an unoccupied orbital, then excited states are

calculated, and the method is known as excitation energy STEOM-CC (EE-STEOM-

CC).64 66 If the final space is over states with two electrons removed, the method is

the double ionization potential STEOM-CC (DIP-STEOM-CC).66 Finally, if the space is









73
over determinants with two electrons added, the method is the double electron attachment

STEOM-CC (DEA-STEOM-CC).66 Each of these methods provides an important tool for

describing interesting chemistry, and although the rest of the dissertation will focus on

the EE-STEOM-CC variant, everything said can be applied to the other two with simple

modifications.


The STEOM-CC Method

In STEOM-CC theory, the orbital space is divided into occupied and virtual orbitals,

sometimes referred to as holes and particles, based on the orbital's occupancy in the

reference determinant. The orbital space is then subdivided into active and inactive

occupied orbitals and active and inactive virtual orbitals. To distinguish the sets, the

labels i, j, k, and I will refer to generic occupied orbitals, and a, b, c, and d will refer

to generic virtual orbitals. Active occupied orbitals will be m and n while active virtual

will be e andf. Explicitly inactive orbitals will be indicated with a prime, such as i' and

a'. The labels p, q, r, and s can refer to any orbital. All such orbitals are spacial orbitals,

as opposed to the previous chapters where spin orbitals have been used.

Although a STEOM-CC calculation could be based on any coupled-cluster or many-

body perturbation theory treatment of the ground state, so far only CCSD and MBPT(2)

have been used.64 This discussion will focus on STEOM-CC with a CCSD9 ground

state, giving the STEOM-CCSD method.64 In coupled-cluster theory, the ground state

of the molecule is represented as


IWo) = erl),


(4-1)









74
where 10) is a closed shell single determinant, often the Restricted Hartree-Fock (RHF)

determinant. For CCSD the operator T consists of pure excitation operators of the form

T = Ti + T2 = {a + atibtj. (4-2)
i,,
The curly brackets mean that the operator is normal ordered with respect to 10).

The untransformed Hamiltonian, in second quantization, is

H = ho + E fpq pfq + Vqprs p,{rqs. (4-3)
p,q Pq
The constant term, ho, is the energy of the reference determinant. The fpq refer to the

appropriate parts of the Fock operator. The Vqprs are non-antisymmetrized two electron

integrals.

After the first transformation the Hamiltonian becomes
H = e-THe

Sho + hfpq pfq} + pqrs.{prqs}
p,q P'q (4-4)
+36 hpqrstuPlsrqttrtu. + ....


With a CCSD reference state, H will contain up to six-body terms, but the four- and higher

body terms do not appear in the EOM-CCSD or STEOM-CCSD equations. Explicit

equations for these H elements, as used here, can be found elsewhere.146 In terms of

H, the CCSD equations can be written as

(I1H >0) = lh, = 0.,
(4-5)
(jiflH10) = /bi = 0.
The CCSD energy is


ECCSD = (0oHl0) = ho.


(4-6)










75

Thus, solving the CCSD equations is equivalent to setting the one- and two-body pure

excitation parts of H to zero.

In STEOM-CC, a second many-body similarity transformation69 is performed in order

to also set to zero selected remaining terms which increase the net excitation level by

one. But since, in diagrammatic language, the necessary operators have a line on the

bottom, they do not commute and can connect to each other. Thus instead of using

a transformation of the form e as originally done by Stolarczyk and Monkhorst,37 a

normal ordered exponential operator {es} is used. The normal ordered exponential,

introduced by Lindgren,33 excludes terms inside the normal ordering from connecting to

each other, and this simplifies the equations.

The S operator used in the transformation consists of two parts S = S+ + S-, where


S+ = S+ + S+ i< {a'te} + + sa {at j}, (4-7)
a',c a.b
Je,

and

S- = S + S2- : Zs{m?} + 1 bs.i{mbtj}. (4-8)
hi',mj

It is the presence of the active q-annihilation operators which cause the components of

S to, in general, not commute.66

The double similarity transformed Hamiltonian takes the form

G = es H es}

= go + gpq{pq} + pqrs {ptrqs} + .... (4-9)
p,q PIq
rs









76

The similarity transformation preserves the zeros for the one- and two-body pure exci-

tation parts, so that66

9ai = hai = 0,
(4-10)
gabij habij = 0.
Also,

90 = ho = ECCSD. (4-11)

In principle, like the CCSD equations for T above, the S equations could be derived

by setting the appropriate parts of G to zero. These are66

9mi' = (i, jIGIm) = 0,

gae = (a' IGIe) = 0,
j (4-12)


gabej = (a IG) = 0.

Note that the second two equations set to zero the primary coupling terms between singly

and doubly excited determinants. In practice, the S coefficients are not calculated this

way. Doing so would lead to systems of nonlinear equations which in some situations can

be numerically unstable.64 Instead, the equivalence of IP-EOM-CC and EA-EOM-CC to

Fock space coupled-cluster theory'47 is exploited to rewrite S in terms of IP-EOM-CCSD

and EA-EOM-CCSD eigenvectors. These equations are66





m b (4-13)

ba= ba -1(A)r

A









77
The A's stand for the set of IP- and EA-EOM-CCSD eigenvectors. There is one

eigenvector per active orbital. The terms r-1 stand for the coefficients in the inverse of

matrix constructed from the active components of the eigenvectors A.66

After the S coefficients are determined, the problem remains of how to determine

the G coefficients. Since the inverse of the normal ordered exponential is not known,69

Eq. (4-9) cannot be used directly. Instead the equation


{es}G = {es} (4-14)

is used. It has been shown69 that this equation is equivalent to


(IeS G) = (ft eS}) (4-15)

or

G= (l{eS}) (es 1 G), (4-16)

where the c means that the equations must be connected.

In practice, all terms involving S1 are dropped.66 Since S1 does not change the net

excitation level, it is not involved in the decoupling of excitation blocks. Ignoring all

of the S1 terms is equivalent to a simple rotation of the eigenvector, and since the G

matrix is diagonalized over all singly excited determinants this "rotation" is automatically

compensated for in the diagonalization. Therefore, including the Si terms does not change

the final answer, and for convenience they are left out. Leaving out S1 also simplifies

Eq. (4-16). Since only the singles-singles block of G is needed, the second term does

not contribute,66 leaving straightforward, linear equations for G in terms of H and S.

The G equations in orbital form can be found in Ref. 66.









78
Finally, the G matrix is diagonalized over the space of single excitations to get the

excited states. This is actually the fastest step in the process, since the final matrix to

be diagonalized is so small. Afterwards, excited state and transition properties can be

calculated, subject to certain approximations.66


Discussion About STEOM-CC

Formally, the introduction of active and inactive orbitals is not needed in STEOM-

CC. All orbitals could be considered active. It would simply mean that S would have

no one-body component and that an IP-EOM-CCSD eigenvector would need to be found

for each occupied orbital and an EA-EOM-CCSD eigenvector would be needed for each

virtual orbital. In practice, this would never work, except for trivial cases. For example,

there is normally a one-to-one relation between IP-EOM-CCSD states and high lying

occupied orbitals, but as the spectrum moves into the inner valence region, shake-up

states appear.13 At that point, many states appear per orbital, and even the concept

of principal ionization potential becomes unclear. The same will happen for EA-EOM-

CCSD. Eventually, the iterative diagonalization procedure99 will have trouble converging.

Even if it were possible to find an IP-EOM-CCSD or EA-EOM-CCSD eigenvector

for every orbital, it would not be desirable. The cost of each IP-EOM-CCSD state is

occ virt times the number of iterations needed for the diagonalization. The cost of each

EA-EOM-CCSD state is n, cNv i per iteration. If the number of states equals the number

of orbitals, then the cost of calculating S becomes as large as the cost of calculating T.

It is not necessary to include all of the orbitals in the active space, anyway. The second

similarity transformation serves to describe the differential dynamic correlation between

the ground and excited states. But relatively few orbitals actually are involved in the









79

excitation, and it is only those orbitals involved in the excitation that need to be in the

active space.

The relationship of STEOM-CCSD to both Fock space coupled-cluster (FSCCSD)

theory3245 and EOM-CCSD58 59 have been discussed at some length.66 Only the

highlights will be reviewed here. If the active space includes all orbitals, then STEOM-

CCSD and FSCCSD are numerically identical for singly excited states.66 For doubly

excited or higher states, the two remain quite different. If the active space does not

include all orbitals, then the two are different for singly excited states. The first difference

is that STEOM-CCSD can be viewed as a series of successive similarity transformations,

whereas FSCCSD has just one transformation. The second difference is that in STEOM-

CCSD the final diagonalization is over all singly excited determinants, while the final

diagonalization in FSCC is over those determinants where only the active orbitals change

occupation.

To consider the difference between STEOM-CCSD and EOM-CCSD, it is useful to

fully expand out the STEOM-CCSD wavefunction. It is

IIr)=e {es}RI0) eo 0)ro

=eT(l+S+ S2)R1 0) +eT ) 'ro (4-17)

where RI is the single excitation part of the STEOM-CCSD eigenvector and r0 is the

contribution of the reference determinant. This should be compared to the EOM-CCSD

wavefunction, which is


I) = eT(R + R2)10) + eO)ro.


(4-18)









80
Since S1 does not change the excitation level, and is not included in the calculation

anyway, it can be ignored in this discussion. It is possible to directly equate terms in

Eqs. (4-17) and (4-18). They match up as follows:

eTI0)ro eT0)ro

eTRIO0) eTjlZ 10)

eTS2R110) eT210) (4-19)

eTS2R10) <
2
The STEOM-CCSD term is on the left and the EOM-CCSD term is on the right. The

primary differences between the two methods are therefore that the doubles term in

EOM-CCSD has the flexibility to determine its optimum weight, while the doubles term

in STEOM-CCSD is a product of a fixed R1 and a fixed S2. On the other hand, the

only triple excitations included in the EOM-CCSD wavefunction are of the form eTR,

while the STEOM-CCSD also includes a S22RI term. Thus STEOM-CCSD has a more

complete description of triple excitations but gives up some flexibility in describing the

doubles.66

The cost of a STEOM-CC calculation in almost dominated by the cost of calculating

the ground state CCSD solution and the ground state left hand state, lambda.107 Lambda

is needed for transition properties.66 Both of these are iterative n2 4Ni steps. The next

most expensive step is the formation of H, which also scales as ncc Nirt but is done only

once. The costs of the IP-EOM-CCSD and EA-EOM-CCSD calculations were discussed

above. Calculating the S coefficients will scale as occNirNirt,aN c where Nvirt,act is

the number of active virtual orbitals. Forming G will scale as noccN irtNvirt,act. The

final diagonalization is only an iterative N4irt step.









81

What all of these scalings imply is that the cost of the calculation is effectively

independent of the number of excited states desired. In fact, the only practical effect that

the number of states has on the cost is that as higher lying excited states are calculated,

more orbitals will need to be included in the active space. This independence of the cost

of the calculation to the number of excited states is in contrast to almost every other

excited state method in use today, and it makes STEOM-CC well suited to calculating

large numbers of excited states.

However, STEOM-CC does have some drawbacks. The primary one is that the

method is only suitable for states dominated by single excitations. It would be possible

to derive a STEOM-CC variant whereby the final diagonalization takes place over the

space of both singles and doubles. This would allow the description of both doubly

excited states and states which are mixed singles and doubles. The problem is that the

equations become much more complicated. For example, the second term in Eq. (4-16)

would contribute. It would also be at least as expensive as an EOM-CCSD calculation.

STEOM-CC also has strong and poorly understood dependencies on the size and the

nature of the active space. It is known that the component of the excited state vector

which lies within the active space should be at least 95-99%. But it is not clear how

the error in the energy and properties depends on that percentage. This is an area that

will require much more practical experience before it becomes clear how best to choose

the active space.

Finally, the calculation of excited state and transition properties in STEOM-CC is

unsatisfactory. Several approximations66 are currently made in the properties calculations.

There is some hope, at least for the excited state properties, in that the advent of










82

gradients for STEOM-CC will give a way to calculate one electron excited state properties

rigorously.














CHAPTER 5
THE SPECTRUM OF FREE BASE PORPHIN


Introduction

Because of their importance in such biological processes as photosynthesis, electron

transfer, and oxygen absorption and transport, the porphyrins have been extensively

studied.148 As the base molecule for the porphyrins, the electronic spectrum of free

base porphin has received much attention, both with semiemperical (see, for example,

Ref. 149) and ab initio methods.150, 75, 151-153, 65 The interesting part of the spectrum

consists of, in order, two visible peaks known as the Q bands, a very intense peak known

as the B (or Soret) band, a shoulder on the B band, called the N band, and two other

small peaks, the L and M bands. The spectrum can be found in Ref. 154 and Ref. 153.

The traditional interpretation of the spectrum is that, with the molecule in the xy plane

and with the two internal hydrogens along the x axis, the lowest energy band, the Qx band,

comes from exciting to the 1'B3u state. The Qy band then comes from the I B2u, with

the B band assigned to the 2'B3u and 2'B2u states.150 Nakatsuji et al.,153 based on their

SAC-CI (symmetry adapted cluster-configuration interaction) calculations,'55 reassigned

the spectrum. They agreed with the assignments of the Q bands, but they claimed that the

B band should be assigned to only the 21B3u state, with the N band being the 2'B2u state.

However, their results have three significant weaknesses. The first is the method used.

In principle SAC-CI energies (but not oscillator strengths) could be equivalent to EOM-

CCSD (equation-of-motion coupled-cluster singles and doubles) 58 and coupled-cluster

singles and doubles linear response53 energies. In practice, though, some non-linear terms









84
in the underlying ground state coupled-cluster result are always omitted. Also, in these

calculations many double excitations were omitted based on a perturbation selection.153

For complex organic molecules even the untruncated EOM-CCSD may not be sufficient.

In a study of benzene and the azabenzenes, EOM-CCSD had an average error of 0.32

eV for the 7r--+7r* states.156 SAC-CI should do no better than this, unless it has some

fortuitous error cancellations.

The second criticism is the basis set used. Nakatsuji et al.153 discuss the importance

of a rearrangement to the excitation energies, yet their basis set only had 2s type and

2p type functions on the carbons and nitrogens, giving it no flexibility to describe the

2s orbitals. They also dropped some of the occupied and virtual orbitals corresponding

to the 2s orbitals, along with all of the Is orbitals on the carbons and nitrogens. The

result is that the basis set had limited flexibility, no polarization functions, and no diffuse

functions. The lack of diffuse and polarization functions corresponds to the conventional

viewpoint that, for at least the lowest states of free base porphin, polarization and diffuse

functions are not needed.152

Finally, there is a problem with the oscillator strengths. The N band appears as a

shoulder to the B band, but Nakatsuji et al.153 calculate the excitation to the 21B2u state to

have an oscillator strength sixty-eight percent larger than the excitation to the 2'B3u state.

They argue that the N band is actually quite broad, with the B band being a narrow peak

on top of it. But if that is true, then the splitting between the vertical excitation energies

for the B and N bands would be less than the 0.32 eV reported.154 To answer some of

these questions, a series of STEOM-CCSD (similarity transformed equation-of-motion

coupled-cluster singles and doubles)64 excited state calculations have been performed.









85

Computational Details

STEOM-CC

In a STEOM-CCSD64 66 calculation (i.e. a STEOM-CC calculation64 based on a

CCSD9 ground state), two similarity transformations are applied to the second quantized

Hamiltonian, such that the one- and two-body terms in the Hamiltonian that increase

the excitation level are set to zero. This effectively blocks the Hamiltonian matrix by

excitation level, so that the single excitations can be accurately calculated with a diago-

nalization over just the singles-singles block. The higher excitation parts of the excited

state wavefunction are then implicitly included through the similarity transformations. In

this way, STEOM-CC has the same conceptual appeal as CI singles, but now in a fully

correlated structure. The STEOM-CCSD method has been implemented within ACES II.

The first similarity transformation, the same as in EOM-CCSD,58 involves the T

amplitudes from a ground state CCSD9 calculation. The new Hamiltonian is then

H = e-THeT. (5-1)

The next step is to solve for a set of states with one electron removed and one electron

added by means of the IP-EOM-CCSD (ionization potential EOM-CCSD)67 and the EA-

EOM-CCSD (electron attachment EOM-CCSD)68 methods. These involve diagonalizing

H over the space of lh-2hlp and lp-2plh determinants, respectively. These eigenvectors

are used to determine the S coefficients in the second similarity transformation


G = eS1 }. (5-2)

The final double similarity transformed Hamiltonian G is then diagonalized over the space

of single excitations. For a detailed description of the method, see Ref. 66.










86

In the current implementation the final wavefunction in terms of G has only single

excitations, limiting the method to only singly excited states, but in terms of the normal

Hamiltonian, the wavefunction actually includes all possible excited determinants (with

only the singles having optimized coefficients). Specifically, the wavefunction contains

double excitations in the form of S2R1 and T1R1, and triple excitations in the form of

S2R1, T2R1, and T1S2R1, where R1 is the STEOM-CC single excitation eigenvector.

Since the higher excitations and the differential correlation between the ground

and excited states are described via the S operator, it is critical that the excitation be

described within the set of active orbitals chosen for the IP-EOM-CC and EA-EOM-CC

calculations and included in S. In the calculations presented here, the active components

of the excitations are almost always above 99% and in all singlets are above 98%. In

test calculations the energy seems to be converged to within 0.05 eV when the active

component is 98%.64

Basis Set and Geometry

The basis set used for the calculations presented here come from the large ANO

basis set of Widmark, Malmqvist, and Roos.157 This is a very large, generally contracted

basis set. It has 14s9p4d3f primitives for carbon and nitrogen and 8s4p3d primitives for

hydrogen. From this set, the first 3s and 2p contractions were used for C and N and the

first s contraction was selected for H. This basis set, consisting of 230 functions, was the

same as used by Merchin et al.152 It is the same number of contracted functions as that

used by Nooijen and Bartlett,65 but their basis set had a much smaller set of primitives.

This set is also significantly larger than the one used by Nakatsuji et al.153









87
This basis set was then extended in two different ways. First, to gauge the effect of

polarization functions on the excitation energies, the first d contraction from the ANO set

was added on the carbon and nitrogen atoms, while the second s function was added to

the hydrogens. This gave a [3s2pld] set on C and N and a [2s] set on H, with a total of

364 contracted functions. In other calculations, to gauge the effect of diffuse functions

on the excited states, a set of 2s and 2p uncontracted functions were added to the center

of the molecule and to the geometrical center of each ring. The diffuse exponents were

taken from Ref. 158, and have been used for naphthalene158 and biphenyl,159 where they

were placed at the center of the molecule. Since two functions of each type were used,

the exact exponent chosen should not matter significantly. In all calculations the first

24 occupied orbitals, corresponding to the ls orbitals on carbon and nitrogen, were left

uncorrelated.

Two different geometries were used in this study. The first is an idealized x-ray

structure,160 where the molecule, without the two internal hydrogens, would be of D4h

symmetry. The hydrogens make the structure D2h. It was used in several of the previous

studies153, 65 and will be referred to as "D4h". With the [3s2p/ls] basis set, this geometry

had a SCF energy of -982.955034 Hartrees and a CCSD energy of -985.154443 Hartrees.

The [3s2pld/2s] basis set gave a SCF energy of -983.416854 Hartrees and a CCSD

energy of -986.676253 Hartrees. The basis with the diffuse functions gave a SCF energy

of -982.967070 Hartrees and a CCSD energy of -985.182042 Hartrees.

The other geometry is a B3PW91/6-31G* optimized geometry,'61 which will be re-

ferred to as "Opt". It also has D2h symmetry. This geometry gave a SCF energy of

-982.964858 Hartrees and a CCSD energy of -985.163441 Hartrees with the [3s2p/ls]









88
basis set. The diffuse basis gave a SCF energy of -982.977085 Hartrees and a CCSD en-

ergy of -985.191231 Hartrees. The [3s2p d/2s] basis gave a SCF energy of -983.430681

Hartrees and a CCSD energy of -986.688179 Hartrees, putting this geometry 7.5 kcal/mol

below the D4h geometry.


Results

Ionized and Electron Attached States

In all of the current calculations, the occupied part of the active space consisted of

eleven orbitals. Therefore, eleven IP-EOM-CCSD states were calculated. In Table 5-1

the current results, with the three basis sets, are compared with the SAC-CI results,153

with the previous STEOM-CC results,65 and with experiment.162 Because of the low

resolution of the experimental spectrum, it is very difficult to relate the measured peaks

to the calculated states. Therefore, the assignments given should be considered tentative.

While adding the diffuse functions does little to change the ionization potentials,

every other time that the basis set was enlarged, from the previous STEOM-CC results,65

to the DZ (the [3s2p/ls] basis) results, to the polarized (the [3s2pld/2s] basis) results, the

electron became more bound. Adding the polarization functions had the most dramatic

effect. It caused the first two states to switch and the sixth and seventh states to switch.

However, the differences are close to an order of magnitude smaller than the errors in

the calculations, so nothing can be said definitively. Going from the "D4h" geometry

to the "Opt" geometry had effects of less than 0.2 eV. The "Opt" geometry typically

had larger IP's.

Twenty-three virtual orbitals (4 ag, 2 big, 3 b2g, 3 b3g, 3 au, 2 blu, 3 b2u, and 3 b3u)

were included in the active space for the DZ and the polarized basis sets. For the diffuse









89

basis set several more orbitals had to be included. The diffuse basis added forty Rydberg

orbitals, and most of them had orbital eigenvalues very close to zero. Therefore, to

include all of the equivalent orbitals that were included for the DZ basis, the active space

needed to consist of fifty-nine virtual orbitals. But that many orbitals caused convergence

problems. The final set included in the active space consisted of fifty-one virtual orbitals

(10 ag, 4 big, 5 b2g, 5 b3g, 3 au, 8 blu, 5 b2u, and 5 b3u). If an IP-EOM-CC or an

EA-EOM-CC eigenvector is less than 70% singles, the code automatically excludes it

from the second similarity transformation. For the diffuse basis, one of the B3g states

was excluded, leaving fifty states for the similarity transformation.

In each of the calculations two positive electron affinities were predicted. They are

listed in Table 5-2. The other electron attached states, even though they do not correspond

to stable states, are still essential for the calculation of the excited states. They go into

the second similarity transformation and help describe the differential correlation between

the ground and excited states.

Excited States

The singlet excited states of free base porphin are listed in Tables 5-3 to 5-7, and the

triplet excited states are in Tables 5-8 to 5-12. The important states will be discussed

later, but some general comments can be made. The first is that the Rydberg states

start about 4.5 eV, which is right in the region of most interest. The second is that the

addition of the diffuse functions never changes the energy of the valence states by more

than 0.02 eV. This has several consequences. It implies that the DZ basis set has diffuse

enough tails that the diffuse functions are not needed to describe the valence region.

It also means that there is essentially no mixing between Rydberg and valence excited









90

states. Finally, it suggests that adding diffuse functions to the polarized basis would not

significantly change the energetic. Adding diffuse functions might, however, have an

effect on the oscillator strengths. Those tend to change a little with the addition of the

diffuse functions.

In general, going from the "D4h" geometry to the "Opt" geometry has little effect on

the spectrum. The excitation energy typically increases by about 0.1 eV. The important

exceptions are the IB2u state which drops from 5.22 eV to 5.07 eV and the 1B3u state

which drops from 5.33 eV to 5.23 eV.

A note should be made about the oscillator strengths. The oscillator strengths are

calculated as


f = l,< |/],)( [, (5-3)


where w is the excitation energy. The right-hand ground state is the coupled-cluster

wavefunction, the left-hand ground state is the lambda solution from coupled-cluster

theory,107 and the right-hand excited state is the STEOM-CC state. Currently, several

approximations are introduced when calculating the left-hand excited state.64 These

sometimes can cause numerical problems with the properties calculated. When the

oscillator strength is listed as (-), it means that the approximations are too severe,

and the calculated oscillator strength is unreliable. Since the energies are calculated

with the right hand excited state wavefunction, they are still correct. For the polarized

basis results, another approximation is made. Because the cost of calculating the ground

state lambda vector is prohibitive for that large a basis set, the left-hand ground state is

estimated with the same approximations as the excited state left-hand wavefunction. In









91

a test calculation with the DZ basis, this makes a difference in the oscillator strengths

of less than 20% in every case.

The "D4h" geometry, polarized basis results for the singlet valence states from Tables

5-3 to 5-7 are listed in Table 5-13, along with CASPT2,152 SAC-CI, 153 previous

STEOM-CC,65 and EOM-CCSD results. The EOM-CCSD calculation is at the "D4h"

geometry with the DZ basis. The valence states are numbered for convenience, but these

numbers are only accurate through 3'Big. After that the Rydberg states should enter

into the numbering.

The energy of the first triplet state has been measured in a solvent/ethyl iodide mixture

at 77K.163 The phosphorescence peak is at 1.58 eV, well above the STEOM-CC polarized

basis result of 1.26 eV. The CASPT2 result is between the two at 1.37 eV.152


Discussion

The recent controversy over the assignment of the spectrum centers around how to

assign the N band. The argument of Nakatsuji et al.,153 essentially, is that they calculated

no other optically allowed states in the 3-4 eV range, and therefore, by default, the N

band must be 21B2u. Table 5-13 shows the poor quality of their results. Had they used a

decent basis set and had they not made the approximations that are always used in their

SAC-CI calculations,153 their calculated energies would have approached the EOM-CCSD

excitation energies, meaning that they all would have increased. At that point, it becomes

impossible to draw any meaningful conclusions about the assignment of the spectrum.

The CASPT2152 results for the B band are too low. This could have been caused by

the size of their active space. It has been shown several times164, 151, 153, 65 that for the

21B3u and 2'B2u states, Gouterman's four orbital model165-167 is not sufficient; the 4blu









92

orbital must also be included. Merchan et al.152 did not include it in their active space.

More CASPT2 calculations, with a larger active space and with d functions, would be

very informative.

The Qx and Qy bands belong to the 11B3u and 1 'B2u states. The current results show

fortuitously good agreement with the Qy band. For the Qx band all of the methods predict

too low an excitation energy. A fluorescence spectrum of free base porphin taken in a

supersonic jet expansion168 placed the 0-0 transition for the Qx band at 2.0234 eV and

the 0-0 transition for the Qy band at 2.4653 eV, slightly higher than the vertical excitation

energies reported by Edwards et al.154 The addition of polarization functions increased

the Qx excitation energy by 0.05 eV, so the use of much larger basis sets may further

increase the calculated excitation energy, moving it towards the experimental number.

The current results strongly support the original assignment for the B band being

both the 21B3u and 21B2u states. The N band is then assigned to the 3'B3u state. If

the assignments of Nakatsuji et al. were correct, it would mean that the polarized basis

STEOM-CC energies would have to be 0.19 eV too low for the first L peak and 0.32

eV too low for the second L peak. It is highly unlikely that the STEOM-CC would

consistently underestimate the excitation energy like that. Instead, this assignment puts

21B3u 0.14 eV above the B peak and 2'B2u 0.29 eV above the B peak. The 3'B3u state

is 0.41 eV above the N peak. This error is uncomfortably large. However, the energies

of all of these states dropped substantially when the polarizations functions were added.

Larger basis set calculations should further decrease these gaps.

The intensities still present a problem. The experimental oscillator strength of the

N band is less than 0.1, and the calculated oscillator strength is 0.93. Edwards et al.154









93

mention that the shape of the B and N bands would be consistent with the two states of

the B bands being split by 1500 cm-1 and an intense N band donating intensity into the

B band. These calculations give a splitting of 1200 cm-1 between the 21B3u and 21B2u

states. Thus these calculations agree with the possible interpretation given by Edwards

et al.154 A B band splitting of 240 cm-1 was measured in a low temperature crystal

spectrum and assigned to the energy difference between the two electronic states,169 but

this assignment has been disputed.170

These assignments then leave the two peaks of the L band assigned to the 3'B2u

and 4'B2u states. The diffuse basis energy for the 3'B2u state agrees well with the first

L peak, but the 41B2u state is 0.33 eV higher than the second L peak. This difference

should also be reduced with larger basis sets. The experimental oscillator strength for the

two states combined is about 0.1. The calculated oscillator strengths are well above that.

The 1'Blu state sits under the L peak, but its intensity is so small that it is not visible.

The one remaining problem with this assignment is the M peak. The calculated

excitation energy for the 4'B3u state is 5.17 eV, 0.33 eV lower than the M band maximum.

The intensities match well, though, and there is no other assignment that makes sense.

To say that 4'B3u is part of the L band would require that the excitation energy to drop

by 0.5 eV, and it would cause problems assigning the other states. Still, it is not clear

why the calculated excitation energy would be so low.

Finally, a point should be made about the cost of these calculations. For the polarized

basis calculation there are 57 occupied and 283 virtual functions. This gave 2479 T, and

34 170 895 T2 amplitudes. The time taken for each step on a Cray C90 is listed in Table

5-14. It took only 36 seconds per excited state to calculate the energies and properties.










94

Most of that time was used calculating the properties. The energies took 1.4 seconds

each. This is a very important advantage of the STEOM-CC method. The total time for

the calculation is effectively independent of the number of excited states calculated. That

makes it possible to study entire spectra instead of just a few selected states.


Conclusions

STEOM-CCSD64 calculation on free base porphin using a [3s2pld/2s] basis set are

presented. These are the first reported calculations on excited states above the lowest

that use polarization functions. For the important optically allowed excited states, the

polarization functions have a significant effect. These calculations strongly support the

traditional interpretation of the spectrum, in that the intense B band is assigned to both

the 2'B3u and 21B2u states.

This study is possible only because the second similarity transformation in STEOM-

CCSD means that the excited states can be represented in terms of single excitations

with respect to a highly modified Hamiltonian. Thus it is possible to efficiently calculate

many excited states at once.




Full Text
xml version 1.0 encoding UTF-8
REPORT xmlns http:www.fcla.edudlsmddaitss xmlns:xsi http:www.w3.org2001XMLSchema-instance xsi:schemaLocation http:www.fcla.edudlsmddaitssdaitssReport.xsd
INGEST IEID E8WDYC52Q_IH29RS INGEST_TIME 2013-10-10T02:53:26Z PACKAGE AA00017654_00001
AGREEMENT_INFO ACCOUNT UF PROJECT UFDC
FILES