From dd619a953f045e5e858da2b31babf5dbcb93a4a7 Mon Sep 17 00:00:00 2001 From: Mahmoud Bentriou <mahmoud.bentriou@centralesupelec.fr> Date: Sun, 21 Feb 2021 12:57:16 +0100 Subject: [PATCH] Fix of the segfault generated by the euclidean automaton test. Julia shouldn't crash but rather raise an error about the existence of a function generated by metaprogramming. I didn't manage to isolate the segfault withtout the package. To overcome the issue, I add another level of multiple dispatch/abstract type for synchronized models. Test of the euclidean distance automaton works. --- automata/euclidean_distance_automaton.jl | 21 +++--- core/common.jl | 18 +++-- core/lha.jl | 8 +-- core/model.jl | 86 +++++++++++++----------- core/trajectory.jl | 2 +- 5 files changed, 74 insertions(+), 61 deletions(-) diff --git a/automata/euclidean_distance_automaton.jl b/automata/euclidean_distance_automaton.jl index f4c09aa..720a01b 100644 --- a/automata/euclidean_distance_automaton.jl +++ b/automata/euclidean_distance_automaton.jl @@ -1,17 +1,22 @@ +# Creation of the automaton types +lha_name = :EuclideanDistanceAutomaton +edge_type = :EdgeEuclideanDistanceAutomaton +@everywhere @eval abstract type $(edge_type) <: Edge end +@everywhere @eval $(MarkovProcesses.generate_code_lha_type_def(lha_name, edge_type)) + 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." @assert length(timeline) == length(observations) "Timeline and observations vectors don't have the same length" nbr_observations = length(observations) - # Creation of the automaton types + # Automaton types and functions + model_name = Symbol(typeof(m)) lha_name = :EuclideanDistanceAutomaton edge_type = :EdgeEuclideanDistanceAutomaton check_constraints = Symbol("check_constraints_$(lha_name)") update_state! = Symbol("update_state_$(lha_name)!") - @everywhere @eval abstract type $(edge_type) <: Edge end - @everywhere @eval $(MarkovProcesses.generate_code_lha_type_def(lha_name, edge_type)) # Locations locations = [:l0, :l1, :l2] @@ -37,7 +42,6 @@ function create_euclidean_distance_automaton(m::ContinuousTimeModel, timeline::A to_idx(var::Symbol) = map_var_automaton_idx[var] id = MarkovProcesses.newid() - model_name = Symbol(typeof(m)) basename_func = "$(model_name)_$(id)" edge_name(from_loc::Location, to_loc::Location, edge_number::Int) = Symbol("Edge_$(lha_name)_$(basename_func)_$(from_loc)$(to_loc)_$(edge_number)") @@ -117,13 +121,12 @@ function create_euclidean_distance_automaton(m::ContinuousTimeModel, timeline::A end # Updating next_state! + @everywhere @eval $(MarkovProcesses.generate_code_synchronized_model_type_def(model_name, lha_name)) @everywhere @eval $(MarkovProcesses.generate_code_next_state(lha_name, edge_type, check_constraints, update_state!)) - @everywhere @eval $(MarkovProcesses.generate_code_synchronized_simulation(lha_name, edge_type, m.f!, m.isabsorbing)) + @everywhere @eval $(MarkovProcesses.generate_code_synchronized_simulation(model_name, lha_name, edge_type, m.f!, m.isabsorbing)) - @eval begin - A = $(lha_name)($(m.transitions), $(locations), $(Λ_F), $(locations_init), $(locations_final), - $(map_var_automaton_idx), $(flow), $(map_edges), $(constants), $(m.map_var_idx)) - end + A = EuclideanDistanceAutomaton(m.transitions, locations, Λ_F, locations_init, locations_final, + map_var_automaton_idx, flow, map_edges, constants, m.map_var_idx) return A end diff --git a/core/common.jl b/core/common.jl index 407383a..4a501df 100644 --- a/core/common.jl +++ b/core/common.jl @@ -2,6 +2,7 @@ abstract type Model end abstract type ContinuousTimeModel <: Model end +abstract type SynchronizedModel <: Model end abstract type AbstractTrajectory end abstract type LHA end @@ -67,9 +68,17 @@ mutable struct StateLHA time::Float64 end -mutable struct SynchronizedModel <: Model - m::ContinuousTimeModel - automaton::LHA +function generate_code_synchronized_model_type_def(model_name::Symbol, lha_name::Symbol) + synchronized_model_name = Symbol("$(model_name)SynchronizedWith$(lha_name)") + return quote + mutable struct $(synchronized_model_name) <: SynchronizedModel + m::$(model_name) + automaton::$(lha_name) + end + + Base.:*(m::$(model_name), A::$(lha_name)) = $(synchronized_model_name)(m, A) + Base.:*(A::$(lha_name), m::$(model_name)) = $(synchronized_model_name)(m, A) + end end struct SynchronizedTrajectory <: AbstractTrajectory @@ -122,9 +131,6 @@ LHA(A::LHA, map_var::Dict{VariableModel,Int}) = getfield(Main, Symbol(typeof(A)))(A.transitions, A.locations, A.Λ, A.locations_init, A.locations_final, A.map_var_automaton_idx, A.flow, A.map_edges, A.constants, map_var) -Base.:*(m::ContinuousTimeModel, A::LHA) = SynchronizedModel(m, A) -Base.:*(A::LHA, m::ContinuousTimeModel) = SynchronizedModel(m, A) - function ParametricModel(am::Model, priors::Tuple{ParameterModel,UnivariateDistribution}...) m = get_proba_model(am) params = ParameterModel[] diff --git a/core/lha.jl b/core/lha.jl index acf19da..cd7dc2a 100644 --- a/core/lha.jl +++ b/core/lha.jl @@ -66,7 +66,7 @@ function Base.copyto!(Sdest::StateLHA, Ssrc::StateLHA) Sdest.A = Ssrc.A Sdest.loc = Ssrc.loc for i = eachindex(Sdest.values) - @inbounds Sdest.values[i] = Ssrc.values[i] + Sdest.values[i] = Ssrc.values[i] end Sdest.time = Ssrc.time end @@ -96,7 +96,7 @@ function generate_code_next_state(lha_name::Symbol, edge_type::Symbol, # A push! method implementend by myself because of preallocation of edge_candidates function _push_edge!(edge_candidates::Vector{<:$(edge_type)}, edge::$(edge_type), nbr_candidates::Int) if nbr_candidates < length(edge_candidates) - @inbounds edge_candidates[nbr_candidates+1] = edge + edge_candidates[nbr_candidates+1] = edge else push!(edge_candidates, edge) end @@ -218,9 +218,9 @@ function generate_code_next_state(lha_name::Symbol, edge_type::Symbol, end # Now time flies according to the flow for i in eachindex(values_state) - @inbounds coeff_deriv = flow[ptr_loc_state[1]][i] + coeff_deriv = flow[ptr_loc_state[1]][i] if coeff_deriv > 0 - @inbounds values_state[i] += coeff_deriv*(tnplus1 - ptr_time_state[1]) + values_state[i] += coeff_deriv*(tnplus1 - ptr_time_state[1]) end end ptr_time_state[1] = tnplus1 diff --git a/core/model.jl b/core/model.jl index e56f9d5..6fad681 100644 --- a/core/model.jl +++ b/core/model.jl @@ -4,30 +4,30 @@ import Distributions: insupport, pdf function _resize_trajectory!(values::Vector{Vector{Int}}, times::Vector{Float64}, transitions::Vector{Transition}, size::Int) - for i = eachindex(values) resize!(@inbounds(values[i]), size) end + for i = eachindex(values) resize!((values[i]), size) end resize!(times, size) resize!(transitions, size) end function _finish_bounded_trajectory!(values::Vector{Vector{Int}}, times::Vector{Float64}, transitions::Vector{Transition}, time_bound::Float64) - for i = eachindex(values) push!(@inbounds(values[i]), @inbounds(values[i][end])) end + for i = eachindex(values) push!((values[i]), (values[i][end])) end push!(times, time_bound) push!(transitions, nothing) end function _update_values!(values::Vector{Vector{Int}}, times::Vector{Float64}, transitions::Vector{Transition}, xn::Vector{Int}, tn::Float64, tr_n::Transition, idx::Int) - for k = eachindex(values) @inbounds(values[k][idx] = xn[k]) end - @inbounds(times[idx] = tn) - @inbounds(transitions[idx] = tr_n) + for k = eachindex(values) values[k][idx] = xn[k] end + (times[idx] = tn) + (transitions[idx] = tr_n) end function generate_code_simulation(model_name::Symbol, f!::Symbol, isabsorbing::Symbol) return quote import MarkovProcesses: simulate - + """ `simulate(m)` @@ -132,34 +132,13 @@ function generate_code_simulation(model_name::Symbol, f!::Symbol, isabsorbing::S end end -function simulate(product::SynchronizedModel; - p::Union{Nothing,AbstractVector{Float64}} = nothing, verbose::Bool = false) - m = getfield(product, :m) - A = getfield(product, :automaton) - p_sim = getfield(m, :p) - if p != nothing - p_sim = p - end - return simulate(m, A, product, p_sim, verbose) -end - -function volatile_simulate(product::SynchronizedModel; - p::Union{Nothing,AbstractVector{Float64}} = nothing, verbose::Bool = false) - m = product.m - A = product.automaton - p_sim = getfield(m, :p) - if p != nothing - p_sim = p - end - return volatile_simulate(m, A, p_sim, verbose) -end - -function generate_code_synchronized_simulation(lha_name::Symbol, edge_type::Symbol, f!::Symbol, isabsorbing::Symbol) +function generate_code_synchronized_simulation(model_name::Symbol, lha_name::Symbol, + edge_type::Symbol, f!::Symbol, isabsorbing::Symbol) return quote import MarkovProcesses: simulate, volatile_simulate - - function simulate(m::ContinuousTimeModel, A::$(lha_name), product::SynchronizedModel, + + function simulate(m::$(model_name), A::$(lha_name), product::SynchronizedModel, p_sim::AbstractVector{Float64}, verbose::Bool) x0 = getfield(m, :x0) t0 = getfield(m, :t0) @@ -293,14 +272,7 @@ function generate_code_synchronized_simulation(lha_name::Symbol, edge_type::Symb return SynchronizedTrajectory(S, product, values, times, transitions) end - """ - `volatile_simulate(sm::SynchronizedModel; p, verbose)` - - Simulates a model synchronized with an automaton but does not store the values of the simulation - in order to improve performance. - It returns the last state of the simulation `S::StateLHA` not a trajectory `σ::SynchronizedTrajectory`. - """ - function volatile_simulate(m::ContinuousTimeModel, A::$(lha_name), p_sim::AbstractVector{Float64}, verbose::Bool) + function volatile_simulate(m::$(model_name), A::$(lha_name), p_sim::AbstractVector{Float64}, verbose::Bool) x0 = getfield(m, :x0) t0 = getfield(m, :t0) time_bound = getfield(m, :time_bound) @@ -360,6 +332,38 @@ function generate_code_synchronized_simulation(lha_name::Symbol, edge_type::Symb end end +""" +`volatile_simulate(sm::SynchronizedModel; p, verbose)` + +Simulates a model synchronized with an automaton but does not store the values of the simulation +in order to improve performance. +It returns the last state of the simulation `S::StateLHA` not a trajectory `σ::SynchronizedTrajectory`. +""" +function volatile_simulate(product::SynchronizedModel; + p::Union{Nothing,AbstractVector{Float64}} = nothing, verbose::Bool = false) + m = product.m + A = product.automaton + p_sim = getfield(m, :p) + if p != nothing + p_sim = p + end + S = volatile_simulate(m, A, p_sim, verbose) + return S +end + +function simulate(product::SynchronizedModel; + p::Union{Nothing,AbstractVector{Float64}} = nothing, verbose::Bool = false) + m = getfield(product, :m) + A = getfield(product, :automaton) + p_sim = getfield(m, :p) + if p != nothing + p_sim = p + end + σ = simulate(m, A, product, p_sim, verbose) + return σ +end + + """ `simulate(pm::ParametricModel, p_prior::AbstractVector{Float64}) @@ -384,7 +388,7 @@ function volatile_simulate(pm::ParametricModel, p_prior::AbstractVector{Float64} epsilon::Float64) @assert typeof(pm.m) <: SynchronizedModel # ABC related automata - if pm.m.name in ["ABC euclidean distance"] + if typeof(pm.m.A) in <: EuclideanDistanceABCAutomaton nothing end full_p = copy(get_proba_model(pm).p) @@ -446,7 +450,7 @@ end number_simulations_smc_chernoff(approx::Float64, conf::Float64) = log(2/(1-conf)) / (2*approx^2) function smc_chernoff(sm::SynchronizedModel; approx::Float64 = 0.01, confidence::Float64 = 0.99) - @assert sm.automaton.name in ["F property", "G property", "G and F property"] + @assert typeof(sm.automaton) <: Union{AutomatonF,AutomatonG,AutomatonGandF} nbr_sim = number_simulations_smc_chernoff(approx, confidence) nbr_sim = convert(Int, trunc(nbr_sim)+1) return probability_var_value_lha(sm, nbr_sim) diff --git a/core/trajectory.jl b/core/trajectory.jl index 8a55b46..da97f76 100644 --- a/core/trajectory.jl +++ b/core/trajectory.jl @@ -210,7 +210,7 @@ function Base.show(io::IO, σ::SynchronizedTrajectory) print(io, "End LHA state:\n") print(io, σ.state_lha_end) print(io, "\n") - print(io, "- Model name: $(σ.m.name) \n") + print(io, "- Model: $(typeof(σ.m)) \n") print(io, "- Variable trajectories:\n") for obs_var in σ.m.g print(io, "* $obs_var: $(σ[obs_var])\n") -- GitLab