md.nlist

Overview

md.nlist.cell Cell list based neighbor list
md.nlist.stencil Cell list based neighbor list using stencils
md.nlist.tree Bounding volume hierarchy based neighbor list.

Details

Neighbor list acceleration structures.

Neighbor lists accelerate pair force calculation by maintaining a list of particles within a cutoff radius. Multiple pair forces can utilize the same neighbor list. Neighbors are included using a pairwise cutoff \(r_\mathrm{cut}(i,j)\) that is the maximum of all \(r_\mathrm{cut}(i,j)\) set for the pair forces attached to the list.

Multiple neighbor lists can be created to accelerate simulations where there is significant disparity in \(r_\mathrm{cut}(i,j)\) between pair potentials. If one pair force has a cutoff radius much smaller than another pair force, the pair force calculation for the short cutoff will be slowed down considerably because many particles in the neighbor list will have to be read and skipped because they lie outside the shorter cutoff.

The simplest way to build a neighbor list is \(O(N^2)\): each particle loops over all other particles and only includes those within the neighbor list cutoff. This algorithm is no longer implemented in HOOMD-blue because it is slow and inefficient. Instead, three accelerated algorithms based on cell lists and bounding volume hierarchy trees are implemented. The cell list implementation is usually fastest when the cutoff radius is similar between all pair forces (smaller than 2:1 ratio). The stencil implementation is a different variant of the cell list, and its main use is when a cell list would be faster than a tree but memory demands are too big. The tree implementation is faster when there is large size disparity, but its performance has been improved to be only slightly slower than the cell list for many use cases. Because the performance of these algorithms depends on your system and hardware, you should carefully test which option is fastest for your simulation.

Particles can be excluded from the neighbor list based on certain criteria. Setting \(r_\mathrm{cut}(i,j) \le 0\) will exclude this cross interaction from the neighbor list on build time. Particles can also be excluded by topology or for belonging to the same rigid body (see nlist.reset_exclusions()). To support molecular structures, the body flag can also be used to exclude particles that are not part of a rigid structure. All particles with positive values of the body flag are considered part of a rigid body (see hoomd.md.constrain.rigid), while the default value of -1 indicates that a particle is free. Any other negative value of the body flag indicates that the particles are part of a floppy body; such particles are integrated separately, but are automatically excluded from the neighbor list as well.

Examples:

nl_c = nlist.cell(check_period=1)
nl_t = nlist.tree(r_buff = 0.8)
lj1 = pair.lj(r_cut = 3.0, nlist=nl_c)
lj2 = pair.lj(r_cut = 10.0, nlist=nl_t)
class hoomd.md.nlist.cell(r_buff=0.4, check_period=1, d_max=None, dist_check=True, name=None, deterministic=False)

Cell list based neighbor list

Parameters:
  • r_buff (float) – Buffer width.
  • check_period (int) – How often to attempt to rebuild the neighbor list.
  • d_max (float) – The maximum diameter a particle will achieve, only used in conjunction with slj diameter shifting.
  • dist_check (bool) – Flag to enable / disable distance checking.
  • name (str) – Optional name for this neighbor list instance.
  • deterministic (bool) – When True, enable deterministic runs on the GPU by sorting the cell list.

cell creates a cell list based neighbor list object to which pair potentials can be attached for computing non-bonded pairwise interactions. Cell listing allows for O(N) construction of the neighbor list. Particles are first spatially sorted into cells based on the largest pairwise cutoff radius attached to this instance of the neighbor list. Particles then query their adjacent cells, and neighbors are included based on pairwise cutoffs. This method is very efficient for systems with nearly monodisperse cutoffs, but performance degrades for large cutoff radius asymmetries due to the significantly increased number of particles per cell.

Use base class methods to change parameters (set_params), reset the exclusion list (reset_exclusions) or tune r_buff (tune).

Examples:

nl_c = nlist.cell(check_period = 1)
nl_c.set_params(r_buff=0.5)
nl_c.reset_exclusions([]);
nl_c.tune()

