How to make approximation of rotational partition function of diatomic linear molecules?

Question

Using the rigid rotor approximation to the level energies, and such other appropriate assumptions, we can approximate rotational partition function, $Q_{\mathrm{rot}}$, of linear molecules as follows: \begin{eqnarray} Q_{\mathrm{rot}} & = & \sum_{J=0}^{\infty}{(2J+1)\exp\left(-\frac{E_J}{kT}\right)}\\ & \simeq & \sum_{J=0}^{\infty}{(2J+1)\exp\left(-\frac{hB_0J(J+1)}{kT}\right)} \\ \end{eqnarray} In the book, Microwave Molecular Spectra (Gordy & Cook, 1984), Chapter 3, Equation 3.64, the authors made approximation as: $$Q_{rot}\simeq\frac{kT}{hB_0}+\frac{1}{3}+\frac{1}{15}\left(\frac{hB_0}{kT}\right)+\frac{4}{315}\left(\frac{hB_0}{kT}\right)^2$$ Also, there is another approximation made by McDowell, R. S. 1988, Journal of Chemical Physics, 88, 356: $$Q_{rot}\simeq\frac{kT}{hB_0}\exp\left(\frac{hB_0}{3kT}\right)$$

And we can see from the picture that $Q_{\mathrm{rot}}\simeq\frac{kT}{hB_0}+\frac{1}{3}$ is also a good approximation. Then I can guess how to derive McDowell's approximation, but I cannot derive Gordy & Cook's approximation.

What I have got are as follows:

Since $$n\exp(-x(n+1)n)+\int_{(n+1)n}^{(m+1)m}{\exp(-xy)}\,\mathrm dy+(m+1)\exp(-x(m+1)m) \\ \leq\sum_{J=n}^{m}{(2J+1)\exp(-x(J+1)J)} \\ \leq\sum_{i=n^2}^{m^2+2m}{\exp(-xi)}$$ hence, $$\int_{0}^{\infty}{\exp(-xy)}\,\mathrm dy=\frac{1}{x} \\ \leq\sum_{J=0}^{\infty}{(2J+1)\exp(-x(J+1)J)} \\ \leq\sum_{i=0}^{\infty}{\exp(-xi)}=\frac{1}{1-\exp(-x)}$$

i.e., $$\frac{kT}{hB_0}<Q_{\mathrm{rot}}<\frac{1}{1-\exp\left(-\frac{hB_0}{kT}\right)}$$

While the latter one approximates $\frac{1}{1-\exp\left(-\frac{hB_0}{kT}\right)}\approx \frac{kT}{hB_0}+\frac{1}{2}+\frac{1}{12}\frac{hB_0}{kT}-\frac{1}{720}\left(\frac{hB_0}{kT}\right)^3$ when $T$ is rather high.

Is there anyone knowing how to derive Gordy & Cook's approximation? Please tell me and I appreciate it very much.

Besides, another little question, are such approximations still useful as nowadays computers are greatly powerful?

Answer 1

Sorry for the late answer. I just discovered this question.

Preludium

Firstly, I think you might have an error or some non-standard notation in your formula for $Q_{\mathrm{rot}}$. The rotational energy levels are given by

\begin{align} E_{\mathrm{rot}} = \frac{\hbar^2}{2 I} J ( J + 1 ) \ , \end{align}

where $I$ is the moment of inertia. Then the rotational constant $B_{0} = \frac{\hbar^2}{2 I}$ (in units of energy) is introduced, so

\begin{align} E_{\mathrm{rot}} = B_{0} J ( J + 1 ) \ . \end{align}

Thus

\begin{align} Q_{\mathrm{rot}} = \sum_{J=0}^{\infty}{(2J+1)\exp\left(-\frac{B_0 J(J+1)}{kT}\right)} \ , \end{align}

whereas your exponent contains another $h$ which should not be there. At least that's the standard notation I am used to. For convenience I'll also introduce the rotational temperature $\Theta_{\mathrm{rot}} = \frac{B_{0}}{k}$, so that

\begin{align} Q_{\mathrm{rot}} = \sum_{J=0}^{\infty}{(2J+1)\exp\left(-\frac{\Theta_{\mathrm{rot}} J(J+1)}{T}\right)} \ . \end{align}

Main Answer

Ok, now that that's sorted out, let's start the real work. It is actually a pretty standard approximation that is used here: They use the Euler-Maclaurin formula

\begin{align} \sum_{n = a}^{b} f(n) &= \int_{a}^{b} f(n) \, \mathrm{d} n + \frac{1}{2} \bigl( f(b) + f(a) \bigr) + \sum_{i=1}^{\infty} (-1)^{i} \frac{b_{2i}}{(2i)!} \Bigl( f^{(2i - 1)}(a) - f^{(2i - 1)}(b) \Bigr) \end{align}

where $f^{(k)}(a)$ is the $k$th derivative of $f$ evaluated at $a$. The $b_{2i}$s are the absolute values of the even Bernoulli numbers, $b_{2} = \frac{1}{6}$, $b_{4} = \frac{1}{30}$, $b_{6} = \frac{1}{42}$, $...$. It provides a way to approximate the summation of a function over a discrete variable by an analogous integration over a continuous variable.

In the special case at hand,

\begin{align} f(J) = (2J+1)\exp\left(-\frac{\Theta_{\mathrm{rot}} J(J+1)}{T}\right) \ , \end{align}

where the summation limits are $a=0$ and $b=\infty$ and where $f(\infty) = f^{\prime}(\infty) = \ldots = 0$ this formula boils down to

\begin{align} Q_{\mathrm{rot}} = \sum_{J = 0}^{\infty} f(J) &= \int_{0}^{\infty} f(J) \, \mathrm{d} J + \frac{1}{2} f(0) + \sum_{i=1}^{\infty} (-1)^{i} \frac{b_{2i}}{(2i)!} f^{(2i - 1)}(0) \end{align}

The first term is what is usually called the classical or the high-temperature limit. It can be easily calculated by making the substitution $K = J(J + 1)$ $\Rightarrow$ $\mathrm{d}J = \frac{\mathrm{d}K}{2J + 1}$. So

\begin{align} \int_{0}^{\infty} f(J) \, \mathrm{d} J &= \int_{0}^{\infty} (2J+1)\exp\left(-\frac{\Theta_{\mathrm{rot}} J(J+1)}{T}\right) \, \mathrm{d} J \\ &= \int_{0}^{\infty} \exp\left(-\frac{\Theta_{\mathrm{rot}} K}{T}\right) \, \mathrm{d} K \\ &= \frac{T}{\Theta_{\mathrm{rot}}} \end{align}

Now, only the derivatives $f^{(k)}(0)$ need to be sorted out. Those that contain terms below the third order in $\frac{\Theta_{\mathrm{rot}}}{T}$ can be calculated to be:

\begin{align} f(0) &= 1 \\ f^{(1)}(0) &= 2 - \frac{\Theta_{\mathrm{rot}}}{T} \\ f^{(3)}(0) &= - 12 \frac{\Theta_{\mathrm{rot}}}{T} + 12 \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{2} - \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{3} \\ f^{(5)}(0) &= 120 \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{2} - 180 \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{3} + 30 \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{4} - \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{5} \ . \end{align}

Doing this by hand is very tedious. Best you use something like Mathematica. Plugging all this into the Euler-MacLaurin formula and dropping all terms above second order in $\frac{\Theta_{\mathrm{rot}}}{T}$ you get exactly the result of Gordy & Cook

\begin{align} Q_{\mathrm{rot}} &\approx \frac{T}{\Theta_{\mathrm{rot}}} + \frac{1}{2} \cdot 1 + (-1) \frac{1}{6 \cdot 2!} \left( 2 - \frac{\Theta_{\mathrm{rot}}}{T} \right) \\ &\quad + (-1)^{2} \frac{1}{30 \cdot 4!} \left( - 12 \frac{\Theta_{\mathrm{rot}}}{T} + 12 \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{2} - \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{3} \right) \\ &\quad + (-1)^{3} \frac{1}{42 \cdot 6!} \left(120 \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{2} - 180 \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{3} + 30 \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{4} - \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{5} \right) \\ &= \frac{T}{\Theta_{\mathrm{rot}}} + \frac{1}{2} - \frac{1}{6} + \frac{1}{12} \frac{\Theta_{\mathrm{rot}}}{T} - \frac{1}{60} \frac{\Theta_{\mathrm{rot}}}{T} \\ &\quad + \frac{1}{60} \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{2} - \frac{1}{252} \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{2} + \mathcal{O}\left( \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{3} \right) \\ &= \frac{T}{\Theta_{\mathrm{rot}}} + \frac{1}{3} + \frac{1}{15} \frac{\Theta_{\mathrm{rot}}}{T} + \frac{4}{315} \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{2} + \mathcal{O}\left( \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{3} \right) \end{align}

As a side note: If you have a homonuclear diatomic molecule you need to introduce the so called symmetry number $\sigma$ which represents the number of indistinguishable orientations a molecule can have. It originates from symmetry constrains: for homonuclear diatomic molecules there is a two-fold axis of symmetry perpendicular to the internuclear axis which gives rise to two indistinguishable orientations, so the above formula for $Q_{\mathrm{rot}}$ overcounts the number of available quantum states by a factor of two because of the inherent indistinguishability of the nuclear pair. Thus, $Q_{\mathrm{rot}}$ needs to be divided by $\sigma = 2$. For heteronuclear diatomic molecules you simply have $\sigma = 1$. With the introduction of $\sigma$ the equation for $Q_{\mathrm{rot}}$ becomes:

\begin{align} Q_{\mathrm{rot}} &\approx \frac{1}{\sigma} \left( \frac{T}{\Theta_{\mathrm{rot}}} + \frac{1}{3} + \frac{1}{15} \frac{\Theta_{\mathrm{rot}}}{T} + \frac{4}{315} \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{2} + \mathcal{O}\left( \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{3} \right) \right) \\ &= \frac{T}{\sigma \Theta_{\mathrm{rot}}} \left(1 + \frac{1}{3} \frac{\Theta_{\mathrm{rot}}}{T} + \frac{1}{15} \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{2} + \frac{4}{315} \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{3} + \mathcal{O}\left( \left(\frac{\Theta_{\mathrm{rot}}}{T}\right)^{4} \right) \right) \end{align}