• TABLE OF CONTENTS
HIDE
 Title Page
 Acknowledgement
 Table of Contents
 List of Tables
 List of Figures
 Abstract
 Introduction
 Exact Xa method and muffin-tin...
 Muffin-tin charge density...
 Hellmann-Feynman force calcula...
 Dipole moment calculations
 Conclusions
 Appendix A: Derivation of the muffin-tin...
 Appendix B: The multiple scattering...
 Appendix C: Hellmann-Feynman theorem...
 Appendix D: Proof of the independence...
 Appendix E: Expansion theorems...
 Appendix F
 Appendix G: Evaluation of the sine...
 Bibliography
 Biographical sketch






Title: Hellmann-Feynman forces and dipole moments using multiple-scattering Xa wavefunctions for diatomic molecules /
CITATION PDF VIEWER THUMBNAILS PAGE IMAGE ZOOMABLE
Full Citation
STANDARD VIEW MARC VIEW
Permanent Link: http://ufdc.ufl.edu/UF00098149/00001
 Material Information
Title: Hellmann-Feynman forces and dipole moments using multiple-scattering Xa wavefunctions for diatomic molecules /
Alternate Title: Multiple-scattering Xa wavefunctions for diatomic molecules, Hellmann-Feynman forces and dipole moments using
Physical Description: ix, 156 leaves : ill. ; 28 cm.
Language: English
Creator: Li, Choy Hing, 1946-
Publication Date: 1975
Copyright Date: 1975
 Subjects
Subject: Dipole moments   ( lcsh )
Molecules   ( lcsh )
Physics thesis Ph. D
Dissertations, Academic -- Physics -- UF
Genre: bibliography   ( marcgt )
non-fiction   ( marcgt )
 Notes
Thesis: Thesis--University of Florida.
Bibliography: Bibliography: leaves 153-155.
General Note: Typescript.
General Note: Vita.
Statement of Responsibility: by Choy Hing Li.
 Record Information
Bibliographic ID: UF00098149
Volume ID: VID00001
Source Institution: University of Florida
Holding Location: University of Florida
Rights Management: All rights reserved by the source institution and holding location.
Resource Identifier: alephbibnum - 000171233
oclc - 02954297
notis - AAT7658

Downloads

This item has the following downloads:

PDF ( 3 MBs ) ( PDF )


Table of Contents
    Title Page
        Page i
    Acknowledgement
        Page ii
    Table of Contents
        Page iii
        Page iv
    List of Tables
        Page v
    List of Figures
        Page vi
        Page vii
    Abstract
        Page viii
        Page ix
    Introduction
        Page 1
        Page 2
        Page 3
    Exact Xa method and muffin-tin approximate Xa method
        Page 4
        Page 5
        Page 6
        Page 7
        Page 8
        Page 9
        Page 10
        Page 11
        Page 12
        Page 13
        Page 14
        Page 15
        Page 16
        Page 17
        Page 18
    Muffin-tin charge density approximation
        Page 19
        Page 20
        Page 21
        Page 22
        Page 23
        Page 24
        Page 25
        Page 26
        Page 27
        Page 28
        Page 29
        Page 30
        Page 31
        Page 32
        Page 33
        Page 34
        Page 35
    Hellmann-Feynman force calculations
        Page 36
        Page 37
        Page 38
        Page 39
        Page 40
        Page 41
        Page 42
        Page 43
        Page 44
        Page 45
        Page 46
        Page 47
        Page 48
        Page 49
        Page 50
        Page 51
        Page 52
        Page 53
        Page 54
        Page 55
        Page 56
        Page 57
        Page 58
        Page 59
        Page 60
        Page 61
        Page 62
        Page 63
        Page 64
        Page 65
        Page 66
        Page 67
        Page 68
        Page 69
        Page 70
        Page 71
        Page 72
        Page 73
        Page 74
    Dipole moment calculations
        Page 75
        Page 76
        Page 77
        Page 78
        Page 79
        Page 80
        Page 81
        Page 82
        Page 83
        Page 84
        Page 85
        Page 86
        Page 87
        Page 88
        Page 89
        Page 90
        Page 91
        Page 92
        Page 93
        Page 94
        Page 95
        Page 96
        Page 97
        Page 98
        Page 99
        Page 100
        Page 101
        Page 102
        Page 103
        Page 104
        Page 105
        Page 106
        Page 107
        Page 108
        Page 109
    Conclusions
        Page 110
        Page 111
        Page 112
        Page 113
        Page 114
        Page 115
        Page 116
        Page 117
        Page 118
        Page 119
    Appendix A: Derivation of the muffin-tin Xa one-electron equations
        Page 120
        Page 121
        Page 122
        Page 123
        Page 124
        Page 125
        Page 126
        Page 127
        Page 128
        Page 129
    Appendix B: The multiple scattering method
        Page 130
        Page 131
        Page 132
        Page 133
        Page 134
    Appendix C: Hellmann-Feynman theorem and the Xa wavefunctions
        Page 135
        Page 136
        Page 137
        Page 138
    Appendix D: Proof of the independence of dipole moment on the choice of the origin for a neutral molecule
        Page 139
        Page 140
        Page 141
    Appendix E: Expansion theorems for the Bessel functions and the Neumann functions
        Page 142
        Page 143
        Page 144
        Page 145
    Appendix F
        Page 146
        Page 147
        Page 148
    Appendix G: Evaluation of the sine and cosine integrals
        Page 149
        Page 150
        Page 151
        Page 152
    Bibliography
        Page 153
        Page 154
        Page 155
    Biographical sketch
        Page 156
        Page 157
        Page 158
        Page 159
Full Text












HELLMANN-FEYNMAN FORCES
AND DIPOLE MOMENTS USING
MULTIPLE-SCATTERING Xa WAVEFUNCTIONS
FOR DIATOMIC MOLECULES











By

Choy Hing Li


A DISSERTATION PRESENTED TO THE GRADUATE COUNCIL
OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT OF THE
REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY




UNIVERSITY OF FLORIDA


1975













ACKNOWLEDGEMENTS


The author would like to express his gratitude to

Professor John W. D. Connolly who suggested the topic of

this study and would like to thank him for his guidance,

instruction and criticisms during this research.

The author would also like to acknowledge the other

members of his advisory committee for their suggestions

and criticisms. Also, thanks are due to all my other

colleagues of the Quantum Theory Project at the University

of Florida. The author particularly extends his apprecia-

tion to Dr. B. Danese for his many valuable discussions,

continuous encouragement and especially for providing the

author with his non-muffin-tin total energy program.

Finally, the author would like to thank his wife

for her patience and encouragement.

This research was supported in part by the grant from

the National Science Foundation which is also appreciated.














TABLE OF CONTENTS


Page


ACKNOWLEDGEMENTS..........................................ii

LIST OF TABLES. .............................................v

LIST OF FIGURES............................................ vi

ABSTRACT .................................................vii

CHAPTER I INTRODUCTION................................ 1

CHAPTER II EXACT Xa METHOD AND MUFFIN-TIN
APPROXIMATE Xa METHOD ......................4

2.1 Introduction....................................... 4
2.2 Exact Xa One-Electron Schr6dinger Equations.......4
2.3 Muffin-Tin Approximate Xa One-Electron
Equations and the Solutions....................... 9
2.4 Relation Between the Exact Xa Method and
tne Muffin-tin Approximate Xa Method..............13

CHAPTER III MUFFIN-TIN CHARGE DENSITY APPRUXIMATION....19

3.1 Introduction.....................................19
3.2 Muffin-Tin Xa Total Energy Calculations..........20
3.3 Force Calculations...............................21
3.4 Dipole Moment Calculations.......................30

CHAPTER IV HELLMANN-FEYNMAN FORCE CALCULATIONS........36

4.1 Introduction. ................. ................. 36
4.2 Prolate spheroidal Coordinates...................36
4.3 Hellmann-Feynman Theorem..........................39
4.4 Method of Hellmann-Feynman Force Calculations....41
4.5 Results and Discussion ...........................62

CHAPTER V DIPOLE MOMENT CALCULATIONS .................75

5.1 Introduction......................................75
5.2 Method of Dipole Moment Calculations.............75
5.3 Results and Discussion.............................99

CHAPTER VI CONCLUSIONS ............................... 110












APPENDIX A


APPENDIX B

APPENDIX C


APPENDIX D



APPENDIX E


APPENDIX F

APPENDIX G


Page


DERIVATION OF THE MUFFIN-TIN
Xa ONE-ELECTRON EQUATIONS...............120

THE MULTIPLE-SCATTERING METHOD..........130

HELLMANN-FEYNMAN THEOREM AND THE Xa
WAVEFUNCTIONS........................... 135

PROOF OF THE INDEPENDENCE OF DIPOLE
MOMENT ON THE CHOICE OF THE ORIGIN
FOR A NEUTRAL MOLECULE..................139

EXPANSION THEOREMS FOR THE BESSEL
FUNCTIONS AND THE NEUMANN FUNCTIONS.....142

