Randomized initial position

Dear community members,

My data file contains atoms that represent a rigid molecule (a helix). I need to have multiple copies of the molecule with randomized initial position for each in the simulation box with no overlap. Can you please introduce to me a tool or recipe on how can I do this?

In an ideal world, I not only need to randomize their initial translational position but also their orientations along all the three axes.


With kind regards,
Alireza Fazeli

There are different protocols to do such a thing, but there is no simple way to get this in a perfect fashion automatically and without any thinking and careful observing.

  1. you can just create a system without worrying about overlap and then “unoverlap” it. in LAMMPS this could be done by editing your data file with a single “molecule” into a “molecule file” and then use create_atoms to create the desired number of molecules into a given volume. then you first use the “soft” pair style and fix adapt starting with a zero potential height and then gradually raising the repulsion so that overlapping molecules/atom will move away from each other. once the system has reached the desired configuration of non-overlapping molecules, you can switch to the intended pair style and follow a normal equilibration protocol. How gradual the unoverlap has to happen and when to switch needs to be found empirically, as it is very system dependent.

  2. you can create regular non-randomized starting geometry, also with a molecule file and create_atoms, but on a regular (simple cubic) lattice. If lattice and initial geometry are chosen properly, there will be no overlap. Then you gradually “break up” the regular structure by running for a period at high temperature until the “memory” of original regular structure is lost. How high and for how long is something that has to be empirically found.

  3. you can use a tool like packmol to generate an initial, densely packed and randomized structure without overlaps. Packmol knows no bond topology, so that would need to be generated later, if necessary. that initial structure could be converted into a LAMMPS data file and then an equilibration started.

  4. there are scripting tools that can build initial structures from templates with a varying degree of sophistication and ease of use. see the pre-post processing tools page on the LAMMPS home page for some examples.

In all cases, it is recommended to start with a structure at a lower than the final density to ease the process and then - during equilibration you can use fix deform to gradually transform your system to the desired density. It is NOT recommended to start with fix npt as non-equlibrium structures will have high pressures and thus will drive the system to much lower than a desired equilibrium density and once a system has significantly expanded, it can take an extended period until it will - on its own - relax into the equilibrium density (e.g. due to jamming, but also from simple statistical mechanical arguments, i.e. due to entropy).