md.nlist
Overview
Base class neighbor list. |
|
Neighbor list computed via a cell list. |
|
Cell list based neighbor list using stencils. |
|
Bounding volume hierarchy based neighbor list. |
Details
Neighbor list acceleration structures.
Pair forces (hoomd.md.pair
) use neighbor list data structures to perform
efficient calculations. HOOMD-blue provides a several types of neighbor list
construction algorithms that you can select from. Multiple pair force objects
can share a neighbor list, or use independent neighbor list objects. When
neighbor lists are shared, they find neighbors within the the maximum
\(r_{\mathrm{cut},i,j}\) over the associated pair potentials.
- class hoomd.md.nlist.Cell(buffer, exclusions=('bond',), rebuild_check_delay=1, diameter_shift=False, check_dist=True, max_diameter=1.0, deterministic=False)
Neighbor list computed via a cell list.
- Parameters
buffer (float) – Buffer width \([\mathrm{length}]\).
exclusions (tuple[str]) – Defines which particles to exlclude from the neighbor list, see more details in
NList
.rebuild_check_delay (int) – How often to attempt to rebuild the neighbor list.
diameter_shift (bool) – Flag to enable / disable diameter shifting.
check_dist (bool) – Flag to enable / disable distance checking.
max_diameter (float) – The maximum diameter a particle will achieve \([\mathrm{length}]\).
deterministic (bool) – When
True
, sort neighbors to help provide deterministic simulation runs.
Cell
finds neighboring particles using a fixed width cell list, allowing for O(kN) construction of the neighbor list where k is the number of particles per cell. Cells are sized to the largest \(r_\mathrm{cut}\). 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.Examples:
cell = nlist.Cell()
- class hoomd.md.nlist.NList(buffer, exclusions, rebuild_check_delay, diameter_shift, check_dist, max_diameter)
Base class neighbor list.
Methods and attributes provided by this base class are available to all subclasses.
Attention
Users should instantiate the subclasses, using
NList
directly will result in an error.Buffer distance
Set the
buffer
distance to amortize the cost of the neighbor list build. Whenbuffer > 0
, a neighbor list computed on one step can be reused on subsequent steps until a particle moves a distancebuffer/2
. Whencheck_dist
isTrue
,NList
starts checking how far particles have movedrebuild_check_delay
time steps after the last build and performs a rebuild when any particle has moved a distancebuffer/2
. Whencheck_dist
isFalse
,NList
always rebuilds afterrebuild_check_delay
time steps.Exclusions
Neighbor lists nominally include all particles within the specified cutoff distances. The
exclusions
attribute defines which particles will be excluded from the list, even if they are within the cutoff.exclusions
is a tuple of strings that enable one more more types of exclusions. The valid exclusion types are:bond
: Exclude particles that are directly bonded together.angle
: Exclude the first and third particles in each angle.constraint
: Exclude particles that have a distance constraint applied between them.dihedral
: Exclude the first and fourth particles in each dihedral.special_pair
: Exclude particles that are part of a special pair.body
: Exclude particles that belong to the same rigid body.1-3
: Exclude particles i and k whenever there is a bond (i,j) and a bond (j,k).1-4
: Exclude particles i and m whenever there are bonds (i,j), (j,k), and (k,m).
- exclusions
Defines which particles to exlclude from the neighbor list, see more details above.
- property shortest_rebuild
The shortest period between neighbor list rebuilds.
shortest_rebuild
is the smallest number of time steps between neighbor list rebuilds during the previousSimulation.run
.(
Loggable
: category=”scalar”)- Type
- class hoomd.md.nlist.Stencil(cell_width, buffer, exclusions=('bond',), rebuild_check_delay=1, diameter_shift=False, check_dist=True, max_diameter=1.0, deterministic=False)
Cell list based neighbor list using stencils.
- Parameters
cell_width (float) – The underlying stencil bin width for the cell list \([\mathrm{length}]\).
buffer (float) – Buffer width \([\mathrm{length}]\).
exclusions (tuple[str]) – Defines which particles to exlclude from the neighbor list, see more details in
NList
.rebuild_check_delay (int) – How often to attempt to rebuild the neighbor list.
diameter_shift (bool) – Flag to enable / disable diameter shifting.
check_dist (bool) – Flag to enable / disable distance checking.
max_diameter (float) – The maximum diameter a particle will achieve \([\mathrm{length}]\).
deterministic (bool) – When
True
, sort neighbors to help provide deterministic simulation runs.
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 with the given widthcell_width
.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
in how the adjacent cells are searched for particles. One stencil is computed per particle type based on the value of cell_width set by the user, 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 toCell
when there is size disparity in the cutoff radius. The memory demands ofStencil
can also be lower thanCell
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
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 must be set manually.Examples:
nl_s = nlist.Stencil(cell_width=1.5)
- class hoomd.md.nlist.Tree(buffer, exclusions=('bond',), rebuild_check_delay=1, diameter_shift=False, check_dist=True, max_diameter=1.0)
Bounding volume hierarchy based neighbor list.
- Parameters
buffer (float) – Buffer width \([\mathrm{length}]\).
exclusions (tuple[str]) – Defines which particles to exlclude from the neighbor list, see more details in
NList
.rebuild_check_delay (int) – How often to attempt to rebuild the neighbor list.
diameter_shift (bool) – Flag to enable / disable diameter shifting.
check_dist (bool) – Flag to enable / disable distance checking.
max_diameter (float) – The maximum diameter a particle will achieve \([\mathrm{length}]\).
Tree
creates a neighbor list using a bounding volume hierarchy (BVH) tree traversal. 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 thanCell
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_dist=False)