Mixed boundary conditions

[mastrof]

Closes #828. WIP, there are still a lot of things to fix. I tried to represent space periodicity with a new type Periodic{D} to only be used internally, which would just wrap a SVector with a Bool value for each spatial dimension. When given a Bool e.g. Periodic{D}(true) it is internally converted to a SVector as fill(true, SVector{D}). I don't have a smart way to fall back to "full periodic" and "non periodic" via dispatch yet, so as it is now every time distances have to be evaluated, each dimensions needs to be checked for whether or not it is periodic. In low-dimensional systems (which I assume make up the vast majority of use cases) there is no performance difference; in high-dimensional ones it might start to become noticeable (although I don't know to what extent this is really an issue...)

Low-dimensional space:

julia> # mixed-boundary
       space = GridSpace(ntuple(i->5,2); periodic=true)
       model = ABM(GridAgent{2}, space)
       add_agent!(model)
       add_agent!(model)
       @btime euclidean_distance($(model[1]), $(model[2]), $(model))
  2.374 ns (0 allocations: 0 bytes)

julia> # main
       space = GridSpace(ntuple(i->5,2); periodic=true)
       model = ABM(GridAgent{2}, space)
       add_agent!(model)
       add_agent!(model)
       @btime euclidean_distance($(model[1]), $(model[2]), $(model))
  2.304 ns (0 allocations: 0 bytes)

High-dimensional space:

julia> # mixed-boundary
       space = GridSpace(ntuple(i->5,8); periodic=true)
       model = ABM(GridAgent{8}, space)
       add_agent!(model)
       add_agent!(model)
       @btime euclidean_distance($(model[1]), $(model[2]), $(model))
  8.389 ns (0 allocations: 0 bytes)

julia> # main
       space = GridSpace(ntuple(i->5,8); periodic=true)
       model = ABM(GridAgent{8}, space)
       add_agent!(model)
       add_agent!(model)
       @btime euclidean_distance($(model[1]), $(model[2]), $(model))
  5.168 ns (0 allocations: 0 bytes)

An alternative approach (probably better and possibly requiring less work) is to not have a Periodic type, just let spaces accept whatever is passed to them (be it a Bool, a tuple an svector etc...). If it's found to be a Bool, then they will just behave as usual, but if they are passed an iterable type they will do the checking.

[mastrof]

An alternative approach (probably better and possibly requiring less work) is to not have a Periodic type, just let spaces accept whatever is passed to them (be it a Bool, a tuple an svector etc...). If it's found to be a Bool, then they will just behave as usual, but if they are passed an iterable type they will do the checking.

Ok I think this is the way to go. With last commit there absolutely no changes for everything that uses periodic=true or periodic=false. Now I will just need to have some fun adding tests everywhere for the mixed boundary case.

For ContinuousSpace with a mixed boundary, I also reimplemented the display function to write "mixed-periodicity continuous space...". Maybe it's better to explicitly write which dimensions are periodic? I wonder why the display for continuous spaces is different from the more schematic one adopted for grid spaces (where now the display will show e.g. periodic=(true,false))?

[Tortar]

Great PR! It seems that the approach is very clean!

[mastrof]

Is there anything I'm missing concerning places where periodicity is used to compute things? I tried grepping my way through the package but didn't see anything obvious. Although I don't really know too well what goes on e.g. in the nearby id searches or in the pathfinder; I guess missing features will come up while I write tests, but if you know/notice something missing please point it out :)

[Tortar]

you are right, the nearby searches for GridSpaces will be surely need to be extended, since the filtering which is done for finding the right nearby positions should be done differently

[codecov-commenter]

Codecov Report

Merging #879 (4804c53) into main (f21211b) will increase coverage by 0.57%. Report is 1 commits behind head on main. The diff coverage is 89.13%.

@@            Coverage Diff             @@
##             main     #879      +/-   ##
==========================================
+ Coverage   71.88%   72.45%   +0.57%     
==========================================
  Files          42       42              
  Lines        2746     2832      +86     
==========================================
+ Hits         1974     2052      +78     
- Misses        772      780       +8     

Files Changed	Coverage Δ
src/submodules/pathfinding/astar.jl	90.12% <50.00%> (-3.12%)	:arrow_down:
src/spaces/continuous.jl	92.30% <60.00%> (-0.96%)	:arrow_down:
src/spaces/walk.jl	86.86% <63.63%> (-2.03%)	:arrow_down:
src/spaces/grid_general.jl	98.92% <100.00%> (+0.08%)	:arrow_up:
src/spaces/grid_multi.jl	83.15% <100.00%> (+1.54%)	:arrow_up:
src/spaces/grid_single.jl	95.65% <100.00%> (+1.05%)	:arrow_up:
src/spaces/utilities.jl	100.00% <100.00%> (ø)
src/submodules/pathfinding/astar_continuous.jl	93.87% <100.00%> (+0.26%)	:arrow_up:
src/submodules/pathfinding/pathfinding_utils.jl	100.00% <100.00%> (ø)

... and 1 file with indirect coverage changes

:mega: We’re building smart automated test selection to slash your CI/CD build times. Learn more

[mastrof]

Can I get any help with this? I don't understand what's going on. How does this even work if bound_range is only defined for non-periodic spaces?

##########################################################################################
# nearby_stuff with special access r::Tuple
##########################################################################################
# TODO: We can re-write this to create its own `indices_within_radius_tuple`.
# This would also allow it to work for any metric, not just Chebyshev!

function nearby_ids(pos::ValidPos, model::ABM{<:GridSpace}, r::NTuple{D,Int}) where {D}
    # simply transform `r` to the Vector format expected by the below function
    newr = [(i, -r[i]:r[i]) for i in 1:D]
    nearby_ids(pos, model, newr)
