Model Hamiltonians for the calculation of atomic and molecular spectroscopy


Material Information

Model Hamiltonians for the calculation of atomic and molecular spectroscopy
Physical Description:
vi, 169 leaves : ill. ; 28 cm.
Baker, John David, 1961-
Publication Date:


bibliography   ( marcgt )
theses   ( marcgt )
non-fiction   ( marcgt )


Thesis (Ph. D.)--University of Florida, 1991.
Includes bibliographical references (leaves 161-168).
Statement of Responsibility:
by John David Baker.
General Note:
General Note:

Record Information

Source Institution:
University of Florida
Rights Management:
All applicable rights reserved by the source institution and holding location.
Resource Identifier:
aleph - 001683425
notis - AHZ5396
oclc - 25034737
System ID:

Full Text








I would like to thank first and foremost my wife, Rebecca, for her under-

standing and support that was crucial to the success of this work. Secondly, I

would like to thank Tennessee Eastman for the financial support that accelerated

this project. A special thanks goes to my committee who have struggled with

me to get many recommendations and resumes in the mail. I thank my advisor,

Dr. Zerner, for allowing me the freedom I desire to explore my own ideas while

supporting me in the process. This, I feel, is truly unique in the academic business

and I certainly realize the personal sacrifice he has made to allow his guys this

privilege. And finally, I would like to thank the guys in the clubhouse who put

up with my songs and occasional bouts of sarcasm.


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

ABSTRACT.................................... .. v


1. INTRODUCTION ................... .............. 1
Motivation .... .. .............. ............... 1
Time Dependent Perturbation Theory ........... ....... 5
Molecular Response to Radiation ...................... 7
The Reporting of Theoretically Calculated Spectral Probabilities .. 13

2. GENERAL THEORY ...................... ... 20
The Hartree-Fock Equations . . ... 20
The Roothaan Equations . .. . 27
The Calculation of Electronic Spectra ...... .......... 31
Singlet Excitations with a Finite Basis ............... 41

3. THE INDO/RPA METHOD ......................... .44
Theory ..... ...... ..... .............. .. .. 44
Benzene and Pyridine . .. . .. 51
Naphthalene and Quinoline ....... ............ ... .. 58
The Diazines .................. ................ 61
The RPA for Extended Systems ... ..... ........... 65

PORPHIN .................... ............... 70
Motivation ..... ............ ................ 70
Results ....... .. ....... .................. 73
Discussion .. ..... ........................ ... 76



The Concept of Chemical Shift ......................

Theory . . . . .

Localization and Integral Evaluation ...................

Reduced Linear Equations ........................

Diamagnetic Shielding in Molecules ................ ....

The INDO/RPA Method for Chemical Shift ..............

PROPERTIES ............................ .....

Im petus ... .. ..... .... .

Atomic Parameterization ...............

A New Atomic Parameterization ..........

A analysis ... .. .. ... .

ELECTRONIC SPECTRA ...............

Density Partitioning ..................

Atomic Hamiltonian .................

Molecular Hamiltonian ................

Modified SCF Procedure ...............

Results .. .. . .

8. CONCLUSIONS......................

Spectroscopy............ ..........

Development of a Better Model .... ......

BIBLIOGRAPHY .....................





. . .






. . 161

. .. .. .. 169

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




May, 1991

Chairman: Michael C. Zerner
Major Department: Chemistry

Semiempirical models appropriate for the calculation of atomic and molecular
spectroscopy for large systems were developed. Through equation of motion
(EOM) constraints, equivalent to the random phase approximation (RPA), elec-

tronic transitions and oscillator strengths are calculated directly. The RPA
formalism is demonstrated to give accurate transition energies for electronic
excitations. The EOM constraints guarantee equivalence between calculated
oscillator strengths in differing formalisms. Although this provision is only true
for a complete basis, considerable improvement in predicting oscillator strengths
that compare with experiment, and differing formalisms, is demonstrated for a
considerably smaller basis.

The computational tools associated with the RPA were applied to the cal-

culation of the electronic spectra of free base and magnesium porphin. The

ability of the RPA formalism to predict accurate oscillator strengths and selectively

correlate transition energies has allowed direct corroboration of the experimental

spectra associated with these molecules. The RPA formalism is also used

for the calculation of isotropic nuclear magnetic resonance (NMR) chemical

shift. The method is capable of demonstrating the inductive effects associated

with electron-withdrawing substituents and is capable of differentiating types

of bonding environments. The prediction of paramagnetic substituent effects,

however, was found to be deficient due to the small size of the basis utilized.

Minimal basis set (MBS) formalisms are restricted in their ability to predict

accurate results for properties. In order to maximize the capability of model

Hamiltonians and maintain the efficiency of such MBS formalisms, a complete

charge-dependent intermediate neglect of differential overlap (INDO) Hamiltonian

appropriate for the calculation of molecular spectra was developed and examined.

The new model is capable of accurately describing the local chemical potential

of atoms in molecules, and surpasses the current INDO/S model in the prediction

of carbonyl excitation energies. A procedure for the acquisition of atomic and

molecular parameters was developed that makes use of readily available exper-

imental data. Simplified parameter functions were developed which allow more

information to be incorporated into the model with less effort than previously.



Chemistry is the study of interactions at the atomic and molecular level.

These interactions can be between the atoms that constitute a molecule, between

molecules or between a molecule and its environment. Quantum chemistry dedi-

cates itself to the task of predicting observables associated with these interactions.

One specific area of study is the interaction of matter with long wavelength radi-

ation. Molecules react to short-lived pulses of radiation in many ways depending

on the frequency of the radiation. High-energy radiation may result in dissocia-

tion or ionization. Low-energy radiation may result in orientational changes of a

molecule. In all cases, however, specific molecules respond to radiation of dif-

fering wavelength very uniquely and can be characterized by their specific ability

to interact with radiation of a given frequency.

The range of wavelength associated with visible and ultraviolet light usually

corresponds to electronic transitions in molecules. Light of a set wavelength

may be absorbed by the material resulting in a new electronic state. This

absorption will take place with a given probability which is again a function of

the molecule of interest. In addition, the probability for which a molecule absorbs

radiation of fixed wavelength but opposite polarization may differ. The qualitative

and quantitative prediction of these types of interactions has been prominent in

quantum chemistry from its inception. While quantum mechanics can in principle

yield exact results for all molecular interactions, the difficulty in obtaining these

solutions has led workers to seek approximate solutions. Early Hiickel theory

models are successful in predicting the qualitative nature of molecular absorptions

of visible and ultraviolet light for planar conjugated systems [1]. More exact

models such as that of Pariser, Parr and Pople (PPP) for conjugated systems were

developed and the agreement between theory and experiment was extended to the

accurate prediction of transition energies [2-4]. Hiickel theory was also extended

to nonplanar systems [5].

Modern approximate methods usually involve the use of well-characterized

experimental data for small systems to fit quantum mechanical models that can

be extrapolated to much larger systems. The spectroscopic parameterization of

the Intermediate Neglect of Differential Overlap (INDO/S) model Hamiltonian has

been very successful in the prediction of properties associated with electronic tran-

sitions [6, 7]. The use of the INDO/S Hamiltonian for the prediction of electronic

transition moments and properties has been primarily a state-oriented approach. In

other words, electronic states are calculated directly yielding transition energies as

differences in state energies and transition probabilities as transition amplitudes

between states.

In this work a transition-oriented approach is developed with the INDO/S

model Hamiltonian and compared to the more traditional state-oriented approach.

Through the use of polarization propagator (PP) techniques [8], transition energies

and probability amplitudes can be calculated directly without the need to describe

the individual states. If exact state functions can be found then both the state-

and transition-oriented approaches are equivalent. For almost all cases, however,

approximate state solutions are employed and methods that are capable of giving

the most consistent solution to a given property must be exploited. A consistent

approach to the calculation of molecular response to electromagnetic radiation

will involve an equitable treatment of transition moments as calculated in both

the dipole length and dipole velocity formalisms. This consistency requirement

will be sufficient to derive the equations associated with the Random Phase

Approximation (RPA) [9, 10]. In this way a balanced treatment of transition

amplitudes can be achieved.

Many properties are computationally related to the transition energies and

moments associated with electronic excitations. Such properties are usually

associated with the response of the total energy of a molecule to an external

perturbation such as a constant magnetic field [11, 12]. These properties are often

expressed as a perturbative expansion of the total energy and truncated at some

order. For a small perturbation the expansion can be truncated at second order

which usually involves a formal sum over all transition moments and associated


energies. Computational tools used must not become prohibitively expensive

as the size of the system becomes larger regardless of the Hamiltonian used.

For properties that involve a formal sum over all the states of a system, the

techniques necessary to calculate individual state transition moments become too

difficult. The method of reduced linear equations [13] as recently employed in

the calculation of molecular polarizabilities [14] meets the criteria of expediency

for sum over state properties. This technique is used in the calculation of Nuclear

Magnetic Resonance (NMR) shielding tensors at the RPA level of theory for large

compounds which would not be possible otherwise [15].

Accurate ab initio calculations with the RPA for small systems has yielded

excellent results but ab initio techniques are often too computationally intensive

for large molecules. If quantum chemistry is to be predictive for the majority

of chemical systems model Hamiltonians will continue to be a necessary tool.

It is therefore prudent to develop and characterize the basic structure of model

Hamiltonians and the content of the parameters in such models. To this end

the basic form of the atomic parameters in the INDO model are examined.

In conjunction with recent developments in density functional theory [16] this

examination suggests generalization and improvement of the current Hamiltonian

used for electronic spectra.

Time Dependent Perturbation Theory

If a molecule is subjected to a short duration of long wavelength radiation

the response of the molecule to such radiation can be treated as a time dependent

perturbation. Stationary state solutions (Ji) to the time dependent Schr6dinger


Hgo (1-1)
1 at
for time independent Hamiltonians (Ho) can be expressed as,

I o(q, t) = exp(-iE t/h) i(q) (1-2)

where Ef is the total energy of stationary state i. The Oi are orthonormal functions

of space and spin coordinates (q) only and are eigenfunctions of Ho.

