Rigorous derivation of the mean free path in a gas

The following is a self answer concerning the rigorous derivation of the relative velocity in a Maxwellian gas. It is here for comparison against other answers.

Assumptions are

• Two particle collisions are the most probable ones, and collisions involving more than this, contribute to the mean relative velocity by a negligible amount.
• The mean relative velocity is strictly positive.

I think that this derivation differs from the one by @Thorondor because within it, $\mathbf{v_1}^2+\mathbf{v_2}^2\neq\mathbf{v_r}^2$. Which is to say that the dot product of $\mathbf{v_1}$, and $\mathbf{v_2}$ is not assumed to be zero. This forces us to integrate over all possible angles between the velocities, which may account for the ~0.07 difference. It is also why the derivation is significantly longer.

The mean free path of an atom/molecule in a Maxwellian gas, depends upon the average relative velocity of each particle to one another. In order to obtain it, we first find the magnitude of velocity $\mathbf{v}_1$ relative to all other particles moving with $\mathbf{v}_2$. This expression is then averaged for all values of $\mathbf{v}_2$, from zero to infinity, producing a mean relative speed given $\mathbf{v}_1$ exists.

Following modus ponens, the mean is multiplied by the probability that an atom/molecule has velocity $\mathbf{v}_1$, and then averaged for all $\mathbf{v}_1$ from zero to infinity.

Consider two velocities $\mathbf{v}_1$ and $\mathbf{v}_2$ inclined at an angle $\theta$, there relative velocity will be described by

\begin{align} (v_r)_{12}&=\sqrt{v_1^2+v_2^2-2v_1v_2cos(\theta)}. \end{align}

All directions are equally probable for $\mathbf{v}_2$. To calculate the average of $(v_r)_{12}$, we multiply it by the probability that it lies within some solid angle $d\Omega$.

\begin{align} \langle(v_r)_{12}\rangle&=\int_{\Omega}\frac{d\Omega}{4\pi}(v_r)_{12}\\ &=\int_{\Omega}\frac{2\pi sin(\theta)d\theta}{4\pi}(v_r)_{12}\\ &=\frac{1}{2}\int_{0}^{\pi}sin(\theta)d\theta(v_1^2+v_2^2-2v_1v_2cos(\theta))^\frac{1}{2} \end{align}

Consider a change of basis, where $cos(\theta)=x$, so $sin(\theta)d\theta=-dx$. The limits of integration will change to $cos(0)=1$ and $cos(\pi)=-1$.

\begin{align} \therefore~\langle(v_r)_{12}\rangle&=\frac{1}{2}\int_{1}^{-1}dx(v_1^2+v_2^2-2v_1v_2x)^\frac{1}{2}\\ &=\frac{1}{3}\frac{(v_1^2+v_2^2-2v_1v_2x)^\frac{3}{2}}{2v_1v_2}\Biggr|_{1}^{-1}\\ &=\frac{1}{6v_1v_2}\left((v_1^2+v_2^2+2v_1v_2)^\frac{3}{2}-(v_1^2+v_2^2-2v_1v_2)^\frac{3}{2}\right)\\ &=\frac{1}{6v_1v_2}\left((v_1+v_2)^\frac{3}{2}(v_1+v_2)^\frac{3}{2}+(v_1-v_2)^\frac{3}{2}(v_1-v_2)^\frac{3}{2}\right)\\ &=\frac{1}{6v_1v_2}\left((v_1+v_2)^3+|v_1-v_2|^3\right)\\ \end{align}

We take the magnitude of $v_1-v_2$ because $\langle(v_r)_{12}\rangle$ should always be positive, and therefore $\left((v_1-v_2)^2\right)^\frac{3}{2}$ must also be. This splits the average into two parts.

