(We assume an idealized model where there are no interactions, no collisions, no external potential field, and all the particles are considered to be classical.)

We write the distribution of velocities as:

\[ \Delta n = \rho (v) \Delta v , \]

And then, employing the standard sloppiness, we take:

\[ dn = \rho(v) dv \]

to be a continuous approximation when \[\Delta n, \]and \[\Delta v, \]are "small."

This distribution of velocities has one boundary condition:

\[ \rho(\infty) = 0 : \]

there are no particles with infinite velocities.

And one normalization condition:

\[ \int{\rho(v)dv} = N \ .\]

Where N is the total number of particles.

Note that the above integral is really over all velocities in a three dimensional "velocity space," and could therefore be written more explicitly as:

\[ \int_{- \infty}^{\infty} \int_{- \infty}^{\infty} \int_{- \infty}^{\infty} \rho (v_{x}, v_{y}, v_{z}) dv_{x} dv_{y} dv_{z} \]

in Cartesian coordinates, or as:

\[ \int_{0}^{\pi} \int_{ 0}^{2\pi} \int_{0}^{\infty} \sigma (r, \theta , \phi ) r^2 dr \sin{\theta} d\theta d \phi \]

in polar coordinates, where \(\sigma( r, \theta, \phi)\) is simply the same distribution in terms of polar coordinates.

Now \[ \rho(\vec{v}) \]must depend on the magnitude of \[\vec{v}, \]but not on its direction.

Therefore

\[ \sigma (r, \theta , \phi) = \sigma(r). \]

Furthermore:

\[ dv = dxdydz = r^2 dr \sin{ \theta} d\theta d\phi ,\]

so that:

\[ \rho (\vec{v})dxdydz = \sigma(r) r^2 dr \sin{ \theta} d\theta d\phi \]

implies:

\[ \rho (\vec{v}) = \sigma(r) \]

Symmetry dictates that no direction should be different than any other.

Expressing \(\rho(\vec{v})\) in cartesian coordinates, this translates to the statement that \(\rho(\vec{v})\) factors as:

\[ \rho (\vec{v}) = f(v_x)f(v_y)f(v_z) \]

Putting this together with our previous results we obtain:

\[ \rho (\vec{v}) = f(v_x)f(v_y)f(v_z) = \sigma(r) \]

Setting \[ v_{y} = v_{z} = 0 \]in the above equation we obtain:

\[ \rho (\vec{v}) = f(v_x)f(0)f(0) = \sigma(r) = \sigma\left(\sqrt{v_{x}^2}\right) = \sigma(v_{x}) \]

Now setting \[v_{x} = r : \]

\[ f(r)f(0)f(0) = \sigma(r) \]

So that:

\[ \rho (\vec{v}) = f(v_x)f(v_y)f(v_z) = \sigma(r)= \sigma \left(\sqrt{v_{x}^2 + v_{y}^2 + v_{z}^2}\ \right) = f(r)f(0)f(0) \]

Divide both sides by \[ \left(f(0)\right)^{3} \]

\[ \frac{f(v_x)}{f(o)} \frac{f(v_y)}{f(o)} \frac{f(v_z)}{f(o)} = \frac{f(r)}{f(0)} \]

And then take the logarithm:

\[ \ln{ \frac{f(v_x)}{f(o)}} + \ln{ \frac{f(v_y)}{f(o)}} + \ln{ \frac{f(v_z)}{f(o)}} = \ln{ \frac{f(r)}{f(o)}} \]

Now define \[ h(x^2) \equiv \ln{ \frac{f(x)}{f(o)}}. \]

Then:

\[ h(v_{x}^2) + h(v_{y}^2) + h(v_{z}^2) = h( v_{x}^2 + v_{y}^2 + v_{z}^2 ) \]

This will be true only if \(h(x^2) = ax^2\), where \(a\) is some constant.

Therefore:

\[ \frac{f(r)}{f(0)} = e^{ar^2} \]

and,

\[ \sigma (r) = f(r)f(0)f(0) = e^{ar^2} \left( f(0) \right)^3 = C e^{ar^2} \]

In order to achieve the boundary condition:

\[ \sigma (\infty) = 0 \ , \]

\(a\) must be less than zero. We will take \(\kappa = | a | \), so that:

\[ \sigma (r) = C e^{-\kappa r^2} \]

We have:

\[ dn_v = \rho (\vec{v}) dv = \sigma (r) dv = C e^{-\kappa v^2} dv \]

The probability that a particle will have velocity \(v\) is:

\[ p(v)dv = \frac{dn_v}{N} = \frac{C e^{-\kappa v^2} dv}{N} = \frac{C e^{-\kappa v^2} dv}{\int C e^{-\kappa v^2} dv } = \frac{ e^{-\kappa v^2} dv}{\int e^{-\kappa v^2} dv } \]

So in order to find \(p(v)dv\) we don't need to know \(N\), or \(C\), but we do need to know \(\int e^{-\kappa v^2} dv\).

This is:

\[ \int e^{-\kappa v^2} dv = \int_{- \infty}^{\infty} \int_{- \infty}^{\infty} \int_{- \infty}^{\infty} e^{-\kappa v^2} dxdydz \]

\[ = \int_{0}^{\pi} \int_{ 0}^{2\pi} \int_{0}^{\infty}e^{-\kappa r^2} r^2 dr \sin{\theta} d\theta d \phi= 4 \pi \int_{0}^{\infty}e^{-\kappa r^2} r^2 dr \]

This is either the product of three gaussian integrals, or \(4\pi\) times a beta function. Either way the answer is:

\[ \int e^{-\kappa v^2} dv = \left( \frac{\pi}{\kappa} \right)^{\frac{3}{2}} \]

So that:

\[ p(v)dv = \left( \frac{\kappa}{\pi} \right)^{\frac{3}{2}}e^{-\kappa v^2} dv.\]

This defines a Gaussian, or normal, probability distribution.

Note that at \[v=0: \]

\[ \frac{dn_0}{N} = p(0)dv = \left( \frac{\kappa}{\pi} \right)^{\frac{3}{2}} dv \]

In other words, at zero, the Gaussian probability distribution is uniform.

If you have studied statistical mechanics, then you know that \(\kappa = \frac{m}{2KT}\), where \(m\) is the mass of each particle, \(K\) is Boltzman's constant, and \(T\) is temperature.

The probability distribution:

\[ p(\vec{v}) dv = \left( \frac{m}{2\pi kT} \right)^{\frac{3}{2}}e^{-\frac{-m v^2}{2kT}} dxdydz\]

Is called the Maxwell distribution.

As \(T\) approaches zero, the Maxwell distribution approaches a Dirac delta distribution.