Ideal Gas Simulation and the Maxwell–Boltzmann Distribution

Kinetic Energy (Joules)

Slider energy: unset

idk

Temperature (Kelvin)

idk

Maxwell-Boltzmann Distribution for a microcanonical ensemble with no potential energy

Equation for 3d elastic collissions of spheres from here.

We start with the Boltzmann Distribution. For a "system" lightly coupled to a "bath" of infinite heat capacity, it gives an ensemble average that is equal to the time average under the assumption of ergodicity. We assume position and momentum $\vec{q}, \vec{p}$ and Hamiltonian $H$ . We look at the state of one particle, the "system", in contact with all the other particles, the "bath". Note that here we have a finite number of other particles in the bath, so the assumption of infinite heat capacity is wrong, but it's pretty close, and becomes exact in the limit of infinite particles. The distribution says probability of the particle being in any state $\vec{q}, \vec{p}$ is:

ρ (\vec{q}, \vec{p}) = \frac{e^{- \frac{H (\vec{q}, \vec{p})}{k_{B} T}}}{\int_{- \infty}^{\infty} \int_{V} e^{- \frac{H ({\vec{x}}^{'}, {\vec{p}}^{'})}{k_{B} T}} d {\vec{q}}^{'} d {\vec{p}}^{'}} = \frac{e^{- \frac{H (\vec{q}, \vec{p})}{k_{B} T}}}{Z}

If the Hamiltonian has no potential energy, then for a non-relativistic particle it would simply be:

H (\vec{q}, \vec{p}) = \frac{{\vec{p}}^{2}}{2 m} = \frac{\vec{p} \cdot \vec{p}}{2 m} = \frac{p_{x}^{2} + p_{y}^{2} + p_{z}^{2}}{2 m}

Because $\vec{q}$ isn't in the expression, we focus just on the probability over the momentum:

ρ (\vec{p}) = \frac{e^{- \frac{p_{x}^{2} + p_{y}^{2} + p_{z}^{2}}{2 m k_{B} T}}}{Z}

The partition function $Z$ is

Z = \int_{- \infty}^{\infty} \int_{- \infty}^{\infty} \int_{- \infty}^{\infty} e^{- \frac{p_{x}^{2} + p_{y}^{2} + p_{z}^{2}}{2 m k_{B} T}} d p_{x} d p_{y} d p_{z}

It's is just a 3D gaussan and we know what that integral is.

= (τ m k_{B} T)^{3 / 2}

Then substituting $Z$ into the Boltzmann equation gives.

ρ (\vec{p}) = (τ m k_{B} T)^{- 3 / 2} e^{- \frac{p_{x}^{2} + p_{y}^{2} + p_{z}^{2}}{2 m k_{B} T}} = (τ m k_{B} T)^{- 3 / 2} e^{- \frac{{\vec{p}}^{2}}{2 m k_{B} T}}

The above is almost the Maxwell-Boltzmann distribution but it's in terms of momentum instead of velocity/speed. We first do a change of variable into velocity on the probability distribution. Momentum is a function of velocity where $\vec{p} (\vec{v}) = m \vec{v}$ .

ρ (\vec{v}) = ρ (\vec{p} (\vec{v})) (\frac{d {\vec{p}}_{x}}{d {\vec{v}}_{x}} \frac{d {\vec{p}}_{y}}{d {\vec{v}}_{y}} \frac{d {\vec{p}}_{z}}{d {\vec{v}}_{z}}) = m^{3} ρ (\vec{p} (\vec{v}))

ρ (\vec{v}) = (\frac{m}{τ k_{B} T})^{3 / 2} e^{- \frac{m {\vec{v}}^{2}}{2 k_{B} T}}

Now if we want to convert it to speed instead of velocity we have to do another change of variables. We first bijection velocity into spherical coordinates:

s = \sqrt{v_{x}^{2} + v_{y}^{2} + v_{z}^{2}}

θ = \arccos (\frac{v_{z}}{s})

ϕ = \arctan (\frac{v_{y}}{v_{x}})

The determinate of the Jacobian of the inverse function is:

s^{2} \sin (θ)

So we can do the change of basis on the probability distribution.

ρ (s, θ, ϕ) = ρ (\vec{v} = ⟨ s \sin (θ) \cos (ϕ), s \sin (θ) \sin (ϕ), s \cos (θ) ⟩) s^{2} \sin (θ)

Then to find the probability of a speed you just integrate the marginal over the unwanted variables.

ρ (s) = \int_{0}^{τ} \int_{0}^{τ / 2} ρ (\vec{v} = ⟨ s \sin (θ) \cos (ϕ), s \sin (θ) \sin (ϕ), s \cos (θ) ⟩) s^{2} \sin (θ) d θ d ϕ

Each velocity coordinate is independent so we can split the probability distribuion.