\begin{align} \langle(v_r)_{12}\rangle&=\frac{1}{6v_1v_2}\left((v_1+v_2)^3-(v_1-v_2)^3\right),~when~v_1\geq v_2,\\ &=\frac{1}{6v_1v_2}\left[2v_2^3+6v_1^2v_2\right]\\ &=\frac{3v_1^2+v_2^2}{3v_1}.\\ \langle(v_r)_{12}\rangle&=\frac{1}{6v_1v_2}\left((v_1+v_2)^3+(v_1-v_2)^3\right),~when~v_2>v_1,\\ &=\frac{3v_2^2+v_1^2}{3v_2}. \end{align}

The average relative velocity of an atom/molecule possessing magnitude and direction $\mathbf{v_1}$, with respect to all other particles moving with $\mathbf{v_2}$, lying within speeds of zero to infinity is then

\begin{align} \langle(v_r)_1\rangle&=\int_{0}^{\infty}P(v_2)\langle(v_r)_{12}\rangle dv_2,\\ P(v_2)dv_2&=4\pi\left(\frac{m}{2\pi k_BT}\right)^\frac{3}{2}e^{-\frac{mv_2^2}{2k_BT}}v_2^2dv_2. \end{align}

Note that $\langle(v_r)_1\rangle$ has two distinct forms depending on the difference between speeds. Because of this we break the integral into two parts.

\begin{align} \langle(v_r)_1\rangle=4\pi\left(\frac{m}{2\pi k_BT}\right)^\frac{3}{2}&\left[\int_{0}^{v_1}dv_2e^{-\frac{mv_2^2}{2k_BT}}\left(\frac{3v_1^2+v_2^2}{3v_1}\right)v_2^2\right.\\ &+\left.\int_{v_1}^{\infty}dv_2e^{-\frac{mv_2^2}{2k_BT}}\left(\frac{3v_2^2+v_1^2}{3v_2}\right)v_2^2\right] \end{align}

To arrive at the relative velocity of any atom with another, we form a product of the above, and the probability a particle has velocity $v_1+dv_1$.

\begin{align} \langle v_r \rangle&=\int_{0}^{\infty}P(v_1)\langle(v_r)_1\rangle dv_1 \end{align}

\begin{align} \langle v_r \rangle=\left(4\pi\left[\frac{m}{2\pi k_BT}\right]^\frac{3}{2}\right)^2\int_{0}^{\infty}dv_1e^{-\frac{mv_1^2}{2k_BT}}v_1^2&\left[\int_{0}^{v_1}dv_2e^{-\frac{mv_2^2}{2k_BT}}\left(\frac{3v_1^2+v_2^2}{3v_1}\right)v_2^2\right.\\ &+\left.\int_{v_1}^{\infty}dv_2e^{-\frac{mv_2^2}{2k_BT}}\left(\frac{3v_2^2+v_1^2}{3v_2}\right)v_2^2\right] \end{align}

which may be re-written as

\begin{align} \langle v_r\rangle=\left(\frac{2m}{k_BT}\right)^3&\left[\int_{0}^{\infty}dv_1e^{-\frac{mv_1^2}{2k_BT}}v_1^2\int_{0}^{v_1}dv_2e^{-\frac{mv_2^2}{2k_BT}}\left(\frac{3v_1^2+v_2^2}{3v_1}\right)v_2^2\right.\\ &+\left.\int_{0}^{\infty}dv_1e^{-\frac{mv_1^2}{2k_BT}}v_1^2\int_{v_1}^{\infty}dv_2e^{-\frac{mv_2^2}{2k_BT}}\left(\frac{3v_2^2+v_1^2}{3v_2}\right)v_2^2\right], \end{align}

or

\begin{align} \langle v_r\rangle&=\left(\frac{2m}{k_BT}\right)^3\left[I_1+I_2\right], \end{align}

where $I_1$ stands for the first integral in the square brackets of the above, and $I_2$ the second one.

\begin{align} I_1&=\int_{0}^{\infty}dv_1e^{-\frac{mv_1^2}{2k_BT}}v_1^2\int_{0}^{v_1}dv_2e^{-\frac{mv_2^2}{2k_BT}}\left(\frac{3v_1^2+v_2^2}{3v_1}\right)v_2^2\\ I_2&=\int_{0}^{\infty}dv_1e^{-\frac{mv_1^2}{2k_BT}}v_1^2\int_{v_1}^{\infty}dv_2e^{-\frac{mv_2^2}{2k_BT}}\left(\frac{3v_2^2+v_1^2}{3v_2}\right)v_2^2 \end{align}

We wish to show that these two definite integrals are equivalent. To do so, a substitution, then change in the order of integration following Fubini's theorem is employed.

Consider the piecewise function,

\begin{align} \mathbb{I}_{\{v_2\leq v_1\}}&=\begin{cases} 1&v_1\leq v_2\\ 0&v_1>v_2\\ \end{cases} \end{align}

which is zero outside the limit of the second integral of $I_1$. Using it, we can rewrite equation $I_1$ as

\begin{align} I_1&=\int_{0}^{\infty}dv_1e^{-\frac{mv_1^2}{2k_BT}}v_1^2\int_{0}^{\infty}dv_2\mathbb{I}_{\{v_2\leq v_1\}}e^{-\frac{mv_2^2}{2k_BT}}\left(\frac{3v_1^2+v_2^2}{3v_1}\right)v_2^2.\\ \end{align}

We can combine these two integrals with a substitution $a=v_1$, allowing us to write the following

\begin{align} I_1&=\int_{0}^{\infty}\int_{0}^{\infty}\mathbb{I}_{\{v_2\leq v_1\}}e^{-\frac{ma^2}{2k_BT}}a^2e^{-\frac{mv_2^2}{2k_BT}}\left(\frac{3v_1^2+v_2^2}{3v_1}\right)v_2^2dv_2da.\\ \end{align}

Since the integrals limits are constant, we can apply Fubini's theorem exchange there order.

\begin{align} I_1&=\int_{0}^{\infty}\int_{0}^{\infty}\mathbb{I}_{\{v_2\leq v_1\}}e^{-\frac{ma^2}{2k_BT}}a^2e^{-\frac{mv_2^2}{2k_BT}}\left(\frac{3v_1^2+v_2^2}{3v_1}\right)v_2^2dadv_2\\ &=\int_{0}^{\infty}dv_2e^{-\frac{mv_2^2}{2k_BT}}v_2^2\int_{0}^{\infty}dv_1\mathbb{I}_{\{v_2\leq v_1\}}e^{-\frac{mv_1^2}{2k_BT}}\left(\frac{3v_1^2+v_2^2}{3v_1}\right)v_1^2\\ &=\int_{0}^{\infty}dv_2e^{-\frac{mv_2^2}{2k_BT}}v_2^2\int_{v_2}^{\infty}dv_1e^{-\frac{mv_1^2}{2k_BT}}\left(\frac{3v_1^2+v_2^2}{3v_1}\right)v_1^2\\ \end{align}

Then, because the integral is definite, we can interchange $v_1$ with $v_2$, and obtain an equivalent volume.

\begin{align} I_1&=\int_{0}^{\infty}dv_1e^{-\frac{mv_1^2}{2k_BT}}v_1^2\int_{v_1}^{\infty}dv_2e^{-\frac{mv_2^2}{2k_BT}}\left(\frac{3v_2^2+v_1^2}{3v_2}\right)v_2^2\\ &=I_2 \end{align}

Hence, the average relative velocity can be rewritten as

\begin{align} \langle v_r\rangle&=2\left(\frac{2m}{k_BT}\right)^3I_2. \end{align}

To evaluate $I_2$, first consider the integral

\begin{align} I_0&=\int_{v_1}^{\infty}dv_2e^{-\frac{mv_2^2}{2k_BT}}\left(\frac{3v_2^2+v_1^2}{3v_2}\right)v_2^2. \end{align}

We can break this into parts,

\begin{align} I_0&=\int_{v_1}^{\infty}dv_2e^{-\frac{mv_2^2}{2k_BT}}v_2^3+\frac{v_1^2}{3}\int_{v_1}^{\infty}dv_2e^{-\frac{mv_2^2}{2k_BT}}v_2, \end{align}

