diff --git a/core/MarkovProcesses.jl b/core/MarkovProcesses.jl
index 476b9e55c2bf365b8581b64fb2c91a8c5b269777..08c2b4a5d208fab06321f6c498d554185fefafa3 100644
--- a/core/MarkovProcesses.jl
+++ b/core/MarkovProcesses.jl
@@ -22,14 +22,18 @@ export load_automaton, get_index, get_value, length_var
# Model related methods
export simulate, set_param!, get_param, set_observed_var!
-export is_bounded
+export is_bounded, check_consistency
export load_model, get_module_path
+# Utils
+export get_module_path
+
include("common.jl")
include("trajectory.jl")
include("lha.jl")
include("model.jl")
+include("utils.jl")
end
diff --git a/core/common.jl b/core/common.jl
index 0f6cb04decded7f70e825195f9bf26998e404ccb..2e9decc0d8d0da38fddee616f3f5496a33456441 100644
--- a/core/common.jl
+++ b/core/common.jl
@@ -88,4 +88,6 @@ end
LHA(A::LHA, map_var::Dict{String,Int}) = LHA(A.l_transitions, A.l_loc, A.Λ,
A.l_loc_init, A.l_loc_final, A.map_var_automaton_idx, A.l_flow,
A.map_edges, A.l_ctes, map_var)
+Base.:*(m::ContinuousTimeModel, A::LHA) = SynchronizedModel(m, A)
+Base.:*(A::LHA, m::ContinuousTimeModel) = SynchronizedModel(m, A)
diff --git a/core/lha.jl b/core/lha.jl
index 4c586bde4049299a752a72b7a1a18c22eea07d08..3221688a35ad8f906ad4ce878d18ade16d814ae4 100644
--- a/core/lha.jl
+++ b/core/lha.jl
@@ -19,6 +19,8 @@ function Base.show(io::IO, S::StateLHA)
end
end
+isaccepted(S::StateLHA) = S.loc in (S.A).l_loc_final
+
# Methods for synchronize / read the trajectory
function init_state(A::LHA, x0::SubArray{Int,1}, t0::Float64)
S0 = StateLHA(A, "", zeros(length_var(A)), t0)
diff --git a/core/model.jl b/core/model.jl
index 82505371324f9ec8ac24deccf9a0260e5208168c..6824680f16338fca690c70a13d1b3aef32d06af4 100644
--- a/core/model.jl
+++ b/core/model.jl
@@ -1,4 +1,6 @@
+load_model(name_model::String) = include(get_module_path() * "/models/" * name_model * ".jl")
+
# Simulation
function simulate(m::ContinuousTimeModel)
# trajectory fields
@@ -120,14 +122,18 @@ function set_observed_var!(m::Model,g::Vector{String})
m._map_obs_var_idx = _map_obs_var_idx
end
-is_bounded(m::Model) = m.time_bound < Inf
-function check_consistency(m::Model) end
-function simulate(m::Model, n::Int; bound::Float64 = Inf)::AbstractObservations end
-function set_param!(m::Model, p::Vector{Float64})::Nothing end
-function get_param(m::Model)::Vector{Float64} end
-
-get_module_path() = dirname(dirname(pathof(@__MODULE__)))
-function load_model(name_model::String)
- include(get_module_path() * "/models/" * name_model * ".jl")
+is_bounded(m::ContinuousTimeModel) = m.time_bound < Inf
+function check_consistency(m::ContinuousTimeModel)
+ @assert m.d == length(m.map_var_idx)
+ @assert m.k == length(m.map_param_idx)
+ @assert m.k == length(m.p)
+ @assert length(m.g) <= m.d
+ @assert length(m._g_idx) == length(m.g)
+ @assert m.buffer_size >= 0
+ @assert typeof(m.is_absorbing(m.p, view(reshape(m.x0, 1, m.d), 1, :))) == Bool
+ return true
end
+set_param!(m::ContinuousTimeModel, p::Vector{Float64}) = (m.p = p)
+set_param!(m::ContinuousTimeModel, name_p::String, p_i::Float64) = (m.p[m.map_param_idx[name_p]] = p_i)
+get_param(m::ContinuousTimeModel) = m.p
diff --git a/core/utils.jl b/core/utils.jl
new file mode 100644
index 0000000000000000000000000000000000000000..4128525ddd046e04f908a87b609048dd7f9ce5a3
--- /dev/null
+++ b/core/utils.jl
@@ -0,0 +1,3 @@
+
+get_module_path() = dirname(dirname(pathof(@__MODULE__)))
+
diff --git a/models/ER.jl b/models/ER.jl
index 243e0b8aa01b599851d47551462518c5ac6128b6..b427f08a050838e6dc3df89084b343b85a6b5d30 100644
--- a/models/ER.jl
+++ b/models/ER.jl
@@ -5,10 +5,10 @@ d=4
k=3
dict_var = Dict("E" => 1, "S" => 2, "ES" => 3, "P" => 4)
dict_p = Dict("k1" => 1, "k2" => 2, "k3" => 3)
-l_tr = ["R1","R2","R3"]
-p = [1.0, 1.0, 1.0]
-x0 = [100, 100, 0, 0]
-t0 = 0.0
+l_tr_ER = ["R1","R2","R3"]
+p_ER = [1.0, 1.0, 1.0]
+x0_ER = [100, 100, 0, 0]
+t0_ER = 0.0
function ER_f!(mat_x::Matrix{Int}, l_t::Vector{Float64}, l_tr::Vector{String}, idx::Int,
xn::SubArray{Int,1}, tn::Float64, p::Vector{Float64})
a1 = p[1] * xn[1] * xn[2]
@@ -43,8 +43,8 @@ function ER_f!(mat_x::Matrix{Int}, l_t::Vector{Float64}, l_tr::Vector{String}, i
end
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"]
+g_ER = ["P"]
-ER = ContinuousTimeModel(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_ER,p_ER,x0_ER,t0_ER,ER_f!,is_absorbing_ER; g=g_ER)
export ER
diff --git a/models/SIR.jl b/models/SIR.jl
index b3ab3aa31663a921638595a7af803779d82f2a81..d3fff557a33a5c1c0127ecbabc876b576fce853a 100644
--- a/models/SIR.jl
+++ b/models/SIR.jl
@@ -5,10 +5,10 @@ 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 = [0.0012, 0.05]
-x0 = [95, 5, 0]
-t0 = 0.0
+l_tr_SIR = ["R1","R2"]
+p_SIR = [0.0012, 0.05]
+x0_SIR = [95, 5, 0]
+t0_SIR = 0.0
function SIR_f!(mat_x::Matrix{Int}, l_t::Vector{Float64}, l_tr::Vector{String}, idx::Int,
xn::SubArray{Int,1}, tn::Float64, p::Vector{Float64})
a1 = p[1] * xn[1] * xn[2]
@@ -42,8 +42,8 @@ function SIR_f!(mat_x::Matrix{Int}, l_t::Vector{Float64}, l_tr::Vector{String},
l_tr[idx] = "R$(reaction)"
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"]
+g_SIR = ["I"]
-SIR = ContinuousTimeModel(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_SIR,p_SIR,x0_SIR,t0_SIR,SIR_f!,is_absorbing_SIR; g=g_SIR)
export SIR
diff --git a/tests/automata/read_trajectory_last_state_F.jl b/tests/automata/read_trajectory_last_state_F.jl
index 636080d5abe3e04ccbc4ad9b18f1ab115c47cede..17f02b82e2e99c45d443e54072a4ded26b044117 100644
--- a/tests/automata/read_trajectory_last_state_F.jl
+++ b/tests/automata/read_trajectory_last_state_F.jl
@@ -11,7 +11,7 @@ A_F = create_automaton_F(SIR, x1, x2, t1, t2, "I") # <: LHA
function test_last_state(A::LHA, m::ContinuousTimeModel)
σ = simulate(m)
Send = read_trajectory(A, σ)
- test = (get_state_from_time(σ, Send.time)[1] == Send["n"]) & (Send["d"] == 0)
+ test = (get_state_from_time(σ, Send.time)[1] == Send["n"]) && (Send["d"] == 0)
#=
if !test
@show Send
@@ -25,7 +25,7 @@ test_all = true
nbr_sim = 10000
for i = 1:nbr_sim
test = test_last_state(A_F, SIR)
- global test_all = test_all & test
+ global test_all = test_all && test
end
return test_all
diff --git a/tests/run_unit.jl b/tests/run_unit.jl
index f334c5dd9fdfe5f780f52c3429ecdc5ca587a16d..49dfe1f240d1901d6d3945d4d99d216fe3d8e575 100644
--- a/tests/run_unit.jl
+++ b/tests/run_unit.jl
@@ -17,5 +17,8 @@ using Test
@test include("unit/dist_lp.jl")
@test include("unit/l_dist_lp.jl")
@test include("unit/check_trajectory_consistency.jl")
+ @test include("unit/check_model_consistency.jl")
+ @test include("unit/set_param.jl")
+ @test include("unit/side_effects_1.jl")
end
diff --git a/tests/unit/check_model_consistency.jl b/tests/unit/check_model_consistency.jl
new file mode 100644
index 0000000000000000000000000000000000000000..5050eff7b89f6506c81269a424cc444323385e8d
--- /dev/null
+++ b/tests/unit/check_model_consistency.jl
@@ -0,0 +1,8 @@
+
+using MarkovProcesses
+
+load_model("SIR")
+load_model("ER")
+
+return check_consistency(ER) && check_consistency(SIR)
+
diff --git a/tests/unit/set_param.jl b/tests/unit/set_param.jl
new file mode 100644
index 0000000000000000000000000000000000000000..222c845b61c7a11c89a93816897fc3b4b92f242c
--- /dev/null
+++ b/tests/unit/set_param.jl
@@ -0,0 +1,38 @@
+
+using MarkovProcesses
+
+load_model("SIR") # p = 0.0012, 0.05
+load_model("ER") # p = 1; 1; 1
+
+test_all = true
+
+new_param_ER = [1.0, 2.0, 1.5]
+set_param!(ER, new_param_ER)
+test_all = test_all && (ER.p == new_param_ER)
+
+k2 = 4.0
+set_param!(ER, "k2", k2)
+test_all = test_all && (ER.p == [1.0, 4.0, 1.5])
+
+k1 = 0.5
+set_param!(ER, "k1", k1)
+test_all = test_all && (ER.p == [0.5, 4.0, 1.5])
+
+k3 = 10.0
+set_param!(ER, "k3", 10.0)
+test_all = test_all && (ER.p == [0.5, 4.0, 10.0])
+
+new_param_SIR = [0.0013, 0.08]
+set_param!(SIR, new_param_SIR)
+test_all = test_all && (SIR.p == new_param_SIR)
+
+kr = 0.06
+set_param!(SIR, "kr", kr)
+test_all = test_all && (SIR.p == [0.0013, 0.06])
+
+ki = 0.011
+set_param!(SIR, "ki", ki)
+test_all = test_all && (SIR.p == [0.011, 0.06])
+
+return test_all
+
diff --git a/tests/unit/side_effects_1.jl b/tests/unit/side_effects_1.jl
new file mode 100644
index 0000000000000000000000000000000000000000..3392843aa8b2e1ab203523557ebadb4ebd551e8d
--- /dev/null
+++ b/tests/unit/side_effects_1.jl
@@ -0,0 +1,16 @@
+
+using MarkovProcesses
+
+load_model("SIR")
+load_model("ER")
+
+test = SIR.map_var_idx !== ER.map_var_idx &&
+ SIR.map_param_idx !== ER.map_var_idx &&
+ SIR.l_transitions !== ER.l_transitions &&
+ SIR.g !== ER.g &&
+ SIR._g_idx !== ER._g_idx &&
+ SIR.f! != ER.f! &&
+ SIR.is_absorbing != ER.is_absorbing
+
+return test
+