From cefe75b76b15afb2f288e1c14ef49f0efbb8ea8d Mon Sep 17 00:00:00 2001
From: Mahmoud Bentriou <mahmoud.bentriou@centralesupelec.fr>
Date: Thu, 19 Nov 2020 02:28:26 +0100
Subject: [PATCH] Change of different names as the Continuous Time model
object.
Major fix of type unstabilities which improved performances.
For now on, we observe for pygmalion benchmarks:
+2x speed simulation compared to simulation of models.
100x speed in terms of reading values (no allocations).
---
bench/pygmalion/multiple_long_sim_er.jl | 7 +++
bench/pygmalion/read_trajectory_er.jl | 32 +++++++++-----
core/MarkovProcesses.jl | 4 +-
core/model.jl | 31 +++++++------
core/trajectory.jl | 43 ++++++++-----------
models/ER.jl | 2 +-
models/SIR.jl | 4 +-
models/_bench_perf_test/ER_col.jl | 2 +-
models/_bench_perf_test/ER_col_buffer.jl | 2 +-
models/_bench_perf_test/ER_row_buffer.jl | 2 +-
models/_bench_perf_test/SIR_col.jl | 2 +-
models/_bench_perf_test/SIR_col_buffer.jl | 2 +-
models/_bench_perf_test/SIR_row_buffer.jl | 2 +-
tests/dist_lp/dist_case_2.jl | 4 +-
tests/dist_lp/dist_case_3.jl | 4 +-
tests/dist_lp/dist_l1_case_1.jl | 4 +-
tests/dist_lp/dist_sim_sir.jl | 4 +-
tests/simulation/sim_er.jl | 2 +-
tests/simulation/sim_er_row_buffer_bounded.jl | 2 +-
tests/simulation/sim_sir.jl | 2 +-
tests/simulation/sim_sir_bounded.jl | 2 +-
.../simulation/sim_sir_col_buffer_bounded.jl | 2 +-
.../simulation/sim_sir_row_buffer_bounded.jl | 2 +-
tests/unit/simulate_sir_bounded.jl | 2 +-
24 files changed, 87 insertions(+), 78 deletions(-)
diff --git a/bench/pygmalion/multiple_long_sim_er.jl b/bench/pygmalion/multiple_long_sim_er.jl
index 3e27d6a..6729aba 100644
--- a/bench/pygmalion/multiple_long_sim_er.jl
+++ b/bench/pygmalion/multiple_long_sim_er.jl
@@ -31,6 +31,13 @@ MarkovProcesses.load_model("ER")
ER.time_bound = 20.0
@timev MarkovProcesses.simulate(ER)
+println("Default buffer size=10")
+b1 = @benchmark for i=1:$(nb_sim) MarkovProcesses.simulate($ER) end
+b2 = @benchmark for i=1:$(nb_sim) MarkovProcesses.simulate(ER) end
+@show minimum(b1), mean(b1), maximum(b1)
+@show minimum(b2), mean(b2), maximum(b2)
+println("Buffer size 100")
+ER.buffer_size = 100
b1 = @benchmark for i=1:$(nb_sim) MarkovProcesses.simulate($ER) end
b2 = @benchmark for i=1:$(nb_sim) MarkovProcesses.simulate(ER) end
@show minimum(b1), mean(b1), maximum(b1)
diff --git a/bench/pygmalion/read_trajectory_er.jl b/bench/pygmalion/read_trajectory_er.jl
index 36ec8cc..f40df1d 100644
--- a/bench/pygmalion/read_trajectory_er.jl
+++ b/bench/pygmalion/read_trajectory_er.jl
@@ -16,7 +16,7 @@ u = Control(10.0)
tml = 1:400
g_all = create_observation_function([ObserverModel(str_oml, tml)])
so = pygmalion.simulate(f, g_all, x0, u, p_true; on = nothing, full_timeline = true)
-function read_trajectory_v1(so::SystemObservation)
+function read_trajectory_pyg_v1(so::SystemObservation)
vals = to_vec(so, "P")
n_states = length(vals)
res = 0.0
@@ -25,7 +25,7 @@ function read_trajectory_v1(so::SystemObservation)
end
return res
end
-function read_trajectory_v2(so::SystemObservation)
+function read_trajectory_pyg_v2(so::SystemObservation)
idx_P = om_findfirst("P", so.oml)
n_states = length(so.otll[idx_P])
res = 0.0
@@ -35,11 +35,11 @@ function read_trajectory_v2(so::SystemObservation)
return res
end
# Bench
-@timev read_trajectory_v1(so)
-b1_pyg = @benchmark read_trajectory_v1($so)
+@timev read_trajectory_pyg_v1(so)
+b1_pyg = @benchmark read_trajectory_pyg_v1($so)
@show minimum(b1_pyg), mean(b1_pyg), maximum(b1_pyg)
-@timev read_trajectory_v2(so)
-b2_pyg = @benchmark read_trajectory_v2($so)
+@timev read_trajectory_pyg_v2(so)
+b2_pyg = @benchmark read_trajectory_pyg_v2($so)
@show minimum(b2_pyg), mean(b2_pyg), maximum(b2_pyg)
println("MarkovProcesses:")
@@ -47,16 +47,28 @@ println("MarkovProcesses:")
MarkovProcesses.load_model("ER")
ER.time_bound = 10.0
σ = MarkovProcesses.simulate(ER)
-function read_trajectory(σ::AbstractTrajectory)
+function read_trajectory_v1(σ::AbstractTrajectory)
n_states = length_states(σ)
res = 0
for i = 1:n_states
- res += (σ["P"])[i]
+ res += σ["P",i]
+ end
+ return res
+end
+function read_trajectory_v2(σ::AbstractTrajectory)
+ n_states = length_states(σ)
+ vals = σ["P"]
+ res = 0
+ for i = 1:n_states
+ res += vals[i]
end
return res
end
# Bench
-@timev read_trajectory(σ)
-b1 = @benchmark read_trajectory($σ)
+@timev read_trajectory_v1(σ)
+b1 = @benchmark read_trajectory_v1($σ)
@show minimum(b1), mean(b1), maximum(b1)
+@timev read_trajectory_v2(σ)
+b2 = @benchmark read_trajectory_v2($σ)
+@show minimum(b2), mean(b2), maximum(b2)
diff --git a/core/MarkovProcesses.jl b/core/MarkovProcesses.jl
index 029af35..e9d2614 100644
--- a/core/MarkovProcesses.jl
+++ b/core/MarkovProcesses.jl
@@ -1,6 +1,6 @@
module MarkovProcesses
-import Base: +, -, getfield, getindex
+import Base: +, -, getfield, getindex, ImmutableDict
export Model, ContinuousTimeModel, DiscreteTimeModel
export simulate, set_param!, get_param, set_observed_var!
@@ -10,7 +10,7 @@ include("model.jl")
export Observations, AbstractTrajectory, Trajectory
export +, -, δ, dist_lp
-export get_obs_var, length_states, length_obs_var, is_bounded
+export get_obs_var, length_states, length_obs_var, is_bounded, times, transitions
include("trajectory.jl")
end
diff --git a/core/model.jl b/core/model.jl
index 541a4bb..e2b23b0 100644
--- a/core/model.jl
+++ b/core/model.jl
@@ -2,15 +2,14 @@
import StaticArrays: SVector
abstract type Model end
-abstract type ContinuousTimeModel <: Model end
abstract type DiscreteTimeModel <: Model end
-mutable struct CTMC <: ContinuousTimeModel
+mutable struct ContinuousTimeModel <: Model
d::Int # state space dim
k::Int # parameter space dim
- map_var_idx::Dict # maps str to full state space
- _map_obs_var_idx::Dict # maps str to observed state space
- map_param_idx::Dict # maps str in parameter space
+ map_var_idx::Dict{String,Int} # maps str to full state space
+ _map_obs_var_idx::Dict{String,Int} # maps str to observed state space
+ map_param_idx::Dict{String,Int} # maps str in parameter space
l_name_transitions::Vector{String}
p::Vector{Float64}
x0::Vector{Int}
@@ -23,7 +22,7 @@ mutable struct CTMC <: ContinuousTimeModel
buffer_size::Int
end
-function CTMC(d::Int, k::Int, map_var_idx::Dict, map_param_idx::Dict, l_name_transitions::Vector{String},
+function ContinuousTimeModel(d::Int, k::Int, map_var_idx::Dict, map_param_idx::Dict, l_name_transitions::Vector{String},
p::Vector{Float64}, x0::Vector{Int}, t0::Float64,
f!::Function, is_absorbing::Function;
g::Vector{String} = keys(map_var_idx), time_bound::Float64 = Inf, buffer_size::Int = 10)
@@ -42,7 +41,7 @@ function CTMC(d::Int, k::Int, map_var_idx::Dict, map_param_idx::Dict, l_name_tra
@warn "You have possibly redefined a function Model.is_absorbing used in a previously instantiated model."
end
- return CTMC(d, k, map_var_idx, _map_obs_var_idx, map_param_idx, l_name_transitions, p, x0, t0, f!, g, _g_idx, is_absorbing, time_bound, buffer_size)
+ return ContinuousTimeModel(d, k, map_var_idx, _map_obs_var_idx, map_param_idx, l_name_transitions, p, x0, t0, f!, g, _g_idx, is_absorbing, time_bound, buffer_size)
end
function simulate(m::ContinuousTimeModel)
@@ -53,27 +52,27 @@ function simulate(m::ContinuousTimeModel)
transitions = Union{String,Nothing}[nothing]
# values at time n
n = 0
- xn = @view m.x0[:] # View for type stability
+ xn = view(reshape(m.x0, 1, m.d), 1, :) # View for type stability
tn = m.t0
# at time n+1
mat_x = zeros(Int, m.buffer_size, m.d)
l_t = zeros(Float64, m.buffer_size)
l_tr = Vector{String}(undef, m.buffer_size)
- is_absorbing = m.is_absorbing(m.p,xn)::Bool
+ is_absorbing::Bool = m.is_absorbing(m.p,xn)
while !is_absorbing && (tn <= m.time_bound)
i = 0
while i < m.buffer_size && !is_absorbing && (tn <= m.time_bound)
i += 1
m.f!(mat_x, l_t, l_tr, i, xn, tn, m.p)
- xn = @view mat_x[i,:]
+ xn = view(mat_x, i, :)
tn = l_t[i]
- is_absorbing = m.is_absorbing(m.p,xn)::Bool
+ is_absorbing = m.is_absorbing(m.p,xn)
end
- full_values = vcat(full_values, @view mat_x[1:i,:])
- append!(times, @view l_t[1:i])
- append!(transitions, @view l_tr[1:i])
+ full_values = vcat(full_values, view(mat_x, 1:i, :))
+ append!(times, view(l_t, 1:i))
+ append!(transitions, view(l_tr, 1:i))
n += i
- is_absorbing = m.is_absorbing(m.p,xn)::Bool
+ is_absorbing = m.is_absorbing(m.p,xn)
end
if is_bounded(m)
if times[end] > m.time_bound
@@ -86,7 +85,7 @@ function simulate(m::ContinuousTimeModel)
push!(transitions, nothing)
end
end
- values = @view full_values[:,m._g_idx]
+ values = view(full_values, :, m._g_idx)
return Trajectory(m, values, times, transitions)
end
diff --git a/core/trajectory.jl b/core/trajectory.jl
index 5ef4940..6f8ce4a 100644
--- a/core/trajectory.jl
+++ b/core/trajectory.jl
@@ -74,11 +74,11 @@ function dist_lp(σ1::AbstractTrajectory, σ2::AbstractTrajectory, var::String;
if !is_bounded(σ1) || !is_bounded(σ2)
@warn "Lp distance computed on unbounded trajectories. Result should be wrong"
end
- return dist_lp(σ1[var], σ1["times"], σ2[var], σ2["times"]; verbose = false, p = p)
+ return dist_lp(σ1[var], times(σ1), σ2[var], times(σ2); verbose = false, p = p)
end
# Distance function. Vectorized version
-function dist_lp(x_obs::SubArray{Int,1}, t_x::SubArray{Float64,1}, y_obs::SubArray{Int,1}, t_y::SubArray{Float64,1};
+function dist_lp(x_obs::SubArray{Int,1}, t_x::Vector{Float64}, y_obs::SubArray{Int,1}, t_y::Vector{Float64};
verbose::Bool = false, p::Int = 1)
current_y_obs = y_obs[1]
current_t_y = t_y[2]
@@ -179,32 +179,23 @@ is_bounded(σ::AbstractTrajectory) = σ.transitions[end] == nothing
# Access to trajectory values
get_var_values(σ::AbstractTrajectory, var::String) =
-@view σ.values[:,(σ.m)._map_obs_var_idx[var]]
-get_state(σ::AbstractTrajectory, idx::Int) = @view σ.values[idx,:]
+ view(σ.values, :, (σ.m)._map_obs_var_idx[var])
+get_state(σ::AbstractTrajectory, idx::Int) =
+ view(σ.values, idx, :)
get_value(σ::AbstractTrajectory, var::String, idx::Int) =
-σ.values[idx,(σ.m)._map_obs_var_idx[var]]
-function δ(σ1::AbstractTrajectory,t::Float64) end
-
-# Get var values ["I"]
-function getindex(σ::AbstractTrajectory, var::String)
- if var == "times"
- return @view σ.times[:]
- elseif var == "transitions"
- return @view σ.transitions[:]
- else
- return get_var_values(σ, var)
- end
-end
+ σ.values[idx,(σ.m)._map_obs_var_idx[var]]
+
+
+δ(σ::AbstractTrajectory,t::Float64) = true
+
+states(σ::AbstractTrajectory) = σ.values
+times(σ::AbstractTrajectory) = σ.times
+transitions(σ::AbstractTrajectory) = σ.transitions
+
# Get i-th state [i]
getindex(σ::AbstractTrajectory, idx::Int) = get_state(σ, idx)
# Get i-th value of var ["I", idx]
-function getindex(σ::AbstractTrajectory, var::String, idx::Int)
- if var == "times"
- return σ.times[idx]
- elseif var == "transitions"
- return σ.transitions[idx]
- else
- return get_value(σ, var, idx)
- end
-end
+getindex(σ::AbstractTrajectory, var::String, idx::Int) = get_value(σ, var, idx)
+# Get the path of a variable ["I"]
+getindex(σ::AbstractTrajectory, var::String) = get_var_values(σ, var)
diff --git a/models/ER.jl b/models/ER.jl
index 3c694b9..243e0b8 100644
--- a/models/ER.jl
+++ b/models/ER.jl
@@ -45,6 +45,6 @@ is_absorbing_ER(p::Vector{Float64},xn::SubArray{Int,1}) =
(p[1]*xn[1]*xn[2] + (p[2]+p[3])*xn[3]) === 0.0
g = ["P"]
-ER = CTMC(d,k,dict_var,dict_p,l_tr,p,x0,t0,ER_f!,is_absorbing_ER; g=g)
+ER = ContinuousTimeModel(d,k,dict_var,dict_p,l_tr,p,x0,t0,ER_f!,is_absorbing_ER; g=g)
export ER
diff --git a/models/SIR.jl b/models/SIR.jl
index 215a4a8..b3ab3aa 100644
--- a/models/SIR.jl
+++ b/models/SIR.jl
@@ -34,7 +34,7 @@ function SIR_f!(mat_x::Matrix{Int}, l_t::Vector{Float64}, l_tr::Vector{String},
b_sup += l_a[i+1]
end
- nu = @view l_nu[:,reaction] # macro for avoiding a copy
+ nu = view(l_nu, :, reaction) # macro for avoiding a copy
for i = 1:3
mat_x[idx,i] = xn[i]+nu[i]
end
@@ -44,6 +44,6 @@ end
is_absorbing_SIR(p::Vector{Float64}, xn::SubArray{Int,1}) = (p[1]*xn[1]*xn[2] + p[2]*xn[2]) === 0.0
g = ["I"]
-SIR = CTMC(d,k,dict_var,dict_p,l_tr,p,x0,t0,SIR_f!,is_absorbing_SIR; g=g)
+SIR = ContinuousTimeModel(d,k,dict_var,dict_p,l_tr,p,x0,t0,SIR_f!,is_absorbing_SIR; g=g)
export SIR
diff --git a/models/_bench_perf_test/ER_col.jl b/models/_bench_perf_test/ER_col.jl
index f237a78..2546fd4 100644
--- a/models/_bench_perf_test/ER_col.jl
+++ b/models/_bench_perf_test/ER_col.jl
@@ -46,6 +46,6 @@ is_absorbing_ER_col(p::Vector{Float64},xn::AbstractVector{Int}) =
(p[1]*xn[1]*xn[2] + (p[2]+p[3])*xn[3]) === 0.0
g = ["P"]
-ER_col = CTMC(d,k,dict_var,dict_p,l_tr,p,x0,t0,ER_col_f!,is_absorbing_ER_col; g=g)
+ER_col = ContinuousTimeModel(d,k,dict_var,dict_p,l_tr,p,x0,t0,ER_col_f!,is_absorbing_ER_col; g=g)
export ER_col
diff --git a/models/_bench_perf_test/ER_col_buffer.jl b/models/_bench_perf_test/ER_col_buffer.jl
index c1c640a..c7315c4 100644
--- a/models/_bench_perf_test/ER_col_buffer.jl
+++ b/models/_bench_perf_test/ER_col_buffer.jl
@@ -45,6 +45,6 @@ is_absorbing_ER_col_buffer(p::Vector{Float64},xn::AbstractVector{Int}) =
(p[1]*xn[1]*xn[2] + (p[2]+p[3])*xn[3]) === 0.0
g = ["P"]
-ER_col_buffer = CTMC(d,k,dict_var,dict_p,l_tr,p,x0,t0,ER_col_buffer_f!,is_absorbing_ER_col_buffer; g=g)
+ER_col_buffer = ContinuousTimeModel(d,k,dict_var,dict_p,l_tr,p,x0,t0,ER_col_buffer_f!,is_absorbing_ER_col_buffer; g=g)
export ER_col_buffer
diff --git a/models/_bench_perf_test/ER_row_buffer.jl b/models/_bench_perf_test/ER_row_buffer.jl
index bc6104d..5c1489d 100644
--- a/models/_bench_perf_test/ER_row_buffer.jl
+++ b/models/_bench_perf_test/ER_row_buffer.jl
@@ -45,6 +45,6 @@ is_absorbing_ER_row_buffer(p::Vector{Float64},xn::AbstractVector{Int}) =
(p[1]*xn[1]*xn[2] + (p[2]+p[3])*xn[3]) === 0.0
g = ["P"]
-ER_row_buffer = CTMC(d,k,dict_var,dict_p,l_tr,p,x0,t0,ER_row_buffer_f!,is_absorbing_ER_row_buffer; g=g)
+ER_row_buffer = ContinuousTimeModel(d,k,dict_var,dict_p,l_tr,p,x0,t0,ER_row_buffer_f!,is_absorbing_ER_row_buffer; g=g)
export ER_row_buffer
diff --git a/models/_bench_perf_test/SIR_col.jl b/models/_bench_perf_test/SIR_col.jl
index 09d2b2a..3052b5e 100644
--- a/models/_bench_perf_test/SIR_col.jl
+++ b/models/_bench_perf_test/SIR_col.jl
@@ -44,6 +44,6 @@ end
is_absorbing_SIR_col(p::Vector{Float64}, xn::AbstractVector{Int}) = (p[1]*xn[1]*xn[2] + p[2]*xn[2]) === 0.0
g = ["I"]
-SIR_col = CTMC(d,k,dict_var,dict_p,l_tr,p,x0,t0,SIR_col_f!,is_absorbing_SIR_col; g=g)
+SIR_col = ContinuousTimeModel(d,k,dict_var,dict_p,l_tr,p,x0,t0,SIR_col_f!,is_absorbing_SIR_col; g=g)
export SIR_col
diff --git a/models/_bench_perf_test/SIR_col_buffer.jl b/models/_bench_perf_test/SIR_col_buffer.jl
index 9f79b3e..ee70393 100644
--- a/models/_bench_perf_test/SIR_col_buffer.jl
+++ b/models/_bench_perf_test/SIR_col_buffer.jl
@@ -44,6 +44,6 @@ end
is_absorbing_SIR_col_buffer(p::Vector{Float64}, xn::AbstractVector{Int}) = (p[1]*xn[1]*xn[2] + p[2]*xn[2]) === 0.0
g = ["I"]
-SIR_col_buffer = CTMC(d,k,dict_var,dict_p,l_tr,p,x0,t0,SIR_col_buffer_f!,is_absorbing_SIR_col_buffer; g=g)
+SIR_col_buffer = ContinuousTimeModel(d,k,dict_var,dict_p,l_tr,p,x0,t0,SIR_col_buffer_f!,is_absorbing_SIR_col_buffer; g=g)
export SIR_col_buffer
diff --git a/models/_bench_perf_test/SIR_row_buffer.jl b/models/_bench_perf_test/SIR_row_buffer.jl
index fcf6a86..4eaed57 100644
--- a/models/_bench_perf_test/SIR_row_buffer.jl
+++ b/models/_bench_perf_test/SIR_row_buffer.jl
@@ -44,6 +44,6 @@ end
is_absorbing_SIR_row_buffer(p::Vector{Float64}, xn::AbstractVector{Int}) = (p[1]*xn[1]*xn[2] + p[2]*xn[2]) === 0.0
g = ["I"]
-SIR_row_buffer = CTMC(d,k,dict_var,dict_p,l_tr,p,x0,t0,SIR_row_buffer_f!,is_absorbing_SIR_row_buffer; g=g)
+SIR_row_buffer = ContinuousTimeModel(d,k,dict_var,dict_p,l_tr,p,x0,t0,SIR_row_buffer_f!,is_absorbing_SIR_row_buffer; g=g)
export SIR_row_buffer
diff --git a/tests/dist_lp/dist_case_2.jl b/tests/dist_lp/dist_case_2.jl
index ce29df4..be7c86d 100644
--- a/tests/dist_lp/dist_case_2.jl
+++ b/tests/dist_lp/dist_case_2.jl
@@ -23,9 +23,9 @@ for p = 1:2
l_tr = Vector{Nothing}(nothing, length(x_obs))
σ1 = Trajectory(SIR, values, t_x, l_tr)
- test_1 = isapprox(dist_lp(@view(x_obs[:]), @view(t_x[:]), @view(y_obs[:]), @view(t_y[:]); p=p), res; atol = 1E-10)
+ test_1 = isapprox(dist_lp(@view(x_obs[:]), t_x, @view(y_obs[:]), t_y; p=p), res; atol = 1E-10)
test_1 = test_1 && isapprox(dist_lp(σ1,σ2; p=p), res; atol = 1E-10)
- test_1_bis = isapprox(dist_lp(@view(y_obs[:]), @view(t_y[:]), @view(x_obs[:]), @view(t_x[:]); p=p), res; atol = 1E-10)
+ test_1_bis = isapprox(dist_lp(@view(y_obs[:]), t_y, @view(x_obs[:]), t_x; p=p), res; atol = 1E-10)
test_1_bis = test_1_bis && isapprox(dist_lp(σ2,σ1; p=p), res; atol = 1E-10)
f_x(t::Real) = MarkovProcesses._f_step(x_obs, t_x, t)
diff --git a/tests/dist_lp/dist_case_3.jl b/tests/dist_lp/dist_case_3.jl
index c7d9558..b271771 100644
--- a/tests/dist_lp/dist_case_3.jl
+++ b/tests/dist_lp/dist_case_3.jl
@@ -26,9 +26,9 @@ for p = 1:2
int_riemann = MarkovProcesses._riemann_sum(diff_f, 0.0, 20.0, 1E-5)
int_riemann = int_riemann^(1/p)
- res1 = dist_lp(@view(x_obs[:]), @view(t_x[:]), @view(y_obs[:]), @view(t_y[:]); p=p)
+ res1 = dist_lp(@view(x_obs[:]), t_x, @view(y_obs[:]), t_y; p=p)
res2 = dist_lp(σ1,σ2; p=p)
- res1_bis = dist_lp(@view(y_obs[:]), @view(t_y[:]), @view(x_obs[:]), @view(t_x[:]); p=p)
+ res1_bis = dist_lp(@view(y_obs[:]), t_y, @view(x_obs[:]), t_x; p=p)
res2_bis = dist_lp(σ2,σ1; p=p)
test_1 = isapprox(res1, int_riemann; atol = 1E-3)
test_1 = test_1 && isapprox(res2, int_riemann; atol = 1E-3)
diff --git a/tests/dist_lp/dist_l1_case_1.jl b/tests/dist_lp/dist_l1_case_1.jl
index d02b7f3..a9152e0 100644
--- a/tests/dist_lp/dist_l1_case_1.jl
+++ b/tests/dist_lp/dist_l1_case_1.jl
@@ -21,7 +21,7 @@ values[:,1] = y_obs
l_tr = Vector{Nothing}(nothing, length(y_obs))
σ2 = Trajectory(SIR, values, t_y, l_tr)
-test_1 = dist_lp(@view(x_obs[:]), @view(t_x[:]), @view(y_obs[:]), @view(t_y[:]); p=1) == 6.4
+test_1 = dist_lp(@view(x_obs[:]), t_x, @view(y_obs[:]), t_y; p=1) == 6.4
test_1 = test_1 && dist_lp(σ1, σ2; p=1) == 6.4
f_x(t::Real) = MarkovProcesses._f_step(x_obs, t_x, t)
@@ -33,7 +33,7 @@ test_2 = isapprox(6.4, int; atol=err)
# Case 1 bis : inverse of case 1
-test_1_bis = dist_lp(@view(y_obs[:]), @view(t_y[:]), @view(x_obs[:]), @view(t_x[:]); p=1) == 6.4
+test_1_bis = dist_lp(@view(y_obs[:]), t_y, @view(x_obs[:]), t_x; p=1) == 6.4
test_1_bis = test_1_bis && dist_lp(σ2, σ1; p=1) == 6.4
return test_1 && test_1_bis && test_2
diff --git a/tests/dist_lp/dist_sim_sir.jl b/tests/dist_lp/dist_sim_sir.jl
index 5236040..13e3889 100644
--- a/tests/dist_lp/dist_sim_sir.jl
+++ b/tests/dist_lp/dist_sim_sir.jl
@@ -14,8 +14,8 @@ for p = 1:2
d = dist_lp(σ1, σ2, "I"; p = p)
d2 = dist_lp(σ1, σ2; p = p)
- f_x(t::Float64) = MarkovProcesses._f_step(σ1["I"], σ1["times"], t)
- f_y(t::Float64) = MarkovProcesses._f_step(σ2["I"], σ2["times"], t)
+ f_x(t::Float64) = MarkovProcesses._f_step(σ1["I"], times(σ1), t)
+ f_y(t::Float64) = MarkovProcesses._f_step(σ2["I"], times(σ2), t)
diff_f(t) = abs(f_x(t) - f_y(t))^p
int_riemann = MarkovProcesses._riemann_sum(diff_f, SIR.t0, SIR.time_bound, 1E-3)
res_int_riemann = int_riemann^(1/p)
diff --git a/tests/simulation/sim_er.jl b/tests/simulation/sim_er.jl
index 09d4575..21e18e2 100644
--- a/tests/simulation/sim_er.jl
+++ b/tests/simulation/sim_er.jl
@@ -6,7 +6,7 @@ load_model("ER")
σ = simulate(ER)
plt.figure()
-plt.step(σ["times"], σ["P"], "ro--", marker="x", where="post", linewidth=1.0)
+plt.step(times(σ), σ["P"], "ro--", marker="x", where="post", linewidth=1.0)
plt.savefig(get_module_path() * "/tests/simulation/res_pics/sim_er.png", dpi=480)
plt.close()
diff --git a/tests/simulation/sim_er_row_buffer_bounded.jl b/tests/simulation/sim_er_row_buffer_bounded.jl
index 6da8376..d684462 100644
--- a/tests/simulation/sim_er_row_buffer_bounded.jl
+++ b/tests/simulation/sim_er_row_buffer_bounded.jl
@@ -8,7 +8,7 @@ ER_row_buffer.time_bound = 10.0
σ = _simulate_row_buffer(ER_row_buffer)
plt.figure()
-plt.step(σ["times"], _get_values_row(σ,"P"), "ro--", marker="x", where="post", linewidth=1.0)
+plt.step(times(σ), _get_values_row(σ,"P"), "ro--", marker="x", where="post", linewidth=1.0)
plt.savefig(get_module_path() * "/tests/simulation/res_pics/sim_er_row_buffer_bounded.png")
plt.close()
diff --git a/tests/simulation/sim_sir.jl b/tests/simulation/sim_sir.jl
index 2a1f00f..74fd911 100644
--- a/tests/simulation/sim_sir.jl
+++ b/tests/simulation/sim_sir.jl
@@ -6,7 +6,7 @@ load_model("SIR")
σ = simulate(SIR)
plt.figure()
-plt.step(σ["times"], σ["I"], "ro--", marker="x", where="post", linewidth=1.0)
+plt.step(times(σ), σ["I"], "ro--", marker="x", where="post", linewidth=1.0)
plt.savefig(get_module_path() * "/tests/simulation/res_pics/sim_sir.png")
plt.close()
diff --git a/tests/simulation/sim_sir_bounded.jl b/tests/simulation/sim_sir_bounded.jl
index 95d6434..469e784 100644
--- a/tests/simulation/sim_sir_bounded.jl
+++ b/tests/simulation/sim_sir_bounded.jl
@@ -7,7 +7,7 @@ SIR.time_bound = 100.0
σ = simulate(SIR)
plt.figure()
-plt.step(σ["times"], σ["I"], "ro--", marker="x", where="post", linewidth=1.0)
+plt.step(times(σ), σ["I"], "ro--", marker="x", where="post", linewidth=1.0)
plt.savefig(get_module_path() * "/tests/simulation/res_pics/sim_sir_bounded.png")
plt.close()
diff --git a/tests/simulation/sim_sir_col_buffer_bounded.jl b/tests/simulation/sim_sir_col_buffer_bounded.jl
index 3343495..524487b 100644
--- a/tests/simulation/sim_sir_col_buffer_bounded.jl
+++ b/tests/simulation/sim_sir_col_buffer_bounded.jl
@@ -8,7 +8,7 @@ SIR_col_buffer.time_bound = 100.0
σ = _simulate_col_buffer(SIR_col_buffer)
plt.figure()
-plt.step(σ["times"], _get_values_col(σ,"I"), "ro--", marker="x", where="post", linewidth=1.0)
+plt.step(times(σ), _get_values_col(σ,"I"), "ro--", marker="x", where="post", linewidth=1.0)
plt.savefig(get_module_path() * "/tests/simulation/res_pics/sim_sir_col_buffer_bounded.png")
plt.close()
diff --git a/tests/simulation/sim_sir_row_buffer_bounded.jl b/tests/simulation/sim_sir_row_buffer_bounded.jl
index b6a3898..b49c589 100644
--- a/tests/simulation/sim_sir_row_buffer_bounded.jl
+++ b/tests/simulation/sim_sir_row_buffer_bounded.jl
@@ -8,7 +8,7 @@ SIR_row_buffer.time_bound = 100.0
σ = _simulate_row_buffer(SIR_row_buffer)
plt.figure()
-plt.step(σ["times"], _get_values_row(σ,"I"), "ro--", marker="x", where="post", linewidth=1.0)
+plt.step(times(σ), _get_values_row(σ,"I"), "ro--", marker="x", where="post", linewidth=1.0)
plt.savefig(get_module_path() * "/tests/simulation/res_pics/sim_sir_row_buffer_bounded.png")
plt.close()
diff --git a/tests/unit/simulate_sir_bounded.jl b/tests/unit/simulate_sir_bounded.jl
index 0a7a10f..c872c02 100644
--- a/tests/unit/simulate_sir_bounded.jl
+++ b/tests/unit/simulate_sir_bounded.jl
@@ -11,7 +11,7 @@ d2 = Dict("I" => 1)
bool_test = SIR.g == ["I"] && SIR._g_idx == [2] &&
SIR.map_var_idx == d1 &&
- SIR._map_obs_var_idx == d2 && σ["times"][end] <= SIR.time_bound
+ SIR._map_obs_var_idx == d2 && times(σ)[end] <= SIR.time_bound
return bool_test
--
GitLab