In this animation we start with all particles at the same energy but with different velocity directions. We will see how the particle energy distribution evolves to become the Maxwell-Boltzmann (MB) distribution. The energies of the majority of the particles must be reduced to achieve this. For the first scatter where the particles are at equal energy, `E`, the energy of particle 1 becomes `E+ delta E` and the energy of particle 2 becomes `E- delta E` so the distribution of particle energy changes has to be symmetrical about zero energy change. After the first scatter, the two particles no longer have equal energies so their energies become `E_1+delta E_(s2)` and `E_2-delta E_(s2)` so the distribution function will generally no longer be symmetrical about zero energy change. In this animation we will examine how this asymmetry leads to greater energy numbers at smaller energies. Starting with all particles at the same energy may seem like a very special energy distribution but, if we solve that problem, we can also solve the evolution for any initial energy distribution. The evolution of the energy distribution of a gas is a very important part of the evolution of the second law of thermodynamics. The second law states that, in a closed system, the motional (kinetic) energy cannot go from a lower energy density region to a higher energy density region. Note that I have avoided the use of the word "temperature" because temperature is a very poorly defined concept. For a simple 3D monoatomic gas like this, the following expression is the Maxwell-Boltzmann energy distribution:

`p(E)=b sqrt(E/bar E)exp(-{3E}/{2bar E})`

where `E` is kinetic energy and `bar E` is the average kinetic energy of all particles.
In this expression `p(E)` is the probability of a given energy and `b` is a normalization constant such that the integral of `p(E)` over all energies becomes 1.

`b int_{0}^{oo}p(E)dE=1`

Explanation of the Animation and Plots

Canvas 2 shows a green Plot of the expected MB distribution. It also, in red, gives a choice of either the histogram of the energy distribution of the particles or the change in energy of the particles from their starting mono-energetic distribution. Soon after starting the animation, note that the change in energy distribution is skewed to lower energies. This skew becomes less noticeable as the animation progresses because the MB energy distribution has already been reached.

Deriving the Equations for Final Velocities After Collision

In the following equations please be aware that bold letters denote vectors which have either (x,y) or (x,y,z) components depending on whether we are considering 2D or 3D calculations. For hard sphere collisions one must consider both conservation of momentum and energy. If the collision was mediated by a potential between the particles, one would only have to consider conservation of kinetic and potential energy. Conservation of momentum requires that the final momentum be the same as the initial momentum. For a collision between particles 1 and 2:

`m_1 bb"v"_{1f}+m_2 bb"v"_{2f}=m_1 bb"v"_{1i}+m_2 bb"v"_{2i} (1)`

where the subscript f stands for final value and the subscript i stands for initial value and the arrowed bold letter `vec bb"v"` stands for the velocity vector. The conservation of kinetic energy equation is as follows:

`m_1 bb"v"_{1f} cdot bb"v"_{1f} +m_2 bb"v"_{2f} cdot bb"v"_{2f}=m_1 bb"v"_{1f} cdot bb"v"_{1i} +m_2 bb"v"_{2f} cdot bb"v"_{2i} (2)`

where the `cdot` between the velocity vectors stands for the inner or dot product of the two vectors. When we combine equations 1 and 2 we get equations that relate the components of the velocity vectors. The most important concept in these calculations is that of the unit vector ` hat bb"u"` between the centers of the two particles that are colliding

`hat bb "u"=[(x_2-x_1)hat bb"x"+(y_2-y_1)hat bb"y"+(z_2-z_1)hat bb"z"]/r_{12}`

where the x, y, and z that have hats are unit vectors along their respective axes and

`r_{12}=sqrt[(x_2-x_1)^2+(y_2-y_1)^2+(z_2-z_1)^2]`

Note that the unit vector points from sphere 1 to sphere 2 so, if `x_2` is greater than `x_1`, the unit vector's x component is positive.
Applying equation 1, we know that the changes of momentum of the spheres are equal and opposite and are along the unit vector `hat bb"u"` so we can write:

`m_1 delta bb"v"_1=M delta v hat bb"u"=-m_2 delta bb"v"_2 (3)`

where `M` has units of mass and `M` and the scalar speed increment `delta v` are still to be determined. Now we may re-write equation 2 using equation 3 and obtain:

