Mixed boundary conditions

CHANGELOG.md

@@ -5,6 +5,7 @@
 - Two new functions `random_id_in_position` and `random_agent_in_position` can be used to select a random id/agent in a position in discrete spaces (even with filtering). 
 - A new argument `alloc` can be used to select a more performant version in relation to the expensiveness of the filtering for all random methods selecting ids/agents/positions.
 - The `sample!` function is up to 2x faster than before.
+- All spaces support boundaries with mixed periodicity, specified by tuples with a `Bool` value for each dimension, e.g. `GridSpace((5,5); periodic=(true,false))` is periodic along the first dimension but not along the second.
 
 ## BREAKING
 - `NTuple` in `ContinuousSpace` is officially deprecated, but backward compatibility is *mostly* maintained. Known breakages include the comparison of agent position and/or velocity with user-defined tuples, e.g., doing `agent.pos == (0.5, 0.5)`. This will always be `false` in v6 as `agent.pos` is an `SVector`. The rest of the functionality should all work without problems, such as moving agents to tuple-based positions etc.

src/spaces/continuous.jl

@@ -13,11 +13,16 @@ Base.eltype(s::ContinuousSpace{D,P,T,F}) where {D,P,T,F} = T
 no_vel_update(a, m) = nothing
 spacesize(space::ContinuousSpace) = space.extent
 function Base.show(io::IO, space::ContinuousSpace{D,P}) where {D,P}
-    s = "$(P ? "periodic" : "") continuous space with $(spacesize(space)) extent"*
+    periodic = get_periodic_type(space)
+    s = "$(periodic)continuous space with $(spacesize(space)) extent"*
     " and spacing=$(space.spacing)"
     space.update_vel! ≠ no_vel_update && (s *= " with velocity updates")
     print(io, s)
 end
+get_periodic_type(space::ContinuousSpace{D,true}) where {D} = "periodic "
+get_periodic_type(space::ContinuousSpace{D,false}) where {D} = ""
+get_periodic_type(space::ContinuousSpace{D,P}) where {D,P} = "mixed-periodicity "
+
 
 """
     ContinuousAgent{D,T} <: AbstractAgent

src/spaces/grid_general.jl

@@ -157,6 +157,26 @@ function nearby_positions(
                 n .+ pos : mod1.(n .+ pos, space_size) for n in nindices)
     end
 end
+function nearby_positions(
+    pos::ValidPos, space::AbstractGridSpace{D,P}, r = 1,
+    get_indices_f = offsets_within_radius_no_0 # NOT PUBLIC API! For `ContinuousSpace`.
+) where {D,P}
+    stored_ids = space.stored_ids
+    nindices = get_indices_f(space, r)
+    space_size = size(space)
+    # check if we are far from the wall to skip bounds checks
+    if all(i -> r < pos[i] <= space_size[i] - r, 1:D)
+        return (n .+ pos for n in nindices)
+    else
+        return (
+            checkbounds(Bool, stored_ids, (n .+ pos)...) ?
+            n .+ pos : mod1.(n .+ pos, space_size)
+            for n in nindices
+            if all(P[i] || checkbounds(Bool, axes(stored_ids,i), n[i]+pos[i]) for i in 1:D)
+        )
+    end
+end
+
 
 function random_nearby_position(pos::ValidPos, model::ABM{<:AbstractGridSpace{D,false}}, r=1; kwargs...) where {D}
     nindices = offsets_within_radius_no_0(abmspace(model), r)

src/spaces/grid_multi.jl

@@ -64,7 +64,9 @@ example for usage of this advanced specification of dimension-dependent distance
 where one dimension is used as a categorical one.
 """
 function GridSpace(
-        d::NTuple{D,Int}; periodic::Bool = true, metric::Symbol = :chebyshev
+        d::NTuple{D,Int};
+        periodic::Union{Bool,NTuple{D,Bool}} = true,
+        metric::Symbol = :chebyshev
     ) where {D}
     stored_ids = Array{Vector{Int},D}(undef, d)
     for i in eachindex(stored_ids)
