GranularMedium simulations made by
Nicholas Guttenberg with his GRNLR particle simulator. The program was orignally used for studying granular jets, wet granular droplet impact, and tipping icebergs.
Ref. "An approximate hard sphere method for densely packed granular flows" (link  283 kb)
Ref.: "Grains and gas flow: Molecular dynamics with hydrodynamic interactions" (link  207 kb)
In
WooDEM we had the problem that particles slowed down due to dissipation (Coefficient of Restitution < 1). So energy was being constantly removed from the simulation due to collisions, to keep the particles going Nicholas used 3 solutions in his program for adding energy:
A. Swimmers: after interaction a force is added to keep the particle's velocity close to 1
F=a*v*(1v^2)
 The first ‘v’ is a vector,
 v^2 is a scalar,
 ’a’ is a scalar that controls how strong this force is.
This turns the particles into active swimmers, and you can create interesting structure by doing this in the presence of a high density of particles (basically they have to move because of the force, but they want to be still because they're near the jamming point, so they end up moving in largescale vortices similar to those in the
Twirls video).

B. Thermalisation: adding a force to keep energy despite dissipation
F=sqrt(T)*eta/sqrt(dt)
 ’T’ is the desired temperature (normally there'd be a measure of dissipation, but your dissipation is due to inelastic collisions),
 ‘dt’ a timestep used by the algorithm (needed for of how noise scales),
 ‘eta’ a random number generated every time you call the function to get the force.
This will cause the particle motions to have some relatively constant level of energy despite dissipation.

C. Rescale timestep: so motion occurs at a constant velocity
For hard grains there is no inherent energy scale, so a system of grains moving at 1mm/year and a system of grains moving at 100 m/s have the same physics, just strewn out over a different time scale.
What you could do is snap a frame at irregular intervals. First at every second, as the grains slow down, then every 2, 4, 8, ... The right interval to use would be based on the square root of the system's average energy (e.g. sum up v^2 for every grain, divide by the number of grains, then take the square root and multiply by a constant to make it look good).
GRNLR  GUI