Electronic Properties of Materials/Quantum Mechanics for Engineers/Many Electron Atoms

This is the tenth chapter of the first section of the book Electronic Properties of Materials.

We now have a solutions for Hydrogen.

$${\begin{aligned}H&=T_{electron}+V_{ion-electric}\\&={-\hbar ^{2} \over 2m}\nabla ^{2}-{|e|^{2} \over r}\\H\psi _{nlmm_{s}}(r)&=E_{n}\psi _{nlmm_{s}}\\E_{n}&={-1 \over n^{2}}\mathbb {R} \\\psi _{nlmm_{s}}&=R_{nl}Y_{lm}\chi _{m_{s}}\end{aligned}}$$

In undergraduate chemistry class, we were shown this solution and told "... and so it goes for the rest of the periodic table..." However, the truth is, it isn't so simple.

Consider Lithium:

<FIGURE> "Diagram of Lithium Ion" (Description)

We can apply center of mass corrections but for simplicity, assume that the nucleus has $M_{noc}=\infty$.

What does hydrogen look like?

$$H=T_{electron}+V_{ion-electron}+V_{electron-electron}$$

$${\begin{aligned}T_{electron}&=T_{e_{1}}+T_{e_{2}}+T_{e_{3}}\\&={-\hbar ^{2} \over 2m}\left({\nabla _{1}}^{2}+{\nabla _{2}}^{2}+{\nabla _{3}}^{2}\right)\\\\V_{ion-electron}&={-3e^{2} \over r_{1}}+{-3e^{2} \over r_{2}}+{-3e^{2} \over r_{3}}\\V_{electron-electron}&={e^{2} \over |r_{1}-r_{2}|}+{e^{2} \over |r_{2}-r_{3}|}+{e^{2} \over |r_{3}-r_{1}|}\end{aligned}}$$

Without the electron-electron interaction the problem would be much simpler.

$${\begin{aligned}H&=T_{1}+V_{ie1}+T_{2}+V_{ie2}+T_{3}+V_{ie3}\\&=h_{1}+h_{2}+h_{3}\end{aligned}}$$

This makes for a separable solution:

${\begin{aligned}H\Psi &=E\Psi \\\Psi &=\psi _{1}\psi _{2}\psi _{3}\end{aligned}}$

Here, each $\psi _{i}$ terms is just the hydrogen wave function with a slight modification for $V={3e^{2} \over r}={ze^{2} \over r}$.

${\begin{aligned}E=\varepsilon _{1}+\varepsilon _{2}+\varepsilon _{3}\\\varepsilon _{i}=-{Z^{2}\mathbb {R} \over \eta ^{2}}\end{aligned}}$

Here $Z$ is the atomic number.

Since electrons are fermions, $\Psi$ must be totally anti-symmetric as seen in 2 Particles in a Box in Chapter 4.

$$\Psi (r_{1},r_{2},r_{3})={1 \over {\sqrt {6}}}\left|{\begin{aligned}\psi _{1}(r_{1})\quad \psi _{2}(r_{1})\quad \psi _{3}(r_{1})\\\psi _{1}(r_{2})\quad \psi _{2}(r_{2})\quad \psi _{3}(r_{2})\\\psi _{1}(r_{3})\quad \psi _{2}(r_{3})\quad \psi _{3}(r_{3})\end{aligned}}\right|$$

This is known as the Slater Determinant.

Approximating Lithium

Unfortunately, we can't just discard electron to electron interactions. The magnitude of ${\textstyle V_{ion-electron}}$ and ${\textstyle V_{electron-electron}}$ are equivalent. So, what do we do? Solve the problem as is? Well, after you solve $H_{Li}$, you could quite possibly earn a Nobel Prize and/or a Fields Medal which makes this a wholly impractical method for this course. Instead we will use approximate methods for quantum mechanical problems more complicated than our earlier Hydrogen solution, and we currently have two primary methods for doing this.

