|
| 1 | +@doc raw""" |
| 2 | + PrecomputedNeighborhoodSearch{NDIMS}(search_radius, n_particles; |
| 3 | + periodic_box_min_corner = nothing, |
| 4 | + periodic_box_max_corner = nothing) |
| 5 | +
|
| 6 | +Neighborhood search with precomputed neighbor lists. A list of all neighbors is computed |
| 7 | +for each particle during initialization and update. |
| 8 | +This neighborhood search maximizes the performance of neighbor loops at the cost of a much |
| 9 | +slower [`update!`](@ref). |
| 10 | +
|
| 11 | +A [`GridNeighborhoodSearch`](@ref) is used internally to compute the neighbor lists during |
| 12 | +initialization and update. |
| 13 | +
|
| 14 | +# Arguments |
| 15 | +- `NDIMS`: Number of dimensions. |
| 16 | +- `search_radius`: The uniform search radius. |
| 17 | +- `n_particles`: Total number of particles. |
| 18 | +
|
| 19 | +# Keywords |
| 20 | +- `periodic_box_min_corner`: In order to use a (rectangular) periodic domain, pass the |
| 21 | + coordinates of the domain corner in negative coordinate |
| 22 | + directions. |
| 23 | +- `periodic_box_max_corner`: In order to use a (rectangular) periodic domain, pass the |
| 24 | + coordinates of the domain corner in positive coordinate |
| 25 | + directions. |
| 26 | +""" |
| 27 | +struct PrecomputedNeighborhoodSearch{NDIMS, NHS, NL, PB} |
| 28 | + neighborhood_search :: NHS |
| 29 | + neighbor_lists :: NL |
| 30 | + periodic_box :: PB |
| 31 | + |
| 32 | + function PrecomputedNeighborhoodSearch{NDIMS}(search_radius, n_particles; |
| 33 | + periodic_box_min_corner = nothing, |
| 34 | + periodic_box_max_corner = nothing) where { |
| 35 | + NDIMS |
| 36 | + } |
| 37 | + nhs = GridNeighborhoodSearch{NDIMS}(search_radius, n_particles, |
| 38 | + periodic_box_min_corner = periodic_box_min_corner, |
| 39 | + periodic_box_max_corner = periodic_box_max_corner) |
| 40 | + |
| 41 | + neighbor_lists = Vector{Vector{Int}}() |
| 42 | + |
| 43 | + new{NDIMS, typeof(nhs), |
| 44 | + typeof(neighbor_lists), |
| 45 | + typeof(nhs.periodic_box)}(nhs, neighbor_lists, nhs.periodic_box) |
| 46 | + end |
| 47 | +end |
| 48 | + |
| 49 | +@inline function Base.ndims(::PrecomputedNeighborhoodSearch{NDIMS}) where {NDIMS} |
| 50 | + return NDIMS |
| 51 | +end |
| 52 | + |
| 53 | +function initialize!(search::PrecomputedNeighborhoodSearch, |
| 54 | + x::AbstractMatrix, y::AbstractMatrix) |
| 55 | + (; neighborhood_search, neighbor_lists) = search |
| 56 | + |
| 57 | + # Initialize grid NHS |
| 58 | + initialize!(neighborhood_search, x, y) |
| 59 | + |
| 60 | + initialize_neighbor_lists!(neighbor_lists, neighborhood_search, x, y) |
| 61 | +end |
| 62 | + |
| 63 | +function update!(search::PrecomputedNeighborhoodSearch, |
| 64 | + x::AbstractMatrix, y::AbstractMatrix; |
| 65 | + particles_moving = (true, true)) |
| 66 | + (; neighborhood_search, neighbor_lists) = search |
| 67 | + |
| 68 | + # Update grid NHS |
| 69 | + update!(neighborhood_search, x, y, particles_moving = particles_moving) |
| 70 | + |
| 71 | + # Skip update if both point sets are static |
| 72 | + if any(particles_moving) |
| 73 | + initialize_neighbor_lists!(neighbor_lists, neighborhood_search, x, y) |
| 74 | + end |
| 75 | +end |
| 76 | + |
| 77 | +function initialize_neighbor_lists!(neighbor_lists, neighborhood_search, x, y) |
| 78 | + # Initialize neighbor lists |
| 79 | + empty!(neighbor_lists) |
| 80 | + resize!(neighbor_lists, size(x, 2)) |
| 81 | + for i in eachindex(neighbor_lists) |
| 82 | + neighbor_lists[i] = Int[] |
| 83 | + end |
| 84 | + |
| 85 | + # Fill neighbor lists |
| 86 | + for_particle_neighbor(x, y, neighborhood_search) do particle, neighbor, _, _ |
| 87 | + push!(neighbor_lists[particle], neighbor) |
| 88 | + end |
| 89 | +end |
| 90 | + |
| 91 | +@inline function foreach_neighbor(f, system_coords, neighbor_system_coords, |
| 92 | + neighborhood_search::PrecomputedNeighborhoodSearch, |
| 93 | + particle; search_radius = nothing) |
| 94 | + (; periodic_box, neighbor_lists) = neighborhood_search |
| 95 | + (; search_radius) = neighborhood_search.neighborhood_search |
| 96 | + |
| 97 | + particle_coords = extract_svector(system_coords, Val(ndims(neighborhood_search)), |
| 98 | + particle) |
| 99 | + for neighbor in neighbor_lists[particle] |
| 100 | + neighbor_coords = extract_svector(neighbor_system_coords, |
| 101 | + Val(ndims(neighborhood_search)), neighbor) |
| 102 | + |
| 103 | + pos_diff = particle_coords - neighbor_coords |
| 104 | + distance2 = dot(pos_diff, pos_diff) |
| 105 | + |
| 106 | + pos_diff, distance2 = compute_periodic_distance(pos_diff, distance2, search_radius, |
| 107 | + periodic_box) |
| 108 | + |
| 109 | + distance = sqrt(distance2) |
| 110 | + |
| 111 | + # Inline to avoid loss of performance |
| 112 | + # compared to not using `for_particle_neighbor`. |
| 113 | + @inline f(particle, neighbor, pos_diff, distance) |
| 114 | + end |
| 115 | +end |
0 commit comments