Skip to content
Snippets Groups Projects
Commit 886e9d98 authored by Bentriou Mahmoud's avatar Bentriou Mahmoud
Browse files

Creation of automata has been change again to improve performance. Now 

check_constraint and update_state functions are created inside
create_automaton_*(), by the evaluation of meta_elementary_functions <:
Expr. check_constraint and update_state don't depend on constants::Dict
anymore: constants are interpolated inside the expression.

Add of tests for euclidean distance automaton.
Add a new function euclidean_distance that computes the euclidean
distance of a Trajectory compared to observations along a discrete
timeline.
parent b1c40b0f
No related branches found
No related tags found
No related merge requests found
Expand all Hide whitespace changes
Inline Side-by-side
Showing
with 827 additions and 748 deletions
automata/automaton_F.jl
+ 99
− 99
View file @ 886e9d98
@everywhere istrue(val::Float64) = convert(Bool, val)
# Invariant predicates functions
@everywhere true_inv_predicate(x::Vector{Int}) = true
# Check constraints and update state functions
# l0 loc : we construct the edges of the form l0 => (..)
# "cc" as check_constraints and "us" as update_state
# l0 => l1
@everywhere cc_aut_F_l0l1_1(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) = true
# us_aut_F_l0l1_1! inside create_automaton_F
# l1 loc
# l1 => l2
@everywhere cc_aut_F_l1l2_1(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =
getfield(S, :time) >= constants[:t1] &&
(constants[:x1] <= S[:n] <= constants[:x2])
@everywhere us_aut_F_l1l2_1!(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =
(S.loc = :l2;
S[:d] = 0)
@everywhere cc_aut_F_l1l2_4(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =
getfield(S, :time) >= constants[:t1] &&
S[:d] == 0
@everywhere us_aut_F_l1l2_4!(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =
(S.loc = :l2)
@everywhere cc_aut_F_l1l2_2(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =
(getfield(S, :time) >= constants[:t2]) &&
(S[:n] < constants[:x1] || S[:n] > constants[:x2])
@everywhere us_aut_F_l1l2_2!(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =
(S.loc = :l2;
S[:d] = min(abs(S[:n] - constants[:x1]), abs(S[:n] - constants[:x2])))
@everywhere cc_aut_F_l1l2_3(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =
istrue(S[:isabs]) && getfield(S, :time) <= constants[:t2]
@everywhere us_aut_F_l1l2_3!(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =
(S.loc = :l2)
# l1 => l3
@everywhere cc_aut_F_l1l3_1(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =
(constants[:x1] <= S[:n] <= constants[:x2])
@everywhere us_aut_F_l1l3_1!(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =
(S.loc = :l3;
S[:d] = 0;)
@everywhere cc_aut_F_l1l3_2(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =
(S[:n] < constants[:x1] || S[:n] > constants[:x2]) &&
(getfield(S, :time) <= constants[:t1])
@everywhere us_aut_F_l1l3_2!(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =
(S.loc = :l3;
S[:d] = min(sqrt((getfield(S, :time) - constants[:t1])^2 + (S[:n] - constants[:x2])^2),
sqrt((getfield(S, :time) - constants[:t1])^2 + (S[:n] - constants[:x1])^2)))
@everywhere cc_aut_F_l1l3_3(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =
(S[:n] < constants[:x1] || S[:n] > constants[:x2]) &&
(constants[:t1] <= getfield(S, :time) <= constants[:t2])
@everywhere us_aut_F_l1l3_3!(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =
(S.loc = :l3;
S[:d] = min(S[:d], min(abs(S[:n] - constants[:x1]), abs(S[:n] - constants[:x2]))))
# l3 loc
# l3 => l1
@everywhere cc_aut_F_l3l1_1(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) = true
# l3 => l2
@everywhere cc_aut_F_l3l2_1(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =
(getfield(S, :time) >= constants[:t2])
@everywhere us_aut_F_l3l2_1!(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =
(S.loc = :l2)
function create_automaton_F(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1::Float64, t2::Float64, sym_obs::VariableModel)
# Requirements for the automaton
@assert sym_obs in m.g "$(sym_obs) is not observed."
......@@ -81,6 +9,7 @@ function create_automaton_F(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1
locations = [:l0, :l1, :l2, :l3]
## Invariant predicates
@everywhere true_inv_predicate(x::Vector{Int}) = true
Λ_F = Dict(:l0 => getfield(Main, :true_inv_predicate), :l1 => getfield(Main, :true_inv_predicate),
:l2 => getfield(Main, :true_inv_predicate), :l3 => getfield(Main, :true_inv_predicate))
......@@ -106,53 +35,124 @@ function create_automaton_F(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1
map_edges[loc] = Dict{Location, Vector{Edge}}()
end
sym_isabs_func = Symbol(m.isabsorbing)
idx_obs_var = getfield(m, :map_var_idx)[sym_obs]
nbr_rand = rand(1:1000)
idx_var_n = map_var_automaton_idx[:n]
idx_var_d = map_var_automaton_idx[:d]
idx_var_isabs = map_var_automaton_idx[:isabs]
nbr_rand = rand(1:100000)
basename_func = "$(replace(m.name, ' '=>'_'))_$(nbr_rand)"
basename_func = replace(basename_func, '-'=>'_')
sym_isabs_func = Symbol(m.isabsorbing)
func_name(type_func::Symbol, from_loc::Location, to_loc::Location, edge_number::Int) =
Symbol("$(type_func)_aut_F_$(basename_func)_$(from_loc)$(to_loc)_$(edge_number)$(type_func == :us ? "!" : "")")
meta_elementary_functions = quote
@everywhere istrue(val::Float64) = convert(Bool, val)
## Check constraints and update state functions
# l0 loc : we construct the edges of the form l0 => (..)
# "cc" as check_constraints and "us" as update_state
# l0 => l1
@everywhere $(func_name(:cc, :l0, :l1, 1))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) = true
@everywhere $(func_name(:us, :l0, :l1, 1))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
(S.loc = :l1;
getfield(S, :values)[$(idx_var_n)] = x[$(idx_obs_var)];
setindex!(getfield(S, :values), Inf, $(idx_var_d));
setindex!(getfield(S, :values), getfield(Main, $(Meta.quot(sym_isabs_func)))(p, x), $(idx_var_isabs)))
# l1 loc
# l1 => l2
@everywhere $(func_name(:cc, :l1, :l2, 1))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
getfield(S, :time) >= $t1 &&
($x1 <= getfield(S, :values)[$(idx_var_n)] <= $x2)
@everywhere $(func_name(:us, :l1, :l2, 1))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
(S.loc = :l2;
setindex!(getfield(S, :values), 0, $(idx_var_d)))
@everywhere $(func_name(:cc, :l1, :l2, 4))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
getfield(S, :time) >= $t1 &&
getfield(S, :values)[$(idx_var_d)] == 0
@everywhere $(func_name(:us, :l1, :l2, 4))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
(S.loc = :l2)
@everywhere $(func_name(:cc, :l1, :l2, 2))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
(getfield(S, :time) >= $t2) &&
(getfield(S, :values)[$(idx_var_n)] < $x1 || getfield(S, :values)[$(idx_var_n)] > $x2)
@everywhere $(func_name(:us, :l1, :l2, 2))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
(S.loc = :l2;
setindex!(getfield(S, :values), min(abs(getfield(S, :values)[$(idx_var_n)] - $x1), abs(getfield(S, :values)[$(idx_var_n)] - $x2)), $(idx_var_d)))
@everywhere $(func_name(:cc, :l1, :l2, 3))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
istrue(getfield(S, :values)[$(idx_var_isabs)]) && getfield(S, :time) <= $t2
@everywhere $(func_name(:us, :l1, :l2, 3))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
(S.loc = :l2)
# l1 => l3
@everywhere $(func_name(:cc, :l1, :l3, 1))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
($x1 <= getfield(S, :values)[$(idx_var_n)] <= $x2)
@everywhere $(func_name(:us, :l1, :l3, 1))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
(S.loc = :l3;
setindex!(getfield(S, :values), 0, $(idx_var_d)))
@everywhere $(func_name(:cc, :l1, :l3, 2))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
(getfield(S, :values)[$(idx_var_n)] < $x1 || getfield(S, :values)[$(idx_var_n)] > $x2) &&
(getfield(S, :time) <= $t1)
@everywhere $(func_name(:us, :l1, :l3, 2))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
(S.loc = :l3;
setindex!(getfield(S, :values), min(sqrt((getfield(S, :time) - $t1)^2 + (getfield(S, :values)[$(idx_var_n)] - $x2)^2),
sqrt((getfield(S, :time) - $t1)^2 + (getfield(S, :values)[$(idx_var_n)] - $x1)^2)), $(idx_var_d)))
@everywhere $(func_name(:cc, :l1, :l3, 3))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
(getfield(S, :values)[$(idx_var_n)] < $x1 || getfield(S, :values)[$(idx_var_n)] > $x2) &&
($t1 <= getfield(S, :time) <= $t2)
@everywhere $(func_name(:us, :l1, :l3, 3))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
(S.loc = :l3;
val_min = min(getfield(S, :values)[$(idx_var_d)],
min(abs(getfield(S, :values)[$(idx_var_n)] - $x1), abs(getfield(S, :values)[$(idx_var_n)] - $x2)));
setindex!(getfield(S, :values), val_min, $(idx_var_d)))
# l3 loc
# l3 => l1
@everywhere $(func_name(:cc, :l3, :l1, 1))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) = true
@everywhere $(func_name(:us, :l3, :l1, 1))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
(S.loc = :l1;
getfield(S, :values)[$(idx_var_n)] = x[$(idx_obs_var)];
setindex!(getfield(S, :values), getfield(Main, $(Meta.quot(sym_isabs_func)))(p, x), $(idx_var_isabs)))
# l3 => l2
@everywhere $(func_name(:cc, :l3, :l2, 1))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
(getfield(S, :time) >= $t2)
@everywhere $(func_name(:us, :l3, :l2, 1))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
(S.loc = :l2)
end
eval(meta_elementary_functions)
# l0 loc
# l0 => l1
sym_func_us_l0l1_1 = Symbol("us_aut_F_$(basename_func)_l0l1_1!")
str_us_l0l1_1 = "
@everywhere $(sym_func_us_l0l1_1)(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) = \n
(S.loc = :l1; \n
S[:n] = x[$(idx_obs_var)];\n
S[:d] = Inf; \n
S[:isabs] = getfield(Main, $(Meta.quot(sym_isabs_func)))(p, x))"
eval(Meta.parse(str_us_l0l1_1))
edge1 = Edge([nothing], getfield(Main, :cc_aut_F_l0l1_1), getfield(Main, sym_func_us_l0l1_1))
edge1 = Edge([nothing], getfield(Main, func_name(:cc, :l0, :l1, 1)), getfield(Main, func_name(:us, :l0, :l1, 1)))
map_edges[:l0][:l1] = [edge1]
# l1 loc
# l1 => l2
edge1 = Edge([nothing], getfield(Main, :cc_aut_F_l1l2_1), getfield(Main, :us_aut_F_l1l2_1!))
edge2 = Edge([nothing], getfield(Main, :cc_aut_F_l1l2_2), getfield(Main, :us_aut_F_l1l2_2!))
edge3 = Edge([nothing], getfield(Main, :cc_aut_F_l1l2_3), getfield(Main, :us_aut_F_l1l2_3!))
edge4 = Edge([nothing], getfield(Main, :cc_aut_F_l1l2_4), getfield(Main, :us_aut_F_l1l2_4!))
edge1 = Edge([nothing], getfield(Main, func_name(:cc, :l1, :l2, 1)), getfield(Main, func_name(:us, :l1, :l2, 1)))
edge2 = Edge([nothing], getfield(Main, func_name(:cc, :l1, :l2, 2)), getfield(Main, func_name(:us, :l1, :l2, 2)))
edge3 = Edge([nothing], getfield(Main, func_name(:cc, :l1, :l2, 3)), getfield(Main, func_name(:us, :l1, :l2, 3)))
edge4 = Edge([nothing], getfield(Main, func_name(:cc, :l1, :l2, 4)), getfield(Main, func_name(:us, :l1, :l2, 4)))
map_edges[:l1][:l2] = [edge1, edge2, edge3, edge4]
# l1 => l3
edge1 = Edge([nothing], getfield(Main, :cc_aut_F_l1l3_1), getfield(Main, :us_aut_F_l1l3_1!))
edge2 = Edge([nothing], getfield(Main, :cc_aut_F_l1l3_2), getfield(Main, :us_aut_F_l1l3_2!))
edge3 = Edge([nothing], getfield(Main, :cc_aut_F_l1l3_3), getfield(Main, :us_aut_F_l1l3_3!))
edge1 = Edge([nothing], getfield(Main, func_name(:cc, :l1, :l3, 1)), getfield(Main, func_name(:us, :l1, :l3, 1)))
edge2 = Edge([nothing], getfield(Main, func_name(:cc, :l1, :l3, 2)), getfield(Main, func_name(:us, :l1, :l3, 2)))
edge3 = Edge([nothing], getfield(Main, func_name(:cc, :l1, :l3, 3)), getfield(Main, func_name(:us, :l1, :l3, 3)))
map_edges[:l1][:l3] = [edge1, edge2, edge3]
# l3 loc
# l3 => l1
sym_func_us_l3l1_1 = Symbol("us_aut_F_$(basename_func)_l0l1_1!")
str_us_l3l1_1 =
"@everywhere $(sym_func_us_l3l1_1)(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) = \n
(S.loc = :l1;\n
S[:n] = x[$(idx_obs_var)];\n
S[:isabs] = getfield(Main, $(Meta.quot(sym_isabs_func)))(p, x))"
eval(Meta.parse(str_us_l3l1_1))
edge1 = Edge([:ALL], getfield(Main, :cc_aut_F_l3l1_1), getfield(Main, sym_func_us_l3l1_1))
edge1 = Edge([:ALL], getfield(Main, func_name(:cc, :l3, :l1, 1)), getfield(Main, func_name(:us, :l3, :l1, 1)))
map_edges[:l3][:l1] = [edge1]
# l3 => l2
edge1 = Edge([nothing], getfield(Main, :cc_aut_F_l3l2_1), getfield(Main, :us_aut_F_l3l2_1!))
edge1 = Edge([nothing], getfield(Main, func_name(:cc, :l3, :l2, 1)), getfield(Main, func_name(:us, :l3, :l2, 1)))
map_edges[:l3][:l2] = [edge1]
## Constants
......
automata/automaton_G.jl
+ 158
− 161
View file @ 886e9d98
This diff is collapsed.
automata/automaton_G_and_F.jl
+ 250
− 260
View file @ 886e9d98
This diff is collapsed.
automata/euclidean_distance_automaton.jl
+ 61
− 49
View file @ 886e9d98
# l0 loc
# l0 => l1
@everywhere cc_eucl_dist_aut_l0l1_1(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) = true
# l1 loc
# l1 => l1
@everywhere cc_eucl_dist_aut_l1l1_1(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =
S[:t] >= constants[Symbol("tml_$(S[:idx])")]
@everywhere us_eucl_dist_aut_l1l1_1!(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =
(S[:d] += (S[:n] - constants[Symbol("y_$(S[:idx])")])^2;
S[:idx] += 1.0)
# l1 => l2
@everywhere cc_eucl_dist_aut_l1l2_1(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =
S[:idx] >= constants[:nbr_obs] + 1
@everywhere us_eucl_dist_aut_l1l2_1!(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =
(S.loc = :l2;
S[:d] = sqrt(S[:d]))
@everywhere cc_eucl_dist_aut_l1l1_2(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) = true
function create_euclidean_distance_automaton(m::ContinuousTimeModel, timeline::AbstractVector{Float64}, observations::AbstractVector{Float64}, sym_obs::VariableModel)
# Requirements for the automaton
@assert sym_obs in m.g "$(sym_obs) is not observed."
......@@ -40,16 +19,7 @@ function create_euclidean_distance_automaton(m::ContinuousTimeModel, timeline::A
map_var_automaton_idx = Dict{VariableAutomaton,Int}(:t => 1, :n => 2,
:d => 3, :idx => 4)
nbr_first_vars = 4
for i = 1:nbr_observations
map_var_automaton_idx[Symbol("tml_$(convert(Float64, i))")] = nbr_first_vars + i
end
for i = 1:nbr_observations
map_var_automaton_idx[Symbol("y_$(convert(Float64, i))")] = nbr_first_vars + nbr_observations + i
end
vector_flow = zeros(nbr_first_vars + 2*nbr_observations)
vector_flow[1] = 1.0
vector_flow = [1.0, 0.0, 0.0, 0.0]
flow = Dict{Location,Vector{Float64}}(:l0 => vector_flow,
:l1 => vector_flow,
:l2 => vector_flow)
......@@ -61,36 +31,78 @@ function create_euclidean_distance_automaton(m::ContinuousTimeModel, timeline::A
end
idx_obs_var = getfield(m, :map_var_idx)[sym_obs]
idx_var_t = map_var_automaton_idx[:t]
idx_var_n = map_var_automaton_idx[:n]
idx_var_d = map_var_automaton_idx[:d]
idx_var_idx = map_var_automaton_idx[:idx]
nbr_rand = rand(1:1000)
basename_func = "$(replace(m.name, ' '=>'_'))_$(nbr_rand)"
basename_func = replace(basename_func, '-'=>'_')
func_name(type_func::Symbol, from_loc::Location, to_loc::Location, edge_number::Int) =
Symbol("$(type_func)_eucl_dist_aut_$(basename_func)_$(from_loc)$(to_loc)_$(edge_number)$(type_func == :us ? "!" : "")")
meta_elementary_functions = quote
# l0 loc
# l0 => l1
@everywhere $(func_name(:cc, :l0, :l1, 1))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) = true
@everywhere $(func_name(:us, :l0, :l1, 1))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
(setfield!(S, :loc, Symbol("l1"));
setindex!(getfield(S, :values), x[$(idx_obs_var)], $(idx_var_n));
setindex!(getfield(S, :values), 0.0, $(idx_var_d));
setindex!(getfield(S, :values), 1.0, $(idx_var_idx)))
# l1 loc
# l1 => l1
# Defined below
@everywhere $(func_name(:cc, :l1, :l1, 1))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
getfield(S, :values)[$(idx_var_t)] >= $(timeline)[convert(Int, getfield(S, :values)[$(idx_var_idx)])]
@everywhere $(func_name(:us, :l1, :l1, 1))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
(setindex!(getfield(S, :values), getfield(S, :values)[$(idx_var_d)] +
(getfield(S, :values)[$(idx_var_n)] - $(observations)[convert(Int, getfield(S, :values)[$(idx_var_idx)])])^2,
$(idx_var_d));
setindex!(getfield(S, :values), getfield(S, :values)[$(idx_var_idx)] + 1.0, $(idx_var_idx)))
@everywhere $(func_name(:cc, :l1, :l1, 2))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) = true
@everywhere $(func_name(:us, :l1, :l1, 2))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
(setindex!(getfield(S, :values), x[$(idx_obs_var)], $(idx_var_n)))
# l1 => l2
@everywhere $(func_name(:cc, :l1, :l2, 1))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
getfield(S, :values)[$(idx_var_idx)] >= ($nbr_observations + 1)
@everywhere $(func_name(:us, :l1, :l2, 1))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
(setfield!(S, :loc, Symbol("l2"));
setindex!(getfield(S, :values), sqrt(getfield(S, :values)[$(idx_var_d)]), $(idx_var_d)))
end
eval(meta_elementary_functions)
# l0 loc
# l0 => l1
sym_func_us_l0l1_1 = Symbol("us_eucl_dist_$(basename_func)_l0l1_1!")
str_us_l0l1_1 = "
@everywhere $(sym_func_us_l0l1_1)(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) = \n
(S.loc = :l1; \n
S[:n] = x[$(idx_obs_var)];\n
S[:d] = 0.0;\n
S[:idx] = 1.0)"
eval(Meta.parse(str_us_l0l1_1))
edge1 = Edge([nothing], getfield(Main, :cc_eucl_dist_aut_l0l1_1), getfield(Main, sym_func_us_l0l1_1))
edge1 = Edge([nothing], getfield(Main, func_name(:cc, :l0, :l1, 1)), getfield(Main, func_name(:us, :l0, :l1, 1)))
map_edges[:l0][:l1] = [edge1]
# l1 loc
# l1 => l1
sym_func_us_l1l2_2 = Symbol("us_eucl_dist_$(basename_func)_l1l2_2!")
str_us_l1l2_2 = "
@everywhere $(sym_func_us_l1l2_2)(S::StateLHA, constants::Dict{Symbol,Float64}, x::Vector{Int}, p::Vector{Float64}) =\n
(S[:n] = x[$(idx_obs_var)])"
eval(Meta.parse(str_us_l1l2_2))
edge1 = Edge([nothing], getfield(Main, :cc_eucl_dist_aut_l1l1_1), getfield(Main, :us_eucl_dist_aut_l1l1_1!))
edge2 = Edge([:ALL], getfield(Main, :cc_eucl_dist_aut_l1l1_2), getfield(Main, sym_func_us_l1l2_2))
#=
edge1 = Edge([:ALL], getfield(Main, func_name(:cc, :l1, :l1, 1)), getfield(Main, func_name(:us, :l1, :l2, 1)))
map_edges[:l1][:l1] = [edge1]
for i = 1:nbr_observations
meta_edge_i = quote
@everywhere $(func_name(:cc, :l1, :l1, 1+i))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
S[:t] >= $(timeline[i])
@everywhere $(func_name(:us, :l1, :l1, 1+i))(S::StateLHA, x::Vector{Int}, p::Vector{Float64}) =
(setindex!(getfield(S, :values), S[:d] + (S[:n] - $(observations[i]))^2, $(idx_var_d));
setindex!(getfield(S, :values), S[:idx] + 1.0, $(idx_var_idx)))
end
eval(meta_edge_i)
push!(map_edges[:l1][:l1], Edge([nothing], getfield(Main, func_name(:cc, :l1, :l1, 1+i)), getfield(Main, func_name(:us, :l1, :l1, 1+i))))
end
=#
edge1 = Edge([nothing], getfield(Main, func_name(:cc, :l1, :l1, 1)), getfield(Main, func_name(:us, :l1, :l1, 1)))
edge2 = Edge([:ALL], getfield(Main, func_name(:cc, :l1, :l1, 2)), getfield(Main, func_name(:us, :l1, :l1, 2)))
map_edges[:l1][:l1] = [edge1, edge2]
# l1 => l2
edge1 = Edge([nothing], getfield(Main, :cc_eucl_dist_aut_l1l2_1), getfield(Main, :us_eucl_dist_aut_l1l2_1!))
edge1 = Edge([nothing], getfield(Main, func_name(:cc, :l1, :l2, 1)), getfield(Main, func_name(:us, :l1, :l2, 1)))
map_edges[:l1][:l2] = [edge1]
## Constants
......
automata/period_automaton.jl
+ 169
− 162
View file @ 886e9d98
This diff is collapsed.
core/MarkovProcesses.jl
+ 1
− 1
View file @ 886e9d98
......@@ -17,7 +17,7 @@ export VariableModel, ParameterModel, Transition
export LHA, StateLHA, Edge, Location, VariableAutomaton
# Trajectory related methods
export +, -, δ, dist_lp
export +, -, δ, dist_lp, euclidean_distance
export get_obs_var, length_states, length_obs_var
export get_state_from_time, get_var_from_time, vectorize
export isbounded, times, transitions
......
core/lha.jl
+ 17
− 14
View file @ 886e9d98
......@@ -104,19 +104,19 @@ end
function _find_edge_candidates!(edge_candidates::Vector{Edge},
edges_from_current_loc::Dict{Location,Vector{Edge}},
Λ::Dict{Location,Function},
Snplus1::StateLHA, constants::Dict{Symbol,Float64},
Snplus1::StateLHA,
x::Vector{Int}, p::Vector{Float64},
only_asynchronous::Bool)
nbr_candidates = 0
for target_loc in keys(edges_from_current_loc)
for edge in edges_from_current_loc[target_loc]
if Λ[target_loc](x) && getfield(edge, :check_constraints)(Snplus1, constants, x, p)
if edge.transitions[1] == nothing
if Λ[target_loc](x) && getfield(edge, :check_constraints)(Snplus1, x, p)
if getfield(edge, :transitions)[1] == nothing
_push_edge!(edge_candidates, edge, nbr_candidates)
nbr_candidates += 1
return nbr_candidates
end
if !only_asynchronous && edge.transitions[1] != nothing
if !only_asynchronous && getfield(edge, :transitions)[1] != nothing
_push_edge!(edge_candidates, edge, nbr_candidates)
nbr_candidates += 1
end
......@@ -133,13 +133,13 @@ function _get_edge_index(edge_candidates::Vector{Edge}, nbr_candidates::Int,
for i = 1:nbr_candidates
edge = edge_candidates[i]
# Asynchronous edge detection: we fire it
if edge.transitions[1] == nothing
if getfield(edge, :transitions)[1] == nothing
return (i, bool_event)
end
# Synchronous detection
if !detected_event && tr_nplus1 != nothing
if (edge.transitions[1] == :ALL) ||
(tr_nplus1 in edge.transitions)
if (getfield(edge, :transitions)[1] == :ALL) ||
(tr_nplus1 in getfield(edge, :transitions))
ind_edge = i
bool_event = true
end
......@@ -163,7 +163,6 @@ function next_state!(Snplus1::StateLHA, A::LHA,
@show Sn
@show Snplus1
end
constants = getfield(A, :constants)
# In terms of values not reference, Snplus1 == Sn
# First, we check the asynchronous transitions
while first_round || length(edge_candidates) > 0
......@@ -172,14 +171,14 @@ function next_state!(Snplus1::StateLHA, A::LHA,
current_loc = getfield(Snplus1, :loc)
edges_from_current_loc = getfield(A, :map_edges)[current_loc]
# Save all edges that satisfies transition predicate (asynchronous ones)
nbr_candidates = _find_edge_candidates!(edge_candidates, edges_from_current_loc, getfield(A, :Λ), Snplus1, constants, xn, p, true)
nbr_candidates = _find_edge_candidates!(edge_candidates, edges_from_current_loc, getfield(A, :Λ), Sn, xn, p, true) #Snplus1=>Sn
# Search the one we must chose, here the event is nothing because
# we're not processing yet the next event
ind_edge, detected_event = _get_edge_index(edge_candidates, nbr_candidates, detected_event, nothing)
# Update the state with the chosen one (if it exists)
# Should be xn here
if ind_edge > 0
getfield(edge_candidates[ind_edge], :update_state!)(Snplus1, constants, xn, p)
getfield(edge_candidates[ind_edge], :update_state!)(Snplus1, xn, p)
end
first_round = false
if verbose
......@@ -193,6 +192,7 @@ function next_state!(Snplus1::StateLHA, A::LHA,
break
end
# For debug
#=
if turns > 100
println("Number of turns in next_state! is suspicious")
@show first_round, detected_event
......@@ -200,10 +200,11 @@ function next_state!(Snplus1::StateLHA, A::LHA,
@show tnplus1, tr_nplus1, xnplus1
@show edge_candidates
for edge in edge_candidates
@show getfield(edge, :check_constraints)(Snplus1, constants, x, p)
@show getfield(edge, :check_constraints)(Snplus1, xn, p)
end
error("Unpredicted behavior automaton")
end
=#
end
if verbose
@show Snplus1
......@@ -227,7 +228,7 @@ function next_state!(Snplus1::StateLHA, A::LHA,
current_loc = getfield(Snplus1, :loc)
edges_from_current_loc = getfield(A, :map_edges)[current_loc]
# Save all edges that satisfies transition predicate (synchronous ones)
nbr_candidates = _find_edge_candidates!(edge_candidates, edges_from_current_loc, getfield(A, :Λ), Snplus1, constants, xnplus1, p, false)
nbr_candidates = _find_edge_candidates!(edge_candidates, edges_from_current_loc, getfield(A, :Λ), Snplus1, xnplus1, p, false)
# Search the one we must chose
ind_edge, detected_event = _get_edge_index(edge_candidates, nbr_candidates, detected_event, tr_nplus1)
# Update the state with the chosen one (if it exists)
......@@ -237,7 +238,7 @@ function next_state!(Snplus1::StateLHA, A::LHA,
@show ind_edge, detected_event, nbr_candidates
end
if ind_edge > 0
getfield(edge_candidates[ind_edge], :update_state!)(Snplus1, constants, xnplus1, p)
getfield(edge_candidates[ind_edge], :update_state!)(Snplus1, xnplus1, p)
end
first_round = false
if verbose
......@@ -249,6 +250,7 @@ function next_state!(Snplus1::StateLHA, A::LHA,
break
end
# For debug
#=
if turns > 100
println("Number of turns in next_state! is suspicious")
@show detected_event
......@@ -256,10 +258,11 @@ function next_state!(Snplus1::StateLHA, A::LHA,
@show tnplus1, tr_nplus1, xnplus1
@show edge_candidates
for edge in edge_candidates
@show getfield(edge, :check_constraints)(Snplus1, constants, x, p)
@show getfield(edge, :check_constraints)(Snplus1, x, p)
end
error("Unpredicted behavior automaton")
end
=#
end
if verbose println("##### End next_state!") end
end
......
core/model.jl
+ 0
− 1
View file @ 886e9d98
......@@ -295,7 +295,6 @@ function volatile_simulate(product::SynchronizedModel;
getfield(m, :f!)(vec_x, l_t, l_tr, xn, tn, p_sim)
if l_t[1] > time_bound
tn = l_t[1]
i -= 1
break
end
if vec_x == xn
......
core/trajectory.jl
+ 11
− 1
View file @ 886e9d98
## About distances
# Top-level Lp distance function
"""
`dist_lp(l_σ1, l_σ2; verbose, p, str_stat_list, str_stat_trajectory)`
......@@ -68,6 +70,7 @@ function dist_lp(x_obs::Vector{Int}, t_x::Vector{Float64}, y_obs::Vector{Int}, t
verbose::Bool = false, p::Int = 1)
current_y_obs = y_obs[1]
current_t_y = t_y[2]
nbr_y_obs = length(y_obs)
idx = 1
res = 0.0
for i = 1:(length(x_obs)-1)
......@@ -78,7 +81,7 @@ function dist_lp(x_obs::Vector{Int}, t_x::Vector{Float64}, y_obs::Vector{Int}, t
@show t_x[i+1]
end
last_t_y = t_x[i]
while current_t_y < t_x[i+1]
while current_t_y < t_x[i+1] && idx <= nbr_y_obs
rect = abs(current_y_obs - x_obs[i])^p * (current_t_y - last_t_y)
res += rect
if verbose
......@@ -172,6 +175,13 @@ function vectorize(σ::AbstractTrajectory, sym_var::Symbol,
return trajectory_points
end
function euclidean_distance(σ::AbstractTrajectory, sym_obs::Symbol,
tml::AbstractVector{Float64}, observations::AbstractVector{Float64})
traj_obs = vectorize(σ, sym_obs, tml)
diff_obs = observations - traj_obs
return sqrt(dot(diff_obs, diff_obs))
end
function check_consistency(σ::AbstractTrajectory)
test_length_var = true
for i = 1:σ.m.dim_obs_state
......
tests/automata/euclidean_distance.jl 0 → 100644
+ 42
− 0
View file @ 886e9d98
using MarkovProcesses
import LinearAlgebra: dot
import Distributions: Uniform
load_automaton("euclidean_distance_automaton")
load_model("SIR")
load_model("ER")
test_all = true
# SIR model
nbr_sim = 10
for i = 1:nbr_sim
set_param!(SIR, [:ki, :kr], [rand(Uniform(5E-5, 3E-3)), rand(Uniform(5E-3, 0.2))])
let tml_obs, y_obs, sync_SIR, σ, test
tml_obs = rand(Uniform(0.0, 5.0)):1.0:rand(Uniform(50.0, 100.0))
y_obs = vectorize(simulate(SIR), :I, tml_obs)
sync_SIR = SIR * create_euclidean_distance_automaton(SIR, tml_obs, y_obs, :I)
σ = simulate(sync_SIR)
test = euclidean_distance(σ, :I, tml_obs, y_obs) == σ.state_lha_end[:d]
#@show test, euclidean_distance(σ, tml_obs, y_obs, :I), σ.state_lha_end[:d]
global test_all = test_all && test
end
end
# ER model
for i = 1:nbr_sim
let tml_obs, y_obs, sync_SIR, σ, test
set_param!(ER, :k3, rand(Uniform(0.0, 100.0)))
tml_obs = rand(Uniform(0.0, 0.2)):1.0:rand(Uniform(0.5,10.0))
y_obs = vectorize(simulate(ER), :P, tml_obs)
sync_ER = ER * create_euclidean_distance_automaton(ER, tml_obs, y_obs, :P)
σ = simulate(sync_ER)
test = euclidean_distance(σ, :P, tml_obs, y_obs) == σ.state_lha_end[:d]
#@show test, euclidean_distance(σ, tml_obs, y_obs, :P), σ.state_lha_end[:d]
global test_all = test_all && test
end
end
return test_all
tests/automata/euclidean_distance_single.jl 0 → 100644
+ 17
− 0
View file @ 886e9d98
using MarkovProcesses
import LinearAlgebra: dot
import Distributions: Uniform
load_automaton("euclidean_distance_automaton")
load_model("SIR")
tml_obs = 0:0.5:200
set_time_bound!(SIR, 200.0)
y_obs = vectorize(simulate(SIR), :I, tml_obs)
sync_SIR = SIR * create_euclidean_distance_automaton(SIR, tml_obs, y_obs, :I)
σ = simulate(sync_SIR)
test = euclidean_distance(σ, tml_obs, y_obs, :I) == σ.state_lha_end[:d]
@show test, euclidean_distance(σ, tml_obs, y_obs, :I), σ.state_lha_end[:d]
return test
tests/run_automata.jl
+ 2
− 0
View file @ 886e9d98
......@@ -4,6 +4,8 @@ using Test
@testset "Automata tests" begin
@test include("automata/absorbing_x0_p.jl")
@test include("automata/accept_R5.jl")
@test include("automata/euclidean_distance.jl")
@test include("automata/euclidean_distance_single.jl")
@test include("automata/read_trajectory_last_state_F.jl")
@test include("automata/read_trajectory_last_state_G.jl")
@test include("automata/sync_simulate_last_state_F.jl")
......
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