1.6: One Electron Moving About a Nucleus

• • Contributed by Jack Simons
• Professor Emeritus and Henry Eyring Scientist (Chemistry) at University of Utah

The Schrödinger equation for a single particle of mass m moving in a central potential (one that depends only on the radial coordinate r) can be written as

$\dfrac{-\hbar^2}{2\mu}\left( \dfrac{\partial^2}{\partial x^2} + \dfrac{\partial^2}{\partial y^2} + \dfrac{\partial^2}{\partial z^2} \right)\psi + V \left( \sqrt{x^2 + y^2 + z^2} \right) \psi = E\psi$

This equation is not separable in cartesian coordinates (x,y,z) because of the way x,y, and z appear together in the square root. However, it is separable in spherical coordinates

$\dfrac{\hbar^2}{2\mu r^2}\left( \dfrac{\partial}{\partial r} \left(r^2 \dfrac{\partial \psi}{\partial r} \right)\right) + \dfrac{1}{r^2 \sin \theta} \dfrac{\partial}{\partial \theta} \left( \sin \theta \dfrac{\partial \psi}{\partial \theta} \right) + \dfrac{1}{r^2 sin^2 \theta} \dfrac{\partial^2 \psi}{\partial \phi^2} + V(r)\psi = E\psi.$

Subtracting V(r)y from both sides of the equation and multiplying by $$\dfrac{-2\mu r^2}{\hbar^2}$$ then moving the derivatives with respect to r to the right-hand side, one obtains

$\dfrac{1}{sin \theta} \dfrac{\partial}{\partial \theta} \left( \sin \theta \dfrac{\partial \psi}{\partial \theta} \right) + \dfrac{1}{sin^2 \theta} \dfrac{\partial^2 \psi}{\partial \phi^2} = \dfrac{-2\mu r^2}{\hbar^2} (E-V(r))\psi - \dfrac{\partial}{\partial r} \left( r^2 \dfrac{\partial \psi}{\partial r} \right).$

Notice that the right-hand side of this equation is a function of r only; it contains no q or f dependence﻿. Let's call the entire right hand side F(r) to emphasize this fact.

To further separate the $$\theta$$ and $$\phi$$ dependence, we multiply by sin$$^2 \theta$$ and subtract the $$\theta$$ derivative terms from both sides to obtain

$\dfrac{\partial^2 \psi}{\partial \phi^2} = \textbf{F}(r)\psi sin^2 \theta - \sin \theta \dfrac{\partial}{\partial \theta} \left( \sin \theta \dfrac{\partial \psi}{\partial \theta} \right)$

Now we have separated the $$\phi$$ dependence from the $$\phi$$ and r dependence. If we now substitute $$\psi = \Phi(f) Q(r,\theta)$$ and divide by $$Psi$$ Q, we obtain

$\dfrac{1}{\Phi} \dfrac{\partial ^2\Phi}{\partial \phi^2} = \dfrac{1}{Q}\left( \textbf{F}(r) sin^2 \theta - \sin \theta \dfrac{\partial}{\partial \theta} \left( \sin \theta \dfrac{\partial Q}{\partial \theta}\right) \right).$

Now all of the $$\phi$$ dependence is isolated on the left hand side; the right hand side contains only r and $$\theta$$ dependence.

Whenever one has isolated the entire dependence on one variable as we have done above for the $$\phi$$ dependence, one can easily see that the left and right hand sides of the equation must equal a constant. For the above example, the left hand side contains no r or $$\theta$$ dependence and the right hand side contains no $$\phi$$ dependence. Because the two sides are equal, they both must actually contain no r, $$\theta$$, or $$\phi$$ dependence; that is, they are constant.

For the above example, we therefore can set both sides equal to a so-called separation constant that we call -m$$^2$$ . It will become clear shortly why we have chosen to express the constant in this form.

The Hydrogenic atom problem forms the basis of much of our thinking about atomic structure. To solve the corresponding Schrödinger equation requires separation of the r, $$\theta$$, and $$\phi$$ variables

The $$\Phi$$ Equation

$\Phi'' + m^2\Psi = 0$

which has as its most general solution