For a time dependent potential energy perturbation (H') Equation (1-1)


(H + H') = --- (1-3)
i at
For a fixed point in time, solutions to Equation (1-3) can be expressed as a linear

combination of %0 since they are eigenfunctions of the Hermitian operator HI

and thus form a complete set. For different points in time the coefficients of the

expansion will in general be different; therefore, solutions for Equation (1-3) of

the form

T = aj(t)1(q, t)= a(t) exp(-iEt/h)jJ(q) (1-4)
j j

must be sought. Substitution of Equation (1-4) into Equation (1-3) yields,
[E aj(t) + H'aj(t)] exp(-iE t/h)Vkj(q) =
h daj(t) (1-5)
Eaj(t) d ] exp(-iEjt/h)n (q) (
Utilizing the orthonormality of Oi one obtains,

aj(t) exp(-iEt/h)(kIi'J h dak(t) exp(-iE t/h) (1-6)
so that Equation (1-3) becomes a series of differential equations in time for the

coefficients in Equation (1-4) without loss of generality [17].
dak(t) i
dt 1 a-j(t)exp[i(E Ef)t/lh](kH' j) (1-7)

For a system in state m in the absence of a time dependent perturbation we

can choose initial conditions such that am ; 1 and anm w 0. For a perturbation

applied at time 0 of sufficiently short duration that the coefficients are close to

their initial values we can approximate Equation (1-7) as

dak(t) i
dat) exp[i(E' Em)t/h](klIH' Im) (1-8)
dt h
and integrate from t=0 to t' to obtain an approximate expression for the coef-

t' to
Jdak(t) /exp[i(E E'm)t/h](k IH'Ilm)dt
0 0 (1-9)

ak(t) ak(O) exp[i(E Em)t/h](k|lH'lm)dt

If the perturbation ceases after t' then the final state function can be expressed as

Equation (1-4) with constant ak(t). This result is a superposition of eigenstates

of the unperturbed Hamiltonian. The probability of observing the system in state

n with energy E after the perturbation is then [18],

Pn = an (t') exp(-iEot/h)|2 = an (t')12 (1-10)

The effect of the perturbation is a transition from the initial state m to final state

n with probability Pn.

Molecular Response to Radiation

Electromagnetic radiation consists of both electric (E) and magnetic (B) fields

that can interact with the electrons within a molecule. The classical potential

energy of interaction between these fields and a system of n identically charged

particles (q) can be written as [19],

Ep= -q Et.r-+B_ -E.d- ./ (1-11)
i \ 2mic

where the position, momentum and mass are represented by r, p and m, respec-

tively, and c is the speed of light. The second equality expresses the interaction

energy in terms of the system dipole and magnetic moments, respectively. If a

classical description of the electromagnetic fields is maintained and we replace

the dipole and magnetic moments with their quantum mechanical operators the

interaction energy becomes for electrons of charge --e

-ielir- x Vi r eL
2mc i=1 2mrnc (1-12)

In the above L is the angular momentum operator.

The space and time variations of the electromagnetic fields of transverse

plane waves in a vacuum can be expressed in a right-handed coordinate system

(e'1, e2, e3) with propagation along ^3 as

E(r, t) = (Ele1 + E2e2) exp(ike3 r iwt)
B(r, t) = e^3 x E
The magnitude of the wave vector is denoted k and has dimensions of reciprocal

wavelength. The variation in coordinate for an electron in a molecule is of the

order of Angstoms while the wavelength of radiation in the ultraviolet region is

of the order of thousands of Angstroms; therefore, the spatial variation of the

electromagnetic field can be neglected to obtain [20]

E(, t)= (Elei + E2e2)exp(-iwt) = Er exp(-iwt)
B(, t) = (Ele2 E2e6i)exp(-iwt) =Br exp(-iwt) .

The terms El and E2 are, in general, complex.

Two special cases are of interest, that of linear and circularly polarized

radiation. Let El be real and E2 purely imaginary and of equal magnitude as


E1. Through the use of the identity exp(iO) = cos 0 + i sin 0 the real component

of the fields can be expressed

EL (, t) = El(e1 cos wt e2 sinwt)
BL (, t) = E(e&2 cos wt F ei1 sinwt)

The upper sign designates left circularly polarized light as the rotation of the

polarization vector is counterclockwise as observed towards the origin from the

positive axis of propagation. The lower sign designates right circular polariza-

tion. The electric and magnetic vectors maintain mutual orthogonality. Linearly

polarized radiation can be expressed as an average of pure left and right circular


Combining Equations (1-11) and (1-15) the interaction Hamiltonian of in-

terest is obtained.

H'L = -E [(ei d + e2* fl)coswt (e2 d -d ) sin wt]
H'L = -EI(HI1 cos wt H2 sin t)


where HI and H2 are independent of time. Defining wkm (E' E )/h Equa-

tion (1-9) becomes
akL (t) =k() + i exp(iwkmit) cos t( klilI m)dt
iE- J exp(iwkmt) sin wt(bklI2 |m)dt
0 t(1-17)

akL (t') ak(0) + (bk|Hli ti m) i- exp(iwkmt)cos wtdt
(bk IH21im)-- / exp(iikmt) sin wtdt
The sine and cosine terms can be written as exponentials
cos wt = -(exp(iwt) + exp(-iwt))
sin wt = -(exp(-iwt) exp(iwt))
and the explicit integration over time performed [19].
ak(t) = ak(0)
SliEl exp[it'(wkm + _)] 1 exp([it'(wkm -)] -
2+( kkm + O km + (1-19)
T fi2m) _El [exp[it'(wkm w)] 1 exp[it'(wkm + w)]- 1
T (0k|H2| m)^ ------------- ^ --
2h wkm O Wkm + W

Equation (1-19) was derived on the premise that t' is small; therefore, the

only cases for which ak becomes of appreciable magnitude is if w approaches

Wkm. Since

exp(itIO) 1
iexp(it = it' (1-20)
a0-o 0

the expression does not become infinite [18]. If E )E0 then the case of w=wk

corresponds to an absorption of radiation to go from state m to state k. For

the opposite case we obtain the probability of stimulated emission from a higher

state to a lower state. The probability for absorption for k 5 m is then given by

Equation (1-10) recalling the initial conditions on ak.

Pk, = akL (t')ak (t') = k (H1,H I2 x Pk(E, t') =
{l(klHil lm)12 + 1(OkljA2[m)j2
Ti[(k|Hl|Cm)( k|H2|m) (kklHllm)*(k IH2ilm)]}
E 2 exp[it'(wkm w)] exp[-it'(wkm w)]
412 km )2
The time dependent portion of the probability can be simplified to yield

( ) E sin [.5(wkm -- (1)t']
Pk L = Pki Hk) ,I x -[Si2[. m -- W-)2 (1-22)
R R Ih (km W)2

For a source of radiation that has equal portions of left and right circularly

polarized light we can define two different types of absorption processes. The

average of left and right polarizations (Pklin = Pk = .5(PkL + PkR)) yields linear

polarization and the absorption process is the total response of the material to the

perturbation. The remaining component corresponds to the difference in absorp-

tion between left and right polarizations (Pkdif PkCD = .5(PkL PkR)). This

is known as circular dichroism [21]. The spatial components of the probability

of absorptions becomes,

Pk(Hil,H2) = I(kHllim) 12

= 1(kklei dl m)|2 + I(^kl2 *l m)12 (1-23)

= dle dkm2 + |e2 Fkml2

Pkco (Hi, IH2)

= [(ObklI3 2 1m)(k)kI2jlmbmm* (1kHllm) *(iPkI72zlm)] (1-24)
= -i[(1i dkm)(e'l km) + (e2 dkm)(e2 k Am)]

Assuming that the energy for the absorption process is supplied by the external

radiation of frequency v then
Wkm W = (E' Em hv)/h = 27r(E~ Eo hv)/h
= 27(vkm v)
Normal radiation is not composed of a single frequency but is better described by

a frequency distribution. The magnitude of the field E1 is related to the radiation

per unit volume (U1) via UI = (E1)2/8r which is also in general a function of the

frequency [22]. Let ul(v)dv be the electromagnetic radiation per unit volume for

frequency v in the differential range v to v + dv. In terms of probability per unit

frequency the field dependent portion of the transition probability becomes [18]

8x sin2[ (km u)t'] d -
Pk(El,t',dv) = -ul() (m dv (1-26)
h- (Vkm V)2

Only for v close to vUk does the frequency dependence of the probability have

appreciable magnitude. The frequency of visible light is on the order of 1015

cycles per second or Hertz (Hz). We can therefore approximate ul(v) to be

constant in the region of Vkm and integrate over all frequencies to obtain, via

f dx sin2(ax)/x2 = ar,

Pk(E, t') = Pk(Vkm, t') = 27r1(km)t (1-27)

The probability per unit time is obtained by dividing by t'.

Thus far we have only been dealing with radiation that is propagating along

a specific axis. For isotropic radiation (ul = u2 = u3 = u/3) we have to sum

the probabilities over all orientations to obtain the total probability of absorption.

The final expressions for the probability of absorption of a single molecule to go

from state m to state k per unit time is

Pk = [ldkmi2 + kml2] -(Vkm)
2w (1-28)
PkcD = -2i[dkm /km] 2U(V(km)

The Reporting of Theoretically Calculated Spectral Probabilities

The total probability per unit time for the absorption of a system of N

molecules is simply N times the single molecule probability. There are two

parts of the transition probability that are characteristic of the molecule and not a


function of time or the radiation field in this treatment. The first is the frequency

or energy of the transition. The second is the proportionality constant that is a

function of the transition matrix elements of the dipole and magnetic moment

operators. The probability for absorption can be formulated as a proportionality

constant times the density of the radiation times the number of molecules. The

proportionality constant then becomes [18]
Bmk [dkml + -kmI 2 i2-
3 (1-29)
BmkcD = -2i[dkm (km]2

The problem that must be addressed theoretically, therefore, is the attainment of

the transition energies and the transition moments themselves.

The transition moments have been expressed as matrix elements between wave

functions of the total time independent Hamiltonian operator.
Si + E q (1-30)
2m; ij
i i,j(i
These state functions therefore include both the nuclear as well as the electronic

components. In this treatment the electronic portion of these transition moments

is of interest. The mass of the proton is approximately 1836 times the mass of

the electron; therefore it is reasonable to assume that the electronic motion can be

treated separately from the nuclear motion. In other words, since the electronic

motion should be very fast with respect to the nuclear motion the nuclei can be

treated as fixed with respect to the electronic motion and treated as a potential

field. This is commonly known as the Born-Oppenheimer approximation [17, 23].

For n electrons and M nuclei the Hamiltonian for fixed nuclei becomes
S2 -"2 ( 1 M M 2
Hei= m + 2 ZE Z Z (1-31)
i mi j(i rij i ,l a r#

Solutions to Equation (1-31) for arbitrary nuclear positions then defines a potential

energy function (V(ri)) from which the nuclear motion Hamiltonian (HT) can be

M 2
Hoi __ V2
2M (1-32)

l= 2m- V + V(ri)

Under the Born-Oppenheimer approximation the Hamiltonian operator is par-

titioned between the electronic and nuclear motion. Therefore the state functions

for the total Hamiltonian can be expressed as a product of nuclear and electronic

functions in which the electronic wave function depends parametrically on the

coordinates of the nuclei.

Ok = 0kel (ri, {r.}) ONk (r.) (1-33)

Transition elements for perturbation operators such as dipole then become

(OklbIdm) = (~ elk/Nk mIdl mEel Nm,)

(= ('Nk, (elk Ide|,)em Nm,) + (melkI(Nk, kdN Nmi)'eim) (1-34)

= (bNk, (Welk delCJel mbN,,), k $ m

The nuclear contribution to the dipole operator vanishes for transitions between

electronic states and will not be considered further.

If for given initial nuclear and electronic state m and specific electronic state

k, (ielk Id.elm,) is treated as an average electronic moment for all excited nuclear

states, k', the electronic integral can be removed and the squared matrix element


Idkm2 = 2 I elk dl elm)l2 X (Nk INm,)( I'N,) (1-35)

For a given electronic transition from state m to k many possible nuclear transitions

are possible. The squared overlap between nuclear state functions is known as the

Franck-Condon factor for the electronic transition. The nuclear states correspond

to various vibrational and rotational states for the molecule. The Frank-Condon

factors then give the relative intensities of electronic transitions as divided among

the various final nuclear states [18]. Considering that the ONm, are eigenfunctions

of a Hermitian operator corresponding to electronic state m, they form a complete

set of orthonormal functions and therefore the total electronic transition probability

over all final nuclear states is

ldkml2 = I()elk dl el m)2 X ZiN NIN )(NmIlNk,)
m' m' (1-36)
= I(elk el X N,,)N,, = I(|< elkl d|I l)I2
Therefore the total probability of an electronic transition under the Born-

Oppenheimer approximation can be expressed at the average electronic state


probability. This average is usually further approximated as the value obtained

between electronic states at a fixed geometry.

With the approximations to the transition moments presented above, the

theoretical problem is reduced to the determination of the state functions and

energies of the electronic Hamiltonian for a molecule. The frequency of a

transition is obtained from the difference in state energies and the probabilities

of interest are found as straightforward integration of the perturbation operators

over the electronic state functions. The total probability for electronic transitions

is proportional to the sum of both dipole and magnetic moments. For dipole

allowed transitions the magnetic moments are usually 10-4 those of the dipole

transition moments and are usually neglected.

Consider a three dimensional harmonic oscillator of the same mass and charge

of an electron initially in its ground vibrational state with vibrational frequency the

same as that of the electronic transition in question. The proportionality constant

for a change of one vibrational quantum number from the ground state is given

by [18]

7re2 we2
B01 = (1-37)
Vkmmh Ekmm

The convention for reporting total electronic absorption probabilities is in the form

of oscillator strengths which is defined as the ratio of the electronic proportionality

constant to that of the pure harmonic oscillator described above. The result is a

dimensionless constant usually designated f,.

S47rVkmm 2 2Ekmm 2
f = IK2 Ikm = e2Idkm
3he2 3h2e2
2Eknm 12 (1-38)

The "r" indicates the use of the dipole length operator to calculate the transition

matrix elements [24].

For circular dichroism the proportionality constant can be rewritten

27r 44r
BmkcD = -2i[dkm km]3h-- = Rmk 2

mk = -ikdbm}) (bklIkm) (1-39)
n 1
2mc1l k m i x i|,k)
i i
Rrk is known as the rotary strength of the electronic transition from state m to k.

Rotary strength has both sign and magnitude since it represents the difference

between the probabilities of absorption for left and right polarized light [21,

25]. For equal absorption of left and right circularly polarized radiation the

rotary strength would be zero even though the total oscillator strength could be

substantial. As with oscillator strength rotary strength is proportional to the total

electronic probability averaged over all final nuclear states.

Oscillator strength is a dimensionless quantity but some system of units must

be used in order to account for the magnitude of rotary strength [26]. The most


commonly used system is cgs units. Rotary strength in cgs units becomes

n n
Rnk(cgs) = 4.45458 x 10-30o(mI I i I'k)(bm i X x i i|k) (1-40)
i i
It is customary in quantum mechanical calculations to use atomic units to calculate

transition energies and moments. In atomic units h = e = m = 1. The atomic unit

of energy is the Hartree and the atomic unit of length is the Bohr. Using atomic

units to calculate transition moments and energies, oscillator strength and rotary

strength becomes

S 2Ekm(au) nm| 2
fr = 3- I(lm|I Sri|k)l2

n n (1-41)
Rmk(cgs) = 235.7262 x 10-40(m I r|k)( m i X Vi 'k)
i i
Common units for reporting transition energies are electron volts and reciprocal

centimeters (1 au = 27.212 eV = 219474.6 cm--1). The conversion constant for

rotary strength from cgs to SI inits is 1.19x106.


The Hartree-Fock Equations

In order to determine the electronic transition moments and energies for a

molecule it is necessary to determine in some fashion the eigenfunctions and

energies of the electronic Hamiltonian. In atomic units the electronic Hamiltonian

n 2 n 1 M
o -+ -- E + Vnuc. (2-1)
ii j(i a
The first term is the electron kinetic energy operator. The second sum contains

the potential energy of repulsion of the electrons and the potential energy of

attraction between the electrons and nuclei. The final term is the nuclear potential

energy and is a constant for fixed geometry regardless of electronic state. Since

we are interested in differences in electronic states this term is dropped. We seek

solutions to Equation (2-1) of the form

H j) = Ejljj) (2-2)

The eigenstates Oj are functions of the electronic coordinates and spin. The

nonrelativistic Hamiltonian used in this treatment is not a function of the spin

coordinates of the electrons. Since electrons are spin 1/2 and therefore fermions

the eigenstates for a system of electrons must be antisymmetric with respect to the

interchange of the coordinates of any two particles [17]. It is necessary therefore

to seek solutions that have the property

j(ql, q2," "- ) = -j(q2, q1,- ) (2-3)

where ql is the space and spin coordinate for electron one.

According to the variation principle, for an arbitrary well behaved normalized

wave function, 14), the expectation value for the Hamiltonian of a system is an

upper bound to the exact ground state energy. In other words

Eo < (I|HIO|4) (2-4)

If 1() contains adjustable parameters they can be varied so that E =(|IHoj) is

a minimum and thus 14) becomes the best possible approximation to the ground

state wave function and E the best possible approximation to the exact ground

state energy. The calculated energy (E) only becomes the exact if 1j)) is the exact

wave function [19].

In general it is not possible to solve Equation (2-2) for more than two

electrons. It becomes necessary to seek expedient approximations that allow

generalized methods of solution. The most common such approximation is

the molecular orbital or Hartree-Fock (HF) approximation [23]. In the HF


approximation Oj is expressed as the best possible determinant of one-electron

spin orbitals. Since a determinant changes sign when two rows or columns are

transposed, which in this case represents a single electron function, the fermion

property of antisymmetry is attained whereas a simple product of spin orbitals

would not have this property. For n electrons

j) 4 ll(qi)02( (q2) m(qm)... n(qn))
n! (2-5)
= (-1)"Pij 1(ql1)2(q) m(qm) n(qn) (2-5)

and is known as a Slater determinant. The symbol Pi is the permutation operator

for the electron labels and the sum is over all n! permutations. The number of

single label interchanges necessary to restore the i'th permutation to the initial

order is denoted pi. The 1/vn! factor is the normalization constant so that

( j lIj) = 1.

Utilizing the variational principle with a Slater determinant of orthonormal

spin orbitals yields

Eo < E = (o0H110|o) (2-6)

The calculated energy becomes a function of the molecular orbitals that make

up the determinant. The remaining discussion is restricted to closed shell deter-

minants. Insertion of Equation (2-1) into (2-6) and integrating over space and

spin yields

E = O q*() Z-i()dq
i=l 27 L i ql)d-q
=1 =1
+ f / O(ql)Oj(q2)- [i(qi)Oj(q2) Oj(qi)Oi(q2)]dqidq2
n n (2-7)
E Oilf--o) + 2 [(E jlij) (Oij0j0\
i=l i,j=l
n n
(iIlhi) + (iij
i=1 i,j=

The operator h is a one-electron operator consisting of kinetic energy and nuclear

electronic attraction. The i'th term is usually denoted hii. The two terms relating to

the electron repulsion are known as the coulomb (Jij) and exchange (Kij) integrals

respectively. The exchange integral is the direct result of using a determinantal

form for the wave function and would be absent if the wave function were a

simple product. In simplified notation Equation (2-7) becomes
n n
E = E(ili) + 1 (ijllij)
i=l i.j=1
n n (2-8)
= hii + 2 [Jij Kij]
i=i i.j=1

In order to obtain the best possible solution to Equation (2-8) the spin orbitals

are varied under the constraint of maintaining orthonormality. Let eij be a set of

Lagrange multipliers so that the Lagrange variational expression becomes

L = E eji[(iIj) 6ij]
ij (2-9)
6L = 6{E eji[(ilj) bij]} = 0
Differentiating the orbitals and rearranging terms yields

J dqis l(1){h(1) i(1) + [dq2O(2) ji(2) ]i(1)
i=1 j=1

-dq2o (2) -Pi2j(2) i(l) + complex conjugate (2-10)
j=lV Jj J J
n n
= i dql6 (1) e jiij(1) + complex conjugate
i=1 j=-
The terms in brackets are known as the coulomb and exchange operators respec-

tively. The variation of the spin orbital is arbitrary so Equation (2-10) must hold

term by term or

1h(1) + E [ j] }ki(1) = f(1)i 1) = Esjij() (2-11)
j=1 j=1
where the first equality defines the one-electron Fock operator. Multiplying by

Ok and integrating one obtains (k [fli) = Iki. The Lagrange multipliers are
matrix elements of the Fock operator in the basis of the spin orbitals. Since

the Fock operator is invariant to a unitary transformation and the matrix of

Lagrange multipliers is Hermitian it is possible to choose the Oi so that the matrix

of Lagrange multipliers is diagonal. That set of Oi and ei are known as the


canonical Hartree-Fock equations. The solution of the HF equations reduces to

the determination of the eigenvalues and eigenfunctions of the Fock operator.

f i) = ii) (2-12)

with the Fock operator defined in Equation (2-11). Since the Fock operator itself

is a function of the spin orbitals the methods for the solution of Equation (2-12)

are in general iterative in nature.

The energy minimization that led to Equation (2-12) pertained only to the

n "occupied" molecular orbitals directly associated with the electrons. Once

these orbitals are determined, the Fock operator becomes a well-defined Hermitian

operator and thus possesses a complete set of solutions. Those orbitals beyond

the active occupied orbitals are known as virtual orbitals. The eigenvalue for an

arbitrary canonical orbital can be expressed

Ei = (iIf Ii) = hii + [Jij Kij] (2-13)

Equation (2-8) represents the total electronic energy expression for n elec-

trons. As a reasonable first guess, the orbital solutions for n electrons can be used

to estimate the energy for a system of (n+1) or (n-1) electrons. The difference

in energy between an n electron state and an (n-1) electron state is known as an

ionization potential (IP). The difference with an n and (n+1) state is known as an

electron affinity (EA). Using the orbitals obtained for an n electron system and

Equation (2-8), one obtains a molecular ionization potential for arbitrary occupied

orbital a and an electron affinity for arbitrary virtual orbital b.

E,(n 1) E =IP = -e
E Eb(n + 1) =EAb = -Eb
Equation (2-14) is an expression of Koopmans' theorem which interprets the

eigenvalues of the occupied orbitals as the negative of the IP's for those orbitals

and the eigenvalues of the virtual orbitals as the negative of the EA's for those

orbitals [23].

The use of the molecular orbital approximation yields an approximate ground

state for a molecular system. Since the solutions for the Fock operator form a

complete set the exact wave function for a system can be expressed in terms of an

expansion of the these solutions. From the HF ground state the occupied orbitals

can be successively replaced by the virtual orbitals and classed in terms of single,

double and higher replacements [27].

o) = [Cio0) + CI0) + c a C ab) +" -] (2-15)
aa aafab

The determinant generated by the replacement of occupied orbital a with virtual

orbital a is denoted I|a) with corresponding expansion coefficient Ca. The coef-

ficients of the determinants become the parameters for the use of the variational

principle. Matrix elements of the Hamiltonian in the basis of determinants must


be formed to solve the problem. The matrix elements between the HF ground

state and a single replacement can be expressed
(oI| -) =haa + (ailai) =
i=1 (2-16)
(0aIf #a) = Ea(0,a1a) = 0
The determinants formed from a single replacements do not directly interact with

the HF solution for the ground state. This is an expression of Brillouin's theorem.

The lowest order correction to the HF ground state is through double replacements

only. The eigenvalue problem in the complete basis of determinants yields not

only the exact ground state energy but also the wave functions and eigenvalues

of the excited states of the system [28].

The Roothaan Equations

The Hartree-Fock equations in the previous section do not suggest a method

by which to obtain the molecular orbitals that define the Fock operator. The

most common way is to assume a basis set that is pertinent to the system of

interest. For molecules this usually involves an atomic-like basis. Each atom in

a molecule contributes a set of functions to the total basis for the solution of the

HF equations. The number and type of functions that an atom can contribute is

limited only to the computational resources of the researcher. The smallest basis

that an atom can contribute is called the minimal basis set (MBS) and is usually

that basis associated with the shell structure of an atom. Each principle quantum

number on an atom contributes a fixed number of basis functions of specific

angular momentum. The molecular orbitals associated with the Fock equations

become linear combinations of atomic orbitals (LCAO-MO).

i = I Cix (2-17)
The expansion coefficients become the variational parameters. The xy are the

atomic basis functions and in general are a function of the electronic coordinates

only. The total number of basis functions is N. The spin of a molecular orbital

is denoted by the existence a bar over the symbol for orbitals that are spin -1/2.

No bar is present for spin 1/2. Alternatively, spin -1/2 is often denoted as /

and spin 1/2 as a. Orbitals of differing spin are orthogonal.

Placing Equation (2-17) into Equation (2-12) and multiplying by the atomic

basis functions and integrating yields

(Xl1 // CEiIXu)= eiCCiIX},
v=l v=l
SCi(xlf xv) = i E Ci(X Xv) (2-18)
v=l v=l
SF rVi = C i Y St Cvri
V=1 V=1
where the elements of the Fock (F) and overlap (S) matrix are defined. In terms

of matrices

FC = SCe



is known as the Roothaan equation [29] and is the most used formalism for the

solution of the HF system of equations. The columns of the coefficient matrix are

the expansion coefficients for the molecular orbitals and e is a diagonal matrix of

orbital eigenvalues. For N basis functions one obtains N molecular spin orbitals,

n of which are occupied and (N-n) are virtual.

A determinant in which each spatial orbital is associated with both spin

a and p functions and both are occupied is known as restricted Hartree-Fock

(RHF) and is the most common type of calculation performed. Introduction of

Equation (2-17) into the expression for the Fock operator and calculating the

matrix elements of the Fock operator in Equation (2-18) for a RHF determinant

n/2 N
F, = h, + [ CC L [2(r a) (Irjov)] (2-20)

The iterative nature of the Fock equations is demonstrated in that the coefficients of

the occupied orbitals must be known in order to define the Fock operator but these

are the same coefficients that are sought in the calculation. The intermediate sum

over the occupied orbitals and coefficients in (2-20) defines the so called density

matrix P = CnCt. Here n represents the diagonal matrix of orbital occupancies

the trace of which is the number of electrons. For an RHF determinant the spatial

orbital occupancies are either 2 or 0.
Par = 2CiCri (2-21)

The elements of the Fock matrix become
FV = h, + Pr[2(rTjv|) (srl 7r)] (2-22)

The matrix elements of the one-electron operator in the atomic basis do not depend

on the density matrix.

For real spatial atomic orbitals the overlap matrix is symmetric and can be

diagonalized and real roots taken of the matrix. Let X be the transformation

matrix that diagonalizes S so that XtSX = s.

S = XsXt = Xsl/2s'/2Xt = [XS1/2Xt][Xsl/2Xt] (2-23)

If there are not any linear dependencies in S then

1 = [Xs-1/2Xt][Xsl/2Xt] = U[Xs/2Xt] (2-24)

defining the matrix U = S-1/2. Consider the transformation of Equation (2-19).

Ut // FIC = FU[Xsl/2Xt]C = [Xsl/2Xt][Xs1/2Xf]Ce = SCe

[UtFU]C' = {Ut[Xsl/2Xt]}C'e (2-25)

F'C' = C'
This procedure by which the atomic basis is orthogonalized is known as symmetric

or L6wdin orthogonalization [30] and the Fock procedure is reduced to the matrix

eigenvalue problem C'tF'C' = C'tC'e = e.

Once the atomic basis is chosen for a problem the solution usually follows

well established lines. The general algorithm is to first evaluate all the integrals in


the atomic basis. The overlap matrix is diagonalized and the transformation matrix

U is formed. An initial guess for the density matrix is formed and the Fock matrix

is built, transformed and diagonalized. The coefficient matrix in the atomic basis

is formed by the back transformation C = UC'. A new density matrix can then

be formed yielding a new Fock matrix and the pattern repeated. The procedure is

considered converged when differences in successive density matrices are within

some set threshold. This technique is known as the self-consistent-field (SCF)


The Calculation of Electronic Spectra

The Hartree-Fock procedure for a given basis provides a convenient set

of orbitals from which Equation (2-15) may be variationally solved for all

the electronic states of a system. The full expansion provided by all orbital

replacements is known as a full configuration interaction (CI). For the special

case of a complete basis set, the full CI is said to be "complete" and represents an

exact solution of the problem represented by the assumed Hamiltonian. For other

than very small systems, however, full CI is not feasible and some truncation is

necessary at some order in the expansion. From Brillouin's theorem the lowest

order correction to the ground state comes at the inclusion of double replacements

only. Since the single replacements do not mix directly with the ground state they

become the basis for the lowest order approximation to singly excited states. It is


the transition moments and energies between such excited states and the ground

state that are of interest here.

The expressions for oscillator and rotary strength in Chapter 1 involve matrix

elements between states for the electron position and angular momentum oper-

ators. These expressions are not unique however. The commutation relation

between the Hamiltonian for the system and the position operator [31] is

[ io = f i (2-26)
i=1 i=1
where Vi is the velocity operator for electron i. The evaluation of matrix elements

between the ground and excited states of the Hamiltonian operator yields
n n n
(001 NolO) = (iol1 rlH0 -H ii|k)
i=l i=l i=l
= (Ek E)(o | i k) (2-27)
n n
(o I Vilk) = EkO(O 1 E i 4l k)
i=1 i=1
Equation (2-27) indicates that the matrix elements for the position operator can

be replaced by those of the velocity operator divided by the transition energy [32,

33]. Equations (1-38) and (1-39) for oscillator and rotary strength in terms of

the velocity operator become

2Ek0 n -,

i=l (2-28)
f_ 2 n-,|

n n
R = CR( O ri kii 0k)(oI >i x Vi|b) =
i=1 i=1 (2-29)
EkO (0 dik)(0 i XVi k)
i=- i=1
For exact eigenfunctions of the Hamiltonian, transition moments in both length

and velocity formalisms are equivalent. Since approximate solutions to ground

and excited states do not necessarily obey Equation (2-27), the equivalence of

oscillator and rotary strength in Equations (2-28) and (2-29) can not be assumed.

Equation (2-27) constitutes a formal constraint that must be incorporated into

the calculation of approximate excited and ground states if equivalence is to be


For practical purposes the HF approximation and the HF plus doubles approx-

imation will be considered for the ground state and singles only for the excited

states [10]. Truncations at higher order are not feasible for use with large sys-

IYo = YKo[I'HF)+{y I C b a# }a#)
aofab (2-30)

Transition moments for a general one-electron operator 0 is expressed in terms

of the approximate states in Equation (2-30) as

1 k Cab*C( o I (2-31)

Matrix elements for one-electron operators are identically zero if the two deter-

minants differ by more that one orbital replacement.

( 1-,o,^|.) = (fl^ c) [6,. + ,, + 6 c + ~ ] (2-32)

The indices of summation in Equation (2-31) are arbitrary so that

(0o I0lk) = ^ KiC ,(rHF IOI)
C; I (2-33)
but (Oa6)aOI) (,l (IO HF) yielding

(0Ck)= KC(HF= Klt )
7c (2-34)
+ { KiC Cc (,16OrHF)}
aa 0 C
For excitation k the following notation is adopted.

xk = K Cc(k), yk, = ZK aC (k)CC(k) (2-35)
This yields,

(0i^O14k) = [xcHFI=OC)
yc (2-36)
+{ YI C(O/6B1HF)}]

The commutation relation between the general operator 6 and the Hamiltonian

operator with the exact state functions [9] is according to Equation (2-27)

(4'oIN ) k)Eko = (Ool[0,HO ]01k)


If the commutation relation itself yields a one-electron operator as in Equation
(2-26) the right hand side of Equation (2-37) becomes for the approximate states
defined by Equation (2-30)

RHS = E [Xc((HF6H|IO ) ( HF 061C))
ye (2-38)

+{Y (( |6f2o HF) (Oc|H6|HF))}
Since the eigenfunctions of the Fock operator form a complete set we can resolve
the identity

1 = IHF)() H + { )( }+ {I )( 1} I )(a}+ (2-39)

between the Hamiltonian and 0. Utilizing Brillouin's theorem and the fact that
matrix elements of the Hamiltonian between two determinants that differ by more
than two orbitals is zero the insertion of Equation (2-39) into (2-38) gives

RHS = X X( (HF I Sd)( i C)
7c bd
S-( ~HI HF HF HF- If6d)(o (IOI))
SHF H(2-40)
+YC((|6 HF)HF|HF) + 4 I11OHF)

Z(4 IHI ) IO HF))}]
bd 6 Q H)
Ad,,c (dIHt (- HFiH|RHF)(-1dc
B6d,,e = (o|H|^HF)

and simplifying according to Equation (2-32) Equation (2-37) takes the form

Eko E d[X HF ) {Y bd 1OHF)}]

7c L 6d

+{Y ((HF 6 ^d )Bbd,7c (4|6Oi|zHF)A de,6d)

The operator 0 can in general be Hermitian as is the position operator or anti-

Hermitian as is the velocity and angular momentum operators. Most calculations

are performed in a basis of real orbitals in which the A and B matrices are

symmetric. For Hermitian operators and real orbitals Equation (2-42) becomes
( lnHVH )[(Xk y,)(A6d, Bsd,7c)
7c 6d
-Eko(Xd + Ykd)76cd] = 0 (2-43)

(Xc -- Yc)(Asd,7C Bsd,7c) = Eko(X d Yd)
and for anti-Hermitian operators

(Xkc + Ykc)(Abd,7c + B6d,yc) = Eko(Xk kd) (2-44)
Equations (2-43) and (2-44) thus define a set of simultaneous linear equations

that must be solved in order to obtain the coefficients necessary for the calculation

of transition moments for both Hermitian and anti-Hermitian operators. In matrix

form they become

(A B)(X- Y) = (X + Y)EK


for Hermitian operators and

(A + B)(X + Y) = (X Y)EK (2-46)

for anti-Hermitian operators. The columns of X and Y are the expansion co-

efficients for the excitations and EK is a diagonal matrix of transition energies.

Equations (2-45) and (2-56) are the RPA Equations [8, 11].

The first case to be considered is for excitations from the Hartree-Fock

approximation to the ground state. This is obtained by setting the matrix Y

identically to zero obtaining

(A B)X = XEK (2-47)


(A + B)X = XEK (2-48)

which are two eigenvalue problems for X and EK. Since A and B do not commute

they do not posses a common set of solutions. Therefore transition moments cal-

culated in different formalisms will not in general be equal. Equation (2-37) can

be solved for Hermitian or anti-Hermitian operators separately, however, by the

diagonalization of (A-B) and (A+B) respectively. If the matrix B is neglected,

Equations (2-47) and (2-48) reduce to the diagonalization of A which is exactly

the matrix obtained from the variational solution of the singles only configura-

tion interaction (CIS) or also known as the Tamm-Dancoff approximation (TDA)

[34]. The TDA is the most commonly used approximation for the calculation of

transition moments and energies but as seen here cannot guarantee equivalence

between transition moments as calculated in different formalisms since Equations

(2-27) and (2-37) do not in general hold. Only in the limit of a full CI can a

variational procedure guarantee this equivalence.

The second case is the direct solution of Equations (2-45) and (2-46) which

were derived on the basis of a ground state consisting Hartree-Fock and double

replacements. Taking the linear sum and difference of Equations (2-45) and

(2-46) yields

BX + AY = -YEK
which is the linear form of the matrix equation

[A B ] =[ X EK (2-50)

Since we have assumed real orbitals the matrix of A and B is real and symmetric

and can be solved as [27],

[Xt Yt] [A B]j ] =[Xt ytj X ]EK = [XX YtY]EK = EK (2-51)

with the special normalization constraint [XtX YtY] = 1. Direct solution of

Equation (2-50) would involve the diagonalization of a matrix of twice the di-

mension of the number of single replacements in order to obtain the coefficients

and transition energies. The problem can be reduced to two successive diago-

nalizations of the same dimension as the number of single replacements. Since

diagonalization scales as the cube of the dimension of the problem, two diagonal-

izations of dimension N are preferable to one diagonalization of dimension 2N.

According to Ullah and Rowe [35] let M = A + B, N = A B, U = X + Y

and V = X Y so that Equations (2-45) and (2-46) become

NV = UEK, MU = VEK (2-52)

The normalization constraint becomes

VtU = UtV = 1 (2-53)

since with real orbitals and real coefficients

vtU = (Xt Yt)(X + Y) = XtX yty + Xty YtX

= xtx yty + [XtY (Xty)t]

= XtX yy = 1 (2-54)

Xty = 4[UU VtV (UtV V'U)]

= [UtU VtV]
If N is positive definite then V = N-'UEK and



Equation (2-55) is similar to the Roothaan equation and an analogous treat-

ment is sought. Let T be a unitary matrix that diagonalizes M so that

TtMT = m, M = Tml/2m1/2Tt. Equation (2-55) can be transformed

ml/2Tt // NTml/2m/2TtU = UEK

[ml/2TtNTml/2][m/2TIU] = [ml/2TtU]E (2-56)

N'U' = U'EK
into a Hermitian eigenvalue problem. To obtain U the reverse transformation

U = Tm-1/2U' is performed. V can be obtained form Equation (2-52).

TTt // MU = MTm-1/2U' = VEK

Tml/2U' = VEK (2-57)

V = Tm/2U'Ei E
In order to obtain X and Y the normalization constraint must be taken into

UtV = [U'tm-1/2Tt][TmI/2U'EK] = 1
UV = U'tU'E1 = 1
The normalization condition is obeyed provided the vectors in U' are scaled so


U'tU = EK (2-59)

With U' properly normalized and U and V obtained as above, X and Y are simply

obtained from

X= (U +V), Y = (U-V) (2-60)
2 2

The use of HF plus doubles as a reference for the ground state and singles

only for the excited states is successful in the maintenance of Equation (2-27)

as a constraint. The resulting solutions are equal to those of the random phase

approximation (RPA) or polarization propagator (PP). The RPA in a complete

HF basis is the lowest level of theory in which a formal equivalence according

to Equations (2-28) and (2-29) can be expected. With the RPA the transition

moments and energies are obtained directly. No formal correlated ground state

or excited states are explicitly formed.

Singlet Excitations with a Finite Basis

The RPA yields formally equivalent expressions for the calculation of optical

properties for a complete HF basis. In practice, however, a finite basis must be

chosen in which Equation (2-39) does not rigorously hold. Calculated optical

properties in differing formalisms become a test for the completeness of the basis

used. Often a full singles space calculation is prohibitive even in a finite basis

and some truncation of the active orbitals is used to calculate transition moments

and energies. If electronic excitation spectra is of primary interest usually only

the lowest excited states are needed. It therefore becomes necessary to study the

stability [36] and equivalence of the calculated transition energies and moments

for the lowest energy transitions as a function of the size and nature of the active


For this study singlet excitations from a closed shell ground state will be

considered. The RPA gives two sets of transition vectors [37]. For Hermitian


(l016KHIk) = (KHF OH| )[Xkc + Yk] (2-61)
and for anti-Hermitian operators

(0olOAlIk) = E(]HFIOA| c)[Xkc yk ] (2-62)
The singlet spin adapted excitation 1Cryc) in terms of spatial orbitals is

1l07c) = 1[114c) 7) I+lcE)] (2-63)

The transition moments in Equations (2-61) and (2-62) become matrix elements

of the one-electron operator in the MO basis.

('lHFlOIlrtc) = [('ObHF6lOl1-rc) + (1CHF1IO61',)]
= +-[(< 610|c) + ( 61c] (2-64)

= 2 (|16|1c)
The operators are usually evaluated in terms of the atomic basis and, after the

SCF coefficients are obtained, transformed into the molecular basis.
(|610|1c) = Ci'r(xi| Ixj)Cjc
Oj (2-65)
Oc(MO) = [CtO(AO)C],

The singlet spin adapted matrix elements of the TDA matrix and B matrix defined

by Equation (2-41) become

Abd,7c = (r6 ec)6766cd + 2(-ybScd) (yc Sd)
B6d,7c = 2(761cd) (7|6dc)
where the two-electron integrals have been suitably transformed from the atomic

basis to the molecular orbital basis and the ei are the SCF orbital eigenvalues.



Use of the LCAO-MO approximation without further simplifying assump-

tions creates a computational barrier. The obstacle is the evaluation of approxi-

mately N4/8 two-electron integrals over the atomic basis set. N is the number of

atomic basis functions. An MBS for benzene would consist of 42 basis functions

and require approximately 400,000 integrals to be evaluated. Accurate calcula-

tions of molecular properties with ab initio techniques, however, often require

extended basis sets. A second set of p functions on benzene increases the num-

ber of integrals to 2,400,000. Commonly used semiempirical Hamiltonians for

electronic structure owe their efficiencies to reducing the number of integrals pro-

cessed to N2. In so doing matrix multiplications and diagonalizations, both N3

steps, dominate, at least at the HF level. In addition most semiempirical methods

utilize a MBS of valence type orbitals only, further reducing N in the more com-

plex N3 steps. These models compensate for the reduction of integrals and small

basis sets by parameterizing directly on atomic phenomena and molecular prop-

erties, and are thus often able to reproduce experiment in a more accurate fashion

than the MBS ab initio calculations on which they are modeled. A valence MBS

for benzene consists of 30 basis functions and only 900 integrals.

The zero differential overlap (ZDO) approximation at the intermediate neglect

of differential overlap (INDO) level is adopted for this work [6, 38-42]. The

ZDO approximation assumes the differential element X,(1l)Xv(1)drl to be zero

unless p=v where it occurs in integral evaluations. A complete description of the

approximations for INDO and other levels of ZDO can be found in Sadlej [43].

The basic equations are summarized below. For INDO the ZDO approximation

is applied to all two-electron integrals in which two or more centers, denoted by

capital letters, are involved and all one-electron integrals associated with the SCF

procedure. Certain two-center integrals will be spherically averaged in order to

maintain rotational invariance. The overlap matrix becomes the unit matrix.

Sp,' = S6p, (3-1)

All one-center two-electron integrals remain and the two-center integrals are

spherically averaged and symbolized by 7AB.

= 7ABSACSBD6gASvo A 0 B (3-2)


Several cases exist for the one-electron matrix elements. Kinetic energy

integrals do not contain the differential element for which the ZDO approximation

is applied so these integrals remain. If both p and v are on center A the matrix

element is partitioned into an atomic part and a two-center part. As with the two-

electron integrals all one-center integrals are kept, but with the use of a Slater basis

the one-electron one-center integrals are diagonal by orbital symmetry. The core

orbitals on one center have a negligible effect on the valence orbitals of a second,

therefore, only the valence electrons are included in a calculation. The interaction

between a core and valence electrons on a single center, however, cannot be

neglected [44]. This energy of interaction scales linearly with the number of

valence electrons and can be accounted for by a core-valence potential function,

EA, that is added to the atomic portion of the one-electron matrix elements.

1 1 ZA
hPAVA = (PA| E + EAIA) = AI Y + EAIVAn)6P
2 C= rlC 2 rlA


The symbol Vc is the potential energy of attraction for an electron to center C.

The atomic one-electron one-center matrix elements are symbolized by UU/p and

are called the core integrals. Since differential overlap has been neglected between

centers, the attraction of an electron on one center to the core of a second should

be similar to a like repulsion. This is accomplished by setting the magnitude of

the two-center nuclear attraction integral to 7AC [38].

hAVA =[Up ZCTAC]&, (3-4)

Two-center one-electron matrix elements are not zero due to the presence

of the kinetic energy operator. The two-center nuclear attraction integrals are of

different type than those that occur on the diagonal of the one-electron Hamiltonian

and can be neglected by the ZDO approximation. This matrix element is called

a resonance integral and is symbolized by /3p,.

hAVB= flAVB, A B (3-5)

Use of the ZDO approximation reduces the Roothaan equation to

FC = Ce (3-6)

which is the exact same form as the symmetrically orthogonalized equation

without the ZDO approximation. Symmetric orthogonalization can be derived

by the orthogonalization of a basis that maintains maximum overlap between

the orthogonalized basis and the original basis [30]. The matrix elements in

INDO can be viewed as if they had been evaluated directly in the orthogonal

basis. In this way the coefficients determined from the INDO Fock matrix can

be back transformed to the original MBS from which it is modeled. These

coefficients, eigenvalues and atomic basis functions can be used to calculate

properties such as transition moments and spectra. With this point of view, the


ZDO approximation then applies only to the SCF procedure and not properties

that are to be calculated that rely on the SCF results. The matrix elements that

remain in the INDO Hamiltonian can be evaluated directly in a basis but this

approach results in a severe approximation to the full MBS calculation which is

usually inadequate for the prediction of properties. Through the judicious use

of experimental information the integrals are parameterized so that calculation

of the same properties yields accurate results. A successful parameterization is

capable of extrapolation to systems and properties for which it was not fit. Integral

approximation schemes appropriate for the calculation of electronic spectra are

of interest here.

Two types of parameters need to be determined. The first are the parameters

that are purely atomic in nature. These are the core integrals and one-center two-

electron integrals. Since all atomic integrals are kept in INDO, atomic parameters

are chosen to reproduce atomic phenomena explicitly. The total energy for an

atom can be written in the same form as Equation (2-8) [45]

ny 1 ny
E=Ec + Uii [Jij- Kij] (3-7)
i=1 i,j=1

Ec is the constant energy of the core electrons and nv is the number of valence

electrons. The exchange integrals, expressed in terms of Condon-Shortly integrals

[46], are determined from differences in atomic state energies. The coulomb

integrals are a linear combination of the spherically symmetric coulomb integral,

denoted Fo or 7AA, and the same Condon-Shortly integrals that make up the

exchange integrals. According to the analysis by Pariser, Fo is set to the ionization

potential of the atom less its electron affinity [47].

/7AA = F = (ssss) = IP EA (3-8)

The core terms remain and are determined as that value for which Koopmans'

theorem would predict the correct orbital ionization potential.

There are two terms that must be addressed in the parameterization of two-

center quantities at the INDO level of approximation. These are the two-center

repulsion integrals and the one-electron bond parameters or resonance integrals.

The generalized form of many two-center repulsion integrals used in semiempirical

theory can be expressed as

7AB(RAB) = (3-9)
[(A ) m +RAB] I/m
K m

In this equation fg is known as the Weiss factor [7] and RAB is the internuclear

separation. -AB(0) is the limit of -AB for which RAB=O and is usually expressed

as an average of the two one-center Fo's. The Weiss factor is set to 1.2 and m is

fixed at 1 according to the prescription of Mataga-Nishimoto [48]. With fg and

m fixed no further parameters are necessary for ^AB.

The INDO two-center one-electron matrix element denoted flv is usually

expressed as the average of one-center bond parameters weighted by the overlap

between two Slater orbitals and an appropriate fixed interaction factor, fi, for the

type of orbitals involved. For example,

S.5(3A + B)S, Y,V= s
PAVB = .5(15A + B)S, fs = s, V = p (3-10)
S5(PA + B)[fa(gaS,, + firgrSr,], /, V = p
The interaction factors for spectroscopic INDO [7] are

fs, = 1.0, f, = 1.267, f, = 0.585 (3-11)

The p-p interactions are divided between pure sigma and pure pi contributions.

The Euler rotation factors are denoted by g.

Properties such as spectral transition moments involve the calculation of

integrals that are not involved in the SCF procedure and have not themselves

been parameterized. In order to evaluate property integrals in the assumed

molecular basis, it becomes necessary to view the parameterized SCF procedure

as connected to an atomic basis set through orthogonalization. The majority

of INDO/S calculations performed to date involve the calculation of oscillator

strength in the dipole length formalism only. The dipole integrals are restricted to

one-center contributions and the SCF coefficients are maintained in the assumed

orthogonalized basis. For this work all property integrals will be calculated

over a Slater basis [49] and subsequently orthogonalized. This will lend to a

more systematic comparison between oscillator strength as calculated in differing

formalisms. Other properties such as NMR shielding involve many different

types of property integrals. Since there is currently no systematic way of

parameterizing or truncating these integrals in a balanced fashion, working in

a fully orthogonalized basis is the most logical point of reference to explore.

Benzene and Pyridine

The INDO/S method was originally parameterized to reproduce excitation

energies for benzene and pyridine with a truncated CIS procedure. In this section

the CIS and RPA methods are compared for these molecules. Oscillator strengths

in both length and velocity formalisms are calculated. The effect the size of the

active space has on oscillator strength and transition energies is also explored.

Benzene is of point group D6h of which the highest occupied molecular

orbital (HOMO) is a doubly degenerate pi orbital of Eig symmetry. The lowest

unoccupied molecular orbital (LUMO) is also a doubly degenerate pi orbital of

E2u symmetry. Single excitations from the HOMO to the LUMO generate excited

states that are of Blu, B2u and Elu symmetry. The only dipole allowed transition

from the ground state is to the degenerate E1u state. The experimental excitation

spectra for benzene consists of three bands. The two lowest bands indicate low

probability of absorption and the third has a large absorption confirming the simple

analysis above. The presence of the nitrogen in pyridine reduces the symmetry

of the molecule to C2v. The excitations from the upper two occupied pi orbitals

into the lower two unoccupied pi orbitals are all now symmetry allowed. The

split degeneracies from that of D6h are small however such that there remains

the same three bands as with benzene. In addition, excitations from the highest

occupied nonbonding orbital on the nitrogen to the lowest unoccupied pi orbitals

are also predicted. The lowest of these is allowed and the second is not. The

SCF orbital eigenvalues for the outer valence region for benzene and pyridine

are presented in Tables 3-1 and 3-2 and compared to the experimental orbital

ionization energies. The agreement with experimental ionization potentials is

reasonable with one noted exception. The nonbonded orbital is predicted to lie

below the pi orbitals but is experimentally found to be nearly degenerate with

the highest pi orbital. This represents the failure of Koopmans' theorem for such

orbitals and not the model Hamiltonian itself. Koopmans' theorem does not take

into account orbital relaxation which can be substantial for such orbitals. The

model Hamiltonian correctly accounts for the experimental ordering when such

effects are considered.

The calculated spectra for all configurations below 70,000cm-1 is presented

in Tables 3-3 and 3-4. Both RPA and CIS give good estimates for the excitation

energies. The RPA excitation energies are consistently lower than the CIS

energies. It is well known that CIS dipole length oscillator strengths are usually

too large and dipole velocity oscillator strengths are too small. From the form of

Equation (2-28) a smaller excitation energy is necessary in order to reduce length

oscillator strength in favor of velocity oscillator strength maintaining constant


Table 3-1 SCF orbital energies for benzene.
MO Energy(AU) EXPa
12 -0.4581 -0.4222
13 -0.4581 -0.4222
14 -0.3291 -0.3399
15 -0.3291 -0.3399
16 0.0304 -
17 0.0304 -
a Ref. [50]

Table 3-2 SCF orbital energies for pyridine.
MO Energy(AU) EXpa
12 -0.4784 -0.4575
13 -0.3709 -0.3554
14 -0.3633 -0.3859
15 -0.3328 -0.3601
16 0.0151 -
17 0.0290 -
a Ref. [51]

Eig (Tr)
Elg (7r)
E2u (r*)
E2u (r*)

A1 (n)
B2 (7r)
A2 (r)
B2 (7T*)
A2 (7r*)

transition moments. The RPA method only guarantees equivalence in a complete

basis. It is demonstrated here however, that the RPA method does a good job in

balancing the two formalisms with a finite basis.


Table 3-3 Electronic excitation spectra for D6h benzene.
Sym. CIS ft fy RPA fr fv EXPa

(cm-1) (osc)
'B2u 37797 37306 38090

1B1u 48806 48305 48972

Elu 54644 1.020 0.222 51566 0.678 0.541 55900
'Eln 54644 1.020 0.222 51566 0.678 0.541 (0.69)
a Ref. [52]

The nature of the length operator is to emphasize the long range region of

the basis functions whereas the velocity operator emphasizes the region closer to

the nuclei where the functions change the most. The dipole velocity oscillator

strength for the lower n-7r* transition for pyridine is larger than the length oscillator

strength for both CIS and RPA. The nonbonding orbital is essentially located on the

nitrogen which is more electronegative than carbon. The vicinity of the nitrogen

atom is expected to gain charge and thus the electron density in the region will

relax. If the assumed atomic basis took this fact into account the atomic velocity

matrix elements would be diminished in favor of the length matrix elements for the

pyridine nonbonding orbital. The INDO method however is based on a minimal

basis formalism which does not allow atoms to become polarized in an asymmetric

Table 3-4 Electronic excitation spectra for C2v pyridine.
Sym. CIS fr fv RPA fr f EXPa

(cm-1) (osc)
'B2(n) 35981 0.009 0.210 35804 0.008 0.230 34771

'A2(n) 44158 44125 -
'B1 38751 0.061 0.010 38120 0.055 0.042 38350

'Al 49991 0.067 0.020 49230 0.099 0.083 49750

1A1 56282 0.732 0.037 53970 0.510 0.220 "55000
1B1 56682 0.887 0.210 54045 0.594 0.449 (1.30)
a Ref. [52, 53]

fashion. The relaxation of all the orbitals on an atom can be accomplished however

with some success and will be the subject of other chapters.

In order to do calculations for large systems some truncation of the active

space must be accomplished. The effect the size of the active space has on the

spectra of small systems can be studied. For benzene and pyridine no qualitative

difference in the ordering of the excitation energies is noted in going from a very

small space of all transitions below 55,000cm-1 to a full singles space consisting

of transitions above 400,000cm-1. The very lowest transitions change very little


Table 3-5 Benzene 'Elu transition as a function of the active space.
fI fV RPA fr fV Space
55444 1.249 0.522 52208 0.829 0.786 55000
54644 1.020 0.222 51566 0.678 0.541 65000
53060 1.041 0.259 51026 0.708 0.567 80000
53850 0.994 0.198 50854 0.677 0.512 90000
53394 0.988 0.199 50407 0.681 0.509 100000
51054 0.770 0.043 48077 0.514 0.375 130000
50848 0.776 0.051 47863 0.520 0.382 200000
50139 0.768 0.056 47097 0.523 0.394 300000
50000 .771 0.061 46929 0.528 0.402 400000a
aFull active space

after all configurations below 65,000cm-1 are included. The quantitative nature of

the other excitation energies, however, systematically diminishes with increasing

active space for both RPA and CIS. This is perhaps caused by excitations to

the higher lying virtual orbitals calculated too low in energy when the overlap

is neglected and mix too readily with the lower lying transitions. The relative

magnitudes of the length versus velocity oscillator strengths remain the same for

both CIS and RPA. Their absolute magnitude becomes smaller with the size of

the space and reduces more rapidly for CIS than RPA. The data for the 'Elu

Table 3-6 Pyridine 'Al transition as a function of the active space.
fr fv RPA fr fv Space
57477 1.082 0.452 54674 0.717 0.646 55000
56282 0.732 0.037 53970 0.510 0.220 70000
55562 0.703 0.023 53247 0.493 0.198 85000
54426 0.645 0.001 52111 0.471 0.140 100000
52887 0.546 0.013 50392 0.382 0.122 150000
52685 0.546 0.010 50171 0.383 0.125 170000
52327 0.551 0.010 49743 0.389 0.138 240000
52074 0.564 0.005 49410 0.404 0.156 310000
51971 0.571 0.003 49267 0.411 0.164 430000a
aFull active space

excitation for benzene and the 1A1 transition for pyridine are presented in Tables

3-5 and 3-6 for illustration The other excitations follow the same trend.

On the basis of the above analysis the active space for the current model

should be truncated in the region 65-80000cm-1 in order to maintain quantitative

transition energies. All subsequent calculations will be truncated in this range

unless otherwise noted.


Table 3-7 SCF orbital energies for Naphthalene and Quinoline.
MO Energy(AU) Energy(AU)

Naph. Quin.
20 -0.44105 r
21 -0.42976 x -0.38831 r
22 -0.36977 7r -0.36693 n
23 -0.31545 7r -0.31520 r
24 HOMO -0.28856 7r -0.30397 7r
25 0.00043 7r* -0.00892 r*
26 0.02225 7r* 0.01954 7r*
27 0.05681 7r* 0.04968 7r*
28 0.08862 7*r 0.08374 r*

Naphthalene and Quinoline

One success of the INDO/CIS method is that it extends to molecules larger

than for which it was parameterized. For a model Hamiltonian this extension is

not guaranteed and becomes a test of the robustness of the method. In this section

the CIS and RPA methods are compared for naphthalene and quinoline which can

be viewed as extending the pi system of benzene and pyridine by one ring. Table

3-7 contains the higher occupied and lower unoccupied orbital eigenvalues for

these molecules.


Table 3-8 Electronic excitation spectra for Naphthalene
Type CIS f fv RPA f: f EXPa

(cm-1) (osc.)
7r->7r* 32138 0.004 0.005 31575 0.003 0.010 32000

7r->r* 37034 0.154 0.008 36059 0.139 0.119 37500

7->7r* 45469 1.844 0.561 43304 1.300 1.131 47500

7->r* 44630 0.000 0.000 44392 0.000 0.000 -
7->r* 46153 0.000 0.000 45594 0.000 0.000 -
7->r* 48551 0.621 0.118 46749 0.416 0.352 49500

a Ref. [54]

The electronic spectra of these molecules is more complicated than benzene

and pyridine due to the presence of excitations from the additional pi orbitals.

The calculated excitation spectra for naphthalene is reported in Table 3-8. The

RPA oscillator strength values are not only nearly equivalent, but are also in line

with experiment. This allows more confidence in the assignment of experimental

bands based on calculation. As with benzene and pyridine the RPA transitions

are systematically lower in energy than the corresponding CIS values and the CIS

dipole length oscillator strength is much larger than dipole velocity.

Table 3-9 Electronic excitation spectra for Quinoline.
Type (ff f RPA *fr fV
(cm-t) __________ (OSC.)
7r->r* 32467 0.027 0.002 31894 0.026 0.012 31900
n->7r* 34312 0.007 0.181 34139 0.007 0.198 to

7r->7r* 37525 0.108 0.002 36788 0.106 0.071 36200
7r->7r* 43277 0.597 0.182 42148 0.724 0.622 44300
n->7r* 45004 0.000 0.000 44968 0.000 0.000 (0.517)
7r->r* 46844 0.502 0.116 45618 0.470 0.387
7r->7r* 46997 0.858 0.278 45770 0.336 0.270
7r->7r* 48963 0.457 0.089 47582 0.157 0.135 (0.743)
a Ref. [55]

The results for quinoline are even more complicated due to the presence of

excitations from the nonbonding orbital. The calculated results are presented in

Table 3-9. The calculated pi excitations indicate a similar structure to naphthalene

with all bands allowed. The oscillator strengths for the nonbonded excitations are

like those of pyridine and a similar explanation would apply here.

The RPA performs better than the CIS method for these molecules in terms of

oscillator strength. The inability of the CIS method to guarantee equivalence even

in a complete basis is not grounds to reject its use for a model Hamiltonian. It is


demonstrated however, that it is not as reliable in the prediction of experimental

oscillator strengths. The RPA performs as well as the CIS method for the

qualitative ordering of the excitation energies and both methods do well for the

quantitative prediction of excitation energies.

The Diazines

The diazines consist of a planar six membered ring containing two nitrogen

atoms and, like pyridine, are isoelectronic with benzene. Two nitrogen atoms will

generate two low lying nonbonding orbitals. The presence of the single nitrogen

did not greatly perturb the band structure for the pi excitations of pyridine as

compared to benzene. The presence of a second nitrogen in the ring perturbs

the pi transition energies somewhat more but the band structure remains intact.

There are now four excitations from the nonbonding orbitals most of which are

not allowed and are not resolved experimentally. The calculated SCF energies

and spectra may be found in Tables 3-10 through 3-13.

The RPA is very capable of reproducing the spectra of the diazines with the

INDO/S Hamiltonian. For benzene, pyridine and the diazines, both CIS and RPA

predict the correct order of the lowest 7r-ir* and n-7r* transition energies with

respect to experiment. Oscillator strengths for n-7r* transitions are systematically

larger for dipole length than dipole velocity as previously observed. The highest

energy n-r* transition is predominately from the lower nonbonding orbital. For

Table 3-10 SCF orbital energies for the Diazines.
MO Energy(AU) Energy(AU) Energy(AU)
1,2-diazine 1,3-diazine 1,4-diazine
12(n) -0.43067 -0.42097 -0.42352
13(7r) -0.38078 -0.38745 -0.40242
14(n) -0.35278 -0.36630 -0.36475
15(r) HOMO -0.34776 -0.34754 -0.33096
16(r*) 0.00299 0.00554 0.00347
17(7r*) 0.01297 0.01524 0.01930
18(7r*) 0.09987 0.09927 0.10012

Table 3-11 Electronic excitation spectra for 1,2-diazine.
Type (cm ) fr f RPA fr fv (osc.)

28531 0.014 0.573 28227 0.013 0.671 26
35271 0.000 0.000 35169 0.000 0.000
46893 0.000 0.000 46826 0.000 0.000

51760 0.014 0.109 51705 0.013 0.117 5 5
39374 0.058 0.001 38741 0.054 0.022
r->r* 49600 0.109 0.039 48664 0.139 0.127 (000)
55363 0.464 0.003 53601 0.403 0.078
57635 0.910 0.260 55079 0.565 0.427
a Ref. [52, 53]

Table 3-12 Electronic excitation spectra for 1,3-diazine.

Type (cms) f f RPA fr fv (.)
(cm ) 1 (osc.)
34282 0.016 0.386 34129 0.015 0.426 31
36640 0.000 0.000 36506 0.000 0.000
45396 0.000 0.000 45313 0.000 0.000

51478 0.010 0.082 51433 0.010 0.087 51143
39803 0.068 0.012 39187 0.062 0.043 40310
50181 0.172 0.048 49139 0.189 0.149 52
Rb 56271
57530 0.652 0.078 55533 0.509 0.167 58500
r->r 57594 0.626 0.011 55700 0.416 0.198 (1.)
a Ref. [52, 53]
b Rydberg transition.

this transition the discrepancy in the calculated oscillator strengths is maintained

but the magnitude of the error is diminished from that of the lowest energy

transition. The two nonbonding orbitals are different in that the lower energy

orbital is less localized than the higher. The hypothesis that more tightly bound

orbitals will generate a larger difference in oscillator strengths holds for these

molecules. In all cases studied thus far the RPA 7-Tr* oscillator strengths are


Table 3-13 Electronic excitation spectra for 1,4-diazine.
Type (cm ) fr fv RPA fr fV ( )
(cm_ m L f IV (osc.)
30409 0.013 0.476 30164 0.013 0.540
n->r* 39056 0.000 0.000 38857 0.000 0.000
40706 0.000 0.000 40652 0.000 0.000
54538 0.000 0.000 54534 0.000 0.000 54000
36671 0.177 0.033 35609 0.146 0.108 38808
48358 0.232 0.024 47030 0.207 0.126
r->ir (0.15)
57694 0.253 0.081 56567 0.203 0.000 6
60653 0.857 0.256 58182 0.552 0.420 60
7t->r* 57900
59305 0.000 0.000 57786 0.000 0.000
Rb 55154
a Ref. [52, 53]
b Rydberg transition.

superior to those of CIS with respect to equivalence and experiment, and in most

cases that obtained from the dipole length operator are better when compared to

experiment than that obtained from the velocity operator.


The RPA for Extended Systems

One utility of modem model Hamiltonians is the ability to handle very large

systems. The INDO/CIS method has proven its relevance in this regard. In

this section the RPA and CIS method are compared for a series of extended

aromatic hydrocarbons. The series of molecules beginning with naphthalene

(two rings) and anthracene (three rings) is extended up to a linear molecule of

twenty rings. The twenty ring molecule consists of 82 carbons and 44 hydrogens

for a total of 372 basis functions. The density of spectroscopic states is large

for all molecules beyond anthracene; however, there are two transitions that

are of specific interest. The first of these is the lowest transition and gives

information about the HOMO-LUMO energy gap. The second of these is the

transition of maximum oscillator strength, which is polarized along the major

axis of the molecule and generates a large transition dipole. Therefore as the

molecule increases in length the transition dipole also increases. Data for these

two transitions for the different molecules are presented in Tables 3-14 and 3-15.

As the size of the molecule increases, the ratio of RPA dipole length to dipole

velocity oscillator strengths for the intense transition is almost constant whereas

those of CIS are varied. This indicates that for this type of transition the RPA

oscillator strengths scale properly with size whereas the CIS method does not.

However the CIS method does improve with size. Both CIS and RPA dipole

Table 3-14 Lowest excitation for extended aromatic system.
# rings CIS(cm-1) f f1 RPA ff fV
2 32145 0.003 0.005 31583 0.003 0.009
3 28505 0.012 0.022 27900 0.010 0.036
4 24761 0.241 0.003 23383 0.160 0.151
5 21482 0.252 0.030 20003 0.155 0.147
6 19258 0.269 0.050 17707 0.158 0.158
10 12272 0.284 0.623 10147 0.123 0.145
20 11045 0.541 0.159 8743 0.226 0.279

Table 3-15 Most intense excitation for extended aromatic system.
# rings CIS(cm-1) fI fV RPA fV fV
2 45434 1.832 0.557 43276 1.293 1.124
3 40131 2.703 0.999 38319 1.917 1.728
4 36509 3.479 1.462 34982 2.503 2.319
5 33882 4.171 1.887 32624 3.073 2.834
6 31988 4.829 2.368 30952 3.641 3.372
10 27656 6.948 4.072 27214 5.735 5.250
20 25179 11.712 8.278 25003 10.998 9.697

length oscillator strengths scale linearly with the number of rings as can be seen

in Figure 3-1. The RPA values are perfectly linear over the whole range of

molecular size whereas the CIS values differ to a small degree for the smaller

molecules. The linear regressions for the oscillator strengths are

f,(RPA) = 0.535n + 0.344 R = 0.9998

fC(CIS) = 0.531n + 1.344 R = 0.9945 (3-12)

n = # rings
The dipole length oscillator strengths scale roughly the same for the two methods

with the CIS values being systematically higher than the RPA values. For the low

energy transition the CIS oscillator strengths indicate no new trends. As with the

other molecules studied the RPA oscillator strengths are much more consistent

than those of CIS.

The excitation energies for both the RPA and CIS procedures follow pre-

dictable trends. Each line was fitted to a general exponential model in order to

extrapolate the transition energies to very large systems. The modeled transition

energies are


E(RPA) = 8612 + 28040 exp(-.0646n1625)

E(CIS) = 10910 + 25610 exp(-.0607n1654)

E(RPA) = 24540 + 48420 exp(-0.588n0691)

E(CIS) = 24710 + 47110 exp(-0.488n0751)
Both modeled energies and transition energies can be found in Figure 3-2.

13.00 .



.04 3.25

1 6 11 16 21
Number of rings

Figure 3-1 Dipole length oscillator strengths for the most intense transition.

The RPA and CIS methods scale similarly with size for a given transition.

The difference between naphthalene and the extrapolated values for very large

systems is roughly 21000cm-1 for both lines and both model calculations. The

RPA excitation energies for the lowest transition are systematically smaller than

the CIS values. The difference between the two levels for the larger molecules

remains constant at 2000cm-1. For the most intense transition the CIS and the

RPA transition energies converge to the same value as the chain length is increased.


47.0I ''

37.2 CIS

-' 27.5 Most intense transition

SExponential Model
17.7 -

Lowest energy transition

8.0 ,- ,,i-1 1 I ,
1 6 11 16 21

Number of rings

Figure 3-2 Modeled and calculated excitation energies



Much experimental and theoretical work has gone into characterizing the

electronic structure and spectroscopy of porphines [56]. These systems and

their excited states are of great interest in a variety of areas, and in particular

are essential in understanding the initial photochemical event in photosynthesis.

The electronic spectra of porphyrin compounds are characterized by three basic

regions. The so called Q bands are relatively weak and occur in the visible

region. The Q bands consist of a degenerate electronic transition for divalent

metalloporphines (e.g. MgP) and two separate electronic transitions for free base

porphines (H2P) separated by ~3300cm-1. The intense Soret or B region occurs

in the near UV and is often accompanied by a closely related N band of lower

intensity. The B band for divalent metalloporphines is assumed to be a degenerate

pair [57]. The higher UV bands in the third region are broad and are often of

near uniform intensity.

The earliest attempts to explain the Q and B bands for H2P were free electron

models and were based on assuming specific resonance structures for porphin that

involve an aromatic macrocycle of 18 7r electrons [58, 59]. This simple model

Figure 4-1 Free Base Porphine

correctly accounts for a weak band (Q) followed by a strong band (B). This 18

membered cyclic polyene can be assigned D18h symmetry. Within this group

excitations between the highest occupied molecular orbital e4u and lowest empty


e5g(7r*) molecular orbital gives rise to states of 1Blu, 'B2u and 1'Eu symmetry.

1B1, and 1B2u lie lowest, and give rise to the observed Q1 and QI peaks.

The Soret band, of 'E, type remains degenerate in this model. In contrast,

metalloporphines can be likened to 16 membered cyclic polyenes with 18 pi

electrons. In D16h, the HOMO-LUMO excitations give rise to states of 'E7u

and 'El symmetry, and both visible and Soret are predicted degenerate. Early

studies determined that the two Q bands for H2P were polarized perpendicular to

each other and that the B band was of mixed polarization [60]. Low temperature

studies of the B band resolved a splitting of 240cm-1 between equally intense

transitions [61]. Hiickel calculations by Gouterman using a four orbital CI model

proved adequate in explaining the splitting of the Q bands and predicted that the

lower Q band would be polarized parallel to the inner H-H axis [62, 63]. Further

polarization studies showed this to be the case and the B band was shown to

have parallel polarization at the low energy side of the band and perpendicular

polarization on the high energy side of the band [64].

By 1972 extensive PPP-CI calculations had failed to reproduce the assumed

B band splitting of 240cm-' for H2P [65]. The B splitting was predicted to be

greater than 1500cm-1 and other 7r--r* transitions beyond the four orbital model

were predicted to be part of the B and N bands now assumed to be coupled and

spread over a region 3300cm-1 wide [66]. Sundbom suggested that n-r* may

also contribute in this region [67]. Ab initio calculations, with scaled transition

energies, seemed to confirm the complexity of the Soret region suggesting as

many as five contributing bands including n-r* [68].

The inability to confirm the nature of the Soret region is complicated by

the fact that predicted oscillator strengths in this region do not conform to the

observed spectrum without making assumptions about vibrational broadening of

the BI relative to BI [67, 68]. In addition, calculated oscillator strengths in both

the length and velocity formalisms can differ by an order of magnitude [66]. The

most commonly used dipole length formalism usually predicts oscillator strengths

too large in the Soret region by a factor of two.

Encouraged by the ability of the INDO/RPA method to predict oscillator

strengths in addition to excitation energies for aromatic hydrocarbons and het-

erocycles the spectroscopy of free base and magnesium porphin is reexamined

in detail.


CIS and RPA calculations have been performed for H2P for several active

space partitions. The number of occupied valence orbitals for H2P are 57. Orbitals

48 and 49 are the nonbonding orbitals located on the nonprotonated pyrrole ring

nitrogens and 50-64 are r. No orbitals outside of this range had significant

contribution to the lowest 30 states. Transition energies and oscillator strengths

Table 4-1 Calculated and Experimental Excitation energies for H2P.
Evb ff fE Ev fE fv Trans.c E4
13.72 0.022 0.0001 11.81 0.020 0.021 QI 15.94
16.58 0.033 0.001 15.92 0.033 0.029 Q_ 19.55
27.18 1.616 0.114 24.08 1.146 1.003 B 26.85
28.42 2.411 0.306 24.44 1.228 1.224 Bj_
28.98 28.70 -
29.52 28.98 -
32.97 1.478 0.162 31.69 0.424 0.343 N1i 29.41
33.58 33.41 -
34.75 34.42 -
35.64 0.366 0.063 35.30 0.169 0.115 1
a Ref. [69]
b Energies are reported in 103cm1.
c Designations II and refer to transitions in the plane of the
molecule polarixed parallel of perpendicular to the inner H-H axis.

for the lowest 10 transitions for an active space of 14 occupied and 14 virtual

orbitals are included in Table 4-1.

The lowest occurring n-r* for the CIS calculation is at 39637cm-' with

a dipole length oscillator strength of 0.019 and for the RPA calculation 39600

and 0.019. In order to ensure that the 14x14 active space is adequate for the

determination of this transition a 1060 configuration calculation was performed

consisting of 53 occupied and 20 virtual orbitals. For CIS the Q, B and N


Table 4-2 Calculated and Experimental Excitation energies for MgP.
Evb f fy Ev fr fv Trans.c Ev
15.77 0.025 0.002 14.64 0.023 0.025 Q
15.77 0.025 0.002 14.64 0.023 0.025 Q
28.62 2.495 0.288 24.95 1.351 1.243 B
28.62 2.495 0.288 24.95 1.351 1.243 B
28.71 28.65 CT
29.58 28.98 7r->r*
30.46 0.032 0.049 30.41 0.030 0.055 CT
31.80 -31.49 r->r*
32.38 32.14 -->
33.55 0.133 0.012 33.19 0.029 0.022 N
33.55 0.133 0.012 33.19 0.029 0.022 N
a Ref. [70]
b Energies are reported in 103cm1.
c CT designates a charge transfer from occupied 7r to the Mg atom.

bands were lowered by less than 600cm-1 and the n-7r* transition was lowered

to 38375cm-1 for CIS and 38254 for RPA.

A small 7x7 configuration calculation was performed to check the sensitivity

of the lower 7r-Tr* to the size of the active space. The results for the 6 lowest

allowed states are qualitatively the same as those in Table 4-1 being generally

400cm-1 higher in energy with differences in oscillator strength less that ten


The results for the low energy portion of MgP are found in Table 4-2. The

major features of the spectrum consist of a weak degenerate Q band followed by

an intense degenerate B band. The Mg is set 0.4 A above the ring as suggested

by X-ray structures [71]. With the addition of the metal atom new features appear

in the spectrum, namely the forbidden or nearly forbidden charge transfer bands.

The positions of these bands are very sensitive to the Mg out-of-plane coordinate

which is easily perturbed. Similarly to H2P, several forbidden 7r->7r* transitions

are calculated to fall between the B and N bands.


Calculations using the RPA procedure yield an explanation of the Soret region

for H2P that conforms to experiment without invoking vibrational coupling to

explain the predicted pattern of oscillator strengths or scaling excitation energies

in a least squares manner. While the RPA formalism only ensures equivalence

between length and velocity oscillator strengths in a complete HF basis the

utility of the procedure to deal with transition moments in a balanced fashion

is again demonstrated. This balance of transition moments afforded by the RPA

is important if properties such as optical rotary strength in addition to oscillator

strength are to be calculated. With the CIS procedure the two moments are very

different as already noted by others [66]. In addition the RPA oscillator strengths

are in better accord with experiment [70]. The experimental B band oscillator

strengths are estimated to be 1.1 as predicted here. The experimental oscillator

strengths for the Q bands are less that 0.1 with QjI)Q1 as predicted by both the

CIS and RPA calculation. The PPP Hamiltonian has to be modified in order to

reproduce this result [65, 67]. The experimental N band occurs on the shoulder of

the high energy side of the intense B band with experimental oscillator strength

much less than that of the B band. This fact is inferred by the RPA results but

not CIS.

Both the CIS and RPA procedures predict the correct ordering of the Q and

B bands for H2P as deduced from polarization experiments [64]. In addition both

procedures predict the splitting of the Q bands in accord with experiment. The

splitting from the CIS procedure is 2860cm-1 and from the RPA 4010cm~1.

For gas phase H2P the experimental splitting is 3610cm-1 [69]. For the B band

the RPA calculation yields a splitting of 360cm-1 which is in excellent accord

with the original low temperature splitting of 240cm-1. The failure to reproduce

this splitting theoretically led to the suggestion that this splitting was due to a

vibrational progression and that the second "B" band actually appeared as "N",

leaving unexplained the reasons for the sharp BII and broad B (N) [66]. Even

for the 1060 configuration active space the CIS procedure predicts a splitting as

large as 1000cm-1.

These calculations therefore support the interpretation that the Soret band

for H2P consists of two nearly degenerate transitions of roughly equal oscillator

strength. The N band consists of 7r-Tr* beyond the four model. No evidence is

found to suspect the presence of n-ir* below 38000cm-1 or that other allowed

7r--T* are present in the Soret region.

No new observations can be attached the MgP spectrum aside from the ability

of the RPA to balance transition moments for oscillator strength. As noted

previously, the frequencies calculated by the RPA method are generally lower

that those calculated from the CIS treatment. Both calculations are in very good

accord with experiment.

The ability of the INDO/RPA method to represent a balanced treatment of

transition moments has been established through these calculations and those in

Chapter 3. The correspondence with experimentally observed transition moments

yields confidence in the method to predict and interpret spectra. The INDO/RPA

method has also been shown to not degenerate with molecular size so that very

large systems can be studied with confidence.


The Concept of Chemical Shift

The ability to quantitate and establish the nature of "chemical shift" as it

relates to magnetic resonance spectroscopy has been, and is, of great interest

in quantum chemistry. The phenomenological Hamiltonian most often used for

this purpose is usually attributed to Ramsey [72, 73] and most commonly takes

the form of a perturbative expansion of the energy of a system with respect to an

external magnetic field and the field of the nuclei of interest [74]. The perturbation

is small with respect to the total energy of most chemical systems; therefore, the

perturbative expansion can be truncated at low order and is usually assumed to

be linear. Many empirical and semiempirical models have been developed to

predict specific chemical shift in terms of the chemical properties of the atoms of

interest and the specific bonding environment [75, 76]. These models and rules

are very useful in interpreting the nature of the local chemical environment in a

molecule. Direct nonempirical calculations of chemical shift have met with only

modest success and normally requires very large basis sets to obtain quantitative

accuracy [77]. The failure of smaller basis sets is usually attributed to their

inability to represent the gauge invariance of the origin of the externally applied

field. As a result, in the last few years methods to minimize the effect of gauge

variance have been developed with much success [78-81].

In this Chapter the local orbital, local origin method of Hansen and Bouman

[81] for the calculation of isotropic shielding of the atoms in a closed shell

molecule is developed for the INDO/S Hamiltonian. This method is chosen

because it utilizes the unique properties of the RPA already discussed, and because

the results obtained are easily analyzed in terms of the chemically useful concept

of localized effects in molecules.

The magnetic energy of interaction of a nucleus in an external magnetic field

can be expressed as [82]

U = -.B (5-1)

where i is the magnetic moment of the nucleus and B is the external magnetic

field. The magnetic moment of the nucleus is quantized. The operator form of

j is

i= 7, I = h[(I 1)]1/2 (5-2)

where P is the angular momentum operator and 7 is the gyromagnetic ratio of the

specific nucleus. For an oriented molecule in an external field of magnitude B

aligned along the z axis and with mi the magnetic quantum number of the nucleus

we can express the interaction energy

Uz = -yhmiB


Differences in magnetic quantum numbers represent differences in interaction

energy. The absolute difference in energy for a change of one quantum number

is then

AUz = hv = I7yhB (5-4)

the frequency of which is

v= B = KB (5-5)

Therefore, in the absence of any other factors, energy differences for unit changes

in magnetic quantum number are linearly dependent on the external magnetic field

and are a function of the specific nuclei only. It is this difference in energy with

applied field that is measured in a magnetic resonance experiment. If no other

factors were involved, a single measurement of for a type of nucleus would

be sufficient to describe all situations dealing with the same nucleus and NMR

would be of very little interest in chemistry. Nuclei are not normally isolated,

however, but in the presence of electrons. The external field induces a current in

the electrons. This current in turn creates a magnetic field that is in opposition to

the applied field. The net field (B) at the nucleus of interest is less than the applied

field (Bo). Different electronic environments will modify the total field differently.

Each chemically unique environment for a given type of atom in a molecule at


fixed external field strength will thus result in a different frequency. The induced

field of the electrons is also linearly dependent on the external field so that
vj = KBo(1 aj) = vo(l aj) =

Vo vj (5-6)
'j = ---
oj is the unitless shielding constant for environment j with respect a bare nucleus

and represents the effect of the electrons in said environment. With respect to an

arbitrary reference environment we define relative chemical shift parameters

jf Vj Vref
S= 0ref 17i = (5-7)

Since |vj Vref ((vo and vo e vrf relative chemical shift is usually tabulated in

terms of parts per million

Sj(ppm) = 106 vi Vref (5-8)

In this formalism positive values of relative chemical shift represent an environ-

ment that is less shielded than the reference environment.


We seek a quantum mechanical description of a which, in general, is depen-

dent on the directional components of the applied field and the components of the

nuclear magnetic moment. Equation (5-1) is modified to read

U=-l-.[i ]&B=- .B+j7.B
= Uo + Ui

U represents the total energy of interaction and U1 the electronic perturbation.

A magnetic field can be expressed in terms of a vector potential that must be

incorporated into the Hamiltonian of the system [17]. The kinetic energy portion

of Equation (2-1) is modified to read

ila 1-. 1-
H=1 + [-iVi + -Ai] [-iVi + -Ai] + V (5-10)

where c is the velocity of light and A is the magnetic vector potential. The vector

potential for a uniform magnetic field is

x [i -R] (5-11)

with arbitrary gauge origin R. All coordinates are relative to the nucleus of

interest. The vector potential associated with the magnetic dipole moment of the

nucleus is expressed [20]

iA = (5-12)

We choose the coulomb gauge so that V A= 0 and expand Equation (5-10)

according to Equations (5-11) and (5-12)

H = Ho + :'

ft + e-x (F-i R-) + 2 (5 -13)

i 12 x(i-R)+ r ] (5-13)

= B- {( R) xVi + 2
i=l i

+ ^1 1xiBx(-) ]2
+ c gx(i-I)+ r ]2

i=1 1
where []2 is a shorthand notation for [- ] [ ]. The last equality in Equation (5-13)

represents the perturbation Hamiltonian of interest here. The terms quadratic in

field are dropped due to the observation that shielding is practically linear with

field according to Equation (5-9). The final expression for the perturbation is

n .
f'= _- I(i- ) x i} + 2 (i3 Vi)
i=l ri
1 [X (i tx)]( [/ i] (5-14)
L.2c2 [r
= H'(1) + I1-I'(2) = HB (1) + H'(1) + IH'(2)

Perturbations to molecular energies due to the presence of a magnetic field

are small so a low order expansion of the perturbation Hamiltonian is adequate.

To second order with a Hartree-Fock reference wavefunction [19]

E' =E' + E2 = HFI + (V)HF I flH k ) k If klHF) (5-15)
where the bk represent excited states of energy Ek. Placing Equation (5-14) into

(5-15) yields
E' 2- {Bi ( HF(?i R) X ViJbHF)

+2 ( FWH}(5 16)
n 1 [I x (i )]. [j7xFi]
+ E 22(HFI 3 IFHF)

= E1(1) + E'(2)

(OHF lII' k)(k HF) =

(O&HFIH'(1) + fH'(2) Ibk)(k IH'(1) + H'(2) IHF)
= ('HF| S'(1)'k)(|'kl'(1)|lHF) + (HHF H'(2)|4'k)(kJ|H'(2)|'HF)

+( HFfl'(1)|lk)( klIH'(2)| IHF) + (OHFIH'(2) Ik)(lkkIH'(1) CHF)
Recalling Equation (5-9) we seek terms that are linear in both external and nuclear

field. It is conventional to interpret experiments in terms of the phenomenological

Hamiltonian of Equation (5-9), therefore is convenient to gear this development

to yield numbers comparable to the interpretation of these experiments. It would

be more systematic, however, to calculate the observed experimental spectrum

directly and would make a logical project for the future. The first term in Equation

(5-16) is only linear in a single field and can be dropped. In Equation (5-17)

all terms involving H'(2) are quadratic or higher in field. Dropping these terms


( HF| k '(1)| k |'(l)| IHF) =

(OHFIH'B(1) + HII(1)[ k )(|k IB(1) + H (1)I) HF) =
(OHFIfB(1)|0k(k|IfIB(1)|bHF) + (kHF IH(1)|k)(OkkIH,(1)|bHF)+

(OHFlH'B(1)|k)(IkIHI(1)|OHF) + (FHFHI (1)|lk)('klIH'(1)|VHF)
Collecting those terms in Equation (5-18) that are of the proper form leaves

1 [I x -(B i iR)] [ x Fl)
U1 = 22 (HF x r3 Ir HF)-
20 P
n x( Xi n _H (5-19)
(oo HF I E .[i] R'k)k E B Rfi ) x Vi]IHF)
EV i=l i=1
c2 k= EHF Ek

The second term in Equation (2-19) has the proper form as suggested by Equation

(5-9). Through vector identities the first term can be rearranged and Equation

(5-19) reads

U1 = L -Cc2 (HF i=l r3IHF

oo (lkHF k k)(k E i[(ij R) XVi]IlHF) (5-20)
i=1 i i=1
c2 EHF Ek

identifying the form of shielding desired. Let e, be a unit vector along the

direction of the nuclear magnetic moment and gB along the applied field.
21 (F [gB x(Fi- J)]. [ xfi]
(pB = 2 (HFl r3 HF)

1 (rrl Z k k B (i ) X Vi 1HF (5-21)
ii= i=1

Equation (5-21) is the basic equation for the calculation of chemical shielding

[73]. The first term can be interpreted as that shielding created by the induced

current around the nuclei of interest and is known as the diamagnetic component

of shielding (ad). The second term is known as the paramagnetic component

of shielding (aP) and can be interpreted as the degree to which the specific

chemical environment quenches the induced current around the nucleus. The

two components are opposite in sign and a delicate balance must be achieved in

order to obtain even qualitative results as compared to experiment.

The transition matrix elements in the paramagnetic portion are over anti-

Hermitian one-electron operators. Being one-electron operators the sum can

be truncated at single excitations only. Transition moments for anti-Hermitian

operators according to the analysis in Chapter 2 for a closed shell reference are

(bHylO = 1(Oic) [X k- Y]yc,k
rc A (5-22)
( k O HF) = V E[X_- Y]d,k dlOA()

The energy in the denominator is the transition energy for the excitation (EK).

Pre and post multiplication of Equation (2-45) yields,
(A B)-1 // (A B)(X Y) = (X + Y)EK // EK1

(X Y)EK = (A B)-1(X + Y) // (Xt Yt) (5-23)

(X Y)EK'(Xt yt) = (A B)-
Incorporation of Equations (5-22) and (5-23) into the paramagnetic shielding

term yields

PB = ( | c)(A B) 7 ( 1dO| 6) (5-24)

Let Ry be the location vector of occupied orbital relative to the nucleus of

interest. The field dependent transition integral in Equation (5-24) becomes
( |dIO )=6 d| (OdI(-- R) xV|0).B

= [(Idl(F- 1r) x V|6) + (|dl(I( iR) x V 6)] E-B (5-25)

= [(dl(d-- Ry) x V16) + (R1 R) x (Od|V|6B)] eB
so that
2 R
aB = (0 ( |lOIc)(A B)d(d(- R) x V16) B
2 (5-26)
+- (I~IA c)[(R, R) x [(A B)-6d(d I I)]] B
According to Equation (2-27) and the properties of the RPA equations in a

complete Hartree Fock basis

(A B)lVd(dVl ) = -(cir-7)


The paramagnetic term becomes

S= ( E (^0,1`c,)(A B)'cd(dl(rd RX V) 6 B
-ycd (5-28)
-4 Z(0,16jI10c)(Th -k x (qciioM)]
C2 ( ,KA R R) X c|W B

where all references to the global gauge origin are contained in the second sum.

For a closed shell reference and local origins the diamagnetic portion of

shielding becomes

Bd -1B x r[ Bx -R)]- -) [E xr
orB= B 2^\^1 -- 3 100r)
+1 [B X (R7 R)] ( x [ r- r0)
C2 r
where again all references to the global origin are contained in the second sum.

Combining the origin dependent sums of Equations (5-28) and (5-29) and

using the vector relation A (B x C) = C- (A x B) yields for the global origin

dependent term
1 B x ( )] )

2 E [Fx V]
-21 r3 )(cfl lx) [E'B x (R, R)]
Utilizing the orbital completeness relation [23]

I,3()(YI + c)(CI 5-31)
7 c (5-31)

c 7

gives for the origin dependent term

1 I [4. x ri
B (R) = E (R- R)] -[I I O-

2 [ [i(x V]
T *

{[WB x (R, R)]. [r- ,)]}

2 [Ix ]
+C IE(l5Iii ) [)(1 [B x (R,- R)]
C E 71 r | I ) ( r ) I X

The second sum can be rearranged

2(1 7{6e [FX V] }

{[1B x (R, R)1. [F1)]}

Ell [F ] [Ex 00]1
= [FB x (R, R)] 2(0{ x V]} r)

-_g)].[ -IP I [ ]- (( { *

=[B x (R1, R)]. (0 [x ], )]


and combined with the first sum

aB(R) = EB x (R ] 7-R) {[F x r [ [f x V],r?}|O1)



= [eB x (R

10 ]

Consider the general relationship

j (F x V)[(ji. r|1")] = -j." F x [V(. r-)thI)]

= r3 x [e1-0) + (E," 1 )]

= (- x rI. (Vi,10)) + (j X rI. [( Re.- r-3 10l')]
(j x f)r (4 0)) = e. (F x V)[(ei -l r )]( (E r)(e x [V )]
= E. (r x V)[(Ei. -F)1)] (ES. r)(.j rF) x [V[t,]

I(ej x rx) l [Ej (Fx V), (-i rf)] ) = 0

the last line of which indicates that the first sum in Equation (5-34) vanishes.
The indices of summation for the remaining term are over all occupied orbitals
and can be rewritten as

e ] 1x,] [
^oIAB (R) = 1 E'l ) B X ([F R)]

+( '[x V]r)(31 [B x (R R)]l (5-36)

1c Z 7 [l'rx I)( I)" [(B x (RI, R,6)]
where all references to an external gauge origin have been eliminated from the
expression for the shielding.

The final working expression for shielding becomes

1B = -' .(-r(, g>B)(T(-- rRi) i-- g (?-- l )(g'B r -)}r-al31-)

1 3x,

+4z; (^I1] KA -q$ )(-R) X V|q$)

property integral and the construction and inversion of the RPA (A-B) matrix


The expression for shielding is a 3x3 matrix that is a function of the orientation

of the molecule and the external field. With the principle exception of solid state

NMR, the orientation of molecules with respect to an external field is random and

shielding becomes an average of the principle components of the shielding matrix.

It is this isotropic value in which we are immediately interested. The following

discussion is for 13C isotropic shielding. 13C is chosen for study because the range

of normal shielding values is 200 ppm and is a well characterized experimental


Localization and Integral Evaluation

The occupied orbitals are localized according to the criteria of Foster and

Boys [83]. A symmetric transformation of the orbitals is chosen such that the

sum of the squares of the pairwise differences in their locations is a maximum.

In other words,

max (R = R Zi 612)
7)8 (5-38)

An iterative process of successive two by two rotations is performed until all

successive pairs are stable within a given threshold [84]. Local orbitals obtained

in this fashion are easily identified as bonds between atoms and nonbonding

orbitals. Double and triple bonds are represented by tau or banana bonds. These

are identified by the program and are selectively delocalized to yield the more

familiar local sigma and pi bonds. The property integrals and RPA matrices in

Equation (5-37) are first evaluated in the canonical molecular orbital basis via

transformation from the appropriate atomic integrals and subsequently transformed

to the local basis.

Three choices of local origins in the calculation of molecular shielding are

evaluated. The first is the coupled HF (CHF) approximation in which all orbitals

are assumed to be located at the nucleus of interest. The second choice assigns all

local orbitals to be at their centroid location. As defined by Hansen and Bouman

[81] this is known as the full local orbital, local origin method (FLORG). The

final choice, intermediate between the first two, assigns all local orbitals that are

bonded to an atom to be at said nucleus and all other orbitals are assumed to be at

their centroid location (LORG). Since in a complete basis the results are invariant

to the choice of origins, the difference between results in a finite basis tests the

sensitivity of the results to the basis type and choice of origins.

There are many atomic integrals necessary for the evaluation of Equation

(5-37) that are not normally associated with those required to perform an SCF

procedure. Fortunately all integrals of the type used here can be put in the form of

a normal overlap integral or a nuclear attraction integral that are easily evaluated.

For atomic basis functions of the form [85]

XA) = NA n,, ny, nz, RA) = NA zA exp(--ArA)) (5-39)

located at center A, dipole integrals, velocity integrals, and angular momentum

integrals can be evaluated with overlap routines. The X component of each is

outlined below. For dipole,

= XB(XAIXB) + NB(XAInx+ l,ny,nz, RB)
and for velocity.

(XAIVxIXB) = NBnx(XAnx l,ny,nz, RB)
-2CBNB(XAInx + 1,ny ,nz, RB)