then consider the substitution $\frac{m}{2k_BT}v_2^2=x$, so $dv_2=\frac{k_BT}{mv_2}$, and the lower limit of integration changes to $\frac{mv_1^2}{2k_BT}$.

\begin{align} I_0&=\int_{\frac{mv_1^2}{2k_BT}}^{\infty}dv_2e^{-\frac{mv_2^2}{2k_BT}}v_2^3+\frac{v_1^2}{3}\int_{\frac{mv_1^2}{2k_BT}}^{\infty}dv_2e^{-\frac{mv_2^2}{2k_BT}}v_2\\ &=2\left(\frac{k_BT}{m}\right)^2\int_{\frac{mv_1^2}{2k_BT}}^{\infty}e^{-x}xdx+\frac{v_1^2}{3}\frac{k_BT}{m}\int_{\frac{mv_1^2}{2k_BT}}^{\infty}e^{-x}dx \end{align}

Since

\begin{align} \int_{a}^{b}xe^{-x}dx&=-(x+1)e^{-x}\Biggr|_{a}^{b}, \end{align}

$I_0$ simplifies to

\begin{align} I_0&=2\left(\frac{k_BT}{m}\right)^2\left[\frac{m}{2k_BT}v_1^2+1\right]e^{\frac{-mv_1^2}{2k_BT}}\\ &=e^{-\frac{mv_1^2}{2k_BT}}\left[\frac{4}{3}\left(\frac{k_BT}{m}\right)v_1^2+2\left(\frac{k_BT}{m}\right)^2\right] \end{align}

Substituting this into $I_2$ yields

\begin{align} I_2=&\int_{0}^{\infty}dv_1e^{-\frac{mv_1^2}{2k_BT}}v_1^2\left[\frac{4}{3}\left(\frac{k_BT}{m}\right)v_1^2+2\left(\frac{k_BT}{m}\right)^2\right] \\=&2\left(\frac{k_BT}{m}\right)^2\int_{0}^{\infty}dv_1e^{-\frac{mv_1^2}{2k_BT}}v_1^2\\ &+\frac{4}{3}\frac{k_BT}{m}\int_{0}^{\infty}dv_1e^{-\frac{mv_1^2}{2k_BT}}v_1^4\\ \end{align}

$I_2$ is actually a standard integral of the form,

\begin{align} \int_{0}^{\infty}e^{-ax^2}x^n&=\frac{1}{2a^{\frac{n+1}{2}}}\Gamma\left(\frac{n+1}{2}\right), \end{align}

where $\Gamma$ is the gamma function, a generalization of the factorial. $I_2$ now simplifies to

\begin{align} I_2&=\left(\frac{k_BT}{m}\right)^{\frac{7}{2}}\left[\Gamma\left(\frac{3}{2}\right)+\frac{2}{3}\Gamma\left(\frac{5}{2}\right)\right]\\ &=\left(\frac{k_BT}{m}\right)^{\frac{7}{2}}\left[\frac{\sqrt{\pi}}{2}+\frac{2}{3}\frac{3\sqrt{\pi}}{4}\right]\\ &=\left(\frac{k_BT}{m}\right)^{\frac{7}{2}}\sqrt{\pi} \end{align}

Hence, the average relative speed of any Maxwellian gas particle with respect to any other is

\begin{align} \langle v_r\rangle&=2\left[4\pi\left(\frac{m}{2\pi k_BT}\right)^\frac{3}{2}\right]^2\left(\frac{k_BT}{m}\right)^\frac{7}{2}\sqrt{\pi},\\ &=\frac{2\sqrt{\pi}\cdot16\pi^2}{8\pi^3}\left(\frac{k_BT}{m}\right)^{-3}\left(\frac{k_BT}{m}\right)^\frac{7}{2},\\ &=\frac{4}{\sqrt{\pi}}\sqrt{\frac{k_BT}{m}},\\ &=\sqrt{\frac{16}{\pi}\frac{k_BT}{m}}. \end{align}

And the ratio of mean velocity $\langle v_r\rangle$ to mean relative velocity $\langle v\rangle$ is