end

function nearby_ids(
    pos::ValidPos,
    model::ABM{<:GridSpace},
    r::Vector{Tuple{Int64, UnitRange{Int64}}},
)
    @assert abmspace(model).metric == :chebyshev
    dims = first.(r)
    vidx = []
    for d in 1:ndims(abmspace(model).stored_ids)
        idx = findall(dim -> dim == d, dims)
        dim_range = isempty(idx) ? Colon() :
            bound_range(pos[d] .+ last(r[only(idx)]), d, abmspace(model))
        push!(vidx, dim_range)
    end
    s = view(abmspace(model).stored_ids, vidx...)
    Iterators.flatten(s)
end

function bound_range(unbound, d, space::GridSpace{D,false}) where {D}
    return range(max(unbound.start, 1), stop = min(unbound.stop, size(space)[d]))
end

EDIT: nvm just realised that indeed this is not working with periodic spaces; I guess I'll just ignore it?

[mastrof]

I think I got a decent implementation for the nearby searches; with 6 dimensions the performance for a mixed-boundary model is pretty close to that of a full periodic model.

julia> bt_pure = @belapsed nearby_positions(
        $(1,5,1,5,1,5),
        $(ABM(GridAgent{6}, GridSpace((5,5,5,5,5,5); metric=:euclidean, periodic=true)))
    )

julia> bt_mix = @belapsed nearby_positions(
        $(1,5,1,5,1,5),
        $(ABM(GridAgent{6}, GridSpace((5,5,5,5,5,5); metric=:euclidean, periodic=(true,false,true,false,true,false))))
    )

pure: 1.483567134268537e-8
mix: 1.7302908726178538e-8

I guess AStar is the only missing thing now

[Tortar]

To me seems good, but will take an in-depth look soon, add a changelog entry anyway :-)

[mastrof]

There's one of the tests in the AStar (L176 in astar_tests.jl) which I marked as wrong, but after some further testing I realized that I don't understand what the expected behavior is supposed to be. I wonder if there is a bug somewhere else? Consider this case, using a periodic=true space (so not affected by anything I implemented here):

julia> space = ContinuousSpace((5., 5.); periodic=true)
periodic continuous space with [5.0, 5.0] extent and spacing=0.25

julia> pf = AStar(space; walkmap=trues(10,10))
A* in 2 dimensions, periodic, diagonal, ϵ=0.0, metric=DirectDistance

julia> model = ABM(ContinuousAgent{2,Float64}, space; properties=(pf=pf,))
StandardABM with 0 agents of type ContinuousAgent
 space: periodic continuous space with [5.0, 5.0] extent and spacing=0.25
 scheduler: fastest
 properties: pf

julia> a = add_agent!((0.,0.), model, (0.,0.))
ContinuousAgent{2, Float64}(1, [0.0, 0.0], [0.0, 0.0])

julia> plan_best_route!(a, [SVector(0.0, 1.0)], model.pf)
2-element SVector{2, Float64} with indices SOneTo(2):
 0.0
 1.0

julia> move_along_route!(a, model, model.pf, 0.1)
ContinuousAgent{2, Float64}(1, [0.0316227766016838, 0.09486832980505139], [0.0, 0.0])

When running the last command, I was expecting the agent to move to [0.0, 0.1], instead it moves also along the x axis, even if its target destination has x=0. So if you keep moving along the route, the agent first strays a bit to positive x and then corrects by going back and eventually gets to the correct destination. But why is it moving like this? Is this really the expected behavior?

If we take a different target destination then indeed the agent is moving straight towards it, without any seemingly arbitrary deroute:

julia> a = add_agent!((0.,0.), model, (0.,0.))
ContinuousAgent{2, Float64}(2, [0.0, 0.0], [0.0, 0.0])

julia> plan_best_route!(a, [SVector(1.0, 1.0)], model.pf)
2-element SVector{2, Float64} with indices SOneTo(2):
 1.0
 1.0

julia> move_along_route!(a, model, model.pf, 0.1)
ContinuousAgent{2, Float64}(2, [0.07071067811865477, 0.07071067811865477], [0.0, 0.0])

[mastrof]

Regarding my last comment, atm I just made the pathfinder test consistent with what the library is currently doing. Whether that is correct or not is a different topic. I implemented all the review comments (thanks for that @Tortar), and I think we have everything now

[Tortar]

I think I got a decent implementation for the nearby searches; with 6 dimensions the performance for a mixed-boundary model is pretty close to that of a full periodic model.

julia> bt_pure = @belapsed nearby_positions(
        $(1,5,1,5,1,5),
        $(ABM(GridAgent{6}, GridSpace((5,5,5,5,5,5); metric=:euclidean, periodic=true)))
    )

julia> bt_mix = @belapsed nearby_positions(
        $(1,5,1,5,1,5),
        $(ABM(GridAgent{6}, GridSpace((5,5,5,5,5,5); metric=:euclidean, periodic=(true,false,true,false,true,false))))
    )

pure: 1.483567134268537e-8
mix: 1.7302908726178538e-8

note that you need to run collect otherwise you just compare the creation of the iterators

[Tortar]

Great work anyway :-)

(But this still need the changelog entry)

[mastrof]

You are right :s then there is some noticeable performance difference (1.5e-7 s and 5e-7 s on my machine with the above benchmark in 6 dimensions). Not as bad as I initially thought it would be, but not sure it can get much better.

Changelog was done btw: https://github.com/JuliaDynamics/Agents.jl/pull/879/commits/b2136fc2aafa8c51ed43de21513c5d057de76ad9

[Tortar]

it's okay since in the end for the majority of points the time should be the same, on the boundary this perf difference seems legit