EXPANSION THEOREM FOR Yo(' )/r2 ........ 146
1 A A
EVALUATION OF THE SINE AND COSINE
INTEGRALS. ...............................149


BIBLIOGRAPHY ............................................153

BIOGRAPHICAL SKETCH .................................... 156













LIST OF TABLES


Table Page

4.1 Force Calculations for H 1( : I 2) 64
2 g g
4.2 Comparison of Ftotal and 3
4.3 Z-Dependence of Hellmann-Feynman Force
for H2 at R=1.4 a.u. 68

4.4 i-Values in the MSXa Calculation for N2(1) 70

4.5 Hellmann-Feynman Forces for N2 at R=2 a.u.
Calculated with MSXa and Hartree-Fock
Wavefunctions 71

5.1 Dipole Moment Components for LiH(1 +)
at R=2.75, 2.836, 3.015, 3.25 a.u. 101

5.2 Dipole Moments and Their Derivatives of
LiH at R=3.015 a.u. 103

5.3 6 and e-Dependence of Dipole Moment for
LiH at R=3.015 a.u. 104

5.4 Dipole Moment for BH(1E+) at R=2.336 a.u. 107

5.5 Dipole Moment for CH( 2) at R=2.124 a.u. 108













LIST OF FIGURES


Figures Page

2.1 Flow Chart of the Conceptual Exact
Xa Method 7

2.2 Flow Chart of the Conceptual Exact
Xa Method and the MSXa Method 15

3.1 Partitioning of the Space of a Homonuclear
Diatomic Molecule in the MSXa Method 23

4.1 Representation of a Point in the Atomic
Coordinates 37

4.2 Relation of rA, R and 8A 49

4.3 Boundary of the Intersphere Region in Prolate
Spheroidal Coordinates with A and B as the Foci 49

4.4 Force Curves for H2 67

5.1 Schematic Representation of A Heteronuclear
Diatomic Molecule in the MSXa Calculation 76

5.2 Boundary of the Intersphere Region in
Prolate Spheroidal Coordinates with A and O
as the Foci 83

5.3 Boundary of the Intersphere Region in
Prolate Spheroidal Coordinates with O
and B as the Foci 89

5.4 Boundary of the Intersphere Region in
Prolate Spheroidal Coordinates with A and B
as the Foci 96

5.5 Dipole Moment v.s. Separation for LiH ( 1+) 102

6.1 Schematic Representation of a Triatomic
Molecule in the MSXa Method 113

6.2 Rotation of Atomic Coordinates to Align with
the internuclear Axis 115











Figures Page

D.1 Transformation of Coordinates 140









Abstract of Dissertation Presented to the Graduate Council
of the University of Florida in Partial Fulfillment of the
Requirements for the Degree of Doctor of Philosophy

HELLMANN-FEYNMAN FORCES
AND DIPOLE MOMENTS USING
MULTIPLE-SCATTERING Xa WAVEFUNCTIONS
FOR DIATOMIC MOLECULES

By

Choy Hing Li

December, 1975

Chairman: John W. D. Connolly
Major Department: Physics

The Hellmann-Feynman forces for H2 and N2, and the

dipole moments for LiH, BH and CH are calculated with the

Multiple-Scattering Xa wavefunctions. First, the muffin-tin

charge density approximation is used to compute these two

quantities. The results show that the internuclear forces

for homonuclear diatomic molecules are always repulsive, and

the dipole moments for LiH, BH and CH are far from the correct

values. Thus, it is necessary to use the non-muffin-tin charge

density in such molecular property calculations. The main

difficulty in such non-muffin-tin calculations is to evaluate

some three-dimensional integrals over the intersphere region.

To achieve this, we express the integrands in Prolate Spher-

oidal Coordinates and reduce the integrals into some one-di-

mensional integrals. The results show that both the force cal-

culations and the dipole moment calculations improve a lot

over the use of the muffin-tin charge density approximation.

For H2, the Hellmann-Feynman force curve is obtained and com-

pared with the curve for the derivative of the total energy,


viii










thus checking whether the Hellmann-Feynman Theorem is valid

for the MSXa wavefunctions. For N2, the force components

contributed by each of the electronic states are calculated

at the experimental equilibrium separation and are compared

with the Hartree-Fock results. For LiH, BH and CH, the dipole

moments thus obtained are respectively 7.29 D, -1.24 D and

-1.958 D at the experimental equilibrium separations. Also,

the derivative of the dipole moment for LiH at the equili-

brium separation is obtained to be 1.27 D/a.u. We conclude

that the non-muffin-tin correction is important in the force

and dipole moment calculations. Any further improvement has

to come from the correction to the MSXa wavefunctions.













CHAPTER I


INTRODUCTION



In the investigation of the electronic structure

and properties of atoms, molecules and solids, we often

encounter the problem of treating the exchange-correlation

effect properly. The Xa theory suggested by Slater provides

us a good way to handle such a problem (see Slater(1972a)).

Together with the muffin-tin approximation, the Xa theory

can be used to investigate the molecular system in a prac-

tical and accurate way. Such a method, called the Multiple-

Scattering Xa method or the MSXa method, was first suggested

by Slater(1965) and then was implemented by Johnson(1966).

Many of the results obtained with the MSXa method

show that it gives better results in the description of

the one-electron features than do the ab initio SCF-LCAO

methods, and is significantly much faster in calculation

(see Connolly and Johnson(1973)). However, in many cases,

we find that the total energy as a function of molecular

conformation is inaccurate. The energy curves of many diatom-

ic molecules, e.g., show that the systems are unbound. It

was pointed out by Connolly (see Connolly and Johnson(1973))










that such a discrepancy may be mainly due to the muffin-tin

approximation, not the Xa approximation, in the calculation

of the total energy. The current version for calculating

the total energy assumes the muffin-tin approximation for

the charge density. The result may be improved by taking

the charge density, instead of its muffin-tin part, gererated

from the converged wavefunctions in the calculation. This

was verified by Danese on the calculations of carbon and

neon molecules (see Danese(1973)).

If the discrepancy in the total energy calculation

is indeed due to the use of the muffin-tin part of the

charge density and not the Xa approximation, and the result

does improve by using the complete form of the charge densi-

ty, then it will mean that the converged MSXa wavefunctions

are not too far from the true Xa wavefunctions and it may

be useful to employ the MSXa wavefunctions in other molecular

property calculations. In this dissertation, we shall inves-

tigate the quality of these MSXa wavefunctions directly in

the dipole moment and force calculations of some diatomic

molecules.

It was proved by Slater (see Slater(1972b)) that the

Hellmann-Feynman Theorem is satisfied by the exact Xa wave-

functions. That is, the quantum mechanical force acting on

a nuclear coordinate is equal to the electrostatic force

exerted by the other nuclei and the electronic charge cloud

which is formed by the exact Xa wavefunctions. Thus, it is









useful to check the quality of the MSXa wavefunctions by

testing if the theorem is also valid for these approximated

Xa wavefunctions. Also, as it is well known that the dipole

moment is more sensitive to the wavefunction than is the

energy, it is useful to check the quality of the MSXa wave-

functions by doing some dipole moment calculations.

The difficulty in all these calculations is that

we have to evaluate some three-dimensional multi-center

integrals over an awkward region. In this dissertation, we

solve this difficulty by a computational technique which

reduces the integrals into one-dimensional integrals whose

integrands involve sine and cosine integrals. The main

point in this method is to employ the Prolate Spheroidal

Coordinates. It is hopeful that this method can be used in

some other molecular property calculations. We shall dis-

cuss it in detail in Chapter IV and Chapter V.

In Chapter II, we shall give a brief discussion of

the exact Xa method and the approximate MSXa method. In

Chapter III, we shall discuss the effect of assuming the

muffin-tin charge density in the calculation of total energy,

force and dipole moment, thus showing the necessity of using

the non-muffin charge density in the calculations. Chapter

IV will be devoted on the force calculations and Chapter V

will be devoted on the dipole moment calculations. Finally,

we shall conclude in Chapter VI.
















CHAPTER II


EXACT Xa METHOD AND MUFFIN-TIN APPROXIMATE Xc METHOD



2.1 Introduction

The discussions in this chapter shall mainly concern

the basic Xa theory and also its muffin-tin approximate form

in the application to molecular systems. First, we shall

briefly discuss the basic Xa theory and its exact one-

electron Schr6dinger equations, in which no specific restrict-

ion is placed on the form of the charge density and the

potential. Secondly, we shall discuss the MSXa method, in

which the charge density and the potential are muffin-tin

averaged. The approximate Xa one-electron Schr6dinger equa-

tions and their solutions are discussed. The quality of

the solutions of the approximate Xa one-electron equations

is one of our main concerns in this dissertation.



2.2 Exact Xa One-Electron Schrddinger Equations

The Xa method was suggested by Slater and has been

discussed in detail by many others (see Slater(1972a),

Slater and Johnson(1972)).

The Xa one-electron Schrodinger equations for the









spin-orbitals ui(r) can be derived by the variational princi-

ple from an expression of spin up and spin down charge den-

sity:

p() = Z niui()ui(7 )
i*


p() = niui()ui(r)
iP

where ni are the occupation numbers and the summation is over

the spin up and spin down orbitals respectively. Also we de-

fine the total charge density as

p(r) = p (r) + p (r)
Then the Xa total energy is





f 2- p 1 2p(rl)p(r2) 1d
+ ; P (rl)drI + r dr drl
I J + t J 12

+ P(rl)Ux (r)dr + o(r5)Uxarl)dr



+ 2 2ZpZq(2.1)
Rpq

Here the nuclear coordinates are Rp and their atomic

numbers are Zp. Rydbergs are used as units of energy, atomic

units as units of distance. The first term is the kinetic

energy. The second term is the nuclear-electron coulomb

energy. The third term is the electron-electron coulomb

interaction term. The fourth and the fifth are the "exchange-

correlation" energy terms of either spin. The last one is









the nuclear-nuclear interaction term.

The above expression is exact except for the exchange-

correlation term..Here we define it as:

Ut (?) = -9a((3/47)p (r))1/3 2.2)
S =/3(2.2)


and similarly for U (r).
xa
The Xa total energy is a function of the occupation

numbers n. and a functional of the spin-orbitals u.. It can
1 1
be shown that if we vary the ui to minimize the total energy

with fixed occupation numbers ni, we can get the Xa one-

electron Schrodinger equations as:

2
{-V1 + Vc(Prl) + Vxa(P,r )} ui(r) = Ei ui(rl) (2.3)

where N -2Z p(2
V(p,l) = d + 2 r (2.4)
P'l r1- 12


is the electrostatic potential at position r1 and

1/3
F -6c((3/47)p,(rl))/ for spin up ui

Vxa (P,rl) =
[ 1/3
-6a((3/47)p (rl)) for spin down ui

(2.5)

is the exchange correlation potential at rl.

Since the potential in the Xa one-electron equations

is a functional of the spin-orbitals ui, the Xa method is a

self-consistent method. We start by assuming a potential

Vc+Vx. Solve the differential equations and find the spin-


























































Fig. 2.1 Flow Chart of the Conceptual Exact Xa Method










orbitals. From these spin-orbitals we form the charge densi-

ty. Then, by equations (2.4) and (2.5), we calculate the

potential V,+V,, .which is then compared with the starting

potential Vc+Vx,, If they are "almost" equal to each other,

the calculation is said to have converged. Otherwise we start

with V%+V, and repeat the cycle again until convergence is

achieved. Then we obtain the spin-orbitals ui, and with

equation (2.1) we calculate the total energy . The

conceptual computational scheme is shown in figure (2.1).

There are some difficulties in this scheme, namely:

(1) There is a parameter a in equation (2.2)..In

atomic systems, the a value can be chosen in several ways,

all of which amount to almost the same total energy (see

Lindgren and Schwarz(1972)). One of them is to set the a

value such that the total energy calculated be equal to

the exact Hartree-Fock total energy. Another is to set the a

value such that the Hartree-Fock total energy calculated

with the Xa orbitals is a minimum. Still another is to set

the a value such that the virial theorem is satisfied.

However, for molecular systems, the space may not be as

homogeneous as that of the atomic systems, and we need

some other prescriptions to determine the a value.

(2) In general, the potential is a very structured

function. For atomic systems, the potential is spherically

symmetric and it is possible to separate the solutions in

spherical coordinates. However, for molecular systems, the










differential equations usually cannot be separated in any

coordinates, unless some approximations are made in the form

of the potential.

(3) Suppose the differential equations are solved

and the spin-orbitals are obtained, we still have to calcu-

late the potential by equations (2.4) and (2.5). and to

evaluate the total energy by equation (2.1). All of these

involve three- and six- dimensional integrals with compli-

cated integrands.

Due to all these difficulties, the MSXa method is

suggested. With the muffin-tin approximation and the multi-

ple scattering formalism, the differential equations can

be solved easily and the expressions for the total energy

and the potential become much more simple. We shall discuss

these in the next section.



2.3 Muffin-Tin Approximate Xa One-Electron Equations
and the Solutions

The Multiple-Scattering Xa method (MSXa) was suggest-

ed by Slater in 1965 and its implementation was begun by

Johnson in 1966. There are three basic characteristics in

the MSXa method: the exchange-correlation approximation, the

muffin-tin approximation for the charge density and the poten-

tial, and the multiple scattering formalism.

First, we divide the coordinate space of a molecule

into three regions: (I) the contiguous spherical regions

surrounding the different nuclei, (II) an intersphere region









which is outside of the spheres of region(I) but inside

region (III), which is the exterior of a sphere enclosing

all the spheres in region (I). The muffin-tin approximation

is such that we approximate the potential and the charge

density in region (I) and region (III) by a spherical aver-

age, and those in region (II) by a volume average, namely:

(r) = +ff(t)dQ if r, I or III

f(S) =

-= f(r)dr if rE II
VII II
where f is the potential or the charge density function,

V is the volume of the intersphere region, and f(r)

is the averaged quantity.

With the muffin-tin approximation, it can be shown

that the difficulties mentioned in the previous section

can be solved, namely: the differential equations for the

spin-orbitals are separable, the potential can be computed

easily, and the calculation of the total energy reduces to

the evaluation of one-dimensional integrals.

Starting with the muffin-tin charge density p(r),

we have an expression for the approximate Xa total energy:

r- C i%* 2 ()d +U + UT
= 2 ni ui(r)-V )u(r)dr + UN +


where ui are the spin-orbitals that minimize

UNN is the nuclear-nuclear interaction energy,

UT is the sum of the nuclear-electron coulomb energy,

electron-electron coulomb energy and the exchange-correlation









energy:


UT4 -2Zp ) (rl)drl
pzI I 1i-pI

