The probability density function for molecular speeds \(v\) in thermal equilibrium is given by:
\[f(v) = 4\pi \left(\frac{m}{2\pi k_B T}\right)^{3/2} v^2 \exp\left(-\frac{mv^2}{2k_B T}\right)\]
Because velocity components are orthogonal and independent, the joint probability distribution factorizes into three independent Gaussians:
\[f(v_x, v_y, v_z) = f(v_x)f(v_y)f(v_z) = \left(\frac{m}{2\pi k_B T}\right)^{3/2} \exp\left(-\frac{m(v_x^2+v_y^2+v_z^2)}{2k_B T}\right)\]
Here is the code used to generate the 3D phase-space simulation above: