This webpage contains executable versions of Fortran programs to generate and analyze two- and three-dimensional hard-particle packings using event-driven molecular dynamics, as explained in this paper:
A. Donev, S. Torquato, and F. H. Stillinger, Neighbor List Collision-Driven Molecular Dynamics for Nonspherical Hard Particles: II. Applications to Ellipses and Ellipsoids, Journal of Computational Physics, 202, 765 (2005).
We have also released C++ source codes for much stripped-down version that only handles
monodisperse packings of spheres, however, in an arbitrary dimensional Euclidean space.
The programs given here are general molecular-dynamics codes. They can be used to prepare disordered jammed packings of hard particles, if used with the appropriate options. The programs do not automatically generate packings, yet alone jammed packings. I would be glad to assist in preparing the input configuration files to perform whatever specific task you have, just e-mail me at
The tar file of all Linux executables and input files needed to play around with event-driven MD and analyze the resulting packings is C++ source codes. For now I have only documented the sphere packing options and hidden the rest, in order to avoid additional complexities for ellipsoids. I have also hidden options related to neighbor lists and other optimizations that require close familiarity with the algorithm and are hard to explain in detail. Contact me if you need more, or see the last section in this webpage (Advanced Usage). If you use this to produce published results, please reference the above paper.
Please read the following instructions carefully before asking questions, and also take a look at the examples I provided!
Producing Packings using MD
The binaries for generating jammed packing using hard-particle MD are “./PackLSD.static.x” for packings which consist of a small number of particle species (monodisperse, bidisperse, tridisperse, etc.), or “./PackLSD.poly.static.x” for truly polydisperse packings where each particle has a different shape/size. The executables only link in standard (OpenGL, libtiff, libjpeg, etc.) Linux libraries, and link openglut statically, which may cause failures on some systems. The programs are compiled at medium optimization, for a generic IA32 processor, on a Red Hat EL4 system. They should run on Red Hat (>=9) and Fedora (>=3) systems just fine, and hopefully other similar distributions as well. If something fails to run, I can try to help prepare a version for you with the correct libraries.
The runtime input file “./PackLSD.input” or “./PackLSD.poly.input” is to be piped into the executable. The main input file is a collection of Fortran 90 namelists as in the examples “./Input/PackLSD.nml” and “./Input/PackLSD.poly.nml”.
The name of the packing can be changed inside “PackLSD.input” and is used as the base of all input/output filenames (as in the filename PackLSD.nml). The program usually prints when files are opened or closed. The executable needs to be run in a directory where there will be sub-directories called Input, Output, Stats, Data, Scenes, Plots, and VRML, and the files are placed or read from these subdirectories. The nml file must be in the Input directory, and the output packing
files will be saved to Data. I have documented the options you will need in the nml file and removed most others (some need to be left in since the defaults are not appropriate for generic uses). Most options are the same for both programs “PackLSD.static.x” and “PackLSD.poly.static.x”, and the minor differences can be seen in the example “./Input/PackLSD.poly.nml”. Please ask me if some option is not clear.
The simulation processes a certain number of MD events (collisions), called a stage. At the end of every stage certain statistics are output, possibly a temporary file saved, the OpenGL view is refreshed, particle velocities rescaled, etc. Many stages can be run until jamming is reached, as defined by the options in the namelist file. The executable can use OpenGL to render to screen the packing process so you can see what is happening. You can use the mouse to zoom in, rotate, etc. There are also options to generate snapshot images that can be used to make animations (I highly recommend the MNG file format!)
The produced packings can later be viewed/analyzed by using the programs “ProcessPacking.static.x” or “ProcessPacking.poly.static.x”. In particular, you can render using OpenGL, make plots, or compute the pair correlation function \(g_2(r)\) or the structure factor \(S(k)\), as well as construct and render NNLs (near-neighbor lists) based on proximity or Delaunay/Voronoi neighbors. The input file is the namelist “./Input/ProcessPacking.nml”
in the Input subdirectory. It contains comments describing the various input fields. In 2D one can generate nice EPS figures to include in papers, and in 3D I provide some VRML models in ./VRML. These are not documented in detail (learn by looking at examples, although knowledge of VRML is required). The (anisotropic, vector k) structure factor in 3D is saved as a VTK file that can be looked at using the VisIt rendering program.
The input/output packing files use a format which should be simple to understand, see the example dat file “./Data/LSD.HS.mono.3D.N=1000.1.dat”. For packings that contain a small number of particle species (mono, bidisperse, etc.), this format is an ASCII text file as follows:
n_dims ! Euclidian space dimensionality (2 or 3)
N dispersity ! Total number of particles and number of species
N_1 N_2 ... ! Number of particles of each specie
D_1 D_2... ! Sphere diameters for each specie
Lambda ! Lattice vectors for unit cell
T T T ! Periodic or hard wall boundaries
x1 y1 z1 ! Coordinates of first sphere center
x2 y2 z2 ! second, etc...
The format is slightly different for truly polydisperse packings, the
particle-specific diameter is after each particle’s coordinates:
N_1 N_2 ... ! One for each specie (each specie is a distribution for polydisperse packings, rather than a fixed size)
Lambda ! Unit cell
T T T ! Periodic or hard wall
x1 y1 z1 D1
x2 y2 z2 D2
Contact/Neighbour networks are written to files in ending with
“.nnl.dat” and have the following format:
M2 N M
j1_1 k1_1 l1_1
j1_2 k1_2 l1_2
j2_1 k2_1 l2_1
Here M2 is the total number of contacts listed in the file. This includes both directions for contacts between particles, and M gives the total number of contacts without double-counting (M=M2/2 with periodic BCs, but it may change if hard-wall neighbors are included). N is the total number of particles. Then comes a list of all N_1 contacts of particle 1, listed one after the other, where j1_1 is the neighbour particle, k1_1 is an image unit cell identifier (explained shortly) and l1_1 is the strength of the
interaction between the two particles (such as force—it can be ignored if not needed, but it has to be there as a floating-point number), etc. The image identifiers is used for periodic boundary conditions to indicate contacts crossing the torus and is explained in the paper cited at the top of this page. It is zero if the contact is
between particles in the base unit cell. After all the contacts for particle 1 have been listed, the same is done for particle 2, etc. Note that if a particle neighbour is negative that indicates a hard wall, with
where k indicates which of the 2d walls it is.
Advanced Usage: Full list of namelist options
Once familiar with the basic operation of the above programs, and one has read the necessary references, one can try to use also more advanced features such as near-neighbor lists (NNLs). For the packing-generation codes, the input file “./PackLSD.full.input” contains all the input options. The namelist input file “./Input/PackLSD.2D.nml” contains all of the input options (it is set up for 2D, with obvious generalization to 3D). For analyzing packings, the input file “./ProcessLSD.full.input contains the full list of input options, and the namelist “./Input/ProcessPacking.2D.nml” contains a listing of all the available options. I recommend adding options gradually, becoming familiar with them first, and also making sure that the version of the executable program you are using is up to date and accepts the given namelist option (since I add new options all the time—if an option is not recognized, just comment it out and ask for an updated executable). It is not trivial to decide which option is compatible with other options or which option is required when some flag is set. Play around and ask me questions if needed.
There is also a version of the packing codes that generates packing of superellipses and superellipsoids. This is implemented by hacking the ellipsoid code to incorporate superellipsoid exponents. Note that this code is numerically sensitive for high exponents, and also that it is much slower than the ellipsoid codes. Note also that near neighbor lists and thus anything that requires them does not work for superellipsoids because the so-called inner overlap potential has not yet been implemented (or theoretically generalized) to superellipsoids. Executable programs that work with superellipsoids are “./PackLSD.super.static.x” and “./ProcessPacking.super.static.x”.
As an example usage, try:
PackLSD.super.static.x < PackLSD.super.input
ProcessPacking.super.static.x < ProcessPacking.full.input