Neighbor lists

Neighbor lists can be computed, returning all pairs of particles that are found within the cutoff, and the corresponding distances.

• Non-periodic systems
• Periodic systems
• In-place computation of neighbor lists
• Options

Non-periodic systems

Without periodic boundary conditions, just provide the coordinates and the cutoff:

using CellListMap
x = [ rand(2) for _ in 1:10_000 ]
neighborlist(positions=x, cutoff=0.01)

15572-element Vector{Tuple{Int64, Int64, Float64}}:
 (2764, 6119, 0.00919483071287996)
 (2764, 6211, 0.00830655317147847)
 (6119, 6211, 0.004276970224300144)
 ⋮
 (3784, 3440, 0.007027402417443701)
 (1914, 4441, 0.007618535761149106)

The returning array contains tuples with the index of the first particle, the second particle, and their distance.

If the neighbor lists between two sets of points are required, use the following notation, where xpositions and ypositions define the input coordinates. In this case using coordinates as arrays of static arrays:

using StaticArrays
x = rand(SVector{3,Float64}, 10^4)
y = rand(SVector{3,Float64}, 10^3)
nblist = neighborlist(xpositions=x, ypositions=y, cutoff=0.01)

46-element Vector{Tuple{Int64, Int64, Float64}}:
 (9759, 803, 0.009712756086369652)
 (8289, 123, 0.007630550668940193)
 (2414, 268, 0.00667155369237361)
 ⋮
 (8368, 543, 0.008357982848126954)
 (4125, 137, 0.009470384866656174)

The returning array contains tuples with the index of the particle in the first vector, the index of the particle in the second vector, and their distance.

Periodic systems

If periodic boundary conditions are used, the unitcell can be provided explicitly as keyword parameters:

x = [ rand(2) for _ in 1:10_000 ]
neighborlist(positions=x, cutoff=0.01, unitcell=[1,1])

15739-element Vector{Tuple{Int64, Int64, Float64}}:
 (3294, 5159, 0.008954090432355115)
 (5671, 6748, 0.006677061239714754)
 (5671, 3011, 0.009814909583849993)
 ⋮
 (9671, 4009, 0.0036407607552604144)
 (2200, 7125, 0.009749184237565761)

In the example above, an Orthorhombic cell was assumed automatically from the fact that a vector of sides was provided. For general periodic boundary conditions, a unit cell matrix must be given:

neighborlist(positions=x, cutoff=0.01, unitcell=[1.0 0.5; 0.5 1.0])

23881-element Vector{Tuple{Int64, Int64, Float64}}:
 (1838, 8936, 0.0015714684699078853)
 (1838, 2112, 0.008118367383449374)
 (7156, 7235, 0.008376098999000004)
 ⋮
 (3206, 6361, 0.008295843067827073)
 (3206, 6481, 0.001991970229297705)

Updating the InPlaceNeighborList object

The update! function provides an interface to modify the coordinates, unitcell, and cutoff, and perform a new neighbor list search. In the example above we have updated the positions, positions, of the system. The full keyword list for updating the system are:

In-place computation of neighbor lists

If neighbor lists are computed within a interactive scenario, it is interesting preallocate all the necessary data and just update the lists at every iteration. This can be achieved by constructing the InPlaceNeighborList object in advance. The performance gain of performing the operations in place might vary and may not be important for single runs, as the allocations do not dominate the computing time.

We can illustrate with an example where sucessive computations of the neighborlists are performed:

julia> using CellListMap

julia> using Chairmarks 

julia> system = InPlaceNeighborList(
           positions=rand(3,100000), 
           cutoff=0.1, 
           unitcell=[1,1,1], 
           parallel=false
       )
InPlaceNeighborList with types: 
CellListMap.CellList{3, Float64}
CellListMap.Box{CellListMap.OrthorhombicCell, 3, Float64, Float64, 9, Float64}
Current list buffer length: 0

julia> @b neighborlist!($system)
262.141 ms (4 allocs: 958.750 MiB, 24.69% gc time, without a warmup)

Note that the initial buffer length has size 0. The first time the neighbor lists are computed, the list will be allocated. Subsequent call of neighborlist! can be completelyl non-allocating, if the number of pairs does not increase:

julia> update!(system; positions=rand(3,100000)) # new positions
InPlaceNeighborList with types: 
CellListMap.CellList{3, Float64}
CellListMap.Box{CellListMap.OrthorhombicCell, 3, Float64, Float64, 9, Float64}
Current list buffer length: 20939808

julia> @b neighborlist!($system)
163.896 ms

For parallel runs, the allocations are minimal, but some small auxiliary data is required for the launching of multiple threads. We illustrate here the convergence of the allocations to the minimum required for multi-threaded calculations.

Options

Input parameters: neighborlist and InPlaceNeighborList

Output

The neighbor list is a Vector{Tuple{Int,Int,T}}, where each element (i, j, d) contains the indices of the two particles and their distance d. The floating-point type T is promoted from the cutoff and unitcell types.

All neighborlist functions return this vector directly.