\begin{align} \frac{\langle v\rangle}{\langle v_r\rangle}&=\frac{\sqrt{\frac{8}{\pi}\frac{k_BT}{m}}}{\sqrt{\frac{16}{\pi}\frac{k_BT}{m}}}\\ &=\frac{1}{\sqrt{2}} \end{align}

And thus the mean free path $\ell$ is

\begin{align} \ell=\frac{1}{n\sigma\sqrt{2}} \end{align}

Edit: This answer is wrong. I made an error that is equivalent to implicitly assuming that the particles' velocities are always perpendicular (see below). The question author has asked me to leave it up as it is interesting that the final result is so nearly correct despite the mistake.

The Maxwell-Boltzmann distribution is

$$B(\mathbf{v}) = \left( \frac{m}{2\pi k T} \right)^{3/2} e^{-\frac{m|\mathbf{v}|^2}{2kT}}$$

This function has a few nice features that make the following calculations easier. First, it's normalized, so integrals such as $\int_{\mathbb{R}^3} B(\mathbf{v}) \, d^3\mathbf{v}$ will always evaluate to 1. Second, since it is a radial function, we'll be able to take advantage of the trick

$$\int_{\mathbb{R}^n} f(|\mathbf{x}|) \, d^n\mathbf{x} = \omega_{n-1} \int_0^{\infty}f(r) r^{n-1} \,dr$$

where $\omega_{n-1}$ is the area of an $n$-sphere of radius 1.

All right, there's no way to fully answer your question without doing a few ugly integrals, so let's dive in. The mean speed of a particle is

\begin{align} \langle v \rangle &= \frac{\int_{\mathbb{R}^3} v \, B(\mathbf{v}) \, d^3\mathbf{v}}{\int_{\mathbb{R}^3} B(\mathbf{v}) \, d^3\mathbf{v}} \\ &= \int_{\mathbb{R}^3} v \, B(\mathbf{v}) \, d^3\mathbf{v} \\ &= \left( \frac{m}{2\pi k T} \right)^{3/2} \int_{\mathbb{R}^3} v \, e^{-\frac{m|\mathbf{v}|^2}{2kT}} \, d^3\mathbf{v} \\ &= \left( \frac{m}{2\pi k T} \right)^{3/2} \int_0^{\infty} v \, e^{-\frac{mv^2}{2kT}} \, 4\pi v^2 dv \\ &= \left( \frac{m}{2kT} \right)^{3/2} \frac{4}{\sqrt{\pi}} \int_0^{\infty} v^3 \, e^{-\frac{mv^2}{2kT}} \, dv \\ &= \left( \frac{m}{2kT} \right)^{3/2} 2\pi \left( \frac{m}{2kT} \right)^{-2} \\ &= \left( \frac{m}{2kT} \right)^{3/2} \frac{2}{\sqrt{\pi}} \left( \frac{m}{2kT} \right)^{-2} \\ &= \frac{2}{\sqrt{\pi}} \left( \frac{2kT}{m} \right)^{1/2} \end{align}