Note

d_max should only be set when slj diameter shifting is required by a pair potential. Currently, slj is the only pair potential requiring this shifting, and setting d_max for other potentials may lead to significantly degraded performance or incorrect results.

class hoomd.md.nlist.nlist

Base class neighbor list.

Methods provided by this base class are available to all subclasses.

add_exclusion(i, j)

Add a specific pair of particles to the exclusion list.

Parameters:
  • i (int) – The tag of the first particle in the pair.
  • j (int) – The tag of the second particle in the pair.

Examples:

nl.add_exclusions(system.particles[0].tag, system.particles[1].tag)
query_update_period()

Query the maximum possible check_period.

query_update_period() examines the counts of nlist rebuilds during the previous hoomd.run(). It returns s-1, where s is the smallest update period experienced during that time. Use it after a medium-length warm up run with check_period=1 to determine what check_period to set for production runs.

Warning

If the previous hoomd.run() was short, insufficient sampling may cause the queried update period to be large enough to result in dangerous builds during longer runs. Unless you use a really long warm up run, subtract an additional 1 from this when you set check_period for additional safety.

reset_exclusions(exclusions=None)

Resets all exclusions in the neighborlist.

Parameters:exclusions (list) – Select which interactions should be excluded from the pair interaction calculation.

By default, the following are excluded from short range pair interactions”

  • Directly bonded particles.
  • Directly constrained particles.
  • Particles that are in the same body (i.e. have the same body flag). Note that these bodies need not be rigid.

reset_exclusions allows the defaults to be overridden to add other exclusions or to remove the exclusion for bonded or constrained particles.

Specify a list of desired types in the exclusions argument (or an empty list to clear all exclusions). All desired exclusions have to be explicitly listed, i.e. ‘1-3’ does not imply ‘1-2’.

Valid types are:

  • bond - Exclude particles that are directly bonded together.
  • constraint - Exclude particles that are directly constrained.
  • angle - Exclude the two outside particles in all defined angles.
  • dihedral - Exclude the two outside particles in all defined dihedrals.
  • pair - Exclude particles in all defined special pairs.
  • body - Exclude particles that belong to the same body.

The following types are determined solely by the bond topology. Every chain of particles in the simulation connected by bonds (1-2-3-4) will be subject to the following exclusions, if enabled, whether or not explicit angles or dihedrals are defined:

  • 1-2 - Same as bond
  • 1-3 - Exclude particles connected with a sequence of two bonds.
  • 1-4 - Exclude particles connected with a sequence of three bonds.

Examples:

nl.reset_exclusions(exclusions = ['1-2'])
nl.reset_exclusions(exclusions = ['1-2', '1-3', '1-4'])
nl.reset_exclusions(exclusions = ['bond', 'angle'])
nl.reset_exclusions(exclusions = ['bond', 'angle','constraint'])
nl.reset_exclusions(exclusions = [])
set_params(r_buff=None, check_period=None, d_max=None, dist_check=True)

Change neighbor list parameters.

Parameters:
  • r_buff (float) – (if set) changes the buffer radius around the cutoff (in distance units)
  • check_period (int) – (if set) changes the period (in time steps) between checks to see if the neighbor list needs updating
  • d_max (float) – (if set) notifies the neighbor list of the maximum diameter that a particle attain over the following run() commands. (in distance units)
  • dist_check (bool) – When set to False, disable the distance checking logic and always regenerate the nlist every check_period steps

set_params() changes one or more parameters of the neighbor list. r_buff and check_period can have a significant effect on performance. As r_buff is made larger, the neighbor list needs to be updated less often, but more particles are included leading to slower force computations. Smaller values of r_buff lead to faster force computation, but more often neighbor list updates, slowing overall performance again. The sweet spot for the best performance needs to be found by experimentation. The default of r_buff = 0.4 works well in practice for Lennard-Jones liquid simulations.

