Building neighborlists¶
A neighbor list is a collision detector for spheres: Given a number of spheres of different radius located at different points, it calculates the pairs of spheres that overlap.
ASE provides two implementations of neighbor lists. The newer
linearlyscaling function
neighbor_list()
and
the older quadraticallyscaling class
PrimitiveNeighborList
. The latter will likely
use the former as a backend in the future for linear scaling.
For flexibility, both implementations provide a “primitive”
interface which accepts arrays as arguments rather than the
more complex Atoms
objects.
Both implementations can be used via the NeighborList
class. It also provides easy access to the two implementations methods and functions:

class
ase.neighborlist.
NeighborList
(cutoffs, skin=0.3, sorted=False, self_interaction=True, bothways=False, primitive=<class 'ase.neighborlist.PrimitiveNeighborList'>)[source]¶ Neighbor list object.
 cutoffs: list of float
 List of cutoff radii  one for each atom. If the spheres (defined by
their cutoff radii) of two atoms overlap, they will be counted as
neighbors. See
natural_cutoffs()
for an example on how to get such a list.  skin: float
 If no atom has moved more than the skindistance since the
last call to the
update()
method, then the neighbor list can be reused. This will save some expensive rebuilds of the list, but extra neighbors outside the cutoff will be returned.  self_interaction: bool
 Should an atom return itself as a neighbor?
 bothways: bool
 Return all neighbors. Default is to return only “half” of the neighbors.
 primitive:
PrimitiveNeighborList
orNewPrimitiveNeighborList
class  Define which implementation to use. Older and quadraticallyscaling
PrimitiveNeighborList
or newer and linearlyscalingNewPrimitiveNeighborList
.
Example:
nl = NeighborList([2.3, 1.7]) nl.update(atoms) indices, offsets = nl.get_neighbors(0)

get_neighbors
(a)[source]¶ See
ase.neighborlist.PrimitiveNeighborList.get_neighbors()
orase.neighborlist.PrimitiveNeighborList.get_neighbors()
.

nneighbors
¶ Get number of neighbors.

npbcneighbors
¶ Get number of pbc neighbors.

nupdates
¶ Get number of updates.

class
ase.neighborlist.
PrimitiveNeighborList
(cutoffs, skin=0.3, sorted=False, self_interaction=True, bothways=False, use_scaled_positions=False)[source]¶ Neighbor list that works without Atoms objects.
This is less fancy, but can be used to avoid conversions between scaled and nonscaled coordinates which may affect cell offsets through rounding errors.

build
(pbc, cell, coordinates)[source]¶ Build the list.
Coordinates are taken to be scaled or not according to self.use_scaled_positions.

get_neighbors
(a)[source]¶ Return neighbors of atom number a.
A list of indices and offsets to neighboring atoms is returned. The positions of the neighbor atoms can be calculated like this:
indices, offsets = nl.get_neighbors(42) for i, offset in zip(indices, offsets): print(atoms.positions[i] + dot(offset, atoms.get_cell()))
Notice that if get_neighbors(a) gives atom b as a neighbor, then get_neighbors(b) will not return a as a neighbor  unless bothways=True was used.


class
ase.neighborlist.
NewPrimitiveNeighborList
(cutoffs, skin=0.3, sorted=False, self_interaction=True, bothways=False, use_scaled_positions=False)[source]¶ Neighbor list object. Wrapper around neighbor_list and first_neighbors.
 cutoffs: list of float
 List of cutoff radii  one for each atom. If the spheres (defined by their cutoff radii) of two atoms overlap, they will be counted as neighbors.
 skin: float
 If no atom has moved more than the skindistance since the
last call to the
update()
method, then the neighbor list can be reused. This will save some expensive rebuilds of the list, but extra neighbors outside the cutoff will be returned.  sorted: bool
 Sort neighbor list.
 self_interaction: bool
 Should an atom return itself as a neighbor?
 bothways: bool
 Return all neighbors. Default is to return only “half” of the neighbors.
Example:
nl = NeighborList([2.3, 1.7]) nl.update(atoms) indices, offsets = nl.get_neighbors(0)

