Starting positions: lattices

The idea

Choosing random coordinates sounds easy, but at liquid density some atoms will land almost on top of each other. The Lennard-Jones repulsion grows as r⁻¹², so one close pair gives an enormous force. On the first step those atoms are shot apart, the energy blows up, and the run crashes.

The standard fix is to start from a crystal lattice. Atoms on a regular grid are evenly spaced, so no two are too close.

• Simple cubic (SC): one atom per cubic cell. Easy, but loosely packed.
• Body-centred cubic (BCC): 2 atoms per cell.
• Face-centred cubic (FCC): 4 atoms per cell (8 corners × ⅛ + 6 faces × ½). This is the crystal LJ atoms actually form, and the closest packing. It is the usual choice, so N = 4n³: 32, 108, 256, 500, 864…

The density fixes the box: for N atoms at number density ρ, the box length is L = (N/ρ)^1/3. Then we give the atoms velocities (Lesson 2.3) and let them run. At a liquid temperature the crystal melts within a few picoseconds, and from then on it has no memory of its start.

For molecules such as water or ions in a solvent, place each molecule on a grid point with a random orientation, or use a packing tool such as Packmol. Then run an energy minimisation to remove any remaining close contacts before MD.

The algorithm

Building an FCC start for N = 4n³ atoms at reduced density ρ:

1. Box length L = (N/ρ)^1/3, and cell size a = L/n.
2. Loop over cells (i, j, k) from 0 to n−1.
3. In each cell, place the 4 basis atoms at (0,0,0), (½,½,0), (½,0,½) and (0,½,½), in units of a.
4. Check: the closest pair distance should be a/√2, well above σ at liquid densities.

basis = [(0,0,0), (.5,.5,0), (.5,0,.5), (0,.5,.5)]
for i, j, k in cells:
    for b in basis:
        pos.append(((i, j, k) + b) * a)

Random with rejection (for mixtures and molecules): propose a random position and accept it only if it is at least d_min from every atom already placed. This gets very slow, or fails, near liquid density, which is why packing tools and minimisation exist.