[Edit: if we make the approximation that the particles' velocities are perpendicular,] the mean relative speed of two particles is

\begin{align} \langle v_r \rangle &\approx \frac{\int_{\mathbb{R}^3} \int_{\mathbb{R}^3} v_r \, B(\mathbf{v})B(\mathbf{v}') \, d^3\mathbf{v} \, d^3\mathbf{v}'}{\int_{\mathbb{R}^3}\int_{\mathbb{R}^3} B(\mathbf{v})B(\mathbf{v}') \, d^3\mathbf{v} \, d^3\mathbf{v}'} \\ &= \int_{\mathbb{R}^3}\int_{\mathbb{R}^3} v_r \, B(\mathbf{v})B(\mathbf{v}') \, d^3\mathbf{v} \, d^3\mathbf{v}' \\ &= \int_{\mathbb{R}^3}\int_{\mathbb{R}^3} v_r \, \left( \frac{m}{2\pi k T} \right)^{3/2} e^{-\frac{m|\mathbf{v}|^2}{2kT}} \left( \frac{m}{2\pi k T} \right)^{3/2} e^{-\frac{m|\mathbf{v'}|^2}{2kT}} \, d^3\mathbf{v} \, d^3\mathbf{v}' \\ &= \left( \frac{m}{2\pi k T} \right)^3 \int_{\mathbb{R}^3}\int_{\mathbb{R}^3} v_r \, e^{-\left(\frac{m|\mathbf{v}|^2}{2kT} + \frac{m|\mathbf{v'}|^2}{2kT}\right)} \, d^3\mathbf{v} \, d^3\mathbf{v}' \\ &= \left( \frac{m}{2\pi k T} \right)^3 \int_{\mathbb{R}^3}\int_{\mathbb{R}^3} v_r \, e^{-\frac{mv_r^2}{2kT}} \, d^3\mathbf{v} \, d^3\mathbf{v}' \\ &= \left( \frac{m}{2\pi k T} \right)^3 \int_0^{\infty} v_r \, e^{-\frac{mv_r^2}{2kT}} \, \pi^3 v_r^5 \, dv_r \\ &= \left( \frac{m}{2 kT} \right)^3 \int_0^{\infty} v_r^6 \, e^{-\frac{mv_r^2}{2kT}} \, dv_r \\ &= \left( \frac{m}{2 kT} \right)^3 \frac{15}{16}\sqrt{\pi} \left( \frac{m}{2kT} \right)^{-7/2} \\ &= \frac{15}{16}\sqrt{\pi} \, \left( \frac{2kT}{m} \right)^{1/2} \end{align}

where on line 6, $\pi^3 v_r^5$ is the surface area of a 6-sphere of radius $v_r$, and the integral on line 8 was calculated using Wolfram Alpha.

Thus we find that $\langle v_r \rangle / \langle v \rangle = 15\pi/32 \approx 1.4726$ which is almost, but not exactly equal to the textbook result $\sqrt{2} \approx 1.4142$.

[Edit: the textbook result $\langle v_r \rangle / \langle v \rangle = \sqrt{2}$ is correct. For the full, correct proof without any simplifying assumptions, please see the other answer.]

Another way:

$$\langle v_r \rangle = \mathcal{N}^2 \iint \vert \vec{v}_1 - \vec{v}_2\vert \ e^{-\frac{m(\vec{v}_1^2+\vec{v}_2^2)}{2k_B T}}\ \mathrm{d}\vec{v}_1 \mathrm{d}\vec{v}_2,$$

wherein the normalized Boltzmann distribution function is $\mathcal{N} e^{-\frac{m \vec{v}^2}{2k_B T}}$, $\ \mathcal{N}$ being the numerical factor that normalizes the Gaussian distribution.

Rotating velocity coordinates by $45^{\circ}$, define new coordinates $$\vec{U}= \frac{\vec{v}_1+\vec{v}_2}{\sqrt{2}},$$ $$\vec{V}= \frac{-\vec{v}_1+\vec{v}_2}{\sqrt{2}}.$$

Since the above transformation is a rotation, the Jacobian is 1.

Rewriting $\langle v_r \rangle$ in terms of new coordinates,

\begin{align} \langle v_r \rangle &= \mathcal{N}^2 \iint \vert \sqrt{2} \vec{V}\vert e^{-\frac{m(\vec{V}^2+\vec{U}^2)}{2k_B T}}\ \mathrm{d}\vec{V} \mathrm{d}\vec{U} \\ &= \sqrt{2} \ \underbrace{\int \vert \vec{V}\vert \ \mathcal{N} e^{-\frac{m\vec{V}^2}{2k_B T}}\ \mathrm{d}\vec{V}}_{\text{equals } \langle v \rangle} \ \ \underbrace{\int \mathcal{N} e^{-\frac{m\vec{U}^2}{2k_B T}} \ \mathrm{d}\vec{U}}_{\text{equals }1} \\ & = \sqrt{2} \ \langle v \rangle. \end{align}

The relation, $\langle v_r \rangle = \sqrt{2} \ \langle v \rangle $, seems to hold independent of the spatial dimensionality.