\int_{0}^{τ} \int_{0}^{τ / 2} ρ (v_{x} = s \sin (θ) \cos (ϕ)) ρ (v_{y} = s \sin (θ) \sin (ϕ)) ρ (v_{z} = s \cos (θ)) s^{2} \sin (θ) d θ d ϕ

And since all orientations of $\vec{v}$ are equally likely, the probability only depends on the magnitude.

\int_{0}^{τ} \int_{0}^{τ / 2} ρ ({\vec{v}}^{2} = s^{2}) s^{2} \sin (θ) d θ d ϕ

We can then move it out of the integral.

ρ ({\vec{v}}^{2} = s^{2}) \int_{0}^{τ} \int_{0}^{τ / 2} s^{2} \sin (θ) d θ d ϕ

And the integral is now just the surface area of a sphere which we know.

ρ (s) = 2 τ s^{2} ρ ({\vec{v}}^{2} = s^{2})

Then put in the probability and we get the Maxwell Boltzmann distribution.

= 2 s^{2} \sqrt{\frac{1}{τ}} (\frac{m}{k_{B} T})^{3 / 2} e^{- \frac{m s^{2}}{2 k_{B} T}}

Finding the temperature

To find the temperature of a microcanonical system (NVE ensemble) of a bunch of particles with only kinetic energy, we look at the manifold of all states where the Hamiltonian equal to a constant energy $H (\vec{q}, \vec{p}) = E$ . Since there is only kinetic energy, the positions of the particles $\vec{q}$ don't matter.

H (\vec{q}, \vec{p}) = \sum_{i} \frac{p_{i}^{2}}{2 m} = E

The definition of temperature in a microcanonical ensemble is

\frac{1}{T} = \frac{\partial S (N, V, E)}{E}

where the entropy is $S = - k_{B} \sum_{i} ρ (x_{i}, p_{i}) l n (ρ (x_{i}, p_{i}) = - k_{B} l n (\frac{1}{W (E)}) = k_{B} l n (W (E))$ and $W (E)$ is the number of accessable states accessable at $E$ . First re-arrange the equation for the energy.

\sum_{i} p_{i}^{2} = 2 m E

If we want to find $W (E)$ then assuming all particles are indistunguishable, means it has a factorial. Then we use Plank's constant as a measure of the phase space density.

W (E) = \frac{1}{h^{3 N} N!} \int_{V} \int_{- \infty}^{\infty} δ (H (\vec{q}, \vec{p}) - E) d^{3 N} \vec{p} d^{3 N} \vec{q}

With no potential energy the position doesn't matter and can be factored out

= \frac{V^{3 N}}{h^{3 N} N!} \int_{- \infty}^{\infty} δ (H (\vec{q}, \vec{p}) - E) d^{3 N} \vec{p} d^{3 N}

First switch to "hyperspherical" coordinates where $s = \sqrt{{\vec{p}}^{2}}$ . All the the angles are just integrated over the surface of the hypersphere, since the radius is fixed, so we can just multiply it by the suface area of a $3 N$ dimensional sphere $\frac{2 (τ / 2)^{3 N / 2}}{Γ (3 N / 2)} s^{3 N - 1}$ .

= \frac{V^{3 N}}{h^{3 N} N!} \frac{2 (τ / 2)^{3 N / 2}}{Γ (3 N / 2)} \int_{0}^{\infty} δ (\frac{s^{2}}{2 m} - E) s^{3 N - 1} d s

And the Dirac measure just chooses the one point.

= \frac{V^{3 N}}{h^{3 N} N!} \frac{2 (τ / 2)^{3 N / 2}}{Γ (3 N / 2)} (2 m E)^{\frac{3 N - 1}{2}}

And we just plug it in to the temperature definition.

\frac{1}{T} = k_{B} \frac{\partial l n (W (E))}{\partial E} = k_{B} \frac{1}{W (E)} \frac{\partial W (E)}{\partial E}

Substituting gives:

= k_{B} \frac{\frac{\partial}{\partial X E} \frac{V^{3 N}}{h^{3 N} N!} \frac{2 (τ / 2)^{3 N / 2}}{Γ (3 N / 2)} (2 m E)^{\frac{3 N - 1}{2}}}{\frac{V^{3 N}}{h^{3 N} N!} \frac{2 (τ / 2)^{3 N / 2}}{Γ (3 N / 2)} (2 m E)^{\frac{3 N - 1}{2}}}

Then almost everythign cancels.

= k_{B} \frac{\frac{3 N - 1}{2} E^{\frac{3 N - 1}{2} - 1}}{E^{\frac{3 N - 1}{2}}}

Simplifying means:

\frac{1}{T} = k_{B} \frac{\frac{3 N}{2} - 1}{E}

this tempeature is:

k_{B} T = \frac{E}{\frac{3 N}{2} - 1}

This equation is consistent with the equipartition theorem as each of the $3 N - 1$ degrees of freedom contribute $\frac{3}{2} k_{B} T$ to the total energy.