$\Phi = Ae^{im\phi} + Be^{-im\phi}.$

We must require the function $$\Phi$$ to be single-valued, which means that

$\Phi(\phi) = \Phi(2\pi + \phi) \: \text{or},$

$Ae^{im\phi}\left(1 - e^{2im\pi}\right) + Be^{-im\phi} \left( 1- e^{-2im\pi}\right) = 0.$

This is satisfied only when the separation constant is equal to an integer m = 0, ±1, ± 2, ... . and provides another example of the rule that quantization comes from the boundary conditions on the wavefunction. Here m is restricted to certain discrete values because the wavefunction must be such that when you rotate through 2$$\pi$$ about the z-axis, you must get back what you started with.

The $$\Theta$$ Equation

Now returning to the equation in which the $$\phi$$ dependence was isolated from the r and $$\theta$$ dependence and rearranging the $$\theta$$ terms to the left-hand side, we have

$\dfrac{1}{sin \theta} \dfrac{\partial}{\partial \theta}\left( \sin \theta \dfrac{Q}{\partial \theta}\right) - \dfrac{-m^2Q}{sin^2 \theta} = \textbf{F}(r)Q.$

In this equation we have separated $$\theta$$ and r variations so we can further decompose the wavefunction by introducing $$Q = \Theta(\theta) R(r)$$ , which yields

$\dfrac{1}{\Theta}\dfrac{1}{sin \theta}\dfrac{\partial }{\partial \theta} \left( \sin \theta \dfrac{\partial \Theta}{\partial \theta} \right) - \dfrac{m^2}{sin^2 \theta} = \dfrac{\textbf{F}(r)R}{R} = -\lambda,$

where a second separation constant, -l, has been introduced once the r and q dependent terms have been separated onto the right and left hand sides, respectively.

We now can write the $$\theta$$ equation as

$\dfrac{1}{sin \theta}\dfrac{\partial}{\partial \theta} \left( \sin \theta \dfrac{\partial \Theta}{\partial \theta} \right) - \dfrac{m^2 \Theta}{sin^2 \theta} = -\lambda \: \Theta,$

where m is the integer introduced earlier. To solve this equation for $$\Theta$$ , we make the substitutions z = cos$$\theta$$ and P(z) = $$\Theta(\theta)$$ , so $$\sqrt{1-z^2} = sin\theta$$ , and

$\dfrac{\partial}{\partial \theta} = \dfrac{\partial z}{\partial \theta}\dfrac{\partial }{\partial z} = -sin\theta \dfrac{\partial}{\partial z}$

The range of values for $$\theta$$ was $$0 \leq \theta < \pi$$ , so the range for z is -1 < z < 1. The equation for $$\Theta$$, when expressed in terms of P and z, becomes

$\dfrac{\text{d}}{\text{dz}}\left((1-z^2)\dfrac{\text{dP}}{\text{dz}}\right) - \dfrac{m^2P}{1-z^2} + \lambda P = 0.$

Now we can look for polynomial solutions for P, because z is restricted to be less than unity in magnitude. If m = 0, we first let

$P = \sum\limits_{k=0}^\infty a_kz^k,$

and substitute into the differential equation to obtain

$\sum\limits_{k=1}^\infty (k+2)(k+1)a_{k+2}z^k - \sum\limits_{k=0}^\infty (k+1)k a_kz^k + \lambda \sum\limits_{k=0}^\infty a_kz^k = 0.$

Equating like powers of z gives

$a_{k+2} = \dfrac{a_k(k(k+1)-\lambda)}{(k+2)(k+1)}.$

Note that for large values of k

$\dfrac{a_{k+2}}{a_k} \rightarrow \dfrac{k^2\left( 1+ \dfrac{1}{k} \right)}{k^2\left(1+\dfrac{2}{k} \right)\left(1+\dfrac{1}{k} \right)} = 1.$

Since the coefficients do not decrease with k for large k, this series will diverge for z = ± 1 unless it truncates at finite order. This truncation only happens if the separation constant $$\lambda$$ obeys $$\lambda$$ = l(l+1), where l is an integer. So, once again, we see that a boundary condition (i.e., that the wavefunction be normalizable in this case) give rise to quantization. In this case, the values of $$\lambda$$ are restricted to l(l+1); before, we saw that m is restricted to 0, ±1, ± 2, ... .

