Solve Schrödinger equation for hydrogen atom?
1 Answer
Talk about a tall order. Obviously (?) this will be an obscenely long answer. More like a tome.
Note that each wave function component given as part of the overall wave function is unnormalized, except for when explicit expressions are presented.
SUMMARY
- Separation of variables into
r and(theta,phi) - Separation of variables of
(theta,phi) intotheta andphi - Solving the
phi part - Solving the
theta part - Putting together the
phi andtheta parts - "Solving" the
r part - Putting together the overall wave function
I would hope you know that you only need to know how to use the wave functions (how to normalize them, how to show they are orthogonal, how to find the radial and angular nodes, etc), and maybe kind of know basic steps of the following derivation.
But see below anyway, because it took me 4+ hours...
The following derivation was adapted from here and from Physical Chemistry: A Molecular Approach by McQuarrie & Simon.
0. INITIAL DEFINITIONS
We begin from the time-independent Schrodinger equation (SE)
hatHpsi = Epsi ,
which for hydrogen atom, has the Hamiltonian
hatH = -ℏ^2/(2mu) nabla^2 - e^2/(4piepsilon_0r) ,where:
ℏ = h//2pi is the reduced Planck's constant.mu = (m_p m_e)/(m_p + m_e) is the reduced mass of the electron and proton.nabla^2 = 1/r^2 (del)/(delr) (r^2 (del)/(delr)) + 1/(r^2 sintheta) (del)/(del theta)(sintheta (del)/(del theta)) + 1/(r^2sin^2theta) (del^2)/(delphi^2) is the Laplacian operator in spherical coordinates.e is the elementary charge,1.602 xx 10^(-19) "C" , e.g. positive for a proton, negative for an electron.epsilon_0 = 8.854187817 xx 10^(-12) "F"cdot"m"^(-1) is the vacuum permittivity.r is the radial position. Under the Born-Oppenheimer approximation, we assume the nucleus is fixed so thatr becomes the radial distance of the electron from the nucleus.
1. SEPARATION OF VARIABLES
We first have to get the SE into standard 2nd order partial differential equation form. Plug in the Laplacian:
-ℏ^2/(2mu)[1/r^2 (del)/(delr) (r^2 (del)/(delr)) + 1/(r^2 sintheta) (del)/(del theta)(sintheta (del)/(del theta)) + 1/(r^2sin^2theta) (del^2)/(delphi^2)]psi - e^2/(4piepsilon_0r) psi = Epsi
Subtract over
1/r^2 (del)/(delr) (r^2 (delpsi)/(delr)) + 1/(r^2 sintheta) (del)/(del theta)(sintheta (delpsi)/(del theta)) + 1/(r^2sin^2theta) (del^2psi)/(delphi^2) + (2mu)/(ℏ^2)(E + e^2/(4piepsilon_0r))psi = 0
We assume separability of the wave function such that
Since
Y/r^2 (del)/(delr) (r^2 (delR)/(delr)) + R/(r^2 sintheta) (del)/(del theta)(sintheta (delY)/(del theta)) + R/(r^2sin^2theta) (del^2Y)/(delphi^2) + (2mu)/(ℏ^2)(E + e^2/(4piepsilon_0r))RY = 0
Multiply through by
overbrace(1/R (d)/(dr) (r^2 (dR)/(dr)))^("Radial Component") + overbrace(1/Y 1/(sintheta) (del)/(del theta)(sintheta (delY)/(del theta)) + 1/Y 1/(sin^2theta) (del^2Y)/(delphi^2))^"Angular Component" + overbrace((2mu r^2)/(ℏ^2) (E + e^2/(4piepsilon_0r)))^"Radial Component" = 0
By subtraction, these two functions are equal to each other:
1/R (d)/(dr) (r^2 (dR)/(dr)) + (2mu r^2)/(ℏ^2) (E + e^2/(4piepsilon_0r)) = -[1/Y 1/(sintheta) (del)/(del theta)(sintheta (delY)/(del theta)) + 1/Y 1/(sin^2theta) (del^2Y)/(delphi^2)]
But since they are functions of different variables, that can only work if they are equal to the same separation constant, which we arbitrarily designate
-[1/R (d)/(dr) (r^2 (dR)/(dr)) + (2mu r^2)/(ℏ^2) (E + e^2/(4piepsilon_0r))] = -lambda
1/Y 1/(sintheta) (del)/(del theta)(sintheta (delY)/(del theta)) + 1/Y 1/(sin^2theta) (del^2Y)/(delphi^2) = -lambda
Now we rearrange these ordinary differential equations separately. Take the radial component first. Switch the signs and then multiply through by
ul((d)/(dr) (r^2 (dR)/(dr)) + (2mu r^2)/(ℏ^2) (E + e^2/(4piepsilon_0r))R - lambdaR = 0)
An analogous process follows for the angular equation and we obtain:
ul(1/(sintheta) (del)/(del theta)(sintheta (delY)/(del theta)) + 1/(sin^2theta) (del^2Y)/(delphi^2) + lambdaY = 0)
2. FURTHER SEPARATING THE ANGULAR PART
The angular part is still a function of
Phi/(sintheta) (del)/(del theta)(sintheta (delTheta)/(del theta)) + Theta/(sin^2theta) (del^2Phi)/(delphi^2) + lambdaThetaPhi = 0
Multiply through by
overbrace((sintheta)/(Theta) (d)/(d theta)(sintheta (dTheta)/(d theta)) + lambdasin^2theta)^(Theta) + overbrace(1/Phi (d^2Phi)/(dphi^2))^(Phi) = 0
As before, these are equal to a separation constant. Call it
ul((sintheta)/(Theta) (d)/(d theta)(sintheta (dTheta)/(d theta)) + lambdasin^2theta - B = 0)
ul(1/Phi (d^2Phi)/(dphi^2) + B = 0)
3. SOLVING THE PHI PART
Obviously, the
(d^2Phi)/(dphi^2) + BPhi = 0
We assume a solution to this ordinary differential equation with constant coefficients by choosing
Phi = e^(im_lphi)
This gives the auxiliary equation:
-m_l^2 cancel(e^(im_lphi))^(ne 0) + Bcancel(e^(im_lphi))^(ne 0) = 0
The solution is then a linear combination of
Phi(phi) = c_1e^(im_lphi) + c_2e^(-im_lphi)
For one full rotation on the
To eliminate variables, currently,
Phi(phi) = (c_1 + c_2)e^(im_lphi) = Ae^(im_lphi)
From Euler's relation,
e^(i cdot m_l cdot 0) = cos(0m_l) + isin(0m_l) = e^(i cdot m_l cdot 2pi) = cos(2pim_l) + isin(2pim_l)
This only holds true for integer values of
m_l = 0, pm1, pm2, . . .
The unnormalized
color(green)(barul(|" "stackrel(" ")(Phi(phi) prop e^(im_lphi))" "|))
Note that we currently do not know that
(That fact would be determined the traditional way at the end of solving the
4. SOLVING THE THETA PART
As the separation constant is consistent throughout the process,
1/(sintheta) (d)/(d theta)(sintheta (dTheta)/(d theta)) + (lambda - m_l^2/(sin^2theta))Theta = 0
The "trick" for this is to use a function of
P(costheta) := Theta(theta)
x := costheta
This transformation changes the differential operator by the chain rule from
d/(dx)(sin^2theta (dP)/(dx)) + (lambda - m_l^2/(sin^2theta))P = 0
Since
d/(dx)((1 - x^2) (dP)/(dx)) + (lambda - m_l^2/(1 - x^2))P = 0
And after the product rule and chain rule to rewrite in terms of just
ul((1 - x^2) (d^2P)/(dx^2) - 2x(dP)/(dx) + (lambda - m_l^2/(1 - x^2))P = 0)
This is something called an associated Legendre-type differential equation. This is a known equation type with solutions called the Associated Legendre Polynomials
I'm not going to bother with the general formula for the solution, but the first few are:
ul(l = 0: ) " "" "" "" "" "" "ul(l = 1: )
P_(0)^(0) = 1 " "" "" "" "" "P_(1)^(0) = costheta
" "" "" "" "" "" "" "" "P_(1)^(1) = sintheta
ul(l = 2: )
P_(2)^(0) = 1/2 (3cos^2theta - 1)
P_(2)^(1) = 3sinthetacostheta
P_(2)^(2) = 3 - 3cos^2theta
What we further get from this is that
5. PUTTING TOGETHER THE ANGULAR PARTS
Now that we have two solutions
color(green)(barul(|Y_(l)^(m_l) prop P_(l)^(|m_l|)(costheta) e^(im_lphi)|))
You can look at the first few normalized spherical harmonics here:
Y_(0)^(0) = 1/(4pi)^(1//2)" "" "" "" "" "Y_(1)^(0) = (3/(4pi))^(1//2) costheta
Y_(1)^(1) = (3/(8pi))^(1//2) sintheta e^(iphi)" "" "Y_(1)^(-1) = (3/(8pi))^(1//2) sinthetae^(-iphi)
and you can find the rest in your textbook.
6. "SOLVING" THE RADIAL PART
Now that we finally have
(d)/(dr) (r^2 (dR)/(dr)) + (2mu r^2)/(ℏ^2) (E + e^2/(4piepsilon_0r))R - l(l+1)R = 0
Apply the product rule on the first half, factor out
ul((d^2R)/(dr^2) + 2/r (dR)/(dr) + [(2mu)/(ℏ^2)(E + e^2/(4piepsilon_0r)) - (l(l+1))/r^2]R = 0)
The solution we care about here is near the nucleus. Most textbooks don't show any more detail than this...
Here,
color(green)(barul(|R_(nl)(r) prop r^l e^(-r//na_0)L_(n+1)^(2l+1)((2r)/(na_0))|)) where
a_0 = "52.9177 pm" is the bohr radius.
The first few associated Laguerre polynomials are shown here:
n = 1, " "l = 0, " "L_(1)^(1)(x) = -1
n = 2, " "l = 0, " "L_(1)^(2)(x) = -2!(2-x)
" "" "" "" "l = 1, " "L_(3)^(3)(x) = -3!
and the rest should be in your textbook.
7. PUTTING TOGETHER THE OVERALL WAVE FUNCTION
Thus (thus?), from the first substitution,
color(blue)(barul(|psi_(nlm_l)(r,theta,phi) = R_(nl)(r)Y_(l)^(m_l)(theta,phi)|))
we know the form of the hydrogen atom wave functions. (I hope you recognize that none of the above green rectangled equations are normalized.)
I will list the normalized wave functions
(n,l,m_l) = (1,0,0): " "" "psi_(100) = 1/sqrtpi (1/a_0)^(3//2) e^(-sigma)
(n,l,m_l) = (2,0,0): " "" "psi_(200) = 1/sqrt(32pi)(1/a_0)^(3//2) (2-sigma)e^(-sigma//2)
(n,l,m_l) = (2,1,0): " "" "psi_(210) = 1/sqrt(32pi)(1/a_0)^(3//2) sigma e^(-sigma//2) costheta
(n,l,m_l) = (2,1,pm1): " "psi_(21pm1) = 1/sqrt(64pi)(1/a_0)^(3//2) sigmae^(-sigma//2) sintheta e^(pmiphi)
And lastly, if you wish to find the energy based on a particular atomic orbital wave function, evaluate the triple integral in allspace:
E = int_(0)^(2pi) int_(0)^(pi) int_(0)^(oo) psi_(nlm_l)^"*"(r,theta,phi) hatH psi_(nlm_l)(r,theta,phi)r^2dr sinthetad thetadphi
In atomic units, the ground-state energy for the hydrogen atom is