Thomas-Fermi Model

The Thomas-Fermi Model is a semi-classical statistical method where we replace the exact potential ${\textstyle V_{ion-electron}}$ and ${\textstyle V_{electron-electron}}$ with an effective potential, ${\textstyle V_{eff}}$, the screened coulomb potential.

Consider an electron gas with so many electrons that you don't have to count them one at a time, but instead just consider a continuous charge distribution. If we put a positive charge in the electron cloud then the positive charge would attract the electrons so the electron distribution would no longer be uniform. That said, the negative charges nearest the positive charge "screen" the electrons further out such that the fringe electrons don't feel as strong of an attraction.

Working out the details of this method requires a discussion of electron gasses, provided later in this text. Meanwhile, you should know that this hypothetical situation turns into a spherically distributed electron gas where the charge density, $\rho$, is slowly varying. Although $\rho (r)$ does change with $r$, at each $r$, we can treat $\rho$ as a uniform gas. Additionally, a point charge $+Q$, placed in a uniform electron gas has a screened potential of:

$$\Phi (r)={Q \over r}e^{-k_{o}r}$$

This is sometimes called a Yukawa Potential. Also, ${\textstyle {1 \over k_{o}}}$ is the Thomas-Fermi screening length and is related to the density of quantum states of the highest filled state, called the density of states at the Fermi Level.

Finally, for charge neutral atoms, where the number of electrons equals the atomic number, with a slowly varying electron gas in a spherically symetric central potential:

$$V(r)\simeq \underbrace {-{Ze^{2} \over r}} _{\text{Attractive Potential}}+\underbrace {1.794\ {Z^{4 \over 3} \over a_{o}}} _{\text{Repulsive Potential}}+\ \cdots$$

<FIGURE> "Electron Gas Potential" (Description)

The truth of the matter is that the Thomas-Fermi Model isn't a great solution for the multi-electron atom. Back in the 1920's it was quite successful and is still useful for simple approximations and the initial input to more advanced methods.

Hartree-Fock Method

A better method is the iterative Hartree-Fock Method, or the Self-Consistent Field Method. For this method, we start by assuming that we can write:

$$\Psi =\prod _{i=1}^{n}\psi _{i}$$

This assumption further implies that the can say that ${\textstyle H=\sum _{i=1}^{N}h_{i}}$. So, what is a reasonable ${\textstyle h_{i}}$?

$$h_{i}={-\hbar ^{2} \over 2m}\nabla _{i}^{2}+V_{o}(r_{i})+{1 \over 2}\sum _{j\neq i}{e^{2} \over |r_{i}-r_{j}|}$$

Now, this gives us:

$$H=\sum _{i=1}^{N}\left[{-\hbar ^{2} \over 2m}\nabla _{i}^{2}+V_{o}(r_{i})+{1 \over 2}\sum _{j\neq i}{e^{2} \over |r_{i}-r_{j}|}\right]$$

${\textstyle V_{o}(r_{i})}$ is all the potential terms not involving electron-electron interaction. For the single atom, ${\textstyle V_{o}(r_{i})={Ze^{2} \over r_{i}}}$, but this method is generally applicable to more complex systems such as molecules. Hartree-Fock is the intellectual forbearer of several modern quantum chemical techniques.

<Source> "Primer on Calculus of Variation and Lagrange Multipliers"

Energy Function

The functional we are interested in is the energy, where ${\textstyle E\left[\Psi (r)\right]=\int \Psi ^{*}H\Psi \ dr}$. We want to find the extrema subject to the constraints ${\textstyle \int \psi _{i}^{*}\psi _{i}\ dr=1}$ for all $i$. We will do this by using Lagrangian multipliers $\varepsilon _{i}$ such that the integral we're interested in is:

$$\int \Psi ^{*}H\Psi \ \operatorname {d} \!r-\sum _{i}\varepsilon _{i}\int \psi _{i}^{*}\psi _{i}\ \operatorname {d} \!r$$

