Commit b41636cb authored by Victor's avatar Victor
Browse files

replaced tuple for position by vector. now faster

parent de690e73
Pipeline #84323 failed with stage
in 17 minutes and 21 seconds
...@@ -10,7 +10,7 @@ $(TYPEDEF) ...@@ -10,7 +10,7 @@ $(TYPEDEF)
""" """
mutable struct Agent{A<:Ancestors,R<:Rates,T<:Tuple,U,V} <: AbstractAgent{A,R} mutable struct Agent{A<:Ancestors,R<:Rates,T<:Tuple,U,V} <: AbstractAgent{A,R}
# history of traits for geotraits # history of traits for geotraits
x_history::Array{T,1} x_history::Array{Any,1}
# birth time of ancestors # birth time of ancestors
t_history::Array{U,1} t_history::Array{U,1}
# death rate # death rate
...@@ -26,12 +26,12 @@ eltype(a::Agent{A,R,T,U,V}) where {A,R,T,U,V} = T ...@@ -26,12 +26,12 @@ eltype(a::Agent{A,R,T,U,V}) where {A,R,T,U,V} = T
function initpos(s::S) where {S<:AbstractSpacesTuple} function initpos(s::S) where {S<:AbstractSpacesTuple}
Eltype = eltype.(s) Eltype = eltype.(s)
Dims = ndims.(s) Dims = ndims.(s)
pos = tuple() pos = []
for i in 1:length(Eltype) for i in 1:length(Eltype)
if Dims[i] > 1 if Dims[i] > 1
pos = (pos...,Eltype[i](ones(Dims[i]))) push!(pos,ones(eltype(Eltype[i]),Dims[i]))
else else
pos = (pos...,one(Eltype[i])) pos = push!(pos,one(Eltype[i]))
end end
end end
Tuple{Eltype...},pos Tuple{Eltype...},pos
...@@ -57,7 +57,7 @@ end ...@@ -57,7 +57,7 @@ end
$(SIGNATURES) $(SIGNATURES)
Initialises agent with `pos` provided Initialises agent with `pos` provided
""" """
function Agent(s::S, pos::P;ancestors=false,rates=false) where {P<:Tuple,S <: AbstractSpacesTuple} function Agent(s::S, pos::P;ancestors=false,rates=false) where {P<:Vector,S <: AbstractSpacesTuple}
T = eltype.(s) T = eltype.(s)
for (i,p) in enumerate(pos) for (i,p) in enumerate(pos)
if typeof(p) !== T[i] if typeof(p) !== T[i]
...@@ -175,14 +175,13 @@ import Base.(+) ...@@ -175,14 +175,13 @@ import Base.(+)
function _get_xinc(a::AbstractAgent,s::AbstractSpacesTuple,p::Dict,t::Number) function _get_xinc(a::AbstractAgent,s::AbstractSpacesTuple,p::Dict,t::Number)
@unpack D,mu = p @unpack D,mu = p
_x = get_x(a) _x = copy(get_x(a))
inc = zero(_x)
for (i,ss) in enumerate(s) for (i,ss) in enumerate(s)
if rand() < mu[i] if rand() < mu[i]
inc[i] = get_inc(_x[i],D[i],ss,t) _x[i] += get_inc(_x[i],D[i],ss,t)
end end
end end
tuple((_x .+ inc)...) _x
end end
## Modifiers ## Modifiers
......
...@@ -14,7 +14,7 @@ AbstractSpacesTuple = Tuple{Vararg{AbstractSpace}} ...@@ -14,7 +14,7 @@ AbstractSpacesTuple = Tuple{Vararg{AbstractSpace}}
import Base:ndims,isfinite,eltype import Base:ndims,isfinite,eltype
Base.ndims(x::AbstractSpace{Dim,T,I}) where {Dim,T,I} = Dim Base.ndims(x::AbstractSpace{Dim,T,I}) where {Dim,T,I} = Dim
Base.isfinite(x::AbstractSpace{Dim,T,IsFinite{t}}) where {Dim,T,t} = t #not sure we need this Base.isfinite(x::AbstractSpace{Dim,T,IsFinite{t}}) where {Dim,T,t} = t #not sure we need this
Base.eltype(::AbstractSpace{Dim,T,I}) where {Dim,T,I} = Dim > 1 ? Tuple{Vararg{T,Dim}} : T Base.eltype(::AbstractSpace{Dim,T,I}) where {Dim,T,I} = Dim > 1 ? Vector{T} : T
Base.ndims(ss::AbstractSpacesTuple) = length(ss) Base.ndims(ss::AbstractSpacesTuple) = length(ss)
Base.eltype(ss::AbstractSpacesTuple) where {Dim,T,I} = Tuple{eltype.(ss)...} Base.eltype(ss::AbstractSpacesTuple) where {Dim,T,I} = Tuple{eltype.(ss)...}
...@@ -71,7 +71,7 @@ get_inc(x,D,s::AbstractSpace{Dim,T,I}) where {Dim,T,I<:IsFinite{false}} = get_in ...@@ -71,7 +71,7 @@ get_inc(x,D,s::AbstractSpace{Dim,T,I}) where {Dim,T,I<:IsFinite{false}} = get_in
function get_inc(D,s::AbstractSpace{Dim,T,I}) where {Dim,T<:AbstractFloat,I<:IsFinite{false}} function get_inc(D,s::AbstractSpace{Dim,T,I}) where {Dim,T<:AbstractFloat,I<:IsFinite{false}}
if Dim > 1 if Dim > 1
return Tuple(D .* randn(T,Dim)) return D .* randn(T,Dim)
else else
return D * randn(T) return D * randn(T)
end end
...@@ -79,7 +79,7 @@ end ...@@ -79,7 +79,7 @@ end
function get_inc(D,s::AbstractSpace{Dim,T,I}) where {Dim,T<:Integer,I<:IsFinite{false}} function get_inc(D,s::AbstractSpace{Dim,T,I}) where {Dim,T<:Integer,I<:IsFinite{false}}
if Dim > 1 if Dim > 1
return Tuple(round.(T,D .*randn(Float32,Dim))) return round.(T,D .*randn(Float32,Dim))
else else
return round(D * randn(Float32)) return round(D * randn(Float32))
end end
......
...@@ -12,7 +12,7 @@ sigma_a = .7; ...@@ -12,7 +12,7 @@ sigma_a = .7;
K0 = 1000; K0 = 1000;
b(X,t) = 1. b(X,t) = 1.
d(X,Y,t) = gaussian(X[1],Y[1],sigma_a)/gaussian(X[1],0.,sigma_K)/K0 d(X,Y,t) = gaussian(X[1],Y[1],sigma_a)/gaussian(X[1],0.,sigma_K)/K0
D = (1e-2,) D = [1e-2]
mu = [.1] mu = [.1]
NMax = 2000 NMax = 2000
tend = 1500 tend = 1500
...@@ -20,11 +20,11 @@ tend = 1500 ...@@ -20,11 +20,11 @@ tend = 1500
dm = d([0],[0],0.);bm = 1. dm = d([0],[0],0.);bm = 1.
p = Dict{String,Any}();@pack! p = d,b,D,mu,NMax,dm,bm p = Dict{String,Any}();@pack! p = d,b,D,mu,NMax,dm,bm
myagents = [Agent(myspace,(1e-2 * randn(),),ancestors=false,rates=false) for i in 1:K0] myagents = [Agent(myspace,[1e-2 * randn()],ancestors=false,rates=false) for i in 1:K0]
w0 = World(myagents,myspace,p,0.) w0 = World(myagents,myspace,p,0.)
w1 = copy(w0) w1 = copy(w0)
@info "Running simulation with CFM algorithm" @info "Running simulation with CFM algorithm"
@time sim = run!(w1,CFM(),tend,dt_saving=10.) @time sim = run!(w1,CFM(),tend,b,d,dt_saving=10.)
# ABMEv.clean!(sim) # ABMEv.clean!(sim)
using Plots using Plots
......
...@@ -12,13 +12,13 @@ sigma_a = .7; ...@@ -12,13 +12,13 @@ sigma_a = .7;
K0 = 1000; K0 = 1000;
b(X,t) = gaussian(X[1],0.,sigma_K) b(X,t) = gaussian(X[1],0.,sigma_K)
d(X,Y,t) = gaussian(X[1],Y[1],sigma_a)/K0 d(X,Y,t) = gaussian(X[1],Y[1],sigma_a)/K0
D = (1e-2,) D = [1e-2]
mu = [.1] mu = [.1]
NMax = 10000 NMax = 10000
tend = 1.5 tend = 1.5
p = Dict{String,Any}();@pack! p = D,mu,NMax p = Dict{String,Any}();@pack! p = D,mu,NMax
myagents = [Agent(myspace,(0,),ancestors=true,rates=true) for i in 1:K0] myagents = [Agent(myspace,[0.],ancestors=true,rates=true) for i in 1:K0]
w0 = World(myagents,myspace,p,0.) w0 = World(myagents,myspace,p,0.)
w1 = copy(w0) w1 = copy(w0)
@info "Running simulation with Gillepsie algorithm" @info "Running simulation with Gillepsie algorithm"
...@@ -41,7 +41,7 @@ w1 = copy(w0) ...@@ -41,7 +41,7 @@ w1 = copy(w0)
@testset "Testing update rates matrix" begin @testset "Testing update rates matrix" begin
@testset "birth event" begin @testset "birth event" begin
myagents = [Agent(myspace,(0,),ancestors=true,rates=true) for i in 1:K0] myagents = [Agent(myspace,[0],ancestors=true,rates=true) for i in 1:K0]
w1 = World(myagents,myspace,p,0.);update_rates!(w1,Gillepsie(),b,d,) w1 = World(myagents,myspace,p,0.);update_rates!(w1,Gillepsie(),b,d,)
mum_idx = 1 mum_idx = 1
updateBirthEvent!(w0,Gillepsie(),1,b,d) updateBirthEvent!(w0,Gillepsie(),1,b,d)
...@@ -54,7 +54,7 @@ w1 = copy(w0) ...@@ -54,7 +54,7 @@ w1 = copy(w0)
end end
@testset "death event" begin @testset "death event" begin
myagents = [Agent(myspace,(0,),ancestors=true,rates=true) for i in 1:K0] myagents = [Agent(myspace,[0],ancestors=true,rates=true) for i in 1:K0]
w1 = World(myagents,myspace,p,0.);update_rates!(w1,Gillepsie(),b,d) w1 = World(myagents,myspace,p,0.);update_rates!(w1,Gillepsie(),b,d)
mum_idx = 1 mum_idx = 1
updateDeathEvent!(w0,Gillepsie(),1,d) updateDeathEvent!(w0,Gillepsie(),1,d)
......
...@@ -9,10 +9,10 @@ g = SimpleGraph(1000,4000) ...@@ -9,10 +9,10 @@ g = SimpleGraph(1000,4000)
myspace4 = (RealSpace{1,Float64}(),GraphSpace(g),) myspace4 = (RealSpace{1,Float64}(),GraphSpace(g),)
K0 = 1000; σ = 1e-1 K0 = 1000; σ = 1e-1
a1 = [Agent(myspace1, (σ,) .* randn() .- .5,ancestors=true) for i in 1:K0] a1 = [Agent(myspace1, σ .* randn(1) .- .5,ancestors=true) for i in 1:K0]
a2 = [Agent(myspace2,tuple((σ, σ) .* randn(2) .- .5...),ancestors=true) for i in 1:K0] a2 = [Agent(myspace2,σ .* randn(2) .- .5,ancestors=true) for i in 1:K0]
a3 = [Agent(myspace3, (rand(Int16.(1:10)), tuple((1e-1.* randn(2) .+ 5.5)... )),ancestors=true) for i in 1:K0] a3 = [Agent(myspace3, [rand(Int16.(1:10)), 1e-1.* randn(2) .+ 5.5 ],ancestors=true) for i in 1:K0]
a4 = [Agent(myspace4, (1.,rand(Int64(1):Int64(1000)),),ancestors=true) for i in 1:K0] a4 = [Agent(myspace4, [1.,rand(Int64(1):Int64(1000))],ancestors=true) for i in 1:K0]
D = (1.,); D = (1.,);
mu = [1.,1.] mu = [1.,1.]
...@@ -80,7 +80,7 @@ w4 = World(a4,myspace4,p4) ...@@ -80,7 +80,7 @@ w4 = World(a4,myspace4,p4)
end end
@testset "Geotrait computation" begin @testset "Geotrait computation" begin
a = Agent(myspace3, (Int16.(1), randn() ),ancestors=true); increment_x!(a,myspace3,p3,2.0); a = Agent(myspace3, [Int16.(1), randn(2) ],ancestors=true); increment_x!(a,myspace3,p3,2.0);
@test get_geo(a,2.0) 2.0 @test get_geo(a,2.0) 2.0
end end
...@@ -89,8 +89,8 @@ end ...@@ -89,8 +89,8 @@ end
## testing for real space of dimension d>1 ## testing for real space of dimension d>1
multispace = (DiscreteSegment{Int8}(1,9),RealSpace{3,Float64}(),) multispace = (DiscreteSegment{Int8}(1,9),RealSpace{3,Float64}(),)
K0 = 10000; K0 = 10000;
multia = [Agent(multispace, (rand(Int8(1):Int8(10)),tuple( randn(3) ...),),ancestors=true) for i in 1:K0] multia = [Agent(multispace, [rand(Int8(1):Int8(10)),randn(3) ],ancestors=true) for i in 1:K0]
D = (Int8(1),tuple(fill(1.,3)...),) D = [Int8(1),fill(1.,3)]
mu = [1] mu = [1]
NMax = 1000 NMax = 1000
multip = Dict{String,Any}();@pack! multip = D,mu,NMax multip = Dict{String,Any}();@pack! multip = D,mu,NMax
......
...@@ -42,7 +42,7 @@ tend = 1.5 ...@@ -42,7 +42,7 @@ tend = 1.5
p = Dict{String,Any}();@pack! p = D,mu,NMax p = Dict{String,Any}();@pack! p = D,mu,NMax
## testing cb ## testing cb
myagents = [Agent(myspace,(0,),ancestors=true,rates=true) for i in 1:K0] myagents = [Agent(myspace,[0.],ancestors=true,rates=true) for i in 1:K0]
w0 = World(myagents,myspace,p,0.) w0 = World(myagents,myspace,p,0.)
w1 = copy(w0) w1 = copy(w0)
cb = (names = ["gamma_div"], agg = Function[w -> var(Float64.(get_x(w,1)))]) cb = (names = ["gamma_div"], agg = Function[w -> var(Float64.(get_x(w,1)))])
...@@ -55,7 +55,7 @@ eltype(cb.agg) ...@@ -55,7 +55,7 @@ eltype(cb.agg)
@test typeof(sim[2]) <: Vector @test typeof(sim[2]) <: Vector
##testing t_saving_cb ##testing t_saving_cb
myagents = [Agent(myspace,(0,),ancestors=true,rates=true) for i in 1:K0] myagents = [Agent(myspace,[0.],ancestors=true,rates=true) for i in 1:K0]
w0 = World(myagents,myspace,p,0.) w0 = World(myagents,myspace,p,0.)
w1 = copy(w0) w1 = copy(w0)
cb = (names = ["gamma_div"], agg = Function[w -> var(Float64.(get_x(w,1)))]) cb = (names = ["gamma_div"], agg = Function[w -> var(Float64.(get_x(w,1)))])
......
...@@ -20,7 +20,7 @@ myspace2 = (mysegment,mycontinuoussegment,real2d) ...@@ -20,7 +20,7 @@ myspace2 = (mysegment,mycontinuoussegment,real2d)
@test ndims(real2d) 2 @test ndims(real2d) 2
@test isfinite(mycontinuoussegment) true @test isfinite(mycontinuoussegment) true
@test typeof(myspace) <: AbstractSpacesTuple @test typeof(myspace) <: AbstractSpacesTuple
@test eltype(myspace2) == Tuple{Int64,Float64,Tuple{Float64,Float64}} @test eltype(myspace2) == Tuple{Int64,Float64,Vector{Float64}}
# increment on infinite spaces # increment on infinite spaces
@test ABMEv.get_inc(0.,myline) (0.) @test ABMEv.get_inc(0.,myline) (0.)
...@@ -30,10 +30,10 @@ myspace2 = (mysegment,mycontinuoussegment,real2d) ...@@ -30,10 +30,10 @@ myspace2 = (mysegment,mycontinuoussegment,real2d)
# @test !(get_inc(1,1.,mydiscreteline) ≈ 0.) # @test !(get_inc(1,1.,mydiscreteline) ≈ 0.)
@test typeof(ABMEv.get_inc([1.,0.],real2d)) == Tuple{Float64,Float64} @test typeof(ABMEv.get_inc([1.,0.],real2d)) == Vector{Float64}
@test typeof(get_inc([1.,0.],[1.,0.],real2d)) == Tuple{Float64,Float64} @test typeof(get_inc([1.,0.],[1.,0.],real2d)) == Vector{Float64}
@test typeof(ABMEv.get_inc([1.],real2d)) == Tuple{Float64,Float64} @test typeof(ABMEv.get_inc([1.],real2d)) == Vector{Float64}
@test typeof(ABMEv.get_inc(1.,real2d)) == Tuple{Float64,Float64} @test typeof(ABMEv.get_inc(1.,real2d)) == Vector{Float64}
# ABMEv._get_inc([1.],real2d) # ABMEv._get_inc([1.],real2d)
# ABMEv.initpos(myspace2) # ABMEv.initpos(myspace2)
...@@ -58,8 +58,8 @@ end ...@@ -58,8 +58,8 @@ end
##### AGENTS ####### ##### AGENTS #######
a1 = Agent(myspace) a1 = Agent(myspace)
a2 = Agent(myspace,ancestors = true) a2 = Agent(myspace,ancestors = true)
a3 = Agent(myspace,(1,1,1.)) a3 = Agent(myspace,[1,1,1.])
a4 = Agent(myspace2,(1,1,(1.,1)),rates=true) a4 = Agent(myspace2,[1,1.,[1.,1]],rates=true)
a5 = Agent(myspace2,ancestors=true) a5 = Agent(myspace2,ancestors=true)
@testset "Agent properties" begin @testset "Agent properties" begin
...@@ -77,4 +77,5 @@ a5 = Agent(myspace2,ancestors=true) ...@@ -77,4 +77,5 @@ a5 = Agent(myspace2,ancestors=true)
@test nancestors(increment_x!(a2,myspace,p_myspace,2.)) > 1 @test nancestors(increment_x!(a2,myspace,p_myspace,2.)) > 1
@test !isnothing(increment_x!(a4,myspace2,p_myspace2,2.)) @test !isnothing(increment_x!(a4,myspace2,p_myspace2,2.))
@test !isnothing(increment_x!(a5,myspace2,p_myspace2,2.)) @test !isnothing(increment_x!(a5,myspace2,p_myspace2,2.))
@test typeof(ABMEv._get_xinc(a2,myspace,p,0.)) == typeof(get_x(a2))
end end
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment