Problem Set 1: Semiclassical quantization of molecular vibrations

As an example combining several basic mathematical operations, we consider the problem of describing a diatomic molecule such as O_{2}, which consists of two nuclei bound together by the electrons that orbit about them. Since the nuclei are much heavier than the electrons we can assume that the latter move fast enough to readjust instantaneously to the changing position of the nuclei (Born-Oppenheimer approximation). The problem is therefore reduced to one in which the motion of the two nuclei is governed by a potential V, depending only upon r, the distance between then:

V (r) = 4V_{O} [(\frac{a}{r})^{12} − (\frac{a}{r})^{6}] ,

The vibrational states of relative motion, with energies E_{n}, are described by the bound state solutions, ψ_{n}(r), of a one-dimensional Schroedinger equation,

[−\frac{h̄^{2}}{2m} \frac{d^{2}}{dr^{2}} + V (r)] ψ_{n} = E_{n}ψ_{n}.

Here, m is the reduced mass of the two nuclei.

Our goal in this example is to find the energies E_{n}, given a particular potential. This can be done exactly by solving the differential eigenvalue equation (1.17). However, the great mass of the nuclei implies that their motion is nearly classical, so that approximate values of the vibrational energies E_{n} can be obtained by considering the classical motion of the nuclei V and then applying “quantization rules” to determine the energies.

Confined classical motion of the internuclear separation in the potential V (r) can occur for energies −V_{0} < E < 0. The distance between the nuclei oscillates periodically (but not necessarily harmonically) between inner and outer turning points, r_{in} and r_{out}. During these oscillations, energy is exchanged between the kinetic energy of relative motion and the potential energy such that the total energy,

E = \frac{p^{2}}{2m} + V (r),

is a constant (p is the relative momentum of the nuclei). We can therefore think of the oscillations at any given energy as defining a closed trajectory in phase space (coordinates r and p) along which Eq. (1.18) is satisfied, as shown in the lower portion of Figure 1.3. An explicit equation for this trajectory can be obtained by solving (1.18) for p:

p(r) = ±[2m(E − V (r))]^{\frac{1}{2}} .

The quantization rules state that, at the allowed energies E_{n}, the action is a half-integral multiple of 2π. Thus, upon using (1.19) and recalling that the oscillation passes through each value of r twice (once with positive p and once with negative p), we have

S(E_{n} = 2 (\frac{2m}{h̄^{2}})^{\frac{1}{2}} ∫_{r_{in}}^{r_{out}} [E_{n} − V (r)]^{\frac{1}{2}} dr = (n + \frac{1}{2})2π

where n is a nonnegative integer. At the limits of the integral, the turning points r_{in} and r_{out}, the integrand vanishes.

It is useful to define the dimensionless quantities

ϵ = \frac{E}{V_{0}} , x = \frac{r}{a}, γ = (\frac{2ma^{2}V_{O}}{h̄} )^{\frac{1}{2}} ,

The quantity γ is a dimensionless measure of the quantum nature of the problem. In the classical limit (h̄ small or m large), γ becomes large. By knowing the moment of inertia of the molecule (from the energies of its rotational motion) and the dissociation energy (energy required to separate the molecule into its two constituent atoms), it is possible to determine from observation the parameters a and V_{0} and hence the quantity γ. For the H_{2} molecule, γ = 21.7, while for the HD molecule, γ = 24.8 (only m, but not V_{0}, changes when one of the protons is replaced by a deuteron), and for the much heavier O_{2} molecule made of two ^{16}O nuclei γ = 150. These rather large values indicate that a semiclassical approximation is a valid description of the vibrational motion.

The goal of this problem set is to find all the (semiclassical) bound iiiistates for a given value of γ.

Exercise 1.7 One of the most important aspects of using a computer as a tool to do physics is knowing when to have confidence that the program is giving the correct answers. In this regard, an essential test is the detailed quantitative comparison of results with what is known in analytically soluble situations. Modify the code fo use a parabolic potential for which the Bohr-Sommerfeld quantization gives the exact eigenvalues of the Schroedinger equation: a series of equally-spaced energies, with the lowest being one-half of the level spacing above the minimum of the potential. For several values of γ, compare the numerical results for the case with what you obtain by solving analytically.

Exercise 1.8 Another important test of a working code is to compare its results with what is expected on the basis of physical intuition. Restore the code to use the Lennard-Jones potential and run it for γ = 50. Note that, as in the case of the purely parabolic potential discussed in the previous exercise, the first excited state is roughly three times as high above the bottom of the well as is the ground state and that the spacings between the few lowest states are roughly constant. This is because the Lennard-Jones potential is roughly parabolic about its minimum. By calculating the second derivative of V at the minimum, find the “spring constant” and show that the frequency of small-amplitude motion is expected to be

\frac{h̄ω}{V_{0}} = \frac{6 × 2^{5/6}}{γ} ≈ \frac{10.691}{γ}

Verify that this is consistent withthe numerical results and explore this agreement for different values of γ. Can you understand why the higher energies are ore densely spaced than the lower ones by comparing the Lennard-Jones potential with its parabolic approximation?

Exercise 1.9 Invariance of results under changes in the numerical algorithms or their parameters can give additional confidence in a calculation. Change the tolerances for the turning point and energy searches (line 120) or the number of Simpson’s rule points (line 130) and observe the effects on the results. Note that because of the way in which the expected number of bound states is calculated (lines 1190-1200), this quantity can change if the energy tolerance is varied.

Exercise 1.10 Replace the searches for the inner and outer turning points by the Newton-Raphson method or the secant method. Replace the Simpson’s rule quadrature for s by a higher-order formula (Eqs. (1.13a) or (1.13b)) and observe the improvement.

Exercise l.ll For the H_{2} molecule, observations show that the depth of the potential is V_{0} = 4.747eV and the location of the potential minimum is r_{min} = 0.74166Å. These two quantities, together with Eq. (1.23), imply a vibrational frequency of

h̄ω = 0.492V_{0} = 2.339eV,

more than four times larger than the experimentally observed energy difference between the ground and first vibrational state, 0.515 eV. The Lennard-Jones shape is therefore not a very good description of the potential of the H_{2} molecule. Another defect is that it predicts 6 bound states, while 15 are known to exist. A better analytic form of the potential, with more parameters, is required to reproduce simultaneously the depth and location of the minimum, the frequency of small amplitude vibrations about it, and the total number of bound states. One such form is the Morse potential

V (r) = V_{0}[(1 − e^{−β(r−r_{mm})})^{2} − 1],

The Morse potential has a minimum at the expected location and th parameter β can be adjusted to fit the curvature of the minimum to the observed excitation energy of the first vibrational state. Find the value of β appropriate for the H_{2} molecule, modify the program above to use the Morse potential, and calculate the spectrum of vibrational states.