Molecular dynamics of hard disks


The algorithm for simulating the dynamics of a system of hard disks can be summarized by the following:

  1. For each particle find the particle that it would collide with next and record the collision time.
  2. Determine the minimum collision time.
  3. Move all the particles along a straight line trajectory for a time equal to the minimum collision time.
  4. Update the collision times.
  5. Compute the post collision velocities of the two particles that collide.
  6. Calculate any quantities of interest and accumulate data.
  7. Update the collision partners for the colliding pair and any other particle that would have collided with either particle of the colliding pair if they had not collided.

Because hard disks move in straight lines at constant speed between collisions and change their velocities instantaneously when a collision occurs, the problem becomes finding the next collision. The algorithm is said to be event driven and the trajectories can be computed exactly in principle; in practice the accuracy is limited only by roundoff errors.


  1. Run the simulation using the default parameters. Does the temperature fluctuate during the simulation? What is the effect of doubling the initial temperature without changing any of the other parameters? Is the temperature a useful quantity for hard disks?
  2. Each point of the plot of PA/NkT is averaged over ten collisions. Does the pressure P fluctuate? Compute PA/NkT as a function of the density for a fixed number of particles. Is there any evidence of a phase transition where the slope of the pressure with respect to density changes abruptly?
  3. The mean free path and the mean collision time are also computed. Show that the mean free path is independent of temperature, and the mean collision time is proportional to the square root of temperature. How does the mean free path and collision time change if you double the linear dimensions of the simulation cell?
  4. *We can gain more insight into the structure of a many body system by looking at how the positions of the particles are correlated with one another due to their interactions. The radial distribution function g(r) is a measure of this correlation. Suppose that N particles are in a region of area A with number density ρ = N/A. Choose one of the particles to be the origin. Then the mean number of other particles between r and r + dr is defined to be ρg(r)dr. If the interparticle interaction is spherically symmetric and the system is a gas or a liquid, then g(r) depends only on the separation r = |r|. The normalization condition for g(r) is

    ρ∫g(r)dr = N - 1 ≅ N,

    where dr = 2πrdr in two dimensions. To determine g(r) for a particular configuration of particles, we first compute n(r,Δr), the number of particles in a circular shell between radius r and r + Δ r, with the center of the shell centered about each particle. We then normalize g(r). What is the qualitative nature of g(r) as a function of r for low densities? How does g(r) change as the density is increased? How does g(r) for a system of hard disks compare to g(r) for a system of Lennard-Jones particles at the same density?

Java Classes

Updated 28 December 2009.