As r_buff is changed, check_period must be changed correspondingly. The neighbor list is updated no sooner than check_period time steps after the last update. If check_period is set too high, the neighbor list may not be updated when it needs to be.

For safety, the default check_period is 1 to ensure that the neighbor list is always updated when it needs to be. Increasing this to an appropriate value for your simulation can lead to performance gains of approximately 2 percent.

check_period should be set so that no particle moves a distance more than r_buff/2.0 during a the check_period. If this occurs, a dangerous build is counted and printed in the neighbor list statistics at the end of a hoomd.run().

When using hoomd.md.pair.slj, d_max MUST be set to the maximum diameter that a particle will attain at any point during the following hoomd.run() commands (see hoomd.md.pair.slj for more information). When using in conjunction, hoomd.md.pair.slj will automatically set d_max for the nlist. This can be overridden (e.g. if multiple potentials using diameters are used) by using set_params() after the hoomd.md.pair.slj class has been initialized.

Caution

When not using hoomd.md.pair.slj, d_max MUST be left at the default value of 1.0 or the simulation will be incorrect if d_max is less than 1.0 and slower than necessary if d_max is greater than 1.0.

Examples:

nl.set_params(r_buff = 0.9)
nl.set_params(check_period = 11)
nl.set_params(r_buff = 0.7, check_period = 4)
nl.set_params(d_max = 3.0)
tune(warmup=200000, r_min=0.05, r_max=1.0, jumps=20, steps=5000, set_max_check_period=False, quiet=False)

Make a series of short runs to determine the fastest performing r_buff setting.

Parameters:
  • warmup (int) – Number of time steps to run() to warm up the benchmark
  • r_min (float) – Smallest value of r_buff to test
  • r_max (float) – Largest value of r_buff to test
  • jumps (int) – Number of different r_buff values to test
  • steps (int) – Number of time steps to run() at each point
  • set_max_check_period (bool) – Set to True to enable automatic setting of the maximum nlist check_period
  • quiet (bool) – Quiet the individual run() calls.

tune() executes warmup time steps. Then it sets the nlist r_buff value to r_min and runs for steps time steps. The TPS value is recorded, and the benchmark moves on to the next r_buff value completing at r_max in jumps jumps. Status information is printed out to the screen, and the optimal r_buff value is left set for further hoomd.run() calls to continue at optimal settings.

Each benchmark is repeated 3 times and the median value chosen. Then, warmup time steps are run again at the optimal r_buff in order to determine the maximum value of check_period. In total, (2*warmup + 3*jump*steps) time steps are run.

Note

By default, the maximum check_period is not set for safety. If you wish to have it set when the call completes, call with the parameter set_max_check_period=True.

Returns:(optimal_r_buff, maximum check_period)
class hoomd.md.nlist.stencil(r_buff=0.4, check_period=1, d_max=None, dist_check=True, cell_width=None, name=None, deterministic=False)

Cell list based neighbor list using stencils

Parameters:
  • r_buff (float) – Buffer width.
  • check_period (int) – How often to attempt to rebuild the neighbor list.
  • d_max (float) – The maximum diameter a particle will achieve, only used in conjunction with slj diameter shifting.
  • dist_check (bool) – Flag to enable / disable distance checking.
  • cell_width (float) – The underlying stencil bin width for the cell list
  • name (str) – Optional name for this neighbor list instance.
  • deterministic (bool) – When True, enable deterministic runs on the GPU by sorting the cell list.

stencil creates a cell list based neighbor list object to which pair potentials can be attached for computing non-bonded pairwise interactions. Cell listing allows for O(N) construction of the neighbor list. Particles are first spatially sorted into cells based on the largest pairwise cutoff radius attached to this instance of the neighbor list.

M.P. Howard et al. 2016 describes this neighbor list implementation in HOOMD-blue. Cite it if you utilize this neighbor list style in your work.

