This webpage contains source codes of C++ programs to generate hard-particle packings in a Euclidean space of arbitrary dimension. If you use these codes, please reference this paper:

"Packing hard spheres in high dimensional Euclidean spaces", by M. Skoge, A. Donev, F. H. Stillinger, and S. Torquato, Phys. Rev. E, Vol. 74: 041127 (2006) [ibid 75: 029901 (2007)], [cond-mat/0608362]

Disclaimer: The programs given here are general molecular-dynamics codes. They can be used to prepape disordered jammed packings of hard particles, if used with the appropriate options. The programs do not automatically generate packings, yet alone jammed packings. The above paper should provide background on the algorithm and appropriate parameters to use. Additional more sophisticated executable codes can be found at, including codes to analyze packings.

Please read the following instructions carefully before asking questions, and contact the original author Monica Skoge (mskoge at domain and/or Aleksandar Donev (aleks.donev at domain for help.

Source Codes

The simplest (but not maintained!) version of code, suitable if you want to quickly modify it or learn how the algorithm works, is found in the zip archive, with periodic boundaries. The usage is documented below. Reflecting (hard wall) boundaries are implemented in box_reflecting.C, which replaces box.C. Packings inside cylindrical containers in 3D can be generated using some modifications by Todd Weisgraber (LLNL), which you can find in the archive

The most advanced and latest code is one that can handle bidisperse packings, found in the archive (the input file has the radius and mass ratios and the fraction of each component). These codes in principle support any distribution of sphere radii or masses, and support both periodic and reflecting boundaries, but you will need to customize for your own needs. You can get monodisperse packings by setting the bidispersity ratio to zero.

Instructions for Monodisperse Codes

The files needed to compile are neighbor.C, spheres.C, box.C, sphere.C, event.C, heap.C, and read_ input.C along with corresponding header (.h) files. The main program is spheres.C and calls the core computational functions in box.C. Simply compile with
	g++ -O neighbor.C spheres.C box.C sphere.C event.C heap.C read_input.C
NOTE: With some compilers you may need to rename our own vector.h into something like my_vector.h and change the C++ sources accordingly because it may confuse it with the STL header! Additionally, newer C++ compilers may not include the standard C library headers by default, in which case you need to "#include <cstdlib>" in spheres.C.

The codes use
The dimension of the packing is defined through the constant DIM in vector.h and the input file containing the parameters and input/output file names is in input. In this file, the user can either specify to read in an existing configuration of hard spheres from a file by specifying the appropriate file name in readfile or the user can create a new configuration by setting readfile to new. In the latter case, a configuration of N spheres is created at density initialpf  by random sequential addition.
If the temperature temp is set to zero, then all spheres are given zero initial velocities; otherwise velocities are initially chosen from a Maxwell-Boltzman distribution. The expansion rate growthrate governs the rate at which the sphere diameter increases in units of the simulation box size divided by the sphere velocity. Due to expansion, the kinetic energy increases after collisions. Therefore velocities are rescaled to 1 after every "cycle". The user specifies the number of events per sphere in one cycle with the parameter eventspercyle. Statistics, such as pressure, are calculated after every cycle and outputted to datafile. The simulation runs until reaching either the maximum allowed packing fraction maxpf or the maximum allowed pressure maxpressure. Before terminating, the positions of all spheres are outputted to writefile.

The program produces two output files:
  1. writefile contains the final positions of all the spheres. Its format is explained below.
  2. datafile contains the statistics calculated after each cycle. The statistics include (as columns in this order) the current packing fraction, pressure, total change in kinetic energy before rescaling and the total number of MD events processed.

The ASCII file formats used for storing sphere packing configurations is as follows (in three dimensions):
n_dims ! Euclidian space dimensionality d
N n_species ! Total number of particles and number of species
N1 [N2...] ! Number of particles in each species
D1 [D2...] ! Sphere diameters for each species
Lambda ! Lattice vectors for unit cell, here a unit cube
T T T ! Periodic boundary conditions along all dimensions
x1 y1 z1 ! Coordinates of first sphere center
x2 y2 z2 ! second, etc...
The format may seem weird, but it is chosen to conform to the format described at