@@ -190,6 +192,16 @@ combine_positions(pos, origin, ::GridSpaceIdIterator{false}) = pos .+ origin
     # checking for bounds is less expensive than calling mod1
     checkbounds(Bool, stored_ids, npos...) ? npos : mod1.(npos, space_size)
 end
+@inline function combine_positions(pos, origin, iter::GridSpaceIdIterator{P,D}) where {P,D}
+    npos = pos .+ origin
+    space_size = iter.space_size
+    stored_ids = iter.stored_ids
+    ntuple(i ->
+        (!P[i] || checkbounds(Bool, axes(stored_ids, i), npos[i])) ?
+        npos[i] : mod1(npos[i], space_size[i]),
+        D
+    )
+end
 combine_positions_nocheck(pos, origin, ::GridSpaceIdIterator) = pos .+ origin
 
 # Must return `true` if the access is invalid
@@ -198,6 +210,10 @@ combine_positions_nocheck(pos, origin, ::GridSpaceIdIterator) = pos .+ origin
     return !valid_bounds || @inbounds isempty(iter.stored_ids[pos_index...])
 end
 invalid_access(pos_index, iter::GridSpaceIdIterator{true}) = @inbounds isempty(iter.stored_ids[pos_index...])
+@inline function invalid_access(pos_index, iter::GridSpaceIdIterator{P,D}) where {P,D}
+    valid_bounds = all(P[i] || checkbounds(Bool, axes(iter.stored_ids, i), pos_index[i]) for i in 1:D)
+    return !valid_bounds || @inbounds isempty(iter.stored_ids[pos_index...])
+end
 invalid_access_nocheck(pos_index, iter::GridSpaceIdIterator) = @inbounds isempty(iter.stored_ids[pos_index...])
 
 ##########################################################################################

src/spaces/grid_single.jl

@@ -27,7 +27,10 @@ with non-zero IDs, either positive or negative. This is not checked internally.
 
 All arguments and keywords behave exactly as in [`GridSpace`](@ref).
 """
-function GridSpaceSingle(d::NTuple{D,Int}; periodic = true, metric = :chebyshev) where {D}
+function GridSpaceSingle(d::NTuple{D,Int};
+        periodic::Union{Bool,NTuple{D,Bool}} = true,
+        metric = :chebyshev
+    ) where {D}
     s = zeros(Int, d)
     return GridSpaceSingle{D,periodic}(s, metric,
         Dict{Int,Vector{NTuple{D,Int}}}(),
@@ -117,6 +120,29 @@ function nearby_ids(pos::NTuple{D, Int}, model::ABM{<:GridSpaceSingle{D,false}},
     return ids_iterator
 end
 
+function nearby_ids(pos::NTuple{D, Int}, model::ABM{<:GridSpaceSingle{D,P}}, r = 1,
+        get_offset_indices = offsets_within_radius # internal, see last function
+    ) where {D,P}
+    nindices = get_offset_indices(model, r)
+    stored_ids = abmspace(model).stored_ids
+    space_size = size(stored_ids)
+    position_iterator = (pos .+ β for β in nindices)
+    # check if we are far from the wall to skip bounds checks
+    if all(i -> r < pos[i] <= space_size[i] - r, 1:D)
+        ids_iterator = (stored_ids[p...] for p in position_iterator
+                        if stored_ids[p...] != 0)
+    else
+        ids_iterator = (
+            checkbounds(Bool, stored_ids, p...) ?
+            stored_ids[p...] : stored_ids[mod1.(p, space_size)...]
+            for p in position_iterator
+            if stored_ids[mod1.(p, space_size)...] != 0 &&
+            all(P[i] || checkbounds(Bool, axes(stored_ids, i), p[i]) for i in 1:D)
+        )
+    end
+    return ids_iterator
+end
+
 # Contrary to `GridSpace`, we also extend here `nearby_ids(a::Agent)`.
 # Because, we know that one agent exists per position, and hence we can skip the
 # call to `filter(id -> id ≠ a.id, ...)` that happens in core/space_interaction_API.jl.

src/spaces/utilities.jl

@@ -43,6 +43,31 @@ function euclidean_distance(
     sqrt(sum(min.(direct, spacesize(space) .- direct).^2))
 end
 
+function euclidean_distance(
+    p1::ValidPos,
+    p2::ValidPos,
+    space::Union{ContinuousSpace{D,P},AbstractGridSpace{D,P}}
+) where {D,P}
+    s = spacesize(space)
+    distance_squared = zero(eltype(p1))
+    for i in eachindex(p1)
+        if P[i]
+            distance_squared += euclidean_distance_periodic(p1[i], p2[i], s[i])^2
+        else
+            distance_squared += euclidean_distance_direct(p1[i], p2[i])^2
+        end
+    end
+    return sqrt(distance_squared)
+end
+
+function euclidean_distance_direct(x1::Real, x2::Real)
+    abs(x1 - x2)
+end
+function euclidean_distance_periodic(x1::Real, x2::Real, l::Real)
+    direct = abs(x1 - x2)
+    min(direct, l - direct)
+end
+
 """
     manhattan_distance(a, b, model::ABM)
 