+ CaJ (+(r )4/3 + ( /)4/3d

1 r r + 2( (
+ 2drlp(rl)Jdd 212
2 J r12

9 1/3
Ca = -9a(3/4T)1/3

With this expression for the total energy, it is

shown in Appendix A that the necessary and sufficient

condition for be a minimum is that the hi satisfy

the following differential equations:

2 + r e n, %%
-V + V(p,r) + Vex() ui(r) = iui(r) (2.6)

where V (p, ) and Vex(p,$) are the muffin-tin average of

V (p,r) and Vex(p,r) as defined in equations (2.4) and (2.5).
We want to emphasize that the muffin-tin approximate

one-electron equations (2.6) are derived from the fact that

the charge density is of muffin-tin form. There is only

one approximation (besides the Xa exchange-correlation approx-

imation) involved in this derivation, and that is the muffin-

tin charge density approximation.
To solve the muffin-tin approximate Xa one-electron

equations:


{ -2 + Vxa(p,)} u (r) = Eiui()









where we write V (p,r)+Ve () as V (p,r), we may use the

partial wave expansions of multiple scattering theory in

the different regions of the molecule. Within an atomic

sphere S (including the outer sphere which we call the 0th

sphere), the potential is spherically symmetric. Thus, we

may expand ui for state i into a linear combination of

products of a radial function and a spherical harmonic

(here we use the real spherical harmonics for computational

convenience):


"U m
ur)= CmR(rr, s)Ymr)

where C are expansion coefficients to be determined, and

the radial functions Rj(rg,Ei) have to satisfy the follow-

ing differential equations:

1 d 2d (z+l) -
r dr (ridr V a(p,r)} R (rE i)
rB dr5 B dr5 3S B x

= EiR(rg,Ei)


In the intersphere region, the potential is a

constant, VII, and the differential equation for ui in this

region is:

2 + (V i)} u(r) = 0


which is just the equation of a free electron. Hence we may

express ui as a multicenter partial wave expansion:










r Z (ir)Y(r8) + A0 j(ir0)Y (ro)

