diff --git a/automata/automaton_F.jl b/automata/automaton_F.jl
index b464a7a98c73d78728867ac7d22cb81c1e48d794..19183884f5aea1133895a1b16f6c88c2057518d4 100644
--- a/automata/automaton_F.jl
+++ b/automata/automaton_F.jl
@@ -35,43 +35,62 @@ function create_automaton_F(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1
tuple = ("l0", "l1")
cc_aut_F_l0l1_1(A::LHA, S::StateLHA) = true
us_aut_F_l0l1_1!(A::LHA, S::StateLHA, x::Vector{Int}) =
- (S.loc = "l1"; S["d"] = Inf; S["n"] = get_value(A, x, str_obs))
+ (S.loc = "l1";
+ S["d"] = Inf;
+ S["n"] = get_value(A, x, str_obs);
+ S["isabs"] = m.isabsorbing(m.p,x))
edge1 = Edge([nothing], cc_aut_F_l0l1_1, us_aut_F_l0l1_1!)
map_edges[tuple] = [edge1]
# l1 loc : we construct the edges of the form l1 => (..)
tuple = ("l1", "l2")
cc_aut_F_l1l2_1(A::LHA, S::StateLHA) =
- (A.l_ctes["x1"] <= S["n"] <= A.l_ctes["x2"]) && (A.l_ctes["t1"] <= S.time <= A.l_ctes["t2"])
- us_aut_F_l1l2_1!(A::LHA, S::StateLHA, x::Vector{Int}) = (S["d"] = 0; S.loc = "l2")
+ (A.l_ctes["x1"] <= S["n"] <= A.l_ctes["x2"]) &&
+ (A.l_ctes["t1"] <= S.time <= A.l_ctes["t2"])
+ us_aut_F_l1l2_1!(A::LHA, S::StateLHA, x::Vector{Int}) =
+ (S.loc = "l2";
+ S["d"] = 0)
edge1 = Edge([nothing], cc_aut_F_l1l2_1, us_aut_F_l1l2_1!)
- cc_aut_F_l1l2_2(A::LHA, S::StateLHA) = (S["d"] > 0 && (S.time > A.l_ctes["t2"] ||Â istrue(S["isabs"])))
- us_aut_F_l1l2_2!(A::LHA, S::StateLHA, x::Vector{Int}) = (S.loc = "l2")
+ cc_aut_F_l1l2_2(A::LHA, S::StateLHA) =
+ S["d"] > 0 &&
+ (S.time > A.l_ctes["t2"] ||Â istrue(S["isabs"]))
+ us_aut_F_l1l2_2!(A::LHA, S::StateLHA, x::Vector{Int}) =
+ (S.loc = "l2")
edge2 = Edge([nothing], cc_aut_F_l1l2_2, us_aut_F_l1l2_2!)
- cc_aut_F_l1l2_3(A::LHA, S::StateLHA) = (S["d"] == 0 && S.time >= A.l_ctes["t1"])
- us_aut_F_l1l2_3!(A::LHA, S::StateLHA, x::Vector{Int}) = (S.loc = "l2")
+ cc_aut_F_l1l2_3(A::LHA, S::StateLHA) =
+ S["d"] == 0 &&
+ S.time >= A.l_ctes["t1"]
+ us_aut_F_l1l2_3!(A::LHA, S::StateLHA, x::Vector{Int}) =
+ (S.loc = "l2")
edge3 = Edge([nothing], cc_aut_F_l1l2_3, us_aut_F_l1l2_3!)
map_edges[tuple] = [edge1, edge2, edge3]
tuple = ("l1", "l3")
- cc_aut_F_l1l3_1(A::LHA, S::StateLHA) = (A.l_ctes["x1"] <= S["n"] <= A.l_ctes["x2"])
- us_aut_F_l1l3_1!(A::LHA, S::StateLHA, x::Vector{Int}) = (S["d"] = 0; S.loc = "l3")
+ cc_aut_F_l1l3_1(A::LHA, S::StateLHA) =
+ (A.l_ctes["x1"] <= S["n"] <= A.l_ctes["x2"])
+ us_aut_F_l1l3_1!(A::LHA, S::StateLHA, x::Vector{Int}) =
+ (S.loc = "l3";
+ S["d"] = 0;)
edge1 = Edge([nothing], cc_aut_F_l1l3_1, us_aut_F_l1l3_1!)
cc_aut_F_l1l3_2(A::LHA, S::StateLHA) =
- (A.l_ctes["x1"] > S["n"] || S["n"] > A.l_ctes["x2"]) && (S.time < A.l_ctes["t1"])
+ (S["n"] < A.l_ctes["x1"] || S["n"] > A.l_ctes["x2"]) &&
+ (S.time <= A.l_ctes["t1"])
us_aut_F_l1l3_2!(A::LHA, S::StateLHA, x::Vector{Int}) =
- (S["d"] = min(sqrt((S.time - A.l_ctes["t1"])^2 + (S["n"] - A.l_ctes["x2"])^2),
- sqrt((S.time - A.l_ctes["t1"])^2 + (S["n"] - A.l_ctes["x1"])^2)); S.loc = "l3")
+ (S.loc = "l3";
+ S["d"] = min(sqrt((S.time - A.l_ctes["t1"])^2 + (S["n"] - A.l_ctes["x2"])^2),
+ sqrt((S.time - A.l_ctes["t1"])^2 + (S["n"] - A.l_ctes["x1"])^2)))
edge2 = Edge([nothing], cc_aut_F_l1l3_2, us_aut_F_l1l3_2!)
cc_aut_F_l1l3_3(A::LHA, S::StateLHA) =
- (A.l_ctes["x1"] > S["n"] || S["n"] > A.l_ctes["x2"]) && (A.l_ctes["t1"] <= S.time <= A.l_ctes["t2"])
+ (S["n"] < A.l_ctes["x1"] || S["n"] > A.l_ctes["x2"]) &&
+ (A.l_ctes["t1"] <= S.time <= A.l_ctes["t2"])
us_aut_F_l1l3_3!(A::LHA, S::StateLHA, x::Vector{Int}) =
- (S["d"] = min(S["d"], min(abs(S["n"] - A.l_ctes["x1"]), abs(S["n"] - A.l_ctes["x2"]))); S.loc = "l3")
+ (S.loc = "l3";
+ S["d"] = min(S["d"], min(abs(S["n"] - A.l_ctes["x1"]), abs(S["n"] - A.l_ctes["x2"]))))
edge3 = Edge([nothing], cc_aut_F_l1l3_3, us_aut_F_l1l3_3!)
map_edges[tuple] = [edge1, edge2, edge3]
@@ -79,11 +98,15 @@ function create_automaton_F(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1
tuple = ("l3", "l1")
cc_aut_F_l3l1_1(A::LHA, S::StateLHA) = true
us_aut_F_l3l1_1!(A::LHA, S::StateLHA, x::Vector{Int}) =
- (S["n"] = get_value(A, x, str_obs); S.loc = "l1"; S["isabs"] = m.isabsorbing(m.p,x))
+ (S.loc = "l1";
+ S["n"] = get_value(A, x, str_obs);
+ S["isabs"] = m.isabsorbing(m.p,x))
edge1 = Edge(["ALL"], cc_aut_F_l3l1_1, us_aut_F_l3l1_1!)
tuple = ("l3", "l2")
- cc_aut_F_l3l2_1(A::LHA, S::StateLHA) = (S.time >= A.l_ctes["t2"])
- us_aut_F_l3l2_1!(A::LHA, S::StateLHA, x::Vector{Int}) = (S["n"] = get_value(A, x, str_obs); S.loc = "l2")
+ cc_aut_F_l3l2_1(A::LHA, S::StateLHA) =
+ (S.time >= A.l_ctes["t2"])
+ us_aut_F_l3l2_1!(A::LHA, S::StateLHA, x::Vector{Int}) =
+ (S.loc = "l2")
edge2 = Edge([nothing], cc_aut_F_l3l2_1, us_aut_F_l3l2_1!)
map_edges[tuple] = [edge1, edge2]
diff --git a/automata/automaton_G.jl b/automata/automaton_G.jl
index 9c982c94fd5d64a87f62d40af85f5b6eb0d16711..c2210988432fe6a47de48740dcc204f01ad97c42 100644
--- a/automata/automaton_G.jl
+++ b/automata/automaton_G.jl
@@ -15,7 +15,7 @@ function create_automaton_G(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1
l_loc_final = ["l2"]
## Map of automaton variables
- map_var_automaton_idx = Dict{VariableAutomaton,Int}("t'" => 1, "in" => 2,
+ map_var_automaton_idx = Dict{VariableAutomaton,Int}("tprime" => 1, "in" => 2,
"n" => 3, "d" => 4, "isabs" => 5)
## Flow of variables
@@ -34,15 +34,19 @@ function create_automaton_G(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1
tuple = ("l0", "l1")
cc_aut_G_l0l1_1(A::LHA, S::StateLHA) = true
us_aut_G_l0l1_1!(A::LHA, S::StateLHA, x::Vector{Int}) =
- (S.loc = "l1"; S["d"] = 0; S["n"] = get_value(A, x, str_obs); S["in"] = true; S["isabs"] = m.isabsorbing(m.p,x))
+ (S.loc = "l1";
+ S["d"] = 0;
+ S["n"] = get_value(A, x, str_obs);
+ S["in"] = true;
+ S["isabs"] = m.isabsorbing(m.p,x))
edge1 = Edge([nothing], cc_aut_G_l0l1_1, us_aut_G_l0l1_1!)
map_edges[tuple] = [edge1]
# l1 loc
tuple = ("l1", "l3")
cc_aut_G_l1l3_1(A::LHA, S::StateLHA) =
- (S.time < A.l_ctes["t1"] &&
- (S["n"] < A.l_ctes["x1"] || S["n"] > A.l_ctes["x2"]))
+ S.time <= A.l_ctes["t1"] &&
+ S["n"] < A.l_ctes["x1"] || S["n"] > A.l_ctes["x2"]
us_aut_G_l1l3_1!(A::LHA, S::StateLHA, x::Vector{Int}) =
(S.loc = "l3";
S["d"] = min(abs(A.l_ctes["x1"] - S["n"]), abs(A.l_ctes["x2"] - S["n"]));
@@ -56,7 +60,7 @@ function create_automaton_G(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1
us_aut_G_l1l3_3!(A::LHA, S::StateLHA, x::Vector{Int}) =
(S.loc = "l3";
S["d"] = S["d"] * (S.time - A.l_ctes["t1"]);
- S["t'"] = 0.0)
+ S["tprime"] = 0.0)
edge3 = Edge([nothing], cc_aut_G_l1l3_3, us_aut_G_l1l3_3!)
cc_aut_G_l1l3_2(A::LHA, S::StateLHA) =
@@ -74,7 +78,7 @@ function create_automaton_G(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1
(A.l_ctes["x1"] <= S["n"] <= A.l_ctes["x2"])
us_aut_G_l1l3_4!(A::LHA, S::StateLHA, x::Vector{Int}) =
(S.loc = "l3";
- S["t'"] = 0.0)
+ S["tprime"] = 0.0)
edge4 = Edge([nothing], cc_aut_G_l1l3_4, us_aut_G_l1l3_4!)
map_edges[tuple] = [edge1, edge2, edge3, edge4]
@@ -111,7 +115,23 @@ function create_automaton_G(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1
(S.loc = "l2";
S["d"] = S["d"] * (A.l_ctes["t2"] - A.l_ctes["t1"]))
edge2 = Edge([nothing], cc_aut_G_l1l2_2, us_aut_G_l1l2_2!)
- map_edges[tuple] = [edge1, edge2]
+ cc_aut_G_l1l2_3(A::LHA, S::StateLHA) =
+ istrue(S["isabs"]) &&
+ S.time <= A.l_ctes["t1"]
+ us_aut_G_l1l2_3!(A::LHA, S::StateLHA, x::Vector{Int}) =
+ (S.loc = "l2";
+ S["d"] = S["d"] * (A.l_ctes["t2"] - A.l_ctes["t1"]))
+ edge3 = Edge([nothing], cc_aut_G_l1l2_3, us_aut_G_l1l2_3!)
+ cc_aut_G_l1l2_4(A::LHA, S::StateLHA) =
+ istrue(S["isabs"]) &&
+ (A.l_ctes["t1"] <= S.time <= A.l_ctes["t2"])
+ us_aut_G_l1l2_4!(A::LHA, S::StateLHA, x::Vector{Int}) =
+ (S.loc = "l2";
+ S["d"] += (A.l_ctes["t2"] - S.time) *
+ min(abs(A.l_ctes["x1"] - S["n"]), abs(A.l_ctes["x1"] - S["n"])))
+ edge4 = Edge([nothing], cc_aut_G_l1l2_4, us_aut_G_l1l2_4!)
+
+ map_edges[tuple] = [edge1, edge2, edge3, edge4]
# l3 loc
tuple = ("l3", "l1")
@@ -124,44 +144,29 @@ function create_automaton_G(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1
map_edges[tuple] = [edge1]
tuple = ("l3", "l2")
- cc_aut_G_l3l2_1(A::LHA, S::StateLHA) =
- istrue(S["in"]) &&
- (S.time >= A.l_ctes["t2"] || istrue(S["isabs"]))
- us_aut_G_l3l2_1!(A::LHA, S::StateLHA, x::Vector{Int}) =
- (S.loc = "l2")
- edge1 = Edge([nothing], cc_aut_G_l3l2_1, us_aut_G_l3l2_1!)
cc_aut_G_l3l2_2(A::LHA, S::StateLHA) =
- !istrue(S["in"]) &&
- (S.time >= A.l_ctes["t2"] || istrue(S["isabs"]))
+ istrue(S["in"]) &&
+ S.time >= A.l_ctes["t2"]
us_aut_G_l3l2_2!(A::LHA, S::StateLHA, x::Vector{Int}) =
(S.loc = "l2";
S["d"] = S["d"] * (A.l_ctes["t2"] - A.l_ctes["t1"]))
edge2 = Edge([nothing], cc_aut_G_l3l2_2, us_aut_G_l3l2_2!)
- cc_aut_G_l3l2_3(A::LHA, S::StateLHA) =
- istrue(S["isabs"]) &&
- S.time <= A.l_ctes["t1"]
- us_aut_G_l3l2_3!(A::LHA, S::StateLHA, x::Vector{Int}) =
- (S.loc = "l2";
- S["d"] = S["d"] * (A.l_ctes["t2"] - A.l_ctes["t1"]))
- edge3 = Edge([nothing], cc_aut_G_l3l2_3, us_aut_G_l3l2_3!)
- cc_aut_G_l3l2_4(A::LHA, S::StateLHA) =
- istrue(S["isabs"]) &&
- S.time >= A.l_ctes["t1"]
- us_aut_G_l3l2_4!(A::LHA, S::StateLHA, x::Vector{Int}) =
- (S.loc = "l2";
- S["d"] += (A.l_ctes["t2"] - A.l_ctes["t1"]) *
- min(abs(S["n"] - A.l_ctes["x1"]), abs(S["n"] - A.l_ctes["x2"])))
- edge4 = Edge([nothing], cc_aut_G_l3l2_4, us_aut_G_l3l2_4!)
+ cc_aut_G_l3l2_1(A::LHA, S::StateLHA) =
+ !istrue(S["in"]) &&
+ S.time >= A.l_ctes["t2"]
+ us_aut_G_l3l2_1!(A::LHA, S::StateLHA, x::Vector{Int}) =
+ (S.loc = "l2")
+ edge1 = Edge([nothing], cc_aut_G_l3l2_1, us_aut_G_l3l2_1!)
- map_edges[tuple] = [edge1, edge2, edge3, edge4]
+ map_edges[tuple] = [edge1, edge2]
# l4 loc
tuple = ("l4", "l1")
cc_aut_G_l4l1_1(A::LHA, S::StateLHA) = true
us_aut_G_l4l1_1!(A::LHA, S::StateLHA, x::Vector{Int}) =
(S.loc = "l1";
- S["d"] += S["t'"] * min(abs(A.l_ctes["x1"] - S["n"]), abs(A.l_ctes["x2"] - S["n"]));
- S["t'"] = 0.0;
+ S["d"] += S["tprime"] * min(abs(A.l_ctes["x1"] - S["n"]), abs(A.l_ctes["x2"] - S["n"]));
+ S["tprime"] = 0.0;
S["n"] = get_value(A, x, str_obs);
S["in"] = true;
S["isabs"] = m.isabsorbing(m.p,x))
@@ -170,29 +175,14 @@ function create_automaton_G(m::ContinuousTimeModel, x1::Float64, x2::Float64, t1
tuple = ("l4", "l2")
cc_aut_G_l4l2_1(A::LHA, S::StateLHA) =
- (S.time >= A.l_ctes["t2"] || istrue(S["isabs"]))
+ (S.time >= A.l_ctes["t2"])
us_aut_G_l4l2_1!(A::LHA, S::StateLHA, x::Vector{Int}) =
(S.loc = "l2";
- S["d"] += S["t'"] * min(abs(A.l_ctes["x1"] - S["n"]), abs(A.l_ctes["x2"] - S["n"]));
- S["t'"] = 0.0)
+ S["d"] += S["tprime"] * min(abs(A.l_ctes["x1"] - S["n"]), abs(A.l_ctes["x2"] - S["n"]));
+ S["tprime"] = 0.0)
edge1 = Edge([nothing], cc_aut_G_l4l2_1, us_aut_G_l4l2_1!)
- cc_aut_G_l4l2_2(A::LHA, S::StateLHA) =
- istrue(S["isabs"]) &&
- S.time <= A.l_ctes["t1"]
- us_aut_G_l4l2_2!(A::LHA, S::StateLHA, x::Vector{Int}) =
- (S.loc = "l2";
- S["d"] = S["d"] * (A.l_ctes["t2"] - A.l_ctes["t1"]))
- edge2 = Edge([nothing], cc_aut_G_l4l2_2, us_aut_G_l4l2_2!)
- cc_aut_G_l4l2_3(A::LHA, S::StateLHA) =
- istrue(S["isabs"]) &&
- S.time >= A.l_ctes["t1"]
- us_aut_G_l4l2_3!(A::LHA, S::StateLHA, x::Vector{Int}) =
- (S.loc = "l2";
- S["d"] += (A.l_ctes["t2"] - A.l_ctes["t1"]) *
- min(abs(S["n"] - A.l_ctes["x1"]), abs(S["n"] - A.l_ctes["x2"])))
- edge3 = Edge([nothing], cc_aut_G_l4l2_3, us_aut_G_l4l2_3!)
- map_edges[tuple] = [edge1, edge2, edge3]
+ map_edges[tuple] = [edge1]
## Constants
l_ctes = Dict{String,Float64}("x1" => x1, "x2" => x2, "t1" => t1, "t2" => t2)
diff --git a/core/lha.jl b/core/lha.jl
index fcc3fce8baa2965b832c3608590070a4e96fc9d0..af93550c86add7222ce30f06ee556bf24651bfd7 100644
--- a/core/lha.jl
+++ b/core/lha.jl
@@ -43,34 +43,39 @@ function init_state(A::LHA, x0::Vector{Int}, t0::Float64)
end
function _find_edge_candidates!(edge_candidates::Vector{Edge}, current_loc::Location,
- A::LHA, Snplus1::StateLHA)
-
+ A::LHA, Snplus1::StateLHA, only_asynchronous::Bool)
for loc in A.l_loc
tuple_edges = (current_loc, loc)
if haskey(A.map_edges, tuple_edges)
for edge in A.map_edges[tuple_edges]
if edge.check_constraints(A, Snplus1)
- push!(edge_candidates, edge)
+ if edge.transitions[1] == nothing
+ pushfirst!(edge_candidates, edge)
+ end
+ if !only_asynchronous && edge.transitions[1] != nothing
+ push!(edge_candidates, edge)
+ end
end
end
end
end
end
-function _get_edge_index(edge_candidates::Vector{Edge}, first_round::Bool,
+function _get_edge_index(edge_candidates::Vector{Edge},
detected_event::Bool, tr_nplus1::Transition)
ind_edge = 0
bool_event = detected_event
for i in eachindex(edge_candidates)
edge = edge_candidates[i]
+ # Asynchronous detection
if edge.transitions[1] == nothing
return (i, bool_event)
end
- if first_round || !detected_event
- if (length(edge.transitions) == 1 && tr_nplus1 != nothing && edge.transitions[1] == "ALL") ||
- (tr_nplus1 in edge.transitions)
- ind_edge = i
- bool_event = true
+ # Synchronous detection
+ if !detected_event && tr_nplus1 != nothing
+ if (length(edge.transitions) == 1 && edge.transitions[1] == "ALL") ||
+ (tr_nplus1 in edge.transitions)
+ return (i, true)
end
end
end
@@ -80,41 +85,45 @@ end
function next_state!(Snplus1::StateLHA, A::LHA,
xnplus1::Vector{Int}, tnplus1::Float64, tr_nplus1::Transition,
Sn::StateLHA; verbose::Bool = false)
- for i in eachindex(Snplus1.l_var)
- Snplus1.l_var[i] += (A.l_flow[Sn.loc])[i]*(tnplus1 - Sn.time)
- end
- Snplus1.time = tnplus1
-
# En fait d'apres observation de Cosmos, après qu'on ait lu la transition on devrait stop.
edge_candidates = Edge[]
first_round::Bool = true
detected_event::Bool = (tr_nplus1 == nothing) ? true : false
- turns = 1
- while first_round || !detected_event || length(edge_candidates) > 0
+ turns = 0
+
+ if verbose
+ println("#####")
+ @show xnplus1, tnplus1, tr_nplus1
+ @show Sn
+ @show Snplus1
+ end
+
+ # In terms of values not reference, Snplus1 == Sn
+ # First, we check the asynchronous transitions
+ while first_round || length(edge_candidates) > 0
+ turns += 1
edge_candidates = empty!(edge_candidates)
current_loc = Snplus1.loc
- if verbose @show turns end
- # Save all edges that satisfies transition predicate (synchronous or autonomous)
- _find_edge_candidates!(edge_candidates, current_loc, A, Snplus1)
- # Search the one we must chose
- ind_edge, detected_event = _get_edge_index(edge_candidates, first_round, detected_event, tr_nplus1)
+ # Save all edges that satisfies transition predicate (asynchronous ones)
+ _find_edge_candidates!(edge_candidates, current_loc, A, Snplus1, true)
+ # 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, detected_event, nothing)
# Update the state with the chosen one (if it exists)
if ind_edge > 0
edge_candidates[ind_edge].update_state!(A, Snplus1, xnplus1)
- # Should add something like if edges_candidates[ind_edge].transition != nohting break end ??
- end
- if (ind_edge == 0)
- #if (ind_edge == 0 ||Â detected_event)
- break
end
+ first_round = false
if verbose
- @show first_round detected_event
- @show tnplus1 tr_nplus1 xnplus1
- @show ind_edge
+ @show turns
@show edge_candidates
- @show tuple_candidates
+ @show ind_edge, detected_event
+ println("After update")
@show Snplus1
end
+ if (ind_edge == 0)
+ break
+ end
# For debug
if turns > 100
println("Number of turns in next_state! is suspicious")
@@ -127,8 +136,61 @@ function next_state!(Snplus1::StateLHA, A::LHA,
end
error("Unpredicted behavior automaton")
end
+ end
+ if verbose
+ @show Snplus1
+ println("Time flies with the flow...")
+ end
+ # Now time flies according to the flow
+ for i in eachindex(Snplus1.l_var)
+ coeff_deriv = (A.l_flow[Snplus1.loc])[i]
+ if coeff_deriv > 0
+ Snplus1.l_var[i] += coeff_deriv*(tnplus1 - Snplus1.time)
+ end
+ end
+ Snplus1.time = tnplus1
+ if verbose
+ @show Snplus1
+ end
+ # Now firing an edge according to the event
+ first_round = true
+ while first_round || length(edge_candidates) > 0
turns += 1
+ edge_candidates = empty!(edge_candidates)
+ current_loc = Snplus1.loc
+ # Save all edges that satisfies transition predicate (synchronous ones)
+ _find_edge_candidates!(edge_candidates, current_loc, A, Snplus1, false)
+ # Search the one we must chose
+ ind_edge, detected_event = _get_edge_index(edge_candidates, detected_event, tr_nplus1)
+ # Update the state with the chosen one (if it exists)
+ if verbose
+ @show turns
+ @show edge_candidates
+ @show ind_edge, detected_event
+ end
+ if ind_edge > 0
+ edge_candidates[ind_edge].update_state!(A, Snplus1, xnplus1)
+ end
first_round = false
+ if verbose
+ println("After update")
+ @show Snplus1
+ end
+ if (ind_edge == 0 ||Â detected_event)
+ break
+ end
+ # For debug
+ if turns > 100
+ println("Number of turns in next_state! is suspicious")
+ @show detected_event
+ @show length(edge_candidates)
+ @show tnplus1, tr_nplus1, xnplus1
+ @show edge_candidates
+ for edge in edge_candidates
+ @show edge.check_constraints(A, Snplus1)
+ end
+ error("Unpredicted behavior automaton")
+ end
end
end
diff --git a/core/model.jl b/core/model.jl
index da8faf3abe5b3aa3abe263d7eb6382cd8ad186d0..c8404b717dac6e1a3e9235d83d6750633373d9b5 100644
--- a/core/model.jl
+++ b/core/model.jl
@@ -125,7 +125,8 @@ function simulate(m::ContinuousTimeModel; p::Union{Nothing,AbstractVector{Float6
return Trajectory(m, values, times, transitions)
end
-function simulate(product::SynchronizedModel; p::Union{Nothing,AbstractVector{Float64}} = nothing)
+function simulate(product::SynchronizedModel; p::Union{Nothing,AbstractVector{Float64}} = nothing,
+ verbose::Bool = false)
m = product.m
A = product.automaton
p_sim = m.p
@@ -149,6 +150,11 @@ function simulate(product::SynchronizedModel; p::Union{Nothing,AbstractVector{Fl
Sn = copy(S0)
isabsorbing::Bool = m.isabsorbing(p_sim,xn)
isacceptedLHA::Bool = isaccepted(Sn)
+ # Alloc of vectors where we stock n+1 values
+ vec_x = zeros(Int, m.d)
+ l_t = Float64[0.0]
+ l_tr = Transition[nothing]
+ Snplus1 = copy(Sn)
# If x0 is absorbing
if isabsorbing || isacceptedLHA
_resize_trajectory!(full_values, times, transitions, 1)
@@ -156,13 +162,12 @@ function simulate(product::SynchronizedModel; p::Union{Nothing,AbstractVector{Fl
if isbounded(m)
_finish_bounded_trajectory!(values, times, transitions, m.time_bound)
end
+ if isabsorbing && !isacceptedLHA
+ next_state!(Snplus1, A, xn, m.time_bound, nothing, Sn; verbose = verbose)
+ copyto!(Sn, Snplus1)
+ end
return SynchronizedTrajectory(Sn, product, values, times, transitions)
end
- # Alloc of vectors where we stock n+1 values
- vec_x = zeros(Int, m.d)
- l_t = Float64[0.0]
- l_tr = Transition[nothing]
- Snplus1 = copy(Sn)
# First we fill the allocated array
for i = 2:m.estim_min_states
m.f!(vec_x, l_t, l_tr, xn, tn, p_sim)
@@ -173,7 +178,7 @@ function simulate(product::SynchronizedModel; p::Union{Nothing,AbstractVector{Fl
n += 1
copyto!(xn, vec_x)
tr_n = l_tr[1]
- next_state!(Snplus1, A, xn, tn, tr_n, Sn)
+ next_state!(Snplus1, A, xn, tn, tr_n, Sn; verbose = verbose)
_update_values!(full_values, times, transitions, xn, tn, tr_n, i)
copyto!(Sn, Snplus1)
isabsorbing = m.isabsorbing(p_sim,xn)
@@ -189,6 +194,10 @@ function simulate(product::SynchronizedModel; p::Union{Nothing,AbstractVector{Fl
if isbounded(m)
_finish_bounded_trajectory!(values, times, transitions, m.time_bound)
end
+ if isabsorbing && !isacceptedLHA
+ next_state!(Snplus1, A, xn, m.time_bound, nothing, Sn; verbose = verbose)
+ copyto!(Sn, Snplus1)
+ end
return SynchronizedTrajectory(Sn, product, values, times, transitions)
end
# Otherwise, buffering system
@@ -207,7 +216,7 @@ function simulate(product::SynchronizedModel; p::Union{Nothing,AbstractVector{Fl
end
copyto!(xn, vec_x)
tr_n = l_tr[1]
- next_state!(Snplus1, A, xn, tn, tr_n, Sn)
+ next_state!(Snplus1, A, xn, tn, tr_n, Sn; verbose = verbose)
_update_values!(full_values, times, transitions,
xn, tn, tr_n, m.estim_min_states+size_tmp+i)
copyto!(Sn, Snplus1)
@@ -230,6 +239,10 @@ function simulate(product::SynchronizedModel; p::Union{Nothing,AbstractVector{Fl
# the trajectory is bounded
_finish_bounded_trajectory!(values, times, transitions, m.time_bound)
end
+ if isabsorbing && !isacceptedLHA
+ next_state!(Snplus1, A, xn, m.time_bound, nothing, Sn; verbose = verbose)
+ copyto!(Sn, Snplus1)
+ end
return SynchronizedTrajectory(Sn, product, values, times, transitions)
end
@@ -249,20 +262,20 @@ function volatile_simulate(product::SynchronizedModel;
Sn = copy(S0)
isabsorbing::Bool = m.isabsorbing(p_sim,xn)
isacceptedLHA::Bool = isaccepted(Sn)
- # If x0 is absorbing
- if isabsorbing || isacceptedLHA
- return Sn
- end
# Alloc of vectors where we stock n+1 values
vec_x = zeros(Int, m.d)
l_t = Float64[0.0]
l_tr = Transition[nothing]
Snplus1 = copy(Sn)
- while !isabsorbing &&Â tn <= m.time_bound && !isacceptedLHA
- if verbose
- @show vec_x
- @show Sn
+ # If x0 is absorbing
+ if isabsorbing || isacceptedLHA
+ if !isacceptedLHA
+ next_state!(Snplus1, A, xn, m.time_bound, nothing, Sn; verbose = verbose)
+ copyto!(Sn, Snplus1)
end
+ return Sn
+ end
+ while !isabsorbing &&Â tn <= m.time_bound && !isacceptedLHA
m.f!(vec_x, l_t, l_tr, xn, tn, p_sim)
tn = l_t[1]
if tn > m.time_bound
@@ -271,12 +284,16 @@ function volatile_simulate(product::SynchronizedModel;
end
copyto!(xn, vec_x)
tr_n = l_tr[1]
- next_state!(Snplus1, A, xn, tn, tr_n, Sn)
+ next_state!(Snplus1, A, xn, tn, tr_n, Sn; verbose = verbose)
copyto!(Sn, Snplus1)
isabsorbing = m.isabsorbing(p_sim,xn)
isacceptedLHA = isaccepted(Snplus1)
n += 1
end
+ if isabsorbing && !isacceptedLHA
+ next_state!(Snplus1, A, xn, m.time_bound, nothing, Sn; verbose = verbose)
+ copyto!(Sn, Snplus1)
+ end
return Sn
end
@@ -313,7 +330,7 @@ function mean_value_lha(sm::SynchronizedModel, str_var::String, nbr_sim::Int)
return sum_val / nbr_sim
end
-function distribute_prob_accept_lha(sm::SynchronizedModel, str_var::String, nbr_sim::Int)
+function distribute_prob_accept_lha(sm::SynchronizedModel, nbr_sim::Int)
sum_val = @distributed (+) for i = 1:nbr_sim Int(isaccepted(volatile_simulate(sm))) end
return sum_val / nbr_sim
end
diff --git a/core/plots.jl b/core/plots.jl
index 9e2ac46fefce11b5bdcc84ab885cde700c5572c2..898267c7cbf167d0fe3760330dc466ee473192c5 100644
--- a/core/plots.jl
+++ b/core/plots.jl
@@ -52,8 +52,7 @@ end
function plot!(A::LHA)
x1, x2, t1, t2 = A.l_ctes["x1"], A.l_ctes["x2"], A.l_ctes["t1"], A.l_ctes["t2"]
- rectangle(w, h, x, y) = Shape(x .+ [0,w,w,0], y .+ [0,0,h,h])
- plot!(rectangle(t2, x2, t1, x1), opacity = 0.5)
+ plot!(Shape([(t1,x1), (t1,x2), (t2,x2), (t2,x1), (t1,x1)]), opacity = 0.5)
end
export plot, plot!
diff --git a/tests/cosmos/distance_F/ER_1D.jl b/tests/cosmos/distance_F/ER_1D.jl
index 80c87964c637c109c3b1a239609126c7f1e376de..3b9176fd7b85e3a69f40830e6fcdc2f5b30fbf29 100644
--- a/tests/cosmos/distance_F/ER_1D.jl
+++ b/tests/cosmos/distance_F/ER_1D.jl
@@ -42,6 +42,7 @@ for exp in l_exp
--const k_3=$(k3),x1=$x1,x2=$x2,t1=$t1,t2=$t2
--level $(level) --width $(width)
--verbose 0`
+ #run(command)
run(pipeline(command, stderr=devnull))
dict_values = cosmos_get_values("Result_dist_F_$(str_model).res")
l_dist_cosmos[i] = dict_values["Estimated value"]
@@ -52,8 +53,11 @@ for exp in l_exp
set_param!(ER, "k3", convert(Float64, k3))
sync_ER = ER*A_F
l_dist_pkg[i] = distribute_mean_value_lha(sync_ER, "d", nb_sim)
+ nb_accepts_pkg = distribute_prob_accept_lha(sync_ER, nb_sim)
+ #@info "About accepts" nb_sim nb_accepted nb_accepts_pkg
test = isapprox(l_dist_cosmos[i], l_dist_pkg[i]; atol = width*1.01)
- global test_all = test_all && test
+ test2 = nb_accepts_pkg == (nb_sim / nb_accepted)
+ global test_all = test_all && test && test2
if !test
@show l_dist_pkg[i], l_dist_cosmos[i]
@show exp
diff --git a/tests/cosmos/distance_F/dist_F_ER.lha b/tests/cosmos/distance_F/dist_F_ER.lha
index baa9f54eb7d36ed1364653a087e4bddb77a3f050..47cd3fb1d983dfb4b06bcaa8d70682ae7297dc88 100644
--- a/tests/cosmos/distance_F/dist_F_ER.lha
+++ b/tests/cosmos/distance_F/dist_F_ER.lha
@@ -1,4 +1,4 @@
-NbVariables = 3;
+NbVariables = 4;
NbLocations = 4;
const int x1 = 50;
@@ -7,7 +7,7 @@ const int x2 = 75;
const double t1 = 0.025;
const double t2 = 0.05;
-VariablesList = {t,d,n};
+VariablesList = {t,d,n,test_abs};
LocationsList = {l0, l1, l2, l3};
@@ -24,17 +24,22 @@ Locations={
};
Edges={
-((l0,l1), #, t>=0, {n=P});
+((l0,l1), #, t>=0, {n=P,test_abs=k_1*(E*S)+k_2*ES});
+
+((l1,l2),#, n>=x1 & n<=x2 & t>=t1 & t<=t2 ,{d=0});
+((l1,l2),#, t>=t2 , #);
+((l1,l2),#, test_abs=0 , #);
+((l1,l2),#, d=0 & t>=t1, #);
+
+((l1,l3), #, n>=x1 & n<=x2, {d=0});
((l1,l3), #, t<=t1 & n<=x1-1, {d=min(((t-t1)^2+(n-x2)^2)^0.5,((t-t1)^2+(n-x1)^2)^0.5)});
((l1,l3), #, t<=t1 & n>=x2+1, {d=min(((t-t1)^2+(n-x2)^2)^0.5,((t-t1)^2+(n-x1)^2)^0.5)});
-((l1,l3), #, n>=x1 & n<=x2, {d=0});
-((l1,l3), #, n<=x1-1 & t >= t1 & t<=t2, {d=min(d ,min( ((n-x1)^2)^0.5, ((n-x2)^2)^0.5 ) )});
+((l1,l3), #, n<=x1-1 & t >= t1 & t<=t2, {d=min(d,min(((n-x1)^2)^0.5,((n-x2)^2)^0.5))});
((l1,l3), #, n>=x2+1 & t >= t1 & t<=t2, {d=min(d,min( ((n-x1)^2)^0.5, ((n-x2)^2)^0.5 ) ) });
-((l3,l1), ALL, t>=0, {n=P});
+
+((l3,l1), ALL, t>=0, {n=P,test_abs=k_1*(E*S)+k_2*ES});
((l3,l2), #, t>=t2 ,#);
-((l1,l2),#, n>=x1 & n<=x2 & t>=t1 & t<=t2 ,{d=0});
-((l1,l2),#, t>=t2 , #);
-((l1,l2),#, d=0 & t>=t1, #);
};
+
diff --git a/tests/cosmos/distance_G/ER_R5.jl b/tests/cosmos/distance_G/ER_R5.jl
index b211a297e922dfc75cd5657bad9353855d0f563b..c62185c9bdcea3caccb704c2c8d6559736799131 100644
--- a/tests/cosmos/distance_G/ER_R5.jl
+++ b/tests/cosmos/distance_G/ER_R5.jl
@@ -1,6 +1,7 @@
@everywhere begin
using MarkovProcesses
+ import Distributed: nworkers
absolute_path = get_module_path() * "/tests/cosmos/"
# Values x1, x2 t1, t2
str_model = "ER"
@@ -32,7 +33,7 @@ for i = 1:nb_k1
k1 = l_k1[i]
k2 = l_k2[j]
command = `Cosmos $(absolute_path * "models/" * str_model * ".gspn")
- $(absolute_path * "distance_G/dist_G_" * str_model * ".lha") --njob $(ENV["JULIA_NUM_THREADS"])
+ $(absolute_path * "distance_G/dist_G_" * str_model * ".lha") --njob $(nworkers())
--const k_1=$(k1),k_2=$(k2),x1=$x1,x2=$x2,t1=$t1,t2=$t2
--level $(level) --width $(width)
--verbose 0`
@@ -47,16 +48,19 @@ for i = 1:nb_k1
set_param!(ER, "k2", convert(Float64, k2))
sync_ER = ER*A_G
mat_dist_pkg[i,j] = distribute_mean_value_lha(sync_ER, "d", nb_sim)
+ nb_accepts_pkg = distribute_prob_accept_lha(sync_ER, nb_sim)
+ #@info "About accepts" nb_sim nb_accepted nb_accepts_pkg
test = isapprox(mat_dist_cosmos[i,j], mat_dist_pkg[i,j]; atol = width*1.01)
+ test2 = nb_accepts_pkg == (nb_sim / nb_accepted)
if !test
@info "Distances too different" (k1,k2) mat_dist_pkg[i,j] mat_dist_cosmos[i,j]
end
- global test_all = test_all && test
+ global test_all = test_all && test && test2
end
end
-@info mat_dist_pkg
-@info mat_dist_cosmos
+@info "Distances R5 pkg" mat_dist_pkg
+@info "Distances R5 Cosmos" mat_dist_cosmos
rm("Result_dist_G_$(str_model).res")
rm("Result.res")
diff --git a/tests/cosmos/distance_G/dist_G_ER.lha b/tests/cosmos/distance_G/dist_G_ER.lha
index 72eeb5bf91af100900f60105190fb7504cb7590d..8b45dcc278014cb4a82ee6b38601eb05e807f52c 100644
--- a/tests/cosmos/distance_G/dist_G_ER.lha
+++ b/tests/cosmos/distance_G/dist_G_ER.lha
@@ -47,7 +47,7 @@ Edges={
((l1,l2), #, in=1 & t>=t2, #);
((l1,l2), #, in=0 & t>=t2, {d=d*(t2-t1)});
((l1,l2), #, test_abs=0 & t<=t1, {d=d*(t2-t1)});
-((l1,l2), #, test_abs=0 & t>=t1, {d=d+(t2-t1)*min(((n-x1)^2)^0.5,((n-x2)^2)^0.5)});
+((l1,l2), #, test_abs=0 & t>=t1 & t<=t2, {d=d+(t2-t)*min(((n-x1)^2)^0.5,((n-x2)^2)^0.5)});
((l3,l1), ALL, t>=0, {n=E,test_abs=k_1*(E*S)+k_2*ES});
diff --git a/tests/run_cosmos.jl b/tests/run_cosmos.jl
index 096a7f53738bc03bcea540fe5426810b88b3a6c5..bc2ad792d2af1dcff5769108369f1f5b1be29637 100644
--- a/tests/run_cosmos.jl
+++ b/tests/run_cosmos.jl
@@ -3,5 +3,6 @@ using Test
@testset "Cosmos tests" begin
@test include("cosmos/distance_F/ER_1D.jl")
+ @test include("cosmos/distance_G/ER_R5.jl")
end