`((m_1 bb"v"_1+M delta v hat bb"u") cdot (m_1 bb"v"_1+M delta v hat bb"u"))/(m_1)+ ((m_2 bb"v"_1-M delta vhat bb"u") cdot (m_2 bb"v"_2-M delta vhat bb"u"))/(m_2) =m_1 bb"v"_{1} cdot bb"v"_{1} +m_2 bb"v"_{2} cdot bb"v"_{2}(3)`

After cancellations of the terms on its right side, equation 4 simplifies to:

`M delta v=(2m_1m_2)/(m_1+m_2)( bb"v"_1- bb"v"_2) cdot hat bb"u"`

We can now make the identifications:

`M=(2m_1m_2)/(m_1+m_2)`

and

`delta v=( bb"v"_1- bb"v"_2) cdot hat bb"u"`

where M is known as the "reduced mass" Finally using equation 3 again we can write the expressions for the final velocity vectors:

` bb"v"_{1f}= bb"v"_{1i}+M/m_1 [( bb"v"_1- bb"v"_2) cdot hat bb"u"] hat bb"u"`
` bb"v"_{2f}= bb"v"_{2i}-M/m_2 [( bb"v"_1- bb"v"_2) cdot hat bb"u"] hat bb"u"`

Deriving Equations for the Changes in Energy in a Collision.

Dropping the subcsripts `i` (i.e. `bb "v"_{1i}` becomes just `bb "v"_1`) the changes in energy due to the collision are

`delta E= m_1/2 [ bb"v"_{1}+M/m_1 [( bb"v"_2- bb"v"_1) cdot hat bb"u"]hat bb"u"]cdot [ bb"v"_{1}+M/m_1 [( bb"v"_2- bb"v"_1) cdot hat bb"u"]hat bb"u"] -(m_1v_1^2)/2`

this expression has the form:

`delta E=m_1/2( bb v_1+delta v_1 hat bb"u")cdot ( bb v_1+ delta v_1 hat bb"u")-(m_1v_1^2)/2`

and after the dot product we have:

`delta E=m_1/2[bb v_1 cdot bb v_1+delta v_1^2+2delta v_1 (bb v_1 cdot hat bb"u") ]-c`

So we are left with :

`delta E=m_1/2[delta v_1^2+2 delta v_1 (bb v_1 cdot hat bb"u") ]`

`delta v_1=M/m_1 [( bb"v"_2- bb"v"_1) cdot hat bb"u"]`

so we have

`delta E= m_1/2[{M/m_1 [( bb"v"_2- bb"v"_1) cdot hat bb"u"]}^2 +2{M/m_1 [( bb"v"_2- bb"v"_1) cdot hat bb"u"](bb "v"_1 cdot hat bb"u")}]`

To simplify notation let `v_(u)=bb "v" cdot bb "u"`. After a long algebra session we find that:

`delta E_2=(2m_1m_2(v_(1u)-v_(2u))(m_1v_(1u)+m_2v_(2u))) /(m_1+m_2)^2=M/(m_1+m_2)(v_(1u)-v_(2u))(m_1v_(1u)+m_2v_(2u))=-delta E_1 (6)`

This is a complicated function of the particle velocity directions with respect to the unit vector between the two particles and there is no way that we can expect to compute the energy change per second from it. Without loss of generality, we may have particle 1 velocity along the x axis. That still leaves the magnitude of the velocity of particle 1, the velocity of particle 2 and the direction of the unit vector as variables. We can make it a lot simpler by letting `m_1=m_2=m` and the result is:

`delta E_2=m/2(v_(1u)-v_(2u))(v_(1u)+v_(2u))=m/2(v_(1u)^2-v_(2u)^2)= -delta E_1 \ \ (6)`

Equation 6 is a function of only two variables `m/2v_(1u)^2 and m/2v_(2u)^2`. The range of these variables is from zero to the energy we assign to particles 1 and 2. We can plot this function by assigning a fixed value to one and plotting the other. For `delta E_2` if we fix the value `m/2v_(1u)^2` then we will have a parabola facing downward. Vice-Versa if we fix `m/2v_(2u)^2`.
We will now revert to numerically computing the `delta E` statistics in 3 dimensions.
As always, the learner is encouraged to press F12 (Windows) and select "Sources" to examine the HTML and javascript code that was used to create this page.