Physical interpretation of statistical mechanics What is the physical interpretation of
$$ \langle p_\theta^2\rangle = \dfrac{8\pi^2 I^2}{\beta^2} \tag{1}$$
and
$$\langle p_\varphi^2 \rangle = \dfrac{16\pi^2 I^2}{3\beta^2} \tag{2}$$
where the Hamiltonian (for $1$ particle) that describe the rotation of a diatomic molecule is
$$ H_1 (\theta, \varphi, p_\theta, p_\varphi)=\dfrac{1}{2I} \left(p_\theta^2 +\dfrac{p_\varphi^2}{\sin^2 \theta} \right)
\tag3 $$
where $I$ is the inertia momentum of the molecule ($I=mr_0^2$). Are the units correct? What do the numerical values mean?
Note that
$\beta =1/k_BT$
as usual in statistical mechanics.
For the definition of the expected value in the canonical ensemble, I've used
$$
\langle p_\theta^2 \rangle =\dfrac{1}{h^2}\dfrac{1}{Z(\beta,\alpha,1)}\int_0^{2\pi} \text{d}\varphi \int_0^\pi \text{d}\theta \int_{-\infty}^\infty \text{d}p_\theta \int_{-\infty}^\infty \text{d}p_\varphi \cdot p_\theta^2 \text{ e}^{-\beta H_1},$$
where
$$ Z(\beta,\alpha,1)= \int_0^{2\pi} \text{d}\varphi \int_0^\pi \text{d}\theta \int_{-\infty}^\infty \text{d}p_\theta \int_{-\infty}^\infty \text{d}p_\varphi  \text{ e}^{-\beta H_1}.$$
I don't know if this is right or not.
For example:
[1] Statistical Mechanics, Ryogo Kubo, p. 200, nº 3.
 A: If an particle has moment of inertia $I$ and angular momentum $L=I\omega$ about some axis, its rotational kinetic energy is $K=\frac12 I\omega^2 = \frac{L^2}{2I}$.  So your source seems to be using $p_\theta$ to mean something like "angular momentum associated with motion in the $\hat\theta$ direction," which isn't too surprising.  If you use the usual "physics" spherical polar coordinates where the coordinate $\theta$ is the angle with the $z$-axis and the coordinate $\phi$ is the angle from the $x$-$z$ plane, then a given angular momentum $p_\phi$ about the $z$-axis is associated with more energy when the moment of inertia is reduced by having the rotor near the $z$-axis.
We expect an angular momentum to have units of joule-seconds, which is consistent with your definition of the Hamiltonian:
$$
[p^2] = [IH] = [\rm {kg\,m^2}][J] = [J\cdot s]^2
$$
However, your expectation values for the angular momenta have units equivalent to $[IH]^2$, which is $[p^4]$.  That's because you haven't normalized your expectation value by the "zeroth moment" of your probability distribution.  Instead you should have something more like
$$
\left<p^2\right>
=
\frac{
  \int d\theta\ d\phi\ dp_\theta\ dp_\phi\ p^2 \exp{-\beta H}
}{
  \int d\theta\ d\phi\ dp_\theta\ dp_\phi\ \exp{-\beta H}
}
$$
Remember when you do these sorts of integrals that the differentials (here, $dp$) also carry their own units, which must be accounted for in the result.
My naïve expectation is that you should find
$$
\frac{\left<p^2\right>}{2I} = \frac12 kT
$$
for both angular momenta
since, in thermal equilibrium, the average energy in each quadratic degree of freedom is $\frac12 kT$.  I don't know whether you're doing the right integral yet or not, in that case.
A: By definition we have
$$ \langle a^2 \rangle = \dfrac{\int_0^{2\pi} \text{d}\phi \int_0^\pi \text{d}\theta \int_{-\infty}^\infty \text{d}p_\theta \int_{-\infty}^\infty \text{d}p_\phi \cdot a^2 \text{ e}^{-\beta H}}{\int_0^{2\pi} \text{d}\phi \int_0^\pi \text{d}\theta \int_{-\infty}^\infty \text{d}p_\theta \int_{-\infty}^\infty \text{d}p_\phi  \text{ e}^{-\beta H}} \tag{S.1}$$
but, for $a=p_\theta$, this expression can be simplified
$$\langle p_\theta^2 \rangle = \dfrac{\int_{-\infty}^\infty p_\theta^2 \text{ e}^{-\beta p_\theta^2/(2I)} \text{ d}p_\theta}{\int_{-\infty}^\infty \text{e}^{-\beta p_\theta^2/(2I)}\text{ d}p_\theta}. \tag{S.2}$$
For the Hamiltonian $(3)$ we obtain
$$ Z(\beta,\alpha,1)=8\pi^2 \dfrac{I}{\beta} \tag{S.3}$$
and using this expression
$$\int_\mathbb{R} \text{d}x \text{ e}^{-ax^2} = 2a\int_\mathbb{R} \text{d}x \text{ } x^2 \text{ e}^{-ax^2}=\sqrt{\dfrac{\pi}{a}}\tag{S.4}$$
into $(S.2)$ we have
$$ \boxed{\langle p_\theta^2 \rangle =\dfrac{I}{\beta}= k_B TI}. \tag{S.5}$$
and, for $a=p_\phi$, the expression $(S.1)$ turns into
$$\langle p_\phi^2 \rangle = \dfrac{\int_{0}^\pi \text{d}\theta \int_{-\infty}^\infty \text{d}p_\phi \text{ } p_\phi^2 \text{ e}^{-\beta p_\phi^2/(2I\sin^2 \theta)}}{\int_{0}^\pi \text{d}\theta \int_{-\infty}^\infty \text{d}p_\phi \text{ }\text{e}^{-\beta p_\phi^2/(2I\sin^2 \theta)}}. \tag{S.6}$$
We obtain
$$ \boxed{\langle p_\phi^2\rangle = \dfrac{2}{3}\dfrac{I}{\beta}=\dfrac{2}{3} k_B TI}. \tag{S.7}$$
Equations $(S.5)-(S.7)$ are $\textbf{the dispersion of each (conjugate) angular momentum}$ because
$$ \langle p_\theta \rangle = \langle p_\phi \rangle =0.$$
Therefore,
$$ \sigma_{p_\theta}^2 = \langle p_\theta^2 \rangle= k_B TI > \sigma_{p_\phi}^2 = \langle p_\phi^2 \rangle = \dfrac{2}{3}k_B TI = \dfrac{2}{3}\langle p_\theta^2 \rangle \geq 0.$$
The physical interpretation is that the $\textbf{dispersion}$ of $\theta$-$\textbf{angular momentum}$  $p_\theta$ is bigger than the $\textbf{dispersion}$ of $\phi$-$\textbf{angular momentum} $ $p_\phi$ in a factor of $3/2$:
$$ \boxed{\sigma_{p_\theta}^2 = \dfrac{3}{2}\sigma_{p_\phi}^2}. \tag{S.8}$$
According to @rob the energy of one atom is
$$ E_1 = \langle H_1 \rangle = \dfrac{1}{2I} \left( k_B TI + k_B TI\right)=k_B T \tag{S.9}$$ and for the molecule
$$ E_{1+1}=\langle 2H_1 \rangle = 2k_B T.\tag{S.10}$$
