From 384b4ad75b327fb6c4a26913c6239d72d4e41d4b Mon Sep 17 00:00:00 2001
From: Mahmoud Bentriou <mahmoud.bentriou@centralesupelec.fr>
Date: Sun, 15 Nov 2020 21:06:07 +0100
Subject: [PATCH] Fix of type instability in SIR model: - Real => Float64 for
better performance (Real is an abstract type) - Should not create a type
inside the model because the compiler fails at estimating the type of the
output of the f function (detected with @code_warntype) - Intentionaly ER is
not as well improved as SIR in order to check the differences of performance.
---
core/model.jl | 18 +++++++++---------
core/observations.jl | 6 +++---
models/ER.jl | 8 ++++----
models/SIR.jl | 13 +++++--------
4 files changed, 21 insertions(+), 24 deletions(-)
diff --git a/core/model.jl b/core/model.jl
index 5e87907..7c428ab 100644
--- a/core/model.jl
+++ b/core/model.jl
@@ -12,19 +12,19 @@ mutable struct CTMC <: ContinuousTimeModel
_map_obs_var_idx::Dict # maps str to observed state space
map_param_idx::Dict # maps str in parameter space
l_name_transitions::AbstractVector{String}
- p::AbstractVector{Real}
+ p::AbstractVector{Float64}
x0::AbstractVector{Int}
- t0::Real
+ t0::Float64
f::Function
g::AbstractVector{String} # of dimension dobs
- _g_idx::AbstractVector{Int} # of dimension dobs
+ _g_idx::Vector{Int} # of dimension dobs
is_absorbing::Function
- time_bound::Real
+ time_bound::Float64
end
function CTMC(d::Int, k::Int, map_var_idx::Dict, map_param_idx::Dict, l_name_transitions::AbstractVector{String},
- p::AbstractVector{Real}, x0::AbstractVector{Int}, t0::Real,
- f::Function, is_absorbing::Function; g::AbstractVector{String} = keys(map_var_idx), time_bound::Real = Inf)
+ p::AbstractVector{Float64}, x0::AbstractVector{Int}, t0::Float64,
+ f::Function, is_absorbing::Function; g::AbstractVector{String} = keys(map_var_idx), time_bound::Float64 = Inf)
dobs = length(g)
_map_obs_var_idx = Dict()
_g_idx = Vector{Int}(undef, dobs)
@@ -72,9 +72,9 @@ end
is_bounded(m::Model) = m.time_bound < Inf
function check_consistency(m::Model) end
-function simulate(m::Model, n::Int; bound::Real = Inf)::AbstractObservations end
-function set_param!(m::Model, p::AbstractVector{Real})::Nothing end
-function get_param(m::Model)::AbstractVector{Real} end
+function simulate(m::Model, n::Int; bound::Float64 = Inf)::AbstractObservations end
+function set_param!(m::Model, p::AbstractVector{Float64})::Nothing end
+function get_param(m::Model)::AbstractVector{Float64} end
function load_model(name_model::String)
include(pathof(@__MODULE__) * "/../../models/" * name_model * ".jl")
diff --git a/core/observations.jl b/core/observations.jl
index 7a35b81..5967ea5 100644
--- a/core/observations.jl
+++ b/core/observations.jl
@@ -4,14 +4,14 @@ ContinuousObservations = AbstractVector{AbstractTrajectory}
struct Trajectory <: AbstractTrajectory
m::ContinuousTimeModel
- values::AbstractMatrix{Real}
- times::AbstractVector{Real}
+ values::AbstractMatrix{Float64}
+ times::AbstractVector{Float64}
transitions::AbstractVector{Union{String,Nothing}}
end
function +(σ1::AbstractTrajectory,σ2::AbstractTrajectory) end
function -(σ1::AbstractTrajectory,σ2::AbstractTrajectory) end
-function δ(σ1::AbstractTrajectory,t::Real) end
+function δ(σ1::AbstractTrajectory,t::Float64) end
function get_obs_variables(σ::AbstractTrajectory) end
get_values(σ::AbstractTrajectory, var::String) =
diff --git a/models/ER.jl b/models/ER.jl
index 876b129..03bebee 100644
--- a/models/ER.jl
+++ b/models/ER.jl
@@ -3,7 +3,7 @@
import StaticArrays: SVector, SMatrix, @SMatrix
State = SVector{4, Int}
-Parameters = SVector{3, Real}
+Parameters = SVector{3, Float64}
d=4
k=3
@@ -14,7 +14,7 @@ p = Parameters(1.0, 1.0, 1.0)
x0 = State(100, 100, 0, 0)
t0 = 0.0
-function f(xn::State, tn::Real, p::Parameters)
+function f(xn::State, tn::Float64, p::Parameters)
a1 = p[1] * xn[1] * xn[2]
a2 = p[2] * xn[3]
a3 = p[3] * xn[3]
@@ -45,10 +45,10 @@ function f(xn::State, tn::Real, p::Parameters)
return xnplus1, tnplus1, transition
end
-is_absorbing_sir(p::Parameters,xn::State) =
+is_absorbing_er(p::Parameters,xn::State) =
(p[1]*xn[1]*xn[2] + (p[2]+p[3])*xn[3]) == 0.0
g = SVector("P")
-ER = CTMC(d,k,dict_var,dict_p,l_tr,p,x0,t0,f,is_absorbing_sir; g=g)
+ER = CTMC(d,k,dict_var,dict_p,l_tr,p,x0,t0,f,is_absorbing_er; g=g)
export ER
diff --git a/models/SIR.jl b/models/SIR.jl
index 122889b..451c9c2 100644
--- a/models/SIR.jl
+++ b/models/SIR.jl
@@ -1,18 +1,15 @@
import StaticArrays: SVector, SMatrix, @SMatrix
-State = SVector{3, Int}
-Parameters = SVector{2, Real}
-
d=3
k=2
dict_var = Dict("S" => 1, "I" => 2, "R" => 3)
dict_p = Dict("ki" => 1, "kr" => 2)
l_tr = ["R1","R2"]
-p = Parameters(0.0012, 0.05)
-x0 = State(95, 5, 0)
+p = SVector(0.0012, 0.05)
+x0 = SVector(95, 5, 0)
t0 = 0.0
-function f(xn::State, tn::Real, p::Parameters)
+function f(xn::SVector{3, Int}, tn::Float64, p::SVector{2, Float64})
a1 = p[1] * xn[1] * xn[2]
a2 = p[2] * xn[2]
l_a = SVector(a1, a2)
@@ -37,13 +34,13 @@ function f(xn::State, tn::Real, p::Parameters)
end
nu = l_nu[:,reaction]
- xnplus1 = State(xn[1]+nu[1], xn[2]+nu[2], xn[3]+nu[3])
+ xnplus1 = SVector(xn[1]+nu[1], xn[2]+nu[2], xn[3]+nu[3])
tnplus1 = tn + tau
transition = "R$(reaction)"
return xnplus1, tnplus1, transition
end
-is_absorbing_sir(p::Parameters,xn::State) = (p[1]*xn[1]*xn[2] + p[2]*xn[2]) == 0.0
+is_absorbing_sir(p::SVector{2, Float64}, xn::SVector{3, Int}) = (p[1]*xn[1]*xn[2] + p[2]*xn[2]) == 0.0
g = SVector("I")
SIR = CTMC(d,k,dict_var,dict_p,l_tr,p,x0,t0,f,is_absorbing_sir; g=g)
--
GitLab