get_neighbors
(a)[source]¶ Return neighbors of atom number a.
A list of indices and offsets to neighboring atoms is returned. The positions of the neighbor atoms can be calculated like this:
>>> indices, offsets = nl.get_neighbors(42) >>> for i, offset in zip(indices, offsets): >>> print(atoms.positions[i] + dot(offset, atoms.get_cell()))
Notice that if get_neighbors(a) gives atom b as a neighbor, then get_neighbors(b) will not return a as a neighbor  unless bothways=True was used.

ase.neighborlist.
neighbor_list
(quantities, a, cutoff, self_interaction=False, max_nbins=1000000.0)[source]¶ Compute a neighbor list for an atomic configuration.
Atoms outside periodic boundaries are mapped into the box. Atoms outside nonperiodic boundaries are included in the neighbor list but complexity of neighbor list search for those can become n^2.
The neighbor list is sorted by first atom index ‘i’, but not by second atom index ‘j’.
Parameters:
 quantities: str
Quantities to compute by the neighbor list algorithm. Each character in this string defines a quantity. They are returned in a tuple of the same order. Possible quantities are:
 ‘i’ : first atom index
 ‘j’ : second atom index
 ‘d’ : absolute distance
 ‘D’ : distance vector
 ‘S’ : shift vector (number of cell boundaries crossed by the bond between atom i and j). With the shift vector S, the distances D between atoms can be computed from: D = a.positions[j]a.positions[i]+S.dot(a.cell)
 a:
ase.Atoms
 Atomic configuration.
 cutoff: float or dict
Cutoff for neighbor search. It can be:
 A single float: This is a global cutoff for all elements.
 A dictionary: This specifies cutoff values for element pairs. Specification accepts element numbers of symbols. Example: {(1, 6): 1.1, (1, 1): 1.0, (‘C’, ‘C’): 1.85}
 A list/array with a per atom value: This specifies the radius of
an atomic sphere for each atoms. If spheres overlap, atoms are
within each others neighborhood. See
natural_cutoffs()
for an example on how to get such a list.
 self_interaction: bool
 Return the atom itself as its own neighbor if set to true. Default: False
 max_nbins: int
 Maximum number of bins used in neighbor search. This is used to limit the maximum amount of memory required by the neighbor list.
Returns:
 i, j, …: array
 Tuple with arrays for each quantity specified above. Indices in \(i\) are returned in ascending order 0..len(a), but the order of (i,j) pairs is not guaranteed.
Examples:
Examples assume Atoms object a and numpy imported as np.
Coordination counting:
i = neighbor_list('i', a, 1.85) coord = np.bincount(i)
Coordination counting with different cutoffs for each pair of species:
i = neighbor_list('i', a, {('H', 'H'): 1.1, ('C', 'H'): 1.3, ('C', 'C'): 1.85}) coord = np.bincount(i)
Pair distribution function:
d = neighbor_list('d', a, 10.00) h, bin_edges = np.histogram(d, bins=100) pdf = h/(4*np.pi/3*(bin_edges[1:]**3  bin_edges[:1]**3)) * a.get_volume()/len(a)
Pair potential:
i, j, d, D = neighbor_list('ijdD', a, 5.0) energy = (C/d**6).sum() pair_forces = (6*C/d**5 * (D/d).T).T forces_x = np.bincount(j, weights=pair_forces[:, 0], minlength=len(a))  np.bincount(i, weights=pair_forces[:, 0], minlength=len(a)) forces_y = np.bincount(j, weights=pair_forces[:, 1], minlength=len(a))  np.bincount(i, weights=pair_forces[:, 1], minlength=len(a)) forces_z = np.bincount(j, weights=pair_forces[:, 2], minlength=len(a))  np.bincount(i, weights=pair_forces[:, 2], minlength=len(a))
Dynamical matrix for a pair potential stored in a block sparse format:
from scipy.sparse import bsr_matrix i, j, dr, abs_dr = neighbor_list('ijDd', atoms) energy = (dr.T / abs_dr).T dynmat = (dde * (energy.reshape(1, 3, 1) * energy.reshape(1, 1, 3)).T).T (de / abs_dr * (np.eye(3, dtype=energy.dtype)  (energy.reshape(1, 3, 1) * energy.reshape(1, 1, 3))).T).T dynmat_bsr = bsr_matrix((dynmat, j, first_i), shape=(3*len(a), 3*len(a))) dynmat_diag = np.empty((len(a), 3, 3)) for x in range(3): for y in range(3): dynmat_diag[:, x, y] = np.bincount(i, weights=dynmat[:, x, y]) dynmat_bsr += bsr_matrix((dynmat_diag, np.arange(len(a)), np.arange(len(a) + 1)), shape=(3 * len(a), 3 * len(a)))