This neighbor-list style differs from cell based on how the adjacent cells are searched for particles. The cell list cell_width is set by default using the shortest active cutoff radius in the system. One stencil is computed per particle type based on the largest cutoff radius that type participates in, which defines the bins that the particle must search in. Distances to the bins in the stencil are precomputed so that certain particles can be quickly excluded from the neighbor list, leading to improved performance compared to cell when there is size disparity in the cutoff radius. The memory demands of stencil can also be lower than cell if your system is large and has many small cells in it; however, tree is usually a better choice for these systems.

The performance of the stencil depends strongly on the choice of cell_width. The best performance is obtained when the cutoff radii are multiples of the cell_width, and when the cell_width covers the simulation box with a roughly integer number of cells. The cell_width can be set manually, or be automatically scanning through a range of possible bin widths using tune_cell_width().

Examples:

nl_s = nlist.stencil(check_period = 1)
nl_s.set_params(r_buff=0.5)
nl_s.reset_exclusions([]);
nl_s.tune()
nl_s.tune_cell_width(min_width=1.5, max_width=3.0)

Note

d_max should only be set when slj diameter shifting is required by a pair potential. Currently, slj is the only pair potential requiring this shifting, and setting d_max for other potentials may lead to significantly degraded performance or incorrect results.

set_cell_width(cell_width)

Set the cell width

Parameters:cell_width (float) – New cell width.
tune_cell_width(warmup=200000, min_width=None, max_width=None, jumps=20, steps=5000)

Make a series of short runs to determine the fastest performing bin width.

Parameters:
  • warmup (int) – Number of time steps to run() to warm up the benchmark
  • min_width (float) – Minimum cell bin width to try
  • max_width (float) – Maximum cell bin width to try
  • jumps (int) – Number of different bin width to test
  • steps (int) – Number of time steps to run() at each point

tune_cell_width() executes warmup time steps. Then it sets the nlist cell_width value to min_width and runs for steps time steps. The TPS value is recorded, and the benchmark moves on to the next cell_width value completing at max_width in jumps jumps. Status information is printed out to the screen, and the optimal cell_width value is left set for further runs() to continue at optimal settings.

Each benchmark is repeated 3 times and the median value chosen. In total, (warmup + 3*jump*steps) time steps are run.

Returns:The optimal cell width.
class hoomd.md.nlist.tree(r_buff=0.4, check_period=1, d_max=None, dist_check=True, name=None)

Bounding volume hierarchy based neighbor list.

Parameters:
  • r_buff (float) – Buffer width.
  • check_period (int) – How often to attempt to rebuild the neighbor list.
  • d_max (float) – The maximum diameter a particle will achieve, only used in conjunction with slj diameter shifting.
  • dist_check (bool) – Flag to enable / disable distance checking.
  • name (str) – Optional name for this neighbor list instance.

tree creates a neighbor list using bounding volume hierarchy (BVH) tree traversal. Pair potentials are attached for computing non-bonded pairwise interactions. A BVH tree of axis-aligned bounding boxes is constructed per particle type, and each particle queries each tree to determine its neighbors. This method of searching leads to significantly improved performance compared to cell lists in systems with moderate size asymmetry, but has slightly poorer performance (10% slower) for monodisperse systems. tree can also be slower than cell if there are multiple types in the system, but the cutoffs between types are identical. (This is because one BVH is created per type.) The user should carefully benchmark neighbor list build times to select the appropriate neighbor list construction type.

M.P. Howard et al. 2016 describes the original implementation of this algorithm for HOOMD-blue. M.P. Howard et al. 2019 describes the improved algorithm that is currently implemented. Cite both if you utilize this neighbor list style in your work.

Examples:

nl_t = nlist.tree(check_period = 1)
nl_t.set_params(r_buff=0.5)
nl_t.reset_exclusions([]);
nl_t.tune()

Note

d_max should only be set when slj diameter shifting is required by a pair potential. Currently, slj is the only pair potential requiring this shifting, and setting d_max for other potentials may lead to significantly degraded performance or incorrect results.