Hey everyone, I’m trying to implement the social force model for pedestrian dynamics, which models pedestrians using physics. I’m not sure what the best course of action is, specifically in terms of which packages I could make use of. I know that it would be possible to stay strictly within the DifferentialEquations.jl package by doing everything with Arrays, but I wanted to ask whether there were some convenient functions or features from other packages in the JuliaDynamics ecosystem ( DynamicalSystems.jl , Agents.jl ) that I could use to make the code more readable and modular.
I did some testing* earlier and it seemed that using structs in DifferentialEquations.jl as opposed to arrays had a time penalty of about 2x, which I’m not sure justifies the nicer interface, but if someone can look at the functions and see something obvious that I’m doing wrong I’d appreciate it.
Edit: removed the link to my repo and made it private by request of my supervisor. I will make the repo public after submitting in March 2022
using BenchmarkTools
f(x) = 2x^2
struct Agent{T}
x::T
y::T
v::T
w::T
end
import Base.similar, Base.copy;
copy(u::Agent{T}) where T = Agent{T}(u.x,u.y,u.v,u.w);
similar(u::Agent{T}) where T = copy(u);
function test!(du::T, u::T, N::Int) where {T<:Array{<:Real}}
for i=1:N
du[1, i] = f(u[1, i])
du[2, i] = f(u[2, i])
du[3, i] = f(u[3, i])
du[4, i] = f(u[4, i])
end
end
function test!(du::T, u::T, N::Int) where {T<:Array{<:Agent}}
for i=1:N
du[i] = Agent(f(u[i].x), f(u[i].y), f(u[i].v), f(u[i].w))
end
end
function wrappedAllocations(f, args...)
f(args...)
res = @timed f(args...)
return res
end
function main()
N = 1000000;
u0 = rand(4, N); # Array{Float64, 2}
U0 = [Agent(u0[:,i]...) for i=1:N]; # Array{Agent{Float64}, 1}
du = similar(u0); # Array{Float64, 2}
dU = similar.(U0); # Array{Agent{Float64}, 1}
print("Arrays: ")
@btime test!($du, $u0, $N)
print("Agents: ")
@btime test!($dU, $U0, $N)
end
main();
yields
Arrays: 3.408 ms (0 allocations: 0 bytes)
Agents: 3.068 ms (0 allocations: 0 bytes)
for me. I believe there a way better optimizations available for immutable structures.
Regarding the approach of @goerch : This is something I have considered a lot in Agents.jl: whether having agents as immutable types was worth it or not. https://github.com/JuliaDynamics/Agents.jl/issues/421 . The conclusion was that even though there were some measurable performance benefits, the huge downsides in convenience and intuitive rationalizing about code were simply too much to consider it useful in practice.
@jacobusmmsmit for your original question. First of all, I wasn’t even aware that DifferentialEquations.jl can work with vectors of arbitrary types. I think you would still need to define addition between vectors of agents, etc. Are you sure this “works”, or is it as easy as you think?
More importantly now. For your social force model, what you need to consider is the importance of accuracy. Is an Euler stepping approach accurate enough? Then you can go ahead with Agents.jl. If you need better accuracy for time evolution, you can consider DifferentialEquations.jl. No matter the level of accuracy, DiffEq will of course be always faster than Agents.jl. Things you should consider having a look at:
For using DynamicalSystems.jl: This package is only useful if you want to do an analysis of differential equations that goes beyond integrating these equations. I am assuming that the entire point of the Social Force Model is to set up some differential equations and then solve them. In this scenario DynamicalSystems.jl offers nothing more besides wrapping DiffEq (and hence, just use DiffEq directly).
Actually, looking the paper in more detail it seems like the continuous space of Agents.jl may be the best fit, or the more intuitive fit to be precise. You can also be inspired by the work done in this PR: https://github.com/JuliaDynamics/Agents.jl/pull/456 on Force Based Motion planning.
I dont think this really connects to to what I said about rationalizing about the code representation of an agent based model though. The in-place version with mutable types is more sensible and represents how one thinks about agent based models. When you move an agent you don’t “create a new one at a new location”, but you rather “change the location of an existing one”. If you drive your mom at the supermarket you don’t create a new mom at the supermarket while killing the one at home. I guess the functional programming community would generally agree with this line of thought.
Buuuut regardless, I think this is becoming too off-topic for this post. To not derail further, feel free to write more in the dynamics-bridged channel on slack or the issue link I posted in the original post!