Since this recursion relation links every other coefficient, we can choose to solve for the even and odd functions separately. Choosing a$$_0$$ and then determining all of the even ak in terms of this a$$_0$$, followed by rescaling all of these a$$_k$$ to make the function normalized generates an even solution. Choosing a$$_1$$ and determining all of the odd a$$_k$$ in like manner, generates an odd solution.

For l= 0, the series truncates after one term and results in $$P_o(z) = 1$$. For l= 1 the same thing applies and P$$_1$$(z) = z. For l= 2, a$$_2 = -6 \dfrac{a_o}{2} = -3a_o$$ , so one obtains P$$_2 = 3z^2-1$$, and so on. These polynomials are called Legendre polynomials.

For the more general case where m ¹ 0, one can proceed as above to generate a polynomial solution for the Q function. Doing so, results in the following solutions:

$P^m_1(z) = (1-z^2)^\dfrac{|m|}{2}\dfrac{d^{|m|}P_1(z)}{dz^{|m|}}.$

These functions are called Associated Legendre polynomials, and they constitute the solutions to the $$\Theta$$ problem for non-zero m values.

The above P and e$$^{im\phi}$$ functions, when re-expressed in terms of $$\theta \:\text{and} \: \phi$$, yield the full angular part of the wavefunction for any centrosymmetric potential. These solutions are usually written as

$Y_{l,m}(\theta,\phi) = P_l^m \dfrac{cos \theta}{\sqrt{2\pi}}e^{im\phi}$

These are called spherical harmonics. They provide the angular solution of the r,$$\theta, \phi$$ Schrödinger equation for any problem in which the potential depends only on the radial coordinate. Such situations include all one-electron atoms and ions (e.g., H, He$$^+$$, Li$$^{++}$$ , etc.), the rotational motion of a diatomic molecule (where the potential depends only on bond length r), the motion of a nucleon in a spherically symmetrical "box" (as occurs in the shell model of nuclei), and the scattering of two atoms (where the potential depends only on interatomic distance).

The $$R$$ Equation

Let us now turn our attention to the radial equation, which is the only place that the explicit form of the potential appears. Using our derived results and specifying $$V(r)$$ to be the coulomb potential appropriate for an electron in the field of a nucleus of charge +Ze, yields:

$\dfrac{1}{r^2} \dfrac{d}{dr} \left( r^2\dfrac{dR}{dr} \right) + \left( \dfrac{2\mu}{\hbar^2} \left( E + \dfrac{Ze^2}{r} \right) - \dfrac{l(l + 1)}{r^2} \right) R = 0.$

We can simplify things considerably if we choose rescaled length and energy units because doing so removes the factors that depend on $$\mu$$, $$\hbar$$, and $$e$$. We introduce a new radial coordinate $$rho \: \text{and a quantity} \: \sigma$$ as follows:

$\rho = \sqrt{ \dfrac{-8\mu E}{\hbar^2} }r$

and

$\sigma^2 = -\dfrac{-\mu Z^2e^4}{2E\hbar^2}.$

Notice that if $$E$$ is negative, as it will be for bound states (i.e., those states with energy below that of a free electron infinitely far from the nucleus and with zero kinetic energy), $$rho$$ is real. On the other hand, if $$E$$ is positive, as it will be for states that lie in the continuum, $$rho$$ will be imaginary. These two cases will give rise to qualitatively different behavior in the solutions of the radial equation developed below.

We now define a function $$S$$ such that

$S(\rho) = R(r)$

and substitute $$S$$ for $$R$$ to obtain:

$\dfrac{1}{\rho^2} \dfrac{d}{d \rho}\left( \rho^2 \dfrac{dS}{d \rho} \right) + \left( - \dfrac{1}{4} - \dfrac{l(l+1)}{\rho^2} + \dfrac{\sigma}{\rho} \right) S = 0.$

The differential operator terms can be recast in several ways using

$\dfrac{1}{\rho^2}\dfrac{d}{d \rho}\left( \rho^2 \dfrac{dS}{d \rho} \right) = \dfrac{d^2S}{d \rho^2} + \dfrac{2}{\rho}\dfrac{dS}{d \rho} = \dfrac{1}{\rho} \dfrac{d^2}{d \rho^2}(\rho S).$

It is useful to keep in mind these three embodiments of the derivatives that enter into the radial kinetic energy; in various contexts it will be useful to employ various of these.

The strategy that we now follow is characteristic of solving second order differential equations. We will examine the equation for S at large and small $$\rho$$ values. Having found solutions at these limits, we will use a power series in $$\rho$$ to "interpolate" between these two limits.

Let us begin by examining the solution of the above equation at small values of $$\rho$$ to see how the radial functions behave at small r. As $$\rho \rightarrow$$0, the second term in the brackets will dominate. Neglecting the other two terms in the brackets, we find that, for small values of $$rho$$ (or r), the solution should behave like $$\rho^L$$ and because the function must be normalizable, we must have $$L\geq 0$$. Since L can be any non-negative integer, this suggests the following more general form for S($$\rho$$) :

$S(\rho) \approx \rho^L e^{-a\rho}.$

This form will insure that the function is normalizable since S$$(\rho) \rightarrow 0 \text{as} r\rightarrow \infty$$ for all L, as long as $$rho$$ is a real quantity. If $$\rho$$ is imaginary, such a form may not be normalized (see below for further consequences).

Turning now to the behavior of S for large $$\rho$$, we make the substitution of $$S(\rho)$$ into the above equation and keep only the terms with the largest power of $$\rho$$ (e.g., first term in brackets). Upon so doing, we obtain the equation

$a^2\rho^Le^{-a\rho} = \dfrac{1}{4} \rho^Le^{-a\rho},$

which leads us to conclude that the exponent in the large-$$rho$$ behavior of S is a = $$\dfrac{1}{2}.$$

Having found the small- and large-$$\rho$$ behaviors of S($$\rho$$), we can take S to have the following form to interpolate between large and small $$rho$$-values:

$S(\rho) = \rho^Le^\dfrac{-\rho}{2}P(\rho),$

where the function L is expanded in an infinite power series in $$\rho$$ as $$P(\rho) = \sum a_k\rho^k$$. Again Substituting this expression for S into the above equation we obtain

$P''\rho + P'(2L_2-\rho) + P(\sigma -L-1) =0,$

and then substituting the power series expansion of P and solving for the ak's we arrive at:

$a_{k+1} = \dfrac{(k- \sigma + L +1) a_k}{(k+1)(k+2L+2)}.$

For large k, the ration of expansion coefficients reaches the limit $$\dfrac{a_k+1}{a_k} = \dfrac{1}{k},$$ which has the same behavior as the power series expansion of e$$^\rho$$. Because the power series expansion of P describes a function that behaves like e$$^\rho$$ for large $$\rho$$, the resulting S($$\rho$$) function would not be normalizable because the e$$^{\dfrac{-\rho}{2}}$$ factor would be overwhelmed by this e$$^\rho$$ dependence. Hence, the series expansion of P must truncate in order to achieve a normalizable S function. Notice that if $$\rho$$ is imaginary, as it will be if E is in the continuum, the argument that the series must truncate to avoid an exponentially diverging function no longer applies. Thus, we see a key difference between bound (with r real) and continuum (with $$\rho$$ imaginary) states. In the former case, the boundary condition of non-divergence arises; in the latter, it does not.

To truncate at a polynomial of order n', we must have n' - $$\sigma$$ + L+ l= 0. This implies that the quantity $$\sigma$$ introduced previously is restricted to $$\sigma$$ = n' + L + l , which is certainly an integer; let us call this integer n. If we label states in order of increasing n = 1,2,3,... , we see that doing so is consistent with specifying a maximum order (n') in the P($$\rho$$) polynomial n' = 0,1,2,... after which the l-value can run from l = 0, in steps of unity up to L = n-1.

Substituting the integer n for $$\sigma$$ , we find that the energy levels are quantized because $$\sigma$$ is quantized (equal to n):

$E = -\dfrac{\mu Z^2e^4}{2\hbar^2 n^2} :\ \text{and} :\ \rho = \dfrac{Zr}{a_on}.$

Here, the length a$$_o$$ is the so called Bohr radius $$\left( a_o = \dfrac{\hbar^2}{\mu e^2} \right)$$ ; it appears once the above E-expression is substituted into the equation for $$\rho$$. Using the recursion equation to solve for the polynomial's coefficients a$$_k$$ for any choice of n and l quantum numbers generates a so-called Laguerre polynomial; P$$_{n-L-1}(\rho)$$. They contain powers of $$\rho$$ from zero through n-l-1.

This energy quantization does not arise for states lying in the continuum because the condition that the expansion of P$$(\rho)$$ terminate does not arise. The solutions of the radial equation appropriate to these scattering states (which relate to the scattering motion of an electron in the field of a nucleus of charge Z) are treated on p. 90 of EWK.

In summary, separation of variables has been used to solve the full r,$$\theta ,\phi$$ Schrödinger equation for one electron moving about a nucleus of charge Z. The $$\theta \: \text{and}\: \phi$$ solutions are the spherical harmonics $$Y_{L,m} (\theta,\phi).$$ The bound-state radial solutions

$R_{n,L}(r) = S(\rho) = \rho^Le^{\dfrac{-\rho}{2}}P_{n-L-1}(\rho)$

depend on the n and l quantum numbers and are given in terms of the Laguerre polynomials (see EWK for tabulations of these polynomials).

Summary

To summarize, the quantum numbers l and m arise through boundary conditions requiring that $$\psi (\theta)$$ be normalizable (i.e., not diverge) and $$\psi (\phi) = \psi(\phi+2\pi).$$ In the texts by Atkins, EWK, and McQuarrie the differential equations obeyed by the $$\theta \: \text{and} \: \phi$$ components of Y$$_{l,m}$$ are solved in more detail and properties of the solutions are discussed. This differential equation involves the three-dimensional Schrödinger equation's angular kinetic energy operator. That is, the angular part of the above Hamiltonian is equal to $$\hbar^2 \dfrac{L^2}{2mr^2}$$, where L$$^2$$ is the square of the total angular momentum for the electron.

The radial equation, which is the only place the potential energy enters, is found to possess both bound-states (i.e., states whose energies lie below the asymptote at which the potential vanishes and the kinetic energy is zero) and continuum states lying energetically above this asymptote. The resulting hydrogenic wavefunctions (angular and radial) and energies are summarized in Appendix B for principal quantum numbers n ranging from 1 to 3 and in Pauling and Wilson for n up to 5.

There are both bound and continuum solutions to the radial Schrödinger equation for the attractive coulomb potential because, at energies below the asymptote the potential confines the particle between r=0 and an outer turning point, whereas at energies above the asymptote, the particle is no longer confined by an outer turning point (see the figure below).

The solutions of this one-electron problem form the qualitative basis for much of atomic and molecular orbital theory. For this reason, the reader is encouraged to use Appendix B to gain a firmer understanding of the nature of the radial and angular parts of these wavefunctions. The orbitals that result are labeled by n, l, and m quantum numbers for the bound states and by l and m quantum numbers and the energy E for the continuum states. Much as the particle-in-a-box orbitals are used to qualitatively describe $$\pi$$ - electrons in conjugated polyenes, these so-called hydrogen-like orbitals provide qualitative descriptions of orbitals of atoms with more than a single electron. By introducing the concept of screening as a way to represent the repulsive interactions among the electrons of an atom, an effective nuclear charge Z$$_{eff}$$ can be used in place of Z in the $$\psi_{n,l,m}$$ and E$$_{n,l}$$ to generate approximate atomic orbitals to be filled by electrons in a many-electron atom. For example, in the crudest approximation of a carbon atom, the two 1s electrons experience the full nuclear attraction so Z$$_{eff}$$=6 for them, whereas the 2s and 2p electrons are screened by the two 1s electrons, so Z$$_{eff}$$= 4 for them. Within this approximation, one then occupies two 1s orbitals with Z=6, two 2s orbitals with Z=4 and two 2p orbitals with Z=4 in forming the full six-electron wavefunction of the lowest-energy state of carbon.