From b27aa2b1aa4bc91586e27b8e0c5577631f66f368 Mon Sep 17 00:00:00 2001 From: Erik Faulhaber <44124897+efaulhaber@users.noreply.github.com> Date: Wed, 7 Feb 2024 16:27:41 +0100 Subject: [PATCH] Update docs for NHS (#368) * Update docs for NHS * Rename kwargs --------- Co-authored-by: Niklas Neher <73897120+LasNikas@users.noreply.github.com> --- src/general/semidiscretization.jl | 21 ++++---- src/neighborhood_search/grid_nhs.jl | 68 ++++++++++++++++++++------ src/neighborhood_search/trivial_nhs.jl | 53 +++++++++++++++++--- test/neighborhood_search/grid_nhs.jl | 14 +++--- 4 files changed, 119 insertions(+), 37 deletions(-) diff --git a/src/general/semidiscretization.jl b/src/general/semidiscretization.jl index 85ed6cf2f..57b6e4bc5 100644 --- a/src/general/semidiscretization.jl +++ b/src/general/semidiscretization.jl @@ -15,10 +15,12 @@ the keyword argument `neighborhood_search`. A value of `nothing` means no neighb By default, the [`GridNeighborhoodSearch`](@ref) is used. Use [`TrivialNeighborhoodSearch`](@ref) to loop over all particles (no neighborhood search). -- `periodic_box_min_corner`: In order to use a periodic domain, pass the coordinates - of the domain corner in negative coordinate directions. -- `periodic_box_max_corner`: In order to use a periodic domain, pass the coordinates - of the domain corner in positive coordinate directions. +- `periodic_box_min_corner`: In order to use a (rectangular) periodic domain, pass the + coordinates of the domain corner in negative coordinate + directions. +- `periodic_box_max_corner`: In order to use a (rectangular) periodic domain, pass the + coordinates of the domain corner in positive coordinate + directions. # Examples ```julia @@ -95,18 +97,19 @@ function Base.show(io::IO, ::MIME"text/plain", semi::Semidiscretization) end function create_neighborhood_search(system, neighbor, ::Val{nothing}, - min_corner, max_corner) + periodic_box_min_corner, periodic_box_max_corner) radius = compact_support(system, neighbor) TrivialNeighborhoodSearch{ndims(system)}(radius, eachparticle(neighbor), - min_corner=min_corner, max_corner=max_corner) + periodic_box_min_corner=periodic_box_min_corner, + periodic_box_max_corner=periodic_box_max_corner) end function create_neighborhood_search(system, neighbor, ::Val{GridNeighborhoodSearch}, - min_corner, max_corner) + periodic_box_min_corner, periodic_box_max_corner) radius = compact_support(system, neighbor) search = GridNeighborhoodSearch{ndims(system)}(radius, nparticles(neighbor), - min_corner=min_corner, - max_corner=max_corner) + periodic_box_min_corner=periodic_box_min_corner, + periodic_box_max_corner=periodic_box_max_corner) # Initialize neighborhood search initialize!(search, initial_coordinates(neighbor)) diff --git a/src/neighborhood_search/grid_nhs.jl b/src/neighborhood_search/grid_nhs.jl index c1e30403d..2b9470cfb 100644 --- a/src/neighborhood_search/grid_nhs.jl +++ b/src/neighborhood_search/grid_nhs.jl @@ -1,7 +1,9 @@ @doc raw""" - GridNeighborhoodSearch{NDIMS}(search_radius, n_particles) + GridNeighborhoodSearch{NDIMS}(search_radius, n_particles; + periodic_box_min_corner=nothing, + periodic_box_max_corner=nothing) -Simple grid-based neighborhood search with uniform search radius ``d``. +Simple grid-based neighborhood search with uniform search radius. The domain is divided into a regular grid. For each (non-empty) grid cell, a list of particles in this cell is stored. Instead of representing a finite domain by an array of cells, a potentially infinite domain @@ -11,15 +13,49 @@ indexed by the cell index tuple \left( \left\lfloor \frac{x}{d} \right\rfloor, \left\lfloor \frac{y}{d} \right\rfloor \right) \quad \text{or} \quad \left( \left\lfloor \frac{x}{d} \right\rfloor, \left\lfloor \frac{y}{d} \right\rfloor, \left\lfloor \frac{z}{d} \right\rfloor \right), ``` -where ``x, y, z`` are the space coordinates. +where ``x, y, z`` are the space coordinates and ``d`` is the search radius. -To find particles within a radius around a point, only particles in the neighboring cells -are considered. +To find particles within the search radius around a point, only particles in the neighboring +cells are considered. See also (Chalela et al., 2021), (Ihmsen et al. 2011, Section 4.4). As opposed to (Ihmsen et al. 2011), we do not sort the particles in any way, -since that makes our implementation a lot faster (although less parallelizable). +since not sorting makes our implementation a lot faster (although less parallelizable). + +# Arguments +- `NDIMS`: Number of dimensions. +- `search_radius`: The uniform search radius. +- `n_particles`: Total number of particles. + +# Keywords +- `periodic_box_min_corner`: In order to use a (rectangular) periodic domain, pass the + coordinates of the domain corner in negative coordinate + directions. +- `periodic_box_max_corner`: In order to use a (rectangular) periodic domain, pass the + coordinates of the domain corner in positive coordinate + directions. + +!!! warning "Internal use only" + Please note that this constructor is intended for internal use only. It is *not* part of + the public API of TrixiParticles.jl, and it thus can altered (or be removed) at any time + without it being considered a breaking change. + + To run a simulation with this neighborhood search, just pass the type to the constructor + of [`Semidiscretization`](@ref): + ```julia + semi = Semidiscretization(system1, system2, system3, + neighborhood_search=GridNeighborhoodSearch) + ``` + The keyword arguments `periodic_box_min_corner` and `periodic_box_max_corner` explained + above can also be passed to the [`Semidiscretization`](@ref) and will internally be + forwarded to the neighborhood search: + ```julia + semi = Semidiscretization(system1, system2, + neighborhood_search=GridNeighborhoodSearch, + periodic_box_min_corner=[0.0, -0.25], + periodic_box_max_corner=[1.0, 0.75])) + ``` ## References: - M. Chalela, E. Sillero, L. Pereyra, M.A. Garcia, J.B. Cabral, M. Lares, M. Merchán. @@ -42,8 +78,8 @@ struct GridNeighborhoodSearch{NDIMS, ELTYPE, PB} cell_size :: NTuple{NDIMS, ELTYPE} function GridNeighborhoodSearch{NDIMS}(search_radius, n_particles; - min_corner=nothing, - max_corner=nothing) where {NDIMS} + periodic_box_min_corner=nothing, + periodic_box_max_corner=nothing) where {NDIMS} ELTYPE = typeof(search_radius) hashtable = Dict{NTuple{NDIMS, Int}, Vector{Int}}() @@ -51,13 +87,15 @@ struct GridNeighborhoodSearch{NDIMS, ELTYPE, PB} cell_buffer = Array{NTuple{NDIMS, Int}, 2}(undef, n_particles, Threads.nthreads()) cell_buffer_indices = zeros(Int, Threads.nthreads()) - if (min_corner === nothing && max_corner === nothing) || search_radius < eps() + if search_radius < eps() || + (periodic_box_min_corner === nothing && periodic_box_max_corner === nothing) + # No periodicity periodic_box = nothing n_cells = ntuple(_ -> -1, Val(NDIMS)) cell_size = ntuple(_ -> search_radius, Val(NDIMS)) - elseif min_corner !== nothing && max_corner !== nothing - periodic_box = PeriodicBox(min_corner, max_corner) + elseif periodic_box_min_corner !== nothing && periodic_box_max_corner !== nothing + periodic_box = PeriodicBox(periodic_box_min_corner, periodic_box_max_corner) # Round up search radius so that the grid fits exactly into the domain without # splitting any cells. This might impact performance slightly, since larger @@ -73,8 +111,8 @@ struct GridNeighborhoodSearch{NDIMS, ELTYPE, PB} "Please use no NHS for very small problems.")) end else - throw(ArgumentError("`min_corner` and `max_corner` must either be " * - "both `nothing` or both an array or tuple")) + throw(ArgumentError("`periodic_box_min_corner` and `periodic_box_max_corner` " * + "must either be both `nothing` or both an array or tuple")) end new{NDIMS, ELTYPE, @@ -333,8 +371,8 @@ function copy_neighborhood_search(nhs::GridNeighborhoodSearch, search_radius, u) search = GridNeighborhoodSearch{ndims(nhs)}(search_radius, nparticles(nhs)) else search = GridNeighborhoodSearch{ndims(nhs)}(search_radius, nparticles(nhs), - min_corner=nhs.periodic_box.min_corner, - max_corner=nhs.periodic_box.max_corner) + periodic_box_min_corner=nhs.periodic_box.min_corner, + periodic_box_max_corner=nhs.periodic_box.max_corner) end # Initialize neighborhood search diff --git a/src/neighborhood_search/trivial_nhs.jl b/src/neighborhood_search/trivial_nhs.jl index c915cfd56..ffa5aa03c 100644 --- a/src/neighborhood_search/trivial_nhs.jl +++ b/src/neighborhood_search/trivial_nhs.jl @@ -2,6 +2,43 @@ TrivialNeighborhoodSearch{NDIMS}(search_radius, eachparticle) Trivial neighborhood search that simply loops over all particles. +The search radius still needs to be passed in order to sort out particles outside the +search radius in the internal function `for_particle_neighbor`, but it's not used in the +internal function `eachneighbor`. + +# Arguments +- `NDIMS`: Number of dimensions. +- `search_radius`: The uniform search radius. +- `eachparticle`: `UnitRange` of all particle indices. Usually just `1:n_particles`. + +# Keywords +- `periodic_box_min_corner`: In order to use a (rectangular) periodic domain, pass the + coordinates of the domain corner in negative coordinate + directions. +- `periodic_box_max_corner`: In order to use a (rectangular) periodic domain, pass the + coordinates of the domain corner in positive coordinate + directions. + +!!! warning "Internal use only" + Please note that this constructor is intended for internal use only. It is *not* part of + the public API of TrixiParticles.jl, and it thus can altered (or be removed) at any time + without it being considered a breaking change. + + To run a simulation with this neighborhood search, just pass the type to the constructor + of [`Semidiscretization`](@ref): + ```julia + semi = Semidiscretization(system1, system2, system3, + neighborhood_search=TrivialNeighborhoodSearch) + ``` + The keyword arguments `periodic_box_min_corner` and `periodic_box_max_corner` explained + above can also be passed to the [`Semidiscretization`](@ref) and will internally be + forwarded to the neighborhood search: + ```julia + semi = Semidiscretization(system1, system2, + neighborhood_search=TrivialNeighborhoodSearch, + periodic_box_min_corner=[0.0, -0.25], + periodic_box_max_corner=[1.0, 0.75])) + ``` """ struct TrivialNeighborhoodSearch{NDIMS, ELTYPE, EP, PB} search_radius :: ELTYPE @@ -9,16 +46,18 @@ struct TrivialNeighborhoodSearch{NDIMS, ELTYPE, EP, PB} periodic_box :: PB function TrivialNeighborhoodSearch{NDIMS}(search_radius, eachparticle; - min_corner=nothing, - max_corner=nothing) where {NDIMS} - if (min_corner === nothing && max_corner === nothing) || search_radius < eps() + periodic_box_min_corner=nothing, + periodic_box_max_corner=nothing) where {NDIMS} + if search_radius < eps() || + (periodic_box_min_corner === nothing && periodic_box_max_corner === nothing) + # No periodicity periodic_box = nothing - elseif min_corner !== nothing && max_corner !== nothing - periodic_box = PeriodicBox(min_corner, max_corner) + elseif periodic_box_min_corner !== nothing && periodic_box_max_corner !== nothing + periodic_box = PeriodicBox(periodic_box_min_corner, periodic_box_max_corner) else - throw(ArgumentError("`min_corner` and `max_corner` must either be " * - "both `nothing` or both an array or tuple")) + throw(ArgumentError("`periodic_box_min_corner` and `periodic_box_max_corner` " * + "must either be both `nothing` or both an array or tuple")) end new{NDIMS, typeof(search_radius), diff --git a/test/neighborhood_search/grid_nhs.jl b/test/neighborhood_search/grid_nhs.jl index d0326d23c..d686efc5b 100644 --- a/test/neighborhood_search/grid_nhs.jl +++ b/test/neighborhood_search/grid_nhs.jl @@ -135,7 +135,8 @@ # 3 x 6 cells nhs = GridNeighborhoodSearch{2}(0.1, size(coords, 2), - min_corner=[-0.1, -0.2], max_corner=[0.2, 0.4]) + periodic_box_min_corner=[-0.1, -0.2], + periodic_box_max_corner=[0.2, 0.4]) TrixiParticles.initialize!(nhs, coords) @@ -173,8 +174,8 @@ # 3 x 6 cells nhs = GridNeighborhoodSearch{2}(0.1, size(coords, 2), - min_corner=[-0.1, -0.2], - max_corner=[0.205, 0.43]) + periodic_box_min_corner=[-0.1, -0.2], + periodic_box_max_corner=[0.205, 0.43]) TrixiParticles.initialize!(nhs, coords) @@ -220,7 +221,8 @@ # 5 x 1 cells nhs = GridNeighborhoodSearch{2}(1.0, size(coords, 2), - min_corner=[-1.5, 0.0], max_corner=[2.5, 3.0]) + periodic_box_min_corner=[-1.5, 0.0], + periodic_box_max_corner=[2.5, 3.0]) TrixiParticles.initialize!(nhs, coords) @@ -239,8 +241,8 @@ # 3 x 6 x 3 cells nhs = GridNeighborhoodSearch{3}(0.1, size(coords, 2), - min_corner=[-0.1, -0.2, 0.05], - max_corner=[0.2, 0.4, 0.35]) + periodic_box_min_corner=[-0.1, -0.2, 0.05], + periodic_box_max_corner=[0.2, 0.4, 0.35]) TrixiParticles.initialize!(nhs, coords)