Trouble trying to implement DifferentialEquations.jl integrator in an agent_step() for a model in Agents.jl?

Specific Domains Modelling & Simulations
,
1

I have created a continuous ABM model using Agents.jl, however, I wish to incorporate an Integrator from the DifferentialEquations.jl package to be used in the agent_step!() function. Though I am having trouble in the implementation.

[Package Versions]
Agents v6.2.5
DifferentialEquations v7.16.0

The following is a brief overview of the code, and finally I share the errors that I encounter.

using Agents, Random, DifferentialEquations
### Agent struct
@agent struct Pedestrian(ContinuousAgent{2,Float64})
    uᵢ::Vector{Float64} # desired_velocity : custom property
end

# Example simple_move!
# This works [only written to show, that the agents side of the code works]
function simple_move(agent::Pedestrian, model)
    p = agent.pos + model.dt.*[1,0]
    v = agent.vel + model.dt.*rand(2)
    return v, p
end

# This agent step requires which solver(or stepping method to use)
# either `simple_move` or `ode_step`
function agent_step!(agent::Pedestrian, model, num_solver)
    agent.vel, p = num_solver(agent,model)
    p = normalize_position(p,model)
    move_agent!(agent, p, model)
end

# ABM Model
function initialize(number_of_peds, x_len, y_len, num_solver)
    space = ContinuousSpace((x_len, y_len); periodic = true)
    properties = Dict(:λ => 2, :A => 5, :B => 0.3, :dt => 0.01)
    rng = Random.seed!(42)
    # num_solver = euler_step

    model = StandardABM(
        Pedestrian,
        space;
        container = Vector,
        agent_step! = (agent, model) -> agent_step!(agent::Pedestrian, model, num_solver),
        properties,
        rng,
        scheduler = Schedulers.Randomly()
    )

    for i in 1:number_of_peds
        add_agent!(model;
            pos = [rand()*x_len, rand()*y_len],  # Initial position
            # pos = [i,0],
            vel = [0,0], # Initial velocity
            # uᵢ = i-1 < (number_of_peds ÷ 2) ? [1,0] : [-1,0]
            uᵢ = [1,0]
        )
    end

    return model
end

model = initialize(4,11,5, simple_move)
step!(model) # works
allagents(model)

This much of the code works, however, I wish to replace simple_move to ode_step! as defined below, and I can’t pinpoint where problem really is

### DiffEq Side

function ped_model!(du,u,p,t)
    agent, model = p
    vel, pos = u[1:2], u[3:4]
    dist = (a1,a2) -> euclidean_distance(a1,a2,model)
    potU = (pos1,pos2) -> - ((pos1 - pos2)/dist(pos1,pos2))*model.A*exp(-dist(pos1,pos2)/model.B)
    du[1:2] = model.λ(agent.uᵢ - vel) - sum( [(i.id != agent.id) ? potU(pos, i.pos) : [0,0] for i in allagents(model)] )
    du[3:4] = vel
end

function ode_step(agent, model)
    u0 = [agent.pos; agent.vel]
    println(u0)
    tspan = (0.0, Inf)
    p = [agent, model]
    prob = ODEProblem{false}(ped_model!, u0, tspan, p)
    
    # Use the model's integrator
    integrator = DifferentialEquations.init(prob, Euler(), dt=model.dt)
    
    # Solve and update agent state
    step!(integrator, model.dt)
    pos = integrator.u[1:2]
    vel = integrator.u[3:4]
    # pos = normalize_position(pos,model)
    return vel, pos
end

model = initialize(4,11,5, ode_step)
step!(model) # Error

The error in question.

MethodError: no method matching ph_ped_model!(::SVector{4, Float64}, ::Vector{Any}, ::Float64)
The function `ped_model!` exists, but no method is defined for this combination of argument types.

MethodErrors are usually easy to solve, however this time I don’t really know how to resolve this. I have my doubts about how I have written ode_function or ped_model, since I am new to DifferentialEquations.jl package and julia itself. Any help or suggestion would be greatly appreciated!

2

I fixed it, somehow,
the followig changes are done inside ode_step(agent,model)

function ode_step!(agent,model)
    # must be given as an array of size 4
    u0 = [agent.vel[1], agent.vel[2], agent.pos[1], agent.pos[2]]
    ...
    # set the differential equation to true, i.e. IIP form so it takes du,u,p,t
    # (as opposed to OOP form, which only takes u,p,t as arguments)
    prob = ODEProblem{true}(ped_model!, u0, tspan, p) 
    ...
    # needs to be converted to static vector, 
    # as this datatype is required by 
    # normalize_position() inside agent_step!()
    return vel, SVector{2}(pos)