ui() = if i >VII
Ai k (Kr )Y(r ) + A0 i ( rOYM r
paI Y, m ^U (ki-ro 0 z ro
if ci

where 1i (IViI- i) '

n is the ordinary Neumann function,
j. is the ordinary Bessel function,

k. is the modified spherical Hankel function of the

first kind,

i. is the modified spherical Bessel function,

and Aim are the expansion coefficients to be determined.

With the requirements that the spin-orbitals ui and

their derivatives have to be continuous on the boundaries of

the spheres, we obtain a secular equation from which the

eigenvalues Ei and the expansion coefficients { CZ' I and
i8
{ A m can be determined. Appendix B gives the details

of the procedure to determine the eigenvalues and the expan-

sion coefficients.


2.4 Relation Between the Exact Xa Method
and the Muffin-Tin Approximate Xa Method

In this section, we shall discuss the relation

between the exact Xa spin-orbitals ui and the spin-orbitals

ui derived from the muffin-tin approximate one-electron

equation (henceforth called MSXa spin-orbitals). First of
























Fig. 2,2 Flow Chart of the Conceptual Exact
Xa Method and the M6Xa Method.

------ Exact Xa
0 MSXc
----> Non-Muffin-Tin Calculation

















olve exac
Xc eq. to
obtain u --

I

Calc.
er Form p
Fp ),p converg..






'/


Generate

Vxa(






no










all, we shall look at the differences between the two self-

consistent cycles for the conceptual exact Xa method and

the muffin-tin approximate Xa method (figure 2.2).

Both SCF cycles start with the superposition of

atomic potentials. In the exact Xa method we try to solve

the Schrodinger equations (2.3) to get ui directly from this

potential whereas in the MSXa method we solve the Schrodinger

equations (2.6) to get ui after the potential being muffin-

tin averaged. Then the charge density p is formed from

either ui or ui. In the exact Xa method the potential is

formed from p while in the MSXa method the potential is

formed from the muffin-tin part of the charge density,

namely T. Then we repeat the cycle again until the conver-

gence is achieved. Then the total energy is computed with

p or p as the charge density respectively in the exact Xa

method and the MSXa method.

As we have mentioned in the previous sections, the

basic difference between the exact Xa spin-orbitals ui and

the MSXa spin-orbitals u. is that they minimize two differ-

ent total energy functionals, and . In

other words, we have two different hamiltonians: H, the

hamiltonian for the exact Xa method, and H, the hamiltonian

for the MSXa method. The ui and ui are the eigenfunctions

of H and H respectively.

Since we had employed the approximate hamiltonian

in our calculations, an important question should arise:









how is the quality of the MSXa spin-orbtials ui, or how

good are they compared to the exact Xa spin-orbitals ui?

The calculations of the non-muffin-tin total energy by using

the MSXa spin-orbitals ui done by Danese (see Danese(1973))

have given part of the answer to this question (we shall

give a more detailed discussion about this calculation in the

next chapter). In this dissertation, we want to look at the

quality of ui in two more aspects: the Hellmann-Feynman

force and the dipole moment calculations.

It was proved by Slater that the exact Xa spin-

orbitals (with the restriction that the a parameter be a

constant) satisfy the Hellmann-Feynman Theorem. That is,

the electrostatic force calculated with ui should be equal

to the derivative of the quantum mechanical expectation

value of the hamiltonian H with respect to the coordinate of

the nucleus upon which the force is exerted. Now, instead

of ui, we have ui; and instead of the exact Xa total energy

, we have the non-muffin-tin total energy .

We can check if the ui satisfy the Hellmann-Feynman Theorem

and test their quality.

Also, dipole moment calculations are performed on

some diatomic molecules by using ui as the wavefunctions. As

we know, the dipole moment is much more sensitive to the

wavefunction than the energy is. Therefore, the dipole moment

calculations should be a good check on the quality of the

spin-orbitals too.










Since we do have some reasonable results for the

total energy by assuming the muffin-tin charge density p, we

wonder what would be the result if we assume the same kind

of approximation in the force and dipole moment calculations.

So, before we discuss the force and dipole moment calcula-

tions by using the MSXa spin-orbitals ui, we shall first

investigate these two kinds of calculations by assuming the

muffin-tin charge c(nsity p in the next chapter. We shall

see that the results are not satisfactory and hence it is

necessary to use the p formed from ui, instead of the crude

approximation p, in such calculations. The distinction be-

tween these two different approaches is also shown in figure

(2.2).














CHAPTER III


MUFFIN-TIN CHARGE DENSITY APPROXIMATION



3.1 Introduction

In the last chapter, we discussed the exact Xa

one-electron equation. Instead of solving this exact equa-

tion, we solve the muffin-tin approximate Xa one-electron

equation which is obtained by varying the Xa total energy

functional by assuming a muffin-tin charge density approxi-

mation. Since we had assumed this approximation, we would

like to see how good it is. One way to check the validity

of the muffin-tin charge density approximation is to look

at some of the molecular property calculations using the

muffin-tin averaged charge density. In this chapter, we

shall discuss three such calculations: the total energy

calculations, the electrostatic force calculations and

the dipole moment calculations. With the results obtained,

we conclude that in most cases it is not satisfactory to

use the muffin-tin charge density in the molecular calcula-

tions. Also, we shall discuss the non-muffin-tin total

energy calculations done by Danese. With the success of

his calculations, it seems that we should use the charge

density generated by the spin-orbitals ui, rather than the.










muffin-tin charge density, in the molecular property calcu-

lations.



3.2 Muffin-Tin Xa Total Energy Calculations

As we have shown in the last chapter, we calculate

the Xa total energy by using the muffin-tin charge density.

In solid state systems, especially the close-packed system,

the use of the muffin-tin approximation is usually quite

successful. However, for molecular systems, it is not uni-

formly good. Although the absolute values of the total ener-

gy are quite close to the experimental results, we have some

bad results when we consider the dissociation energy and

the conformation of the systems. The calculation by Connolly

and Sabin(see Connolly and Sabin(1972)) on the water mole-

cule led to a minimum total energy for a linear molecule

rather than for the bent form known to be correct experimen-

tally. The work in many diatomic molecules indicates no

binding too.

It is noticed that the results are worse in loose-

packed systems than in close-packed systems. It was suggest-

ed that this is because the intersphere region is larger in

the previous case than in the latter case, and the muffin-

tin effect is correspondingly more severe. This was verified

by the work done by Danese and Connolly, Danese performed

the calculations for the non-muffin-tin Xa total energy of

carbon and neon molecule by using the MSXa spin-orbitals









rather than the muffin-tin charge density (see Danese

(1973)). The result is dramatically improved for C2. One

finds binding very much closer to the experimental value

than in Hartree-Fock calculation, and somewhat closer than

the result of the configuration interaction calculation.

For Ne2, one also finds binding.

With the calculations on the total energy, one may

conclude that the approximation of charge density has quite

a large effect on the molecular properties. In the later

sections, we shall show such muffin-tin charge density effect

on the electrostatic force and the dipole moments.



3.3 Force Calculations

In the last section, we have shown that the use of

MSXa spin-orbitals ?.causes a big improvement in the total

energy calculations. In this section, we shall look at ui

in a more elaborate way.

As we know, the error in the total energy calcula-

tion is of second order in the error in the wavefunction.

We would like to look at some molecular properties which

are more sensitive to the wavefunction. The force exerted

on a nucleus by the electronic charge cloud and the other

nuclei is a good test of the quality of the wavefunctions.

The Hellmann-Feynman Theorem, which we shall discuss in

more detail in the next chapter, provides a way to calculate

this force. It was proved (see Slater(1972b)) that this









theorem holds for exact Xa spin-orbitals u. (with the re-
1
striction that the a parameter be a constant). However, it

is well known that for approximate wavefunctions the theorem

often gives rather unreliable results. Thus, it should be

a good test of the MSXa spin-orbitals by employing them in

a calculation of the force exerted on a nucleus of a mole-

cule by means of the Hellmann-Feynman Theorem.

We shall discuss such a calculation in Chapter IV.

Before that, we may wonder what the result would be if we

use the muffin-tin charge density in such a calculation,

just as in the total energy calculation. We shall show in

the following that for homonuclear diatomic molecules, the

force obtained with such a muffin-tin approximation would

be always be repulsive. This cannot be physically possible

because that would mean there is no binding for such a

system.

In the MSXa calculation of a homonuclear diatomic

molecule, the space is partitioned into three regions: two

atomic spheres surrounding the two nuclei A and B, each with

nuclear charge Z, an outer sphere enclosing the two atomic

spheres, and the intersphere region II, as shown in figure

(3.1).

By symmetry, the force on nucleus A is along the z-

axis and is given by:

F = 2Z P() cosOA d3r 2 (3.1)
r A
A AB



















































Fig. 3.1 Partitioning of the Space of a Homonuclear
Diatomic Molecule in the MSXa Method.









where p is the charge density in atomic unit,

eA is the azimuthal angle measured with nucleus A at

the center,

rA is the distance from nucleus A,

RAB is the separation of the two nuclei,

FA is in Rydberg/a.u. units.

Now, if the muffin-tin charge density approximation

is assumed, the quantity p in equation (3.1) will be replaced

by p and the calculation for FA becomes much more simple.

First, we decompose the integral into four different parts,

one over each region:

Sp(r) cosOA 3
r dr = IA + IB + out +II


where I p(-r ) coseA d3
where IA = d2 r
A JA rA

p(r) cosA 3
IB r 2 d3r


I p(r) coseA 3
out d r

Sp(r) coseA 3
III = --- r- d3r
1 i rA

We shall discuss each of the above integrals in the follow-

ing:

Integral Over Sphere A: IA

Since p(r) is spherically symmetric inside sphere A,








we have:

Sp(r) cos9A 3
IA 2= r d3r
A

= T(r) rA drA cos sin O dSA d

where RA is the radius of the sphere A and is equal to one

half of the separation of the bond distance R.

Since the integral over the coordinate 0A is equal

to zero, we have IA = 0.

Integral Over Sphere B: IB

To evaluate the other integrals, we need to expand cosA

in terms of functions centered at nucleus B or at the cen-

ter of the outer sphere. Such an expression can be derived

from an expansion theorem for Neumann functions (see Appen-

ix F). The expression is as follows:
0
S9+1,(RA-R B)
-4t I(t+1,0;1,0;1,0) -- r Yo:rB)
Y RB B Z BAB
Y(r RAB
-- = (if rB < RAB) (3.2)

-i Y (r~+)
4 I(-10;11,0;,0) RAOZ1 Y -O(RA YB +


(if r0 > RA0) (3.3)

where Ym are the real spherical harmonics,

I(. ,m ;Z ,m ;Z ,m ) are the Gaunt coefficients for
1 1 L 2 3 J
real spherical harmonics and are defined by:

Y I(r) Y 2(r) Y (r) d (3.4)
t i Z,2 9. j









and RAB, RAO are the distances between A.and B and-

between A and 0 respectively.

Inside sphere B, we have rB < RAB and we can use

equation (3.2) for the evaluation of IB:


Sp(r) cos9A
I = r2 d3r
B B rA


= p(r) (47/3) Yl(r d3r
J rA


= (4i7/3)' p(rB)(-4Ir)" I(t+1,0;1,0;, 0)
-+ --
Y+1(RA-RB) o 3
rx +2 B Y(rB) d r
RAB

Y 01 A-AB)
=(4/3) (-4~) I(t+l,0;1,0;,0) eR+2
x-- RAB

R +2 O
x' (rB) rB drB Y(rB) dQB


0 r
r (4u)Bg) dQB ko) Q


I = (4Tr/3)2 (-4i) I(1,0;1,0;0,0)
B


Since

and


Yl(RA-RB)
R'AB


x (4T) j P(rB) rB drg


I(1,0;1,0;0,0)=(1/4i)
00 ( (3/4)
Y (RA -RB) = Y(0A = 1) = (3/4r)2cos =(3/4a)









1 2 2
therefore I P(r)4r2 dr = Q /R
B R2 B B B B AB
AB

where QB = IP(rB)4T rB drB is the electronic charge
inside sphere B.
Integral Over the Outer Sphere: Iout

In the outer sphere, ro > RAO, so we can use equation
cos0A
(3.3) for the expansion of "- Then we have:
rA

S cosOA 3
Iout = Jo p() d3r
A

= (4T/3)2 P(r )4i I(-1,0;1,0;,0)


-Z1 0 Y-( 0) 3
x RAO Yk-1 RA-R0) t+ d r
r0

= 4n(4/3) I(-l,0;l,0;,0) RA0 Y RA 0


p(r0) 2 0
O 1 r0 dr0 Y(ro) d$0

0 1 )
Since Y(r ) d0 = (4 ) 60

and k j 0 in the summation, therefore Iout = O
Integral Over the Intersphere Region: III

In the intersphere region, the muffin-tin charge

density is a constant, so we can take p out of the integral:










If p( ) cosOA 3 -
S= d r = p
I I r I


We further decompose the integral into three parts:


ScosA 3
J-- d r =
rI A


cosA 3
--- d3r
rA


cose d3
6 A


Just as for IA and IB, we-have:


ScoseA
A -A d3r = 0
A rA


cosA
d3r =


and 1


2
( 4r2 dr )
( B B


2
/RAB = 1TRAB/6


For the first term, we have:


ScosA 3
-7-- d3r
AtB*tII rA

= (4x/3) Y(rA) d3r +
SA


(4T/3)1 {(-4T)
Azoo


JSr+2dr2


0 '
-rY(--A) d3r
, rat R A


y 0(AA-O
I(Z+l,0;l,;Z,)ZA-
KO+2
AO


YO(rO) .dQ


Z-1 0- -
+4,r I(z-l,;lO;Z ,) RAO YZ1l(RA-RO)
P1


cosA d3
rAdr
rA


cos0A 3
---- d3r
rA







29


.x 2 -- -- dr0 Y0(r0)d2O}
AO r0

= 4rRA0/3 = 2TRAB/3


Putting all of them together:


III= =II ( 2nRAB/3 TRAB/6 ) = p 7RAB/2

We can express II in terms of the total charge QI1

in the intersphere region and the volume VII;

3 3 3 3
Vi = (47/3) (RAB RA RB) = TRA


QII QII
Therefore I = V II RAB/2 = 2
II 2 RAB

Combining IA,IB,Iout' II and the nucleus force, we have:


A cosOA 3 2 2
FA = 2Z (r) -- dr 2Z2/RAB


= 2Z (IA + I + It + I ) 2/RB


_B QII 2 2
=2Z (0 +- + 0 2 ) 2Z /RAB
AB AB

2Z
= 2 (% + QII/2 z)
RAB

Since for a homonuclear diatomic molecule, we have:


QII = total 2 QB out









where Qout is the charge in the outer sphere.
out
Also, for a neutral homonuclear diatomic molecule,

Qtotal=2Z. Thus, we finally have:

2Z
F= 2- (Q + Qtota/2 Q Q /2 -Z)
A R B total B out


2Z
AB
= AS <-ut2 + Z Z)

Z out (3.5)

AB

Since Qout is always positive, we have the result

that the force exerted on the nucleus A is always negative,

or away from the other nucleus. Thus we prove that the force

is always repulsive for a homonuclear diatomic molecule if

we assume the muffin-tin charge density.

Due to such a result, we see that it is inadequate

to employ the muffin-tin charge density approximation in the

force calculations. As in the total energy calculations, we

hope that improvements may be achieved if we use the charge

density p generated by the MSXa spin-orbitals, instead of

the muffin-tin charge density. Indeed, we have some improve-

ments as we shall discuss in Chapter IV.



3.4 Dipole Moment Calculations

Another way to check the quality of the wavefunctions

is to look at the dipole moment calculations. As is well

known, the dipole moment is very sensitive to the wavefunc-










tion employed in the calculations. We shall use the MSXa

spin-orbitals for such an investigation. Before that, we

first derive an expression for the dipole moment calculated

from the muffin-tin charge density. We shall see that such a

charge density approximation is also inadequate, and we need

a more elaborate way to evaluate the dipole moment.

Consider a neutral diatomic molecule AB. We partition

the space as we do in the force calculation. Since we have a

heteronuclear diatomic molecule, sphere A and sphere B are

in general of different size.

The dipole moment is defined as:


= p(r)r d r + vN


where P(r) is the muffin-tin charge density,

UN is the nuclear contribution and is defined by:

UN = ZARA + ZBRB.
For a neutral molecule, it can be shown (Appendix D)

that the dipole moment is independent of the choice of the

origin. Thus, for computational convenience, we can choose

the origin at the center of the outer sphere. Furthermore, by

the axial symmetry of the system, we note that the only non-

vanishing component is the z-component. Thus,



S= j -(r)z0 d r + N






32

To evaluate the integral, we decompose it into four
parts again:


S(r)zO d3r


= p(r)z0d r + j (r)Zq d3r + () dr


+ J p(r)zo d3r


We shall consider them one by one:





33
(1) J' p(r)z0 dr


= J p(r) (rAcos0A RAO) d3r


p= (rA)rA drA 1 cosGA dA RAO I (r) d r


= 0 RAQA = RAQA


(2) p (r)zo d3r


= ] P(r) (rBgosOB + RB0) d3r

e 3
S (r )r drB cosoB dQB + RBO p(r) d3r


= 0 + RBOQB = RBOQB









(3) Iot P(?)z0 d3r


= lt o(r0)r0cosO0 d3r


= r p(r0)r dr0J cos00 d00 = 0


(4) J (r)z0 d3r


= PII 11Z d3r


= ( z0 d3r z0 d3r )
0 A ZO d ZO d


= { r" 3 dr0 Jcs0 d0 I (rAcsOA RAO) d3r


I (r cos + R ) d3r
J B B BO


= II ( 0 + RAOVA RBOVB )


where VA and VB are the volume of sphere A and sphere B
respectively.
Combining all the terms, we have:


I p(r)z0 d3r = RAQA + RBOQB + II(VARAO VBRBO


Since RAO = (RA + RB) RA = RB









RBO = (RA + RB) RB = RA


VA = 4r Ri/3


3
VB = 4r RB/3


PI1 = QII/VII = QII/ 4ARARB(RA+RB)


Therefore p(r)z0 d3r


= RBQA +


Also,


QII 4R 3 3RA
RAQB + 4RRARB(RA+RB) W- (ARB RBRA)


= RBQA + RAQB + QII(RA RB)/3



PN = ZBRBO ZARAO = ZgRA -ZAB


Thus, finally we have:


= {- RBQA + RAQB + QII(RA-RB)/3} + ZBRA ZARB


= RB(QA + QII/3 ZA) RA(QB + QII/3 ZB)
(3.6)

We shall discuss the results computed with the above

equation in Chapter V. There we shall find that most of the

results are not satisfactory when compared with the experi-
mental values or with the results obtained by other methods.







35


The reason is simple: the expression (3.5) depends entirely

on the crude approximation of muffin-tin charge density. We

hope that improvements may be achieved by using the charge

density generated by the MSXa spin-orbitals instead of the

crude muffin-tin density. We shall devote Chapter V to such

calculations.














CHAPTER IV


HELLMANN-FEYNMAN FORCE CALCULATIONS



4.1 Introduction

In the last chapter, we saw the necessity of using

the non-muffin-tin charge density generated from the MSXa

spin-orbitals ui in the force and dipole moment calculations.

To carry out these calculations, we need to evaluate a three-

dimensional integral whose integrand is a multi-center func-

tion. The basic difficulty of such calculations is the eval-

uation of the integral over the intersphere region, which is

in general very awkward. In this chapter, we shall develop

an analytical way to perform such integrals for the case of

diatomic molecules, and we shall apply it to the force calcu-

lations.



4.2 Prolate Spheroidal Coordinates

Before we go into the force and dipole moment calcu-

lations, we first introduce the Prolate Spheroidal Coordinates

that are essential to our scheme of evaluating the integral

over the intersphere region.

Let A and B, with separation of R, be the two foci

of our coordinates. Also, let (rA, OA, iA) and (rB, 6B, (B)



















































Fig. 4.1 Representation of a Point in the Atomic Coordinates.









be the two spherical coordinates with centers located at A

and B respectively. OA and OB are both measured from the

major axis in the same direction, as shown in figure (4.1).

The Prolate Spheroidal Coordinates are d: fined .as:



5 = (rA + rB)/R


n = (rA + rB)/R


(4.1)


A = CA = dB
SA B


The reciprocal relations are:


rA = R( + n)/2


cosOA = (1 + +n)/( + n)


{A = 2 )(1 -n2)}
+ n


rB = R( n)/2


and cosO = (1 5n)/(5 n)


= (2- 1)(1 -n2)
sin(B
(4.2)


Also, it can be shown that the volume element is:


d3r = (R/2)3(2 n2) d dn d4


(4.3)


If the integral is over the entire space, the inte-

gration limits would be: 1 $ ( < <, -1 I q $1, and 0 2ir .

However, if we are evaluating an integral over a finite space









only, such as the integral over the intersphere region, the

integration limits have to be modified. We shall use the two

atomic centers or one atomic center together with the center

of the outer sphere as the two foci, depending on the form

of the integrands, in our treatment of the integrals.



4.3 Hellmann-Feynman Theorem

It was proved (see Slater (1972b)) that the exact Xa

spin-orbitals, with the restriction that the a value be a

constant everywhere, rigorously obey the Hellmann-Feynman
9
Theorem. In other words, the force component
aXp
acting on a nuclear coordinate Xp, where is the exact

Xa total energy of the system as calculated quantum-mechani-

cally from the exact Xa spin-orbitals, is the force which

would be calculated by pure classical electrostatics acting

on this nuclear coordinate, resulting from the electronic

charge cloud, formed by the exact Xa spin-orbitals, and all

other nuclei. We can express the theorem in the following

form:


9 ( 2Z ) d3r Z 2ZpZq
-Xp x (r) ( ) d r - ( --T- p )
DX pJx p P @Xp Rp-RqI

(4.4)
th
where Zp is the atomic number of the p nucleus,

is the exact Xa total energy in Rydberg units,

Pxa is the electronic charge density formed by the
exact Xa spin-orbitals,









R R are the nuclear coordinates.

In the above theorem, we want to emphasize that the

spin-orbitals and the total energy are those derived from

the exact Xa method. As we have mentioned in the previous

chapters, the existing MSXa calculations are not exact Xa

calculations. They assume the muffin-tin approximation. With

such an approximation, we can ask the following two questions:

(1) Does the Hellmann-Feynman Theorem hold for the

approximate Xa spin-orbitals ui? If not, then how severe is

the discrepancy? In other words, we want to see how large

is the difference between xa and Fp(pxa), where
Xp
is the Xa total energy obtained with MSXa spin-orbitals ui

and Fp(pxa) is the electrostatic force exerted on the nuclear

coordinate Xp by the other nuclei and the electronic charge

cloud pxa formed from ui.

(2) How well does Fp(pxa) compare with the experimen-

tal value of force acting on the nuclear coordinate Xp?

In order to answer the above questions, we made some

force calculations on homonuclear diatomic molecules. The

reasons that we pick the homonuclear diatomic molecules are

that the a values in these systems are constant throughout

the whole -space which is the restriction of the Hellmann-

Feynman Theorem, and that the computations are easier and

the results can be compared with the experimental values.









4.4 Method of Hellmann-Feynman Force Calculations

To investigate the Hellmann-Feynman Theorem, we have

to calculate two quantities, the total energy and the elec-

trostatic force exerted by the electron cloud as expressed

by the integral on the left hand side of equation (4.4). The

calculation of the total energy can be done by using

Danese's program (see Danese (1973)), as we have discussed

in Section 2 of Chapter III. The main quantity that we wait

to evaluate in this discussion is the electrostatic force

Fe (xC)

2Zp
Fe 'a) = r) ( _--7 ) d3r
3Xp Ir-Rpl

Since p(r) = niui(r)ui(r) where ni is the occupa-

tion number of the ith state, we can express Fe (xa) as:

i p I,2Zp
Fe(xa) = ni uI (i)i() ( TX ) d3r
p r- p1


For a homonuclear diatomic molecule, the expression

for the z-component, which is the only nonvanishing part,

of the force exerted on nucleus A is reduced to:


S1* -t -> cosOA 3
Fe (xa) = 2 ZA ni i(r)ui(r) dr
i rA

To evaluate the above integral, we partition the

region of integration into four parts: one over sphere A,









one over sphere B, one over the outer sphere, and one over

the intersphere region:


S%-* % cosGA 3r
ui() i(r) c d r = F + F + Fot + FII (4.5)


i U COSOAG 3
where F = ()ui() d3r j = A, B, out, II


We shall discuss each of the integrals in the follow-

ing:

Integral Over Sphere A: FA

Inside sphere A, we have


a + iA +k m 9
ui(r) = C R (rA,,i) Y A(rA)

iA
where C.m are the C-coefficients,
km
Ym(rA) are the real spherical harmonics centered at A.
Since we have used real functions as our basis, all
%* t
the C-coefficients are real, and ui = ui. Therefore,


_i f cos6A 3
FA = ui(r)ui(r) -OS- d3r
A rA

= Cmi CiAm ( RyZ(rA)Ym (rA)
Ai f7, JA m 1 i


xR (rA)Y (rA) C }.d3r
r A


S iCA CiA fA R(rA)RZ(rA) r2 drA
R, miPr~rA









SmY r my (rA)COSEA dQA
x Y A 2 A A A


= CA C Am,(4t/3) I((4/2I2,m2;mlO )



IRR (rA)R (rA) drA

where I(L, ,m; 2,m2;i ,m3) are the Gaunt coefficients as

defined in equation (3.4).

Since I(t1,m ;Z2,m2;1,0) = 0 unless the arguments

satisfy the following conditions:


m = m = m
1 2
S+Z +1 is even
IL,- I < 1 < Z + 2


therefore FA = C m C Am (4r/3)t
A M Yi 2 M


x I( ,m;2,m;1,0) R R(rA)R (rA) drA


By the axial symmetry of the diatomic molecule, the

C-coefficients are nonzero only if m=mi for state i (hence-

forth called m). Therefore, we have:


Fi = C iA (4i/3)1
CA 2 , m 2m


x I(,m;Z2,m;l,0) J RZ (rA)R 2(rA)drA

(4.6)








i
Integral Over Sphere B: FB

Inside sphere B, we have:


r C iB y(rg (B)
ui(r) = C R(rB) rB


So FB C Z C J2 R (r )Y (rB)R (r )Y (r)
B Is m B B 2 B B B


cos0A 3
A

2
By equation (3.2), we can expand cosCA/rA in terms
of functions centered at B:

0+
cosegA YI(rA)
r2 = (43/3) I rA
rA rA

= 47 (47/3)2 I( Z+1,0;1,0;l,0 )

0
S+1 (RA-RB)L Y0
R+2 B B
AB

Therefore F = iBC iB(-4)(4/3)


S I(+l,0;l,0;l,0) Y (~1A-+B)
d+2
AB

x R (rB)RL (rB)rg drB


X Y 1(rB)Y (rB )y(rB) dQB
j iB BB











F = C iB CiB (-4T)(4iT/3) {I(+1,0;1,0;Z,0)
B m 2M (-4(4/


x I(e,m;2,m;,0) Y+1(RA-RB)
AB

x R(rB)R2(rB)rB+2drB1 (4.7)
*/


Integral Over the Outer Sphere: Fout

Inside the outer sphere, we have:


ui (r) = C R (r0) Y(r


Thus Fut = C 00 {R (rO)Yi(r)R (r)Y (r)
out 1m 92m LA t xOOkOz


cosA 3
x 2-- d3r
rA

By equation (3.3), we expand cosOA/r around the

center of the outer sphere 0. After integrating over QA'

we have:


Fi = 4(4w/3) C c0i 0 {I(,-1,0;1,;0',o0)
out m '2m
A '.i, 1atPt"


x I(,,m;,m;,) RAO YI-Ar



R r- dg (4.8
Ro r0) ro








i
Integral Over the Intersphere Region: FII
In the intersphere region, ui can be expressed as:

iA m + iBA
iA n (KrA)Y (rA) + A r n (KrB rm
2i1 Az1 A I 2 2 2

iO m .

mui() =

AA km (CrA)Y1 (rA) + TAiB k (rB)Yk (rB)
21 1 2 2 2

io )m
*+ A0 i (Kr0)Y (r0) if Ei < V

For diatomic molecules, the eigenvalues Ei for most
states (except the core.state) are greater than the constant

potential VII in the intersphere region. We shall look at
the positive energy casd (i>VI ) in more detail (the other
case can be treated in a similar way). For core states, the
electronic charge is concentrated around the atomic centers
and the contribution due to the intersphere region can be
neglected in the core states. Thus, we have the following
i
expression for FII:


F + cosOA 3
= ui(r)ui(r) dr
A
'iA AiA i 1 + AiB A iB i
= A m A2m IAA 1Z2m) + m Am BB ( 2m)


+ iO i A A Om I00( 2,m) + 2 A iA iA B i '
Rim ,m m) 2tAgm A2m IAB(-i,2,m)
+ ,.i t.









iA iO i iB iO i
+ 2 AmA2mlAO (1, ,m) + 2 mA BO( ,,m)
(4.9)
i If m
where IA (' m) = in1 (KrA ) (rA)q 2(

cos0A d3r
A






rA
i m m
IAB ( 9M ,m) = An ( I i

COSOA d3
x --- dr
A
i ((KZm) = Ii (Kr )Y mr ()

AO (i A r ,(A ) 22(r 0 2 0

xcosOA d3r
r2



rA

I 0(,)' 2m) 1 2i'n(rB)Y(rB) (Kr0Y(rO)

c osA 3

x d3r
r2
rA
i m m
I0O(1 2,m) = f 1(Kro)Y1(ro)j2(Kro)Y12(r)


cosOA d3r
x -- dr
rA









To evaluate each of the above integrals, we employ

Prolate Spheroidal Coordinates as introduced in Section 2

at the beginning of this chapter. First of all, we want to

see what the limits of integration are for such coordinates.

For the intersphere region, we have the following

restrictions on E and n:


(1) rA + rB > R


(2) 1 rA rB I < R


(3) rA > RA = R/2


(4) rB > RB = R/2


or 5 = (rA + rB)/R > 1


or I n = (rA rB)/Rl 1


or + 1


S- > 1


(5) From figure (4.2), we have:


2 2 2
r0 = rA + (R/2) 2rA(R/2)cos0A


2 2 2
= (R/2) (+n)2 + (R/2) 2(R/2)(+n)/2)(R/2)(+n)/(+n)


= (R/2)2 (E2 + n2 -l)
2 2 2
= (R/2) ( + n -1)


Since r0 H R, we have 2 + (.5


Putting all of the restrictions together, we have





















Fig. 4.2 Relation of rA, R and eA



TI


n-=/5-iZ


n= -/5-


Fig. 4.3 Boundary of the Intersphere Region in
Prolate Spheroidal Coordinates with A
and B as the Foci.









the transformed intersphere region as defined by: (figure

4.3)


-f(S) $ C $ f(S) and 1 4 C J5


where r 1 if 1 5 S 2

f() =

(5 5 ) if 2 C VS 75


Before discussing the evaluation of the integrals,

we note that all the integrands, except for the factor of

the spherical harmonics, are 0-independent, and, because of

the cylindrical symmetry of the intersphere region, all the

integrals can be reduced into two-dimensional integrals

through the 4-integration, such as:


m 2 3
|n (KrA )Y(rA )n2(KrB) (rB)cOSA/rA d r


[ mm 23
-= rJ (Kr )PA (cosS )nZ (Kr )P2(cosOB)cosOA/rAd r/d4


where P'j are the normalized Legendre polynomials.


To evaluate the integrals, we also introduce the fol-
m
lowing expressions for P ri and j :


m = N in)m m -m-2
P (coseA) NZ(sinSA) ) (cos0A






51


i+1 -1-
n (z) = (-1) z {Q(i+-,z)cos(z+-7)- Q'(a+,z)sin(z+inw)}


jk(z) = z-1 {Qo(+,z)sin(z-iir) + Q'(+,z)cos(z-jit)}


m m (2+1l)(-m)!
where N = (-1) { 2 (i+m)!


im =(_ (2i-2v)!
v ( 2(ti-m-2v)!(i-v)!v!


Q(i+L,z) = (-l)k(2+,2k)(2z)-2k (4.10)
k=o
Q'(Z+1,z) = (-1)k(k+J,2k+1)(2z)-2k-1
Q'(i+jz)= Z-h) -2k-1
Z = (-l)k(t+,2k+l)(2z)
Ico

(+k)!
(R+i,k) = k(k)!
k!(9,-k)!


and [ b] = the integer part of the real number b.



With the above expressions and some algebra, we have
i i i i i
the following results for IAA, IBB' IAB' IAO' IBO, and I00

as required by equation (4.9):



(1) IAA: ni(KrA)Y (rA)n (KrA)Y2(rA)cos5A/rA d r


This is a one-center integral. With some straight-

forward algebra, we have:






52



I ( 1 (l) 2 (R/2) d I 2 d,
AA = N 2- 95 Ak IO S.0


{W1 m zm(-_1)rCm(KR) 1+ 2-2m+1-2v+2r C R1+Z-2m+l-2v+2r
V1 V2 r r

(KR)-s s (1 2 L +Z2-2m+1--2v+2r-s



x {4Sgil2 (z+2+ 2-2v+2r-s+4,0)


(2/KR) gl (I1+Z2-2m-2v+2r-s+3,))}} (4.11)
112


where v = v,+ v2


m = m!
Cr = (m-r)!r!


Here the g1 -functions are defined as:


g1 (2(- ) = 5 7 (-l)U+V
k I P 2 3 0


x A(u,v) (u v1) +n+u+v
n 1 2


where the A-coefficients are defined in terms of the (k+,k)
coefficients through the following relations:


S A(uv)(2z)-n = QU(t +,z)QV(e +2,z)
A1 s 0)






53

and P(uv,1)(,0) are defined in terms of the sine and
1 2
cosine integrals as below:


P(u,v,1)(,U) = Icos(Z1+2)r (6u+v, 6u+v,2

KR(t _cost))
+ sin (Z +2 ) 6 }uv,l ( ) t dt
j u+vo1 f-(P)) t


+ {sini( +R2)T (6u+v,2-6u+v,O

rrR(j ( ~ sin t
+ cos((Z +z2 ) 6 sin dt
1 2 u+v,f) t


+ {cos(I -z )7T (6 -++6 )
1 2) ('u+v,0 u+v,2


+ sin(1 2) 6u+v,1(6vO-6v1)


i m m 23
(2) IBB: 1(KrB )Y1(rB)n2 )Y (rB)COSOA/rA dr
(2) IBB: T1/^B ZI B Z2 B Z2 B A A


This is a two-center integral. The result is as

following:


i = N N (-1)z + Z+1(R/2) d "
BB I 2 ki 2o Ro p' d0

Xikm Pm m +2i-2m-2v+3r ,+,-2m-2v+2r
x{Vw 2 r (-) Cs


x (KR)S 2(1_-2) s +2-2m-2v+2r-s {(1+2) g 2(s+)1


(S/KR) g92 2(s,)~) (4.12)
XZ XZ2









where v = v, + v,
and the g2-functions are defined as:


t (1U+V
g 2(s,) = t Y, E (-1)
US= v0o ,=o

x A(uv) p(uv,2)(s+2n+u+v+l,E)
n 2
i
where the A-coefficients are defined as in IAA, and
P1u v,2) (,C) are defined in terms of the sine and cosine

integrals as below:


p ,2 (e,) = .{cosj(i +z )T (6 6 2)
1i2 1 2 u+v,O u+v,2

rkR(Pflq))
+ sinji(1+,)T6 } os t dt
u+v,1 J -KI()) (2KRC-t) t


+{sinI(il+i2)n (6u+v,2-6u+v,)


+ cosj(ZJ+iZ)76 }-t) sin t dt
u+v,1 ) (2KRE-t)2t


+{cosSi(,-,)T (6u+v,o- u+v,2)


+ sin(i1 2 ) (6 vo-6 )6+v,1

fR() rf))
J--KRE 1 dt
J' (tqI) (2





55

i 2
(3) IAB: (Kr )Y (rA (KrA )Y2(rA)OS0A/rA d3r


This is also a two-center integral. The result is:


I N N (-1) ,2(R/2) f d
AB 1 2 .2 av 0 rV

u+v+ -2v2 CCm 2u C -m-2v .1m m
=(-1) C v V1 V2


2 C2-m-2v2 ) +i2-2v -2v2-2u+v-s-r v+r+s
s (KR)


x{2(1- 2) g3 (. _-2v,-2u+v-r+3,Z-2v2-s,E)



Sg3 ,L( 1-2v i-2u+v-r+2,i,-2v2-s,C)}} (4.13)


where the g3-functions are defined as:


gi (Z,(I ) = (-1)u+v+k,+k2 ( 1+1,2k,+u)
111X 2 W=o Vso *=O arO
((u,v,3)
x ( +,2k +v) P1u ,3) (+2k1+u, +2k +v,0


and P2U v3)(r,s,) are defined in terms of the sine and

cosine integrals as below:


(uv 3) ,



= {cos(KRC+-2I2 )(6 +v,o+6u+v
2 u+v,0 u+v,2








rkr(j+kp)
+ sin(KRE+ 9 2 )6u+v'l(6vO-6vl)} f cos t dt
K(J-fCp) (2K R-t)rts

+{sin(KRC+--1 )(6 u+v, 0+u+v,2

+ cos(KR + -) 6u+i,1(6u0o6ul1 r sin r
z -Z2 sin t
+ cos(KRC+-12- u+,v,( uO-6ul)l (2KR -t)rts dt

+
2 rU()-VU U (2cR62-t) t

+{cos(KRS++ 2
2 ) (6u+v,O u+v,2


+ sin(KRE--2- w)6u+v,l (2KR-t)rts dt


(4)1 i I(r r (KO)008r/r 2 d3r
(4) IA :z n(rA )Y(rA)i( 2 (ro)cosOA/rA dr


Here we have a two-center integral. We want to ex-
m +
press j 2(Kr)Y 2(rO) in terms of functions centered at A so
that we can reduce the integral into sums of one-center inte-
grals. Such an expansion is derived in Appendix E:


j2(KrO)Y2(ro) = (-1) 24 7 i


x I(",m" ;!,m;,ii' ,m' ) (-1)'


x j, (KRA)Y (RA l(KrA)Y (rA) (4.14)


for r0 in the intersphere region
With the above equation, we have:










AO A<. '.



0 [ m m M
X j ,(

coseA 3
X Co dr (4.15)
A

i
Thus, we reduce I to an infinite sum of one-center
AO
integrals. Since j decreases rapidly with increasing 2',

we can expect that the sum would converge rapidly.

To evaluate the integral on the RHS of equation

(4.15), we employ the same technique as before and obtain

the following expression:


mI m 2 3
J (Kr)Y1 A) (K A) 2 (rA)coseA/rA d r
N i(KrA)YtN( r A) d({A)YZ(rAIsXA/


r + 1 L v r-
Z ,P 2t 1\ Yt Y' =0D S=O


x(-1)r m 2m Cm C c+2-2m-2v+2r+l (KR) 2 )
-) 2 r S (KR) (1-i )
V 1 V 2 r S

XE1+z,-2m+1-2v+2r-s 2 4E g (s+3,) 4 (s+2,)}
Sz KR (s+2,)}}
(4.16)

where v = v + v2

and the g4-functions are defined as:



SI (-)uA (u, v) (u v,4P )(+2n+u+v,,)
Z 1 P 2 XAL 2






58

where P(u v,4)(,E) are defined in terms of the sine and
e i 2
cosine integrals as below:


S(,,4)(,) = {sini(L1-Z2) (6 +6U )
1i2 u+v,O u+v,2
(KR(fwfP)
+ cos(Z-,2) 6u+v,l(6 -u-ul )} t-dt


+{cosi(1-2)> ( ++v, ++v)



u+v,l ul uO KR[-J(p) t


+{sin(1+ 2) (6U+v,2-6u+v,O)


+ cosi(+a). 6+u2 ) dt

i ( m m 23
(5) BO: (Kr )Y (rB) (Kr )Y2 (r)cosO A/rA d r
(5) .IBOB:2 0 d0rA

m
Expanding j 2(Kr )Y2(r0) in terms of functions cen-

tered at at B, we have,similar to IA0, the following equation:


I1O = (-l)24t 4 (-1)k I(",m;,m;',0)

m + r m -+ .+
x j,(KRB)Y,(RB) n (KTr)Y (r )j,(KrB)Y,,(rB)


x cosOA/r d3r






59

The integral on the RHS of the above equation can be
evaluated as below:


2 3
11 (rB )Y (rB). ,(KrB )"(rB)cosEA/r d3r


= Nm (- ) +1(R/2) d( L
k1 2 I *, =0 1ro :O

+Z-2m-2v+3r 1m ,m m C1+-2m-2v+2r 2(KR)s
x(-l) w Cr s
V1 V2 r s

x.(l-2)s ,+12-2m-2v+2r-s {(i+2 )gS (s+l,) g (s,O)}
12 cKR Z2


where v = v + v
1 2
and the gs-functions are defined as:


gkt (S,) = t Y_ (-I) An (s+2n+u+v+l,C)
g (s, = Az ((unu (uvu) (uv,5) u
1L2 a2l V. ". n

where the P (u,,5)(t,) are defined in terms of sine and
cosine integrals as:


P v,5)( .,) = sin(1,-2)l (6+v,O+6u+v,2)


+ cos(-Z2)r 6 (6-6()}
+ cosi2( ,-z)T 6u+v,1(6uO-6ul


x cos t dt
(2KR(-t) t

+{cosj(1,-z2)T (6u+v,O+6 u+v,2)










u+v,l ul u0
+ sin(t,-2)n 6)u+v,(6ul


x I -sin t dt
.I()5-jt)2KRS-t) 2t

+{sin( 1+2 )T (6u+v,2-6u+v,0)


+ cos--k +k )TT 6dt
+ cosi(+2) u 6u+,1} v(2KR -t)2t


(6) Ii J (KrO)Y (r)Z(KrO)Y (Or)cosOA/r d3r


Expanding j (Kro)Y (r0)about A, we have:

i 0 0 Z ,Z2,m
I100 = E js(

where C ,2 ,m
ZI k3


= (-l) ,+ (47)2


{ (-l) I(P2I",(,m;,,m;R',0)

(m)
x I(,,m;k; ,m;R 3,0) I 1,
Z" I Z.4


and I(m
Z", 4


( i m )y 3
,(rA)Y (rA ) (KrA ) ( )COSA/rA d3r
4 4(


= j I(t",m;4,m;3,0)


0 2 3
3 (rA)cOSA/rA d r


, J"(CKrA)Ji (rA)


n-It
('+-(
V-Ii'-A


Sa~t-Ia


4b O~,4h~2rh









Thus, I00 is reduced to a sum of one-center integrals

which, by employing the Prolate Spheroidal Coordinates,can

be further reduced as below:


0 2 3
j1(KrA )J (
2 L 30 -2v+l r
= (2a3+1)2 .i (R/8) d{ L 30 Cr2v (KR)


x(1 2 r 3-2v+l-r 2 KRf j(t)j0( dt
x(l- 5 ) 5 {26 dt
(<-f))) r+1

-(I/KR) j j (Lit) dt}}
fR(cO-f1))) tr

To conclude this section, we note that:

(i) We have reduced the three-dimensional multi-center inte-

grals into a sum of one-dimensional integrals whose inte-

grands involve the sine and cosine integrals. The sine and

cosine integrals can be evaluated quickly and accurately and

can be called directly from the IBM Scientific Subroutines

Package. We shall discuss such evaluations in Appendix G.

(ii) It looks as though the expressions for I I I,
AA BB' AB'
i i i
I, I and I are very formidable, yet in the actual cal-
AO BO 00
culations they are much more simple. This is because the A-val-

ues used in the basis of the MSXa spin-orbitals are usually

small, and in most cases, do not exceed 3. Thus, for example,

if we use s and p atomic orbitals (. = 0,1) for the a state

(m = 0) in the MSXa calculation, then by equation (4.11),






62


there is only one term in the summation of v,, v2and r since

(i-m)/2 =0 .

(iii) By partitioning the region of integration into differ-

ent parts, we can find the contributions to the total force

by:the electronic charge in each of the regions and for each

molecular state:


Total= Ni


= n(F + F + F + F ut) + FN


where i represents the molecular state i,

FA is the force exerted by the electronic charge
th i i
inside sphere A for the i state. Similarly for F F
B' II
and F
out'
FN is the force exerted by the other nucleus.



4.5 Results and Discussion

The Hellmann-Feynman force calculations have been

carried out for two different systems: the hydrogen molecule

and the nitrogen molecule. The hydrogen molecule is investi-

gated because it is the simplest diatomic molecule and accu-

rate results had been previously obtained (see Kolos and

Wolniewicz (1964)) so that we can make some comparison. The

nitrogen molecule is investigated because its electronic

configuration for the ground state-is a closed shell config-

uration in which case the MSXa method should work better.

Also, it involves the t-states so that we can see how our










method works for a more complicated case.

Hydrogen molecule: H2

(A) The ground state for the hydrogen molecule H2

is z+ and the electronic configuration is log. Forces for

three different bond distances are computed in order to

plot a force curve and to obtain the equilibrium separation

Re and the force constant ke. First, for each bond distance,

a non-spin-polarized MSXa calculation is carried out by

using s and p orbitals in the atomic regions and s and d

.orbitals in the outer sphere. Then, the electrostatic force

is calculated, according to the method described in the pre-

vious section, by using the converged wavefunctions. The re-

sults are listed in table (4.1), and the curve is plotted in

figure (4.4).

(B) In order to check the validity of the Hellmann-

Feynman Theorem of the MSXa spin-orbitals, the non-muffin-tin

total energy of H2 is calculated. From the values
a
of , we obtain the values of -DR at the three

different separations by least squares fit. The results are

listed in table (4.2). Two such curves are plotted, one for

the wavefunction with basis: i=0,1 in the atomicsphere, i=0,2

in the outer sphere; another for the wavefunction with basis:

P=0 in the atomic sphere, Q=0,2 in the outer sphere. Also,
DE
to compare with the experimental result, a curve with -

against R where E is the total energy calculated by Kolos

and Wolniewicz (see Kolos and Wolniewicz (1964)) is plotted.

Such an H2 calculation is considered to be more accurate














TABLE 4.1


1+
Force calculations For H2 ( g



2FAg 2Fag 2Fg 2FCg
A B out II


2
: lg )
g


F Ftota
N total


R=1.4 0.264 0.271 -0.026 0.442 -1.020 -0.068

R=1.44 0.262 0.266 -0.U26 0.417 -U.965 -u.046

R=1.5 0.257 0.258 -0.025 0.403 -0.889 0.004


Basis: L=0,1 in the atomic spheres
=0,2 in the outer sphere.
a=0.77627
R in atomic unit; F in Rydberg/a.u.













TABLE 4.2


3 Comparison of Fto and R for H2
total aR


(a) (b) '(c)
Sa 3
total R dR


(d)
D
3R


R=1.4 -0.068 -0.024 -0.035 -0.004

R=1.44 -0.046 0.004 -0.012 0.032

R=1.5 0.004 0.046 0.024 0.068


Re 1.495 1.43 1.46 1.401

ke 0.72 U.71 0.60 0.73


(a) Hellmann-Feynman Force with basis: =0,1 in the atomic
region; L=0,2 in the outer sphere.
(b) Basis: =0,1 in the atomic region, =0,2 in the outer
sphere.
(c) Basis: =0 in the atomic region; =0,2 in the outer
sphere.
(d) Calculated with the Kolos-Wolniewicz wavefuncqion
Force in Ry/a.u.; R in atomic units;ke in Hy/a.u.
























Fig. 4.4 Force Curves for H2

1. Basis: =0,I in the H sphere
Z=u,2 In the outer sphere
2. Basis: i=0 in the H sphere
X=0,z in the outer sphere.











0.06 X Experimental
Q Exa() 1
aR
SDExa() 2

aR
0.04 Fx (x )





U.02





0





-.0 __





-0.04





-O.U6








-0.08


-0.10


1.5 separation(a.u.)


1.3 1.4













TABLE 4.3


i-D'ependence of Hellmann-Feynman Forces for H2 at R=1.4 a.u.


2FAg 2FBg zF 2F (
A B out II


FN Ftotal


gl 0.265 0.271 -0.026 0.442 -1.020 -0.068

;2 0.268 0.274 -0.026 0.442 -1.020 -0.062


#1: =0,1 in the atomic region; a=0,2 in the outer sphere
#2: =0,1,2 in the atomic region; =0,2 in the outer sphere.
Force in Ry/a/u.









than that experimentally observed. The equilibrium separa-

tions Re and the force constants ke for each of the curves

are obtained. The results are listed in table (4.2) and the

curves are plotted in figure (4.4).

(C) To check the sensitivity of the force with re-

spect to the basis functions, two different calculations, one

with s, p orbitals in the atomic region while the other with

with s, p, d orbitals (both with s, d functions in the outer

sphere), are carried at R=1.4 a.u. The results are in table

(4.3).

Nitrogen Molecule: N2

Bader,Henneker and Cade-have done some force calcula-

tions by using Hartree-Fock wavefunctions as the basis (see

Bader, Henneker and Cade (1967)). In order to compare these

two methods, a force calculation for nitrogen molecule is

also done by using the MSXa wavefunctions.

Both the MSXa and the Hartree-Fock calculations are

carried out at experimental equilibrium separation Re=2.07
1+
a.u. The ground state for N2 is g and the electronic con-

figuration is l1 la 2a 2 21 3a The spin- polarized a value
g u g u u g
for nitrogen, c=0.74522 is used, and the values used for

the different molecular orbitals are listed as in table (4.4).

Here we treat the lag and lau as core states. The results of

the force calculations are as in table(4.5).

Discussion

When we compare Fxa(p), Fx (p) and the experimental

values for the force, we find that we have achieved a big














TABLE 4.4



Z-Vaiues in the MSXa Calculation for N2 ( Zg)


0,1,2

0,2


0,1,2

1,3


0,1,2

0,2


* o0 and lou are both treated as core states.


atomic

outer














TABLE 4.5



Hellmann-Feynman Forces for N2 at R=2 a.u.
Calculated with MSXa and Hartree-Fock Wavefunctions


2 2 2 2 4 2
Io la 2 2 0o ir 30
g u g u u g

fA 0. 1.513 -0.501 1.278 0.462

fB 7. 2.227 0.9b2 2.317 1.685

f 0. U.03 -0.844 0.160 -0.788
out
fI 0. 3.316 -U.0u3 3.274 0.452


fxa 7. 7.059 -0.J95 7.029 1.811

fHar-Fock 7.858 9.387 -1.621 8.512 0.525


Ftotal(Pxa) = -1.996

Ftotal(Hartree-Fock) = -0.263

* Total electronic force
Force in Ry/a.u.










improvement by using the non-muffin-tin charge density to

calculate the force. The force constant k and the equili-
e
brium separation for H2 calculated with this method involve

errors of 1.4% and 7% respectively, while the Fx,(P) indi-

cates no binding at all.
Ex(>) we
However, when we compare Fx(p) and R we

find that there is some difference between them. This means

that the Hellmann-Feynman Theorem is not quite correct with

the approximate MSXa wavefunctions. Some reasons may be:

(1) Since the wavefunctions we used are hot exact

solutions to the Xa Schrodinger equations. It also assumes

the muffin-tin potential approximation and the potential is

generated only from a muffin-tin averaged charge density.

Thus, it is only an approximate Xa wavefunction and we cannot

claim that it should satisfy the Hellman-Feynman Theorem very

well. Actually, due to the polarization effect, the potential

on both sides of the atomic spheres should be different. The

potential in the bonding region facing the other nucleus

should be less shallow than the potential in the outer re-

ion. Thus, there should be more charge in the bonding region

than what is expected by using the muffin-tin averaged poten-

tial. This means that the actual force should be more attrac-

tive than the force we got by assuming the muffin-tin poten-

tial. This polarization effect should be especially severe

for the core state, for which the electronic charge accumu-

lates near the nucleus.









(2) Besides assuming the muffin-tin averaged poten-

tial, we also use only a finite set of basis functions with

small i values. This may not be important in calculating

the total energy but it may be important in calculating the

force, because the force is much more sensitive to the wave-

function than the total energy is. We further note that:

(i) In the force calculation for the hydrogen mole-

cule ( table (4.2) and figure (4.4)), the values of Fxa(P)

by using the basis (H =0,1; out =0,2)are closer to the values

of Exa) calculated with the basis (Hp=0; 9out=0,2) than
3R.
to that obtained with the basis (tH=0,1; zout=0,2). Further-

more, we have the .following relation:


S3Ex(P)
F x() ( =0,1; =0,2) < xa (i =O; R =0,2)
xci H out 'R I out


< 3Exa() ( ; =0,2) < F
aR H=0; out experimental


This can be expected since if the error in Fxa(P) (9H=0,1;

9out=0,2) is of first order of the errror in the wavefunction,

then the error in E x(p) (ZH=0,1; Zout=0.2) is of second

order and that the error in E (p) (. =0; out=0,2) should
xCi H out
be larger than that of Ex (P)(ZH=0,1; 9out=0,2) but smaller

than second order.

(ii) In table (4.3), we see that when we increase

the Z-values, the net force becomes less repulsive and agrees

better with 3Exa(P) The improvement is very slight. However,
3R
if we use only s orbitals in the atomic spheres,then FA=0










and the net force becomes more repulsive. Thus, the p orbital

is very important in the force calculations, although the s

orbital is almost sufficient for the total energy calculation.

In the force calculation for the nitrogen molecule,

we see that the MSXa forces have the same characteristic as

the Hartree-Fock forces. Again, the MSXa forces are too re-

pulsive. However, the relative values of forces contributed

by the different molecular states for both cases are compar-

able, and the ou state in both calculations contributes re-

pulsive force. Of course, we should not expect our result

could be the same as that of the Hartree-Fock calculation.

What we should expect is just a general characteristic of

these two kinds of forces.

In conclusion, the use of MSXa wavefunctions instead

of the muffin-tin averaged charge density to calculate the

force makes a big improvement. It gives a reasonable equili-

brium separation and a good force constant. However, due to

the muffin-tin potential and the truncation of k-values in

the basis set, the existing calculation of MSXa wavefunctions

still has to be modified to satisfy the Hellmann-Feynman

Theorem.













CHAPTER V


DIPOLE MOMENT CALCULATIONS



5.1 Introduction

In Chapter III,we have discussed how to compute the

dipole moment of a diatomic molecule by assuming the muffin-

tin charge density. In this chapter, we shall discuss the

same calculation but using the non-muffin-tin charge density,

i.e., the charge density formed from the MSXa wavefunctions.

We shall first discuss the computational procedure and then

give the results for some diatomic molecules. The results

show that the use of non-muffin-tin charge density indeed

provides big improvements over the use of muffin-tin charge

density.



5.2 Method of Dipole Moment Calculations

To calculate the dipole moment of a diatomic mole-

cule, we need to evaluate the integral p()z d3r contrib-

uted by the electronic cloud. Again, in the MSXa calcula-

tions, we have:


p(r ) = Z ni ui(r) ui(r)

As usual, we partition the coordinate space into three re-

gions: the atomic region, the intersphere region, and the









































II out














Fig. 5.1 Schematic Representation of A Heteronuclear
Diatomic Molecule in the MSXa Calculation.








outer sphere region (figure 5.1). As we have mentioned in

Chapter III, the value of the dipole moment is independent

of the choice of the origin for a neutral molecule, we can

choose the center of the outer sphere as the origin. Then

we have:


3 i i i
p()z dr = n (A B + out + II) (5.1)


i %* 3
where PA = J ui(r)ui(r)z0 d r





outi i *3
B = ui(r)ui (r)z0 d3r






11 = z ui(r)ui(r)zO d r


and the dipole moment u is: P = uN e


where PN = ZBRBO ZARAO is the nuclear contribution,



Pe = J p(r )z0 d3r is the electronic contribution,


ZA, ZB are respectively the atomic numbers of the two

nuclei A and B; p is in atomic units.
i
Integral Over Sphere A: uA

Inside sphere A, we have:









A, iA M
u. (r)= Ci R(r)Y(r )
1 Zm 2 A 9 A


and z0 = A RA0 = (4/3) rAY(rA) RAO

Thus, we have:


1A = ui(r ) ui(r) Z d r


= i CImA R (rAY m (rA)RZ(rA)Yr (rA)


x { (4T/3) rAY (rA) RAO d r


iA iAO
= C ciA c iC I (4r/3)2 I(Z1,m;92,m;l,0)

RA 3
x R (r)R(r)R )rA drA


A 2
RAO 62 RZ(rA)R2(rA)rA drA }


i
Integral Over Sphere B: up

Inside sphere B, we have


ai(?) = : Cis Re(rg) Y (?)
an = Z (r) Ym(rB)


and z0 = zB + RB0 = (4T/3) rB Yl(rB) + RB0


Thus, we have:









i = 3
B =J ui(r) u(r) d r


iB iB
= P CB CiB (47/3)2 I(t,1,m;2,,m;l,0)
= C 9m 92m


x R R(rB)R 2(rB)r drB + RB0o 6 [R (rB)R2(rB)rBdrB}
io

Integral Over the Outer Sphere: out

Inside the outer sphere, we have:


r iO
u.(r) = C R (r ) Y (ro)
1 Lm z 0 z 0

i %* r" 3
Thus = u.(r) u() z dr
out 1 1 0


= iCC. 1 C ou R (r )m(our)R (rOI (r


x (47/3) roY0(r0) d r


= C C m (47/3) I(k1,m;Z2,m;l,0)
ZIm Z2m


x R (ro)R (r0)ro dr0


Integral Over the Intersphere Region: ui

In the intersphere region, we can express ui(r) as

linear combinations of products of Neumann functions and

spherical harmonics as discussed in the force calculations.

In the same way as before, we obtain:









uI = (r) ui(r) z d3r


= Z lAmA 2m AA 1'2,m)' + AmAm IBB(1,2'm)


+ A 0A Om0( 1, 2,m) + 2 AAA AiB I 1 2,m)
Y 2m 00 il.m I2m AB

+2 AiA Ai0 Io( m) + 2 AiB Ai0
+2J, A m 2 m IAO(Il.2'm) + 2 A iBm A Im BO( 2'm)
m A 1 2 2 2m BO(t 1 2
(5.2)
where I AA(, 2 m) = T (rA )Y (r (KrA 2 A )z d 3r
AA 1 2 ti(A rA).(A Z2 A/Z2rA 0

m m ) 3
IBB I 1Z 2m) = (KrB B (r B 2 (KrB rB) zd r


IAB(" Z2,m) = n., (Kr r A )^^2 B9 31B
IAB(,IAm) = (KrA rA)-2(KrB)Y(rB)Zd r
Sm 3
AO ') =2' l (KrA)Ym(r )j (,r )Ym (r )z d r
II(ll 2'm) 21A 2.1 A 2.2 0 9 2 0J0

IO( 2m) = I (Kr)Y ( )j (Kr )Ym ( )z d r
BO J I B i B 2 0 22 0 0

IB0 1' 2' j(Kr )Y B( 22 (Kr )Ymr)d
2m Z 1 (0 kr1 0 Z2 0 2.)2 0 0

To evaluate the above integrals, again we employ the
Prolate Spheroidal.Coordinates as we did in the force calcula-
tions. However, since in this case we have a heteronuclear
diatomic case, so the limits of integration when transformed to
Prolate Spheroidal Coordinates would become more complex.
We shall discuss each of the integrals in the following:









() IAA 1(Kr A Yl rA 2 A Z2 A 0dr


This is a.two-center integral. We shall express the

integrand in terms of Prolate Spheroidal Coordinates with A

and 0 as the two foci (figure(5.1)):


S= (rA + r )/RA


n = (rA r0)/RAO


S= A = 0O


The reciprocal relations are:


r = RAO ( + n) r0 = AO ( )


cosO = (l+5n)/(4+n) and coso = (l-(n)/(5-n)
A 0

(2_)( 2)}2 2
sinOA +n sinG = ( -1)(1-n )
"+n o S-n

(5.3)
and the volume element is:


3 3 2 2
d r = (A RA)3 (E 2) dE dn d4


Since the two atomic spheres are in general of dif-

ferent size, the geometry of the partition is a function of









two parameters, the bond distance R and the ratio of the

radii of any two spheres, in this case, we choose 6 = RB/R.

Then we can express the different quantities in terms of

these two parameters:

RA = (1-)/R: RB = BR; Rout = R; RAO = BR;
RBO = (1-6)/R (5.4)

With the coordinate system as expressed in equation (5.3)

and also by the relations in (5.4), we have the following

restrictions for the intersphere region:



(1) rA > RA



implies ARAO(S+n) 3 RA or c+n > 2(1-B)/B



(2) r0 < Rout


implies 2RAO(S+n) S Rout or E-n < 2/6



(3) rB RB


2 2 2
and r = rA + R 2r Rcose
B A A A

24 242
implies n (1- B)En + E- (B +8-1) ) 0

(4)

(4) E 1


(5) -l < n < 1










(i) B>1/2


(ii) 6<1/2


12 2
r-


Fig. 5.2 Bounoary of the Intersphere Fegion in Prolate
Spheroidal Coordinates with A and O as the Foci.









Combining (1), (2), (3), (4) and (5), we have the

transformed region as defined by the following limits (figure

5.2):


1

0
2 3







1 2


1
if 2 2


if
if 6 <


where f (S) = +{(1-36)C- {(1-B)C + (B2+6-)}1}



2/6 2 C if 5 < 2/8-1

and f (C) =
1
S- 2/ if C > 2/8-1



Expanding the spherical harmonics and the Neumann

functions in terms of 5, q, p as we did in the force calcu-

lation in Chapter IV, we have the following result after

some straightforward algebra:


A= N () 1 t-wi] ((fi~--.)]
IA = k N (-1) z+k2+1(RAO/2)4 dC .


-2m-2v+2r rs-2
((_l)r Wim W2m m CC,+2-2m-2+2r ( RAs-2
vSe r1 2 R)









x 2(1_-2 s E1+Z2-2m-2v+2r-s {2C(1+2 ) hi' (s+l,)
I1
g +1 h *1 (sR) + C/(KRA)2 h2 (-1i}}
RAO 1 2 AO Z012


where v = v + v2

and the h'-functions are defined as:


h' (s,) = Z Z 1)+ A ( )
U12 V:o t rO

x Q(u v,1)(s+2n+u+v,S)
~ 1~2


where the A-coefficients are also defined in Chapter IV and the

Q-functions are defined in terms of sine and cosine integrals,

similarly to the corresponding terms in the force calcula-

tions except for the integration limits:


Q(uv,1)(,') = ( cos(z + )r (6 -6 )
Zi2 1 2 u+v,O u+v,2

I cost
+ sin(Z1+Z)n 2 u)vT1} 6R ,(1) r- dt


+ {sini(1+L,2)7 (Su+v,2-6u+v,O)


+ cosi( ,+e,) 6 u+v,}1 ) intdt
jKg(+tfirp5) tL

+ {sin(1- 2)r 6U+V,l (6vO-vl )


)+ cs(-2) (uv, u ,) ( )) dt
+ cosl(9L-L2)r (6u+v,O u+v,2 1f(5tfp) t
kRllRJf,









(2) IA: J1 nTk(KrA)Y(A) J2(KrO)Y ((0)z0d3r

This is also a two-center integrals, we can use the

same coordinate system as in the evaluation of IAA. The re-

sult is:


tm m
IAA = N NM (-1) (RAO/2)4( 1 2 1

-z J- )u+v m
n,.o -a (so 1


SO v=o =

2m 1 2-2v
V (-1)
2


m 2u ,-m-2v, 9.-m-2V2 1-2v -2u+v-r+L -2v -s
x C C C (KRAO) 1
u v r s AO


x (1-2 )2,-m-2v-r+-2v2-s v+r+s {(1+2



x h 2 (Z -2v -2u+v-r,92-2v2-s, ) S/(KRAO)


x h 2( ,-2v -2u+v-r-l, 2-2v,-s,()}}

J, -tCL-a.)-] t-4'-)J
where h2 2(3,9,) = '
gso VgO kls
1 20 as voc

u+k.+k,
( (-1) 1 +'(Z1+,2k,+u)(2+&,2k+v)


x Q^ Uv,2)(.3+2k,+u,Z,+2k2+v,S) }


and the Q-functions are defined by:









Q(uv,2)(p,q,) = {sin(KRAAo- (Z-_2)I) (6 -6 )
t2 AO u+v,0 u+v,2


+ cos(KRAO-i(i1-2)r) u+v, }

x JK S cost
P -dt
RoI~tf,ct) t(2KRA0,_t)q

+{cos(KRA0-( Z1-Z )) (u+v,2 6u+v,0)


+ sin(KRAO -(ZL1-2)O) u+v,l)


Sjsint dt
J<^o( ,p) tP(2
+{sin(KRA0+( +k2)7) (6u+v,0+6 u+v,2)


+ cos(KRAO+*( +z 2 )T) 6 u,(6u -6 ul)


1
"4o(ftg1) tP(2KRAO-t)q dt

(3) IBB: rliz(rB)Y(rB h(KrB)Y2rB)zO d3r
This is also a two-center integral. We treat this
integral by considering the Prolate Spheroidal Coordinates
with 0 and B as the two foci:


S= (r0+ rB)/RBO = (r rB)/RBO = 0 = B .


With this coordinate system, we have the following restric-









tion for the intersphere region:



(1) rB > RB


implies (RB0/2)(S n) > RB or n 5 28/(1 l )


(2) r0 < Rout


implies (RB0/2)( + n) ( Rout or 5 + n < 2/(1 8)


(3 rA > RA


2 2 2
and rA = rB + R + 2rgR cosOB


2 4 2 4 2
implies n + (1 2) n + + (1 y ) 0
Y Y


where y = 1-6 = RA/R


(4) 5 3 1


(5) -1 $ n ( 1


Combining all the above restrictions, we have the

transformed region as defined by the following (figure 5.3):









(1) Y<1/2 (B>1/2)





1- -
n= +2= ;_=


F t r Ri i\


-1 n=g-()


aii) y> 1/2 (8<1/2)




nn= 2 =

_/ E 2/y~

1 2 13
/ ._-------^ ^
rln=gl (E)
-I"





Fig. 5.3 Boundary of the Intersphere Region in
Prolate Spheroidal Coordinates with O
and B as the Foci.








2/Y -3. if y < (B> )

i >, =
1 if Y (B F )
2 2


C c 2/y + 1


g (C) n $ g2( )


2 2 (1_ 2))}
where g () = {- (1-y/2) + ((-y)2 (-Yy2))


g( 2/y + 2 if 5 < 2/y- 1

g (5) =
2/y 5 if 5 > 2/y- 1


Expanding the spherical harmonics and the Neumann
functions in , n, ), we have the following expression for

IBB:


BEB < < (-1)N+2+1 (RBO/2)4 d"' I
r--o , s o

(_ 1)z+2,-2m+3r-2v 2m 2m
V" 1 V2

xm C sZ+z2-2m+2r-2v 2(Ro ) s-2 E,+,2-2m+2r-2v-s
xr s 2(2BO

2
x (1-s2 ) {2((1+E2) h ) 2(s+1,) 3- +1 h (s,
+ R/ B BO

+ f/(.RB. 2 h 1.2(s,-1,W) }









where v = v, + v2

and the h3-functions are defined as:


h (s, i) = L(-1
2 uo ve a n=O


u+v (u,v) (u,v,3)
An (si+nu+v,C)
n1 2.ij


where the A-coefficients are defined as before and the Q-
functions are defined as below:


= ( cos (Z,+k2)I (Su+v,- u+v,2)


+ sin&(+Z +E2 uZ v sJt dt
u+v,1 t,


+ {sinj(tZ1+2) (6u+v,2- u+v,O)

+ cKRsoQ-Odl))
+ cos(Zi+k2)r au+vT jS }
k,(- 2.'Oap)


sint
-9dt


(4) IBO: f,1


+ icos(-z ( U+v,0+6u+v,2

(*Reo~l -,tp) .dt
+ sin(1^ z) 6u+v, (6v6v)} t



1("KrB)Y'(rB)j2(OKro)Y (r0)z0 d3r


We evaluate this kind of integrals by using the same

coordinate system as defined in the treatment of IBB, i.e.,

with 0 and B as the two foci. The result after some straight-

forward algebra is:


Q(u,v,3)( )
1i2




University of Florida Home Page
© 2004 - 2010 University of Florida George A. Smathers Libraries.
All rights reserved.

Acceptable Use, Copyright, and Disclaimer Statement
Last updated October 10, 2010 - - mvs