CellListMap: Go from list of pairs to list of indices

General Usage Performance
2

I would of course first ask if you really need that conversion. Usually running over that list of tuples is faster to compute properties for the set than running over a list of lists.

Concerning the conversion, I don’t think you can get much faster than something like this:

julia> function convert(x, list)
          out = [ Int[] for _ in x ]
          for (i,j,d) in list
              push!(out[i], j)
              push!(out[j], i)
          end
          return out
       end

Yet, if you need that kind of output, you could implement it using CellListMap lower level interface directly. Something like this, for example, such that you don’t need the conversion at all:

using CellListMap
using StaticArrays

const T = Vector{Vector{Int}}

function add_pair!(i,j,list::T)
    push!(list[i], j)
    push!(list[j], i)
    return list
end

function reduce_lists(list::T, list_threaded::Vector{T})
    for it in eachindex(list_threaded)
        for i in eachindex(list)
            append!(list[i], list_threaded[it][i])
        end
    end
    return list
end

function run()

  sides = [250,250,250]
  cutoff = 10.
  box = Box(sides,cutoff)

  # positions
  x = [ SVector{3,Float64}(sides .* rand(3)) for i in 1:1500 ]

  # Initialize auxiliary linked lists
  cl = CellList(x,box)

  # Initial, empty, list
  list = [ Int[] for _ in x ]

  # Run pairwise computation
  map_pairwise( 
      (x,y,i,j,d2,list) -> add_pair!(i,j,list),
      list,box,cl;
      reduce=reduce_lists
  )
  return list

end

(there are some additional optimizations that can be used if your calling this within an iterative procedure, to preallocate the threaded output and auxiliary arrays for threading).

Finally, if I’m not mistaken, NearestNeighbors.jl default output is the list of lists you want.

The reason CellListMap.jl provides the output as a list of tuples is exactly because building that list is faster, and running over it is faster as well, if one need to run over all pairs. But of course that may not be the best structure for every application.

1 Like