@@ -60,16 +85,40 @@ function manhattan_distance(
     p2::ValidPos,
     space::Union{ContinuousSpace{D,false},AbstractGridSpace{D,false}},
 ) where {D}
-    sum(abs.(p1 .- p2))
+    sum(manhattan_distance_direct.(p1, p2))
 end
 
 function manhattan_distance(
     p1::ValidPos,
     p2::ValidPos,
     space::Union{ContinuousSpace{D,true},AbstractGridSpace{D,true}}
 ) where {D}
-    direct = abs.(p1 .- p2)
-    sum(min.(direct, spacesize(space) .- direct))
+    sum(manhattan_distance_periodic.(p1, p2, spacesize(space)))
+end
+
+function manhattan_distance(
+    p1::ValidPos,
+    p2::ValidPos,
+    space::Union{ContinuousSpace{D,P},AbstractGridSpace{D,P}}
+) where {D,P}
+    s = spacesize(space)
+    distance = zero(eltype(p1))
+    for i in eachindex(p1)
+        if P[i]
+            distance += manhattan_distance_periodic(p1[i], p2[i], s[i])
+        else
+            distance += manhattan_distance_direct(p1[i], p2[i])
+        end
+    end
+    return distance
+end
+
+function manhattan_distance_direct(x1::Real, x2::Real)
+    abs(x1 - x2)
+end
+function manhattan_distance_periodic(x1::Real, x2::Real, s::Real)
+    direct = abs(x1 - x2)
+    min(direct, s - direct)
 end
 
 """
@@ -97,6 +146,21 @@ function get_direction(
     return to .- from
 end
 
+function get_direction(
+    from::ValidPos,
+    to::ValidPos,
+    space::Union{ContinuousSpace{D,P},AbstractGridSpace{D,P}}
+) where {D,P}
+    direct_dir = to .- from
+    inverse_dir = direct_dir .- sign.(direct_dir) .* spacesize(space)
+    return map(
+        i -> P[i] ?
+        (abs(direct_dir[i]) <= abs(inverse_dir[i]) ? direct_dir[i] : inverse_dir[i]) :
+        direct_dir[i],
+        1:D
+    )
+end
+
 #######################################################################################
 # %% Utilities for graph-based spaces (Graph/OpenStreetMap)
 #######################################################################################

src/spaces/walk.jl

@@ -82,15 +82,31 @@ function normalize_position(pos::SVector{D}, space::ContinuousSpace{D,false}) wh
     return clamp.(pos, 0.0, prevfloat.(spacesize(space)))
 end
 
+function normalize_position(pos::SVector{D}, space::ContinuousSpace{D,P}) where {D,P}
+    s = spacesize(space)
+    return SVector{D}(
+        P[i] ? mod(pos[i], s[i]) : clamp(pos[i], 0.0, prevfloat(s[i]))
+        for i in 1:D
+    )
+end
+
 #----
-# for backward compatibility 
+# for backward compatibility
 function normalize_position(pos::NTuple{D}, space::ContinuousSpace{D,true}) where {D}
     return Tuple(mod.(pos, spacesize(space)))
 end
 
 function normalize_position(pos::NTuple{D}, space::ContinuousSpace{D,false}) where {D}
     return Tuple(clamp.(pos, 0.0, prevfloat.(spacesize(space))))
 end
+
+function normalize_position(pos::NTuple{D}, space::ContinuousSpace{D,P}) where {D,P}
+    s = spacesize(space)
+    return ntuple(
+        i -> P[i] ? mod(pos[i], s[i]) : clamp(pos[i], 0.0, prevfloat(s[i])),
+        D
+    )
+end
 #----
 
 function normalize_position(pos::ValidPos, space::AbstractGridSpace{D,true}) where {D}
@@ -101,6 +117,14 @@ function normalize_position(pos::ValidPos, space::AbstractGridSpace{D,false}) wh
     return clamp.(pos, 1, spacesize(space))
 end
 
+function normalize_position(pos::ValidPos, space::AbstractGridSpace{D,P}) where {D,P}
+    s = spacesize(space)
+    return ntuple(
+        i -> P[i] ? mod1(pos[i], s[i]) : clamp(pos[i], 1, s[i]),
+        D
+    )
+end
+
 #######################################################################################
 # %% Random walks
 #######################################################################################

src/submodules/pathfinding/astar.jl

@@ -63,7 +63,7 @@ Utilization of all features of `AStar` occurs in the
 """
 function AStar(
     dims::NTuple{D,T};
-    periodic::Bool = false,
+    periodic::Union{Bool,NTuple{D,Bool}} = false,
     diagonal_movement::Bool = true,
     admissibility::Float64 = 0.0,
     walkmap::BitArray{D} = trues(dims),
@@ -115,13 +115,16 @@ function vonneumann_neighborhood(D)
 end
 
 function Base.show(io::IO, pathfinder::AStar{D,P,M}) where {D,P,M}
-    periodic = P ? "periodic, " : ""
+    periodic = get_periodic_type(pathfinder)
     moore = M ? "diagonal, " : "orthogonal, "
     s =
         "A* in $(D) dimensions, $(periodic)$(moore)ϵ=$(pathfinder.admissibility), " *
         "metric=$(pathfinder.cost_metric)"
     print(io, s)
 end
+get_periodic_type(::AStar{D,false,M}) where {D,M} = ""
+get_periodic_type(::AStar{D,true,M}) where {D,M} = "periodic, "
+get_periodic_type(::AStar{D,P,M}) where {D,P,M} = "mixed periodicity, "
 
 struct GridCell
     f::Int
@@ -195,6 +198,13 @@ end
     (mod1.(cur .+ β.I, size(pathfinder.walkmap)) for β in pathfinder.neighborhood)
 @inline get_neighbors(cur, pathfinder::AStar{D,false}) where {D} =
     (cur .+ β.I for β in pathfinder.neighborhood)
+@inline function get_neighbors(cur, pathfinder::AStar{D,P}) where {D,P}
+    s = size(pathfinder.walkmap)
+    (
+        ntuple(i -> P[i] ? mod1(cur[i] + β[i], s[i]) : cur[i] + β[i], D)
+        for β in pathfinder.neighborhood
+    )
+end
 @inline inbounds(n, pathfinder, closed) =
     all(1 .<= n .<= size(pathfinder.walkmap)) && pathfinder.walkmap[n...] && n ∉ closed
 

src/submodules/pathfinding/astar_continuous.jl

@@ -7,6 +7,14 @@ sqr_distance(from, to, pathfinder::AStar{D,true}) where {D} =
     sum(min.(abs.(from .- to), pathfinder.dims .- abs.(from .- to)) .^ 2)
 sqr_distance(from, to, pathfinder::AStar{D,false}) where {D} =
     sum((from .- to) .^ 2)
+@inline function sqr_distance(from, to, pathfinder::AStar{D,P}) where {D,P}
+    s = pathfinder.dims
+    delta = abs.(from .- to)
+    sum(
+        P[i] ? (min(delta[i], s[i] - delta[i]))^2 : delta[i]^2
+        for i in 1:D
+    )
+end
 
 """
     find_continuous_path(pathfinder, from, to)

src/submodules/pathfinding/pathfinding_utils.jl

@@ -9,6 +9,12 @@ position_delta(pathfinder::GridPathfinder{D,true}, from::Dims{D}, to::Dims{D}) w
 position_delta(pathfinder::GridPathfinder{D,false}, from::Dims{D}, to::Dims{D}) where {D} =
     abs.(to .- from)
 
+@inline function position_delta(pathfinder::GridPathfinder{D,P}, from::Dims{D}, to::Dims{D}) where {D,P}
+    s = size(pathfinder.walkmap)
+    delta = abs.(to .- from)
+    ntuple(i -> P[i] ? min(delta[i], s[i] - delta[i]) : delta[i], D)
+end
+
 """
     Pathfinding.delta_cost(pathfinder::GridPathfinder{D}, metric::M, from, to) where {M<:CostMetric}
 Calculate an approximation for the cost of travelling from `from` to `to` (both of

test/astar_tests.jl

@@ -93,6 +93,18 @@ using Agents.Pathfinding
             a = add_agent!((1, 3), model, 0.)
             @test plan_best_route!(a, [(1, 3), (4, 1)], model.pf) == (1, 3)
             @test isnothing(plan_best_route!(a, [(5, 3), (4, 1)], model.pf))
+
+            sp = GridSpace((5, 5); periodic = (true, false))
+            pf = AStar(sp)
+            model = ABM(Agent3, sp; properties = (pf = pf,))
+            model.pf.walkmap[3, :] .= 0
+            model.pf.walkmap[:, 2] .= 0
+            a = add_agent!((1, 3), model, 0.)
+            # agent cannot move to (1,1) because y is not periodic
+            @test isnothing(plan_best_route!(a, [(1,1)], model.pf))
+            # but it can move to (5,3) because x is periodic
+            # and (5,3) is only one step away while (2,4) is two steps away
+            @test plan_best_route!(a, [(5,3), (2,4)], model.pf) == (5,3)
         end
 
         @testset "ContinuousSpace" begin
@@ -133,6 +145,7 @@ using Agents.Pathfinding
             move_along_route!(a, model, model.pf, 1.0)
             @test all(isapprox.(a.pos, (0.7071, 0.7071); atol))
 
+
             model.pf.walkmap[:, 3] .= 0
             move_agent!(a, SVector(2.5, 2.5), model)
             @test all(plan_best_route!(a, SVector.([(3., 0.3), (2.5, 2.5)]), model.pf) .≈ (2.5, 2.5))
@@ -149,6 +162,17 @@ using Agents.Pathfinding
             @test all(get_spatial_property(random_walkable(model, model.pf), model.pf.walkmap, model) for _ in 1:10)
             rpos = [random_walkable(SVector(2.5, 0.75), model, model.pf, 2.0) for _ in 1:50]
             @test all(get_spatial_property(x, model.pf.walkmap, model) && euclidean_distance(x, SVector(2.5, 0.75), model) <= 2.0 + atol for x in rpos)
+
+            # mixed boundary
+            space = ContinuousSpace((5., 5.); periodic = (false, true))
+            pathfinder = AStar(space; walkmap = trues(10, 10))
+            model = ABM(Agent6, space; properties = (pf = pathfinder,))
+            model.pf.walkmap[5,:] .= 0
+            a = add_agent!(SVector(0., 0.), model, SVector(0., 0.), 0.)
+            @test isnothing(plan_best_route!(a, [SVector(4.0, 0.)], model.pf))
+            @test plan_best_route!(a, [SVector(1.0, 1.0)], model.pf) == SVector(1.0, 1.0)
+            move_along_route!(a, model, model.pf, 1.0)
+            @test isapprox(a.pos, SVector(1/sqrt(2), 1/sqrt(2)); atol)
         end
     end