ase.neighborlist.
primitive_neighbor_list
(quantities, pbc, cell, positions, cutoff, numbers=None, self_interaction=False, use_scaled_positions=False, max_nbins=1000000.0)[source]¶ Compute a neighbor list for an atomic configuration.
Atoms outside periodic boundaries are mapped into the box. Atoms outside nonperiodic boundaries are included in the neighbor list but complexity of neighbor list search for those can become n^2.
The neighbor list is sorted by first atom index ‘i’, but not by second atom index ‘j’.
Parameters:
 quantities: str
Quantities to compute by the neighbor list algorithm. Each character in this string defines a quantity. They are returned in a tuple of the same order. Possible quantities are
 ‘i’ : first atom index
 ‘j’ : second atom index
 ‘d’ : absolute distance
 ‘D’ : distance vector
 ‘S’ : shift vector (number of cell boundaries crossed by the bond between atom i and j). With the shift vector S, the distances D between atoms can be computed from: D = positions[j]positions[i]+S.dot(cell)
 pbc: array_like
 3tuple indicating giving periodic boundaries in the three Cartesian directions.
 cell: 3x3 matrix
 Unit cell vectors.
 positions: list of xyzpositions
 Atomic positions. Anything that can be converted to an ndarray of shape (n, 3) will do: [(x1,y1,z1), (x2,y2,z2), …]. If use_scaled_positions is set to true, this must be scaled positions.
 cutoff: float or dict
Cutoff for neighbor search. It can be:
 A single float: This is a global cutoff for all elements.
 A dictionary: This specifies cutoff values for element pairs. Specification accepts element numbers of symbols. Example: {(1, 6): 1.1, (1, 1): 1.0, (‘C’, ‘C’): 1.85}
 A list/array with a per atom value: This specifies the radius of
an atomic sphere for each atoms. If spheres overlap, atoms are
within each others neighborhood. See
natural_cutoffs()
for an example on how to get such a list.
 self_interaction: bool
 Return the atom itself as its own neighbor if set to true. Default: False
 use_scaled_positions: bool
 If set to true, positions are expected to be scaled positions.
 max_nbins: int
 Maximum number of bins used in neighbor search. This is used to limit the maximum amount of memory required by the neighbor list.
Returns:
 i, j, … : array
 Tuple with arrays for each quantity specified above. Indices in \(i\) are returned in ascending order 0..len(a)1, but the order of (i,j) pairs is not guaranteed.

neighborlist.
get_connectivity_matrix
(sparse=True)¶ Return connectivity matrix for a given NeighborList (dtype=numpy.int8).
A matrix of shape (nAtoms, nAtoms) will be returned. Connected atoms i and j will have matrix[i,j] == 1, unconnected matrix[i,j] == 0. If bothways=True the matrix will be symmetric, otherwise not!
If sparse is True, a scipy csr matrix is returned. If sparse is False, a numpy matrix is returned.
Note that the old and new neighborlists might give different results for periodic systems if bothways=False.
Example:
Determine which molecule in a system atom 1 belongs to.
>>> from ase import neighborlist >>> from ase.build import molecule >>> from ase.utils import natural_cutoffs >>> from scipy import sparse >>> mol = molecule('CH3CH2OH') >>> cutOff = natural_cutoffs(mol) >>> neighborList = neighborlist.NeighborList(cutOff, self_interaction=False, bothways=True) >>> neighborList.update(mol) >>> matrix = neighborList.get_connectivity_matrix() >>> #or: matrix = neighborlist.get_connectivity_matrix(neighborList.nl) >>> n_components, component_list = sparse.csgraph.connected_components(matrix) >>> idx = 1 >>> molIdx = component_list[idx] >>> print("There are {} molecules in the system".format(n_components)) >>> print("Atom {} is part of molecule {}".format(idx, molIdx)) >>> molIdxs = [ i for i in range(len(component_list)) if component_list[i] == molIdx ] >>> print("The following atoms are part of molecule {}: {}".format(molIdx, molIdxs))