Applying Calculus of Variation

We want to use calculus of variation to find the functions that make this function stationary (an extrema). Consider this function:

$$I[\phi _{1}(x)\ \phi _{2}(x)]=\int f(\phi _{1}(x),\ \phi _{2}(x),\ x)\ \operatorname {d} \!x$$

We want to vary the function to look for the stability condition. In other words, we want to change $\phi _{i}(x)$ at an arbitrary constant $x$.

$$\Delta I=I[\phi _{1}(x)+\varepsilon _{1}\eta _{2}(x),\ \phi _{2}(x)+\varepsilon _{2}\eta _{2}(x)]-I[\phi _{1},\ \phi _{2}]$$

Here $\varepsilon _{i}$ is a small constant (we take the limit of $\varepsilon _{i}$ as it goes to zero), and $\eta _{i}$ is an arbitrary function that deviates $\phi _{i}$ (this must be continuous over the limits of integration. This gives us:

$$\Delta I=\int f(\phi _{1}(x)+\varepsilon _{1}\eta _{2}(x),\ \phi _{2}(x)+\varepsilon _{2}\eta _{2}(x),\ x]-\int f(\phi _{1}(x),\ \phi _{2}(x),\ x)$$

The Maclaurin expansion in powers of $\varepsilon _{i}$ is:

$$f(x)=f(0)+f'(0)+{1 \over 2!}f''(0)+{1 \over 3!}f'''(0)+\cdots$$

$${\begin{aligned}\Delta I&=\int f(\phi _{1}(x),\ \phi _{2}(x),\ x)+{\partial f \over \partial \phi _{1}}\ \varepsilon _{1}\eta _{1}(x)+{\partial f \over \partial \phi _{2}}\ \varepsilon _{2}\eta _{2}(x)+{\text{Higher Order Terms}}-f(\phi _{1}(x),\ \phi _{2}(x),x)\\&={\partial f \over \partial \phi _{1}}\ \varepsilon _{1}\eta _{1}(x)+{\partial f \over \partial \phi _{2}}\ \varepsilon _{2}\eta _{2}(x)\\\delta I&=\int {\partial f \over \partial \phi _{1}}\ \delta \phi _{1}+{\partial f \over \partial \phi _{2}}\ \delta \phi _{2}\\&=\int \delta f\qquad {\text{(The variation of the functional}}\ f)\end{aligned}}$$

$\delta \phi _{i}=\varepsilon _{i}\eta _{i}$ in the limit as $\varepsilon _{i}$ goes to zero.

Let's return to our problem:

$I=\int \Psi ^{*}H\Psi -\sum \varepsilon _{j}\psi _{j}^{*}\psi _{j}\operatorname {d} \!r$

Writing this function out in detail gives us:

$$I=\int \underbrace {(\psi _{1}^{*},\ \cdots ,\ \psi _{n}^{*})(h_{1},\ \cdots ,\ h_{n})(\psi _{1},\ \cdots ,\ \psi _{n})-\varepsilon _{1}\psi _{1}^{*}\psi _{1}-\cdots -\varepsilon _{n}\psi _{n}^{*}\psi _{n}} _{f(\psi _{1}^{*},\ \cdots ,\psi _{n}^{*})}\operatorname {d} \!r$$

Note that you could alternatively carry through $\psi$ in this step, but in the end you don't need them. I chose to utilize $\psi ^{*}$ instead of $\psi$ for purely aesthetic reasons. This gives us:

$$\delta I=0=\int \sum _{q=1}^{n}{\partial f \over \partial \psi _{q}^{*}}\delta \psi _{q}^{*}\ \operatorname {d} \!r$$

Looking at the $q^{th}$ term in the summation, ${\textstyle {\partial f \over \partial \psi _{q}^{*}}\ \delta \psi _{q}^{*}}$, we can expand to:

$${\begin{aligned}&\sum _{j=1}^{n}\left[(\psi _{1}^{*},\ \cdots ,\ \delta \psi _{q}^{*},\ \cdots ,\psi _{n}^{*})(h_{j})(\psi _{1},\ \cdots ,\ \psi _{q},\ \cdots ,\psi _{n})\right]-\delta \psi _{q}^{*}\varepsilon _{q}\psi _{q}\\&\sum _{j=1}^{n}\left[(\psi _{1}^{*},\ \cdots ,\ \delta \psi _{q}^{*},\ \cdots ,\psi _{n}^{*})\left[\underbrace {{-\hbar ^{2} \over 2m}\nabla _{j}^{2}} _{T_{j}}+\underbrace {V_{o}(r_{j})} _{V_{j}^{o}}+\underbrace {{1 \over 2}\sum _{k\neq j}^{n}{e^{2} \over |r_{j}-r_{k}|}} _{V_{jk}^{ee}}\right](\psi _{1},\ \cdots ,\ \psi _{q},\ \cdots ,\psi _{n})\right]-\delta \psi _{q}^{*}\varepsilon _{q}\psi _{q}\end{aligned}}$$

For the $T_{j}$ and $V_{j}^{o}$ terms, there exists two conditions:

1. When $j\neq q$, the terms inside the summation become $math$, but integrating this will go to zero because $\delta \psi _{q}^{*}$ and $\psi _{q}$ are orthogonal.
2. The only non-zero terms come about for $j=q$

Thus, dividing through by $\Psi ^{*}\Psi$ gives us $\delta \psi _{q}^{*}(T_{q}+V_{q}^{o})\psi _{q}$. Here, the $V_{jk}^{ee}$ term is again zero for $j\neq q$, because the orthogonal functions come out from $\delta \psi _{q}^{*}\psi _{q}$. The resulting term is:

$$\delta \psi _{q}^{*}\left({1 \over 2}\sum _{k\neq q}\psi _{k}^{*}{e^{2} \over |r_{q}-r_{k}|}\psi _{k}\right)\psi _{q}$$

Because $V_{jk}^{ee}$ depends on both $j_{q}$ and $k$, we retain the wavefunctions $\psi _{k}$ inside the summation.

The functional variation now looks like:

$${\begin{aligned}\delta I&=\int \sum _{q=1}^{n}\delta \psi _{q}^{*}(T_{q}+V_{q}^{o})\psi _{q}+\delta \psi _{q}^{*}\left({1 \over 2}\sum _{k\neq q}\psi _{k}^{*}{e^{2} \over |r_{q}-r_{k}|}\psi _{k}\right)\psi _{q}-\delta \psi _{q}^{*}\varepsilon _{q}\psi _{q}\\&=0\end{aligned}}$$

Which can be solved as $n$ equations.

$$\left[{-\hbar ^{2} \over 2m}\nabla ^{2}+V_{o}(r)+V_{d}(r)\right]\psi _{k}(r)=\varepsilon _{k}\psi _{k}(r)$$

Where $V_{d}$ is the Hartree Term which is: $V_{d}=\sum _{j\neq k}^{n}\int \psi _{j}^{*}(r_{j})\psi _{j}(r_{j}){e^{2} \over |r_{j}-r_{k}|}\ \operatorname {d} \!r$

This looks like a Schrodinger equation, but is it? Yes, sort of... $\psi _{i}$ are Hartree wave functions. Remember that these are single electron $\psi _{i}$ since we started with $H=\sum _{i}h_{i}$. $\varepsilon _{i}$ are Lagrangian multipliers for enforce normalization. However, we can use these $\{\psi _{i}\}$ and $\varepsilon _{i}$ to study the system.

Hartree Equation

The Hartree equation has a physically intuitive form (if anything in quantum mechanics is physically intuitive). THe first term is kinetic energy, the second term is the electron-ion or any external potential energy (such as E or B), and the third term is the electron-electron potential energy defined by the sum:

$$\sum _{j\neq k}\int \underbrace {\psi _{j}^{*}(r_{j})\psi _{j}(r_{j})} _{\rho _{j}}{e^{2} \over |r_{j}-r_{k}|}\ \operatorname {d} \!r$$

Here $\rho _{j}$ is the charge density of the $j^{th}$ electron. So this sum can be thought of as an electron at $r_{k}$ interacting with some charge at $r_{j}$.

What's wrong with this picture? Well, excluding our single particle assumption, we're putting several electrons into this system, which are fermions, but our wavefunction is not totally anti-symmetric. What we need to do is write $\Psi$ using a Slater determinant.

$${\begin{aligned}\Psi &={1 \over {\sqrt {n!}}}\left|{\begin{matrix}\psi _{1}(r_{1})&\cdots &\psi _{n}(r_{1})\\\cdots \\\psi _{1}(r_{n})&\cdots &\psi _{n}(r_{n})\end{matrix}}\right|\\&={1 \over {\sqrt {n!}}}\sum _{\rho }(-1)^{\rho }\psi _{1}(r_{1})\ \cdots \ \psi _{n}(r_{n})\end{aligned}}$$

Put this back into our original energy function and apply calculus of variation methods to the whole expression to find our answer. This was done by Fock (and Slater) in 1930. The resulting set of expressions is exactly the same except for one additional term.

$$\left[{-\hbar ^{2} \over 2m}\nabla ^{2}+V_{o}(r)+V_{d}(r)+V_{\text{exchange}}(r)\right]\psi _{k}(r)=\varepsilon _{k}\psi _{k}(r)$$

The evaluation of $V_{\text{exchange}}$ is quite tricky. The expression for the expectation value is:

$$\langle E_{\text{exchange}}\rangle ={1 \over 2}\sum _{1\neq j}\int \int \psi _{i}^{*}(r_{i})\psi _{j}^{*}(r_{j}){e^{2} \over |r_{i}-r_{j}|}\underbrace {\psi _{i}(r_{j})\psi _{j}(r_{i})} _{\operatorname {d} \!{r_{i}}\operatorname {d} \!{r_{j}}}$$

In many cases it is acceptable to approximate using the so-called "free electron exchange" ${\textstyle V_{\text{exchange}}(r)=-e^{2}\left({3 \over \pi }\right)^{1 \over 3}\rho (r)^{1 \over 3}}$

Self-Consistent Field Loop

These equations are not trivial to solve. Notice that the operator depends on $\rho (r)$, but to get $\rho (r)$ you need $\psi _{i}(r)$ which requires solving the set of equations. This is where we solve by iterative, self-consistent methods. This technique is the basis for all of modern quantum mechanical methods.

Write ${\textstyle (T+V(r))\psi _{k}(r_{k})=\varepsilon _{k}\psi _{k}(r_{k})}$, where ${\textstyle V(r)}$ is all the terms, not $T$. It will depend on $\rho (r)$ or $\psi (r)$ or both.

1. Guess ${\textstyle V^{(1)}(r)}$. You can use known solutions, random guesses, the Thomas-Fermi, or any other method of guessing you wish.
2. Solve the system of equations to find $\{\psi _{i}^{(1)}(r)\}$.
3. Calculate from $\{\psi _{i}^{(1)}(r)\}$.
4. If ${\textstyle V^{(2)}\simeq V^{(1)}}$ end. Otherwise return to step 2 and keep looping until ${\textstyle V^{(n)}\simeq V^{(n+1)}}$ to within a given tolerance.

This is called a self-consistent field loop (S.C.F.).

These type of SCF methods can be applied to all sorts of systems, but since this chapter is about many-electron problems, we will look at the Herman-Skillman atomic data, calculated from the Hartree-Fock method in 1963.

<SOURCE> "Herman-Skillman" atomic data"