==================================================================
=== ===
=== GENESIS Distributed Memory Benchmarks ===
=== ===
=== MD1 ===
=== ===
=== Molecular Dynamics for a Lennard Jones fluid - 3D ver. ===
=== ===
=== Versions: Std F77, PARMACS, Subset HPF ===
=== ===
=== Author: Mark Pinches ===
=== Subset HPF: Bryan Carpenter ===
=== Department of Electronics and Computer Science ===
=== University of Southampton ===
=== Southampton SO9 5NH, U.K. ===
=== fax.:+44-703-593045 e-mail:mrsp@uk.ac.soton.pac ===
=== dbc@uk.ac.soton.ecs ===
=== ===
=== Last update: May 1993; Release: 2.2 ===
=== ===
==================================================================
1. Description
--------------
Molecular dynamics essentially involves solution of the equations of
motion for a system of a large number of interacting particles, i.e.
solving second order differential equations. The problem is solved
numerically by calculating approximate solutions at a large number of
time steps with interval dt, where dt is small in comparison with the
time taken for a molecule with average velocity to travel a typical
interatomic distance. The particles are usually considered to interact
through an effective pair potential which is often adjusted to reproduce
experimental results and is used to model the exact many-body potential.
For monatomic neutral atoms the form of pair potential most often
employed is the Lennard-Jones potential, which is characterised by a
steep repulsive term and a long-range tail. When modelling homogeneous
bulk phases, periodic boundary conditions are usually applied, so the
`box' containing the system to be studied can be considered to be
replicated by an infinite `lattice'.
A typical MD calculation for N particles involves, for each time step:
- calculating the forces on the N particles;
- calculating their new positions;
- calculating properties of interest, such as energies and
radial distribution functions.
In addition, diagnostic output and perhaps particle coordinates (the
latter involving a large amount of output for large systems) may be
required at intervals of tens or hundreds of time steps. A simulation
usually proceeds in two phases:
- Equilibration - the atoms in an initially static lattice are
given random velocities, sampled from a Gaussian distribution
scaled to give the correct temperature, and an equilibrium (i.e.
minimum free energy) configuration is determined by running the
simulation, rescaling the kinetic energies at each timestep to
maintain the required temperature.
- Simulation - the coordinates from the previous step are used
in a simulation, run at constant energy, from which information
on the system being studied is produced. The period of this
simulation must be sufficient to obtain good statistical
averages of the required properties.
For many systems, the number of terms can be limited because the short
range nature of the Lennard-Jones potential means that reasonable
accuracy can be obtained by considering only atom pairs with interatomic
distances within a cutoff range, which is of the order of a few atomic
diameters. Long range corrections are used in the calculation of
thermodynamic properties. For small systems this strategy is sufficient
to ensure that the forces are evaluated efficiently. For large systems,
however, the number of interatomic terms required to evaluate the force
is small compared with the total number of atom pairs, and the time
taken to decide whether the separation of a given pair is less than
r_{c} is itself highly significant. An algorithm has been developed by
Hockney et al, the linked-list algorithm, which minimises the
time taken to find atomic pairs that contribute significantly to the
potential. This algorithm subdivides the system into a number of cells,
which may be 1, 2 or 3 dimensional. Cubic cells are the most appropriate
for large bulk systems. At each timestep each atom is allocated to a
cell, and one atom from each cell is selected as the first atom in
that cell. Provided that the cell size is chosen properly, the
interatomic distances need only be calculated for atoms in the current
and nearest neighbour cells. Note however that for edge and corner
cells, those at the opposite edges and corners are involved because of the
periodic boundary conditions.
Implementation of this algorithm on a parallel machine involves the
allocation of a number of adjacent cells to one processor. Cells at the
edges and corners of the processor must communicate with neighbouring
processors and each cell must communicate with a `driver' which is
responsible for input/output and supervision of the process, which is
synchronised once per timestep.
2. Operating instructions
-------------------------
Changing problem size and numbers of processes:
Most of the parameters are internally fixed for the case of the gas
Krypton. The user has to specify only the number of processors in
each direction and the number of latice cells NC, which will result in
4*NC*NC*NC atoms. These are read from the standard input on channel 5.
Suggested Problem Sizes :
For realistic simulation at least 2000 atoms per node is recommended,
which means that
NC**3 > 500*NP
where NP is the total number of processors (NP = NPX*NPY*NPZ). This is in
fact the limit for the smallest possible problem size. The upper limit
depends on the amount of local memory.
Compiling and running the benchmark:
1) To compile and link the benchmark type: `make' for the distributed
version or `make slave' for the single-node version.
2) To run either sequential or distributed version of the benchmark,
type: md1
3) Input parameters NPX, NPY, NPZ and NC (see above) on standard input.
Output from the benchmark is written to the file "result"
3. Accuracy check
-----------------
After the equilibration phase the total energy should be conserved, i.e.
the first few digits should remain constant and there should certainly
be no net drift in the values. The temperature should not vary by more
than 0.1 if the simulation is to produce sensible simulation data. If
it does, the equilibration time has to be increased by increasing the
number of equilbration steps, NEQUIL, in the host-program.
$Id: ReadMe,v 1.2 1994/04/20 17:41:15 igl Rel igl $
Submitted by Mark Papiani,
last updated on 10 Jan 1995.