From a24fb05571d3a7ace40bd4ad65f10fc60000f7f2 Mon Sep 17 00:00:00 2001 From: Xiaoyan Li Date: Mon, 11 Nov 2024 10:22:57 -0600 Subject: [PATCH 1/3] add statecharts interface, tests and an example SIR networks added StateCharts.jl as weak dependencies of AlgebraicABMs.jl, and also add test files and an example of SIR model with statecharts and networks. --- Project.toml | 3 + docs/Project.toml | 2 + docs/literate/sir_statechart_with_networks.jl | 113 ++++++++ ext/StateChartsInterface.jl | 256 ++++++++++++++++++ src/AlgebraicABMs.jl | 10 +- test/Project.toml | 1 + test/StateChartsInterface.jl | 59 ++++ test/runtests.jl | 4 + 8 files changed, 447 insertions(+), 1 deletion(-) create mode 100644 docs/literate/sir_statechart_with_networks.jl create mode 100644 ext/StateChartsInterface.jl create mode 100644 test/StateChartsInterface.jl diff --git a/Project.toml b/Project.toml index 461ad51..cb9e60f 100644 --- a/Project.toml +++ b/Project.toml @@ -19,13 +19,16 @@ StructEquality = "6ec83bb0-ed9f-11e9-3b4c-2b04cb4e219c" [weakdeps] AlgebraicPetri = "4f99eebe-17bf-4e98-b6a1-2c4f205a959b" Makie = "ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a" +StateCharts = "bfc335c0-d564-425b-98aa-d5c699c8b7ab" [extensions] MakieExt = "Makie" PetriInterface = "AlgebraicPetri" +StateChartsInterface = "StateCharts" [compat] AlgebraicPetri = "0.9" +StateCharts = "0.0.1" AlgebraicRewriting = "^0.4.0" Catlab = "^0.16.16" CompetingClocks = "0.1" diff --git a/docs/Project.toml b/docs/Project.toml index 0e8fff9..c7d8929 100644 --- a/docs/Project.toml +++ b/docs/Project.toml @@ -11,4 +11,6 @@ Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4" Literate = "98b081ad-f1c9-55d3-8b20-4c87d4299306" LiveServer = "16fef848-5104-11e9-1b77-fb7a48bbb589" Makie = "ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a" +Random = "9a3f8284-a2c9-5f02-9a11-845980a1fd5c" +StateCharts = "bfc335c0-d564-425b-98aa-d5c699c8b7ab" StructEquality = "6ec83bb0-ed9f-11e9-3b4c-2b04cb4e219c" diff --git a/docs/literate/sir_statechart_with_networks.jl b/docs/literate/sir_statechart_with_networks.jl new file mode 100644 index 0000000..2ff3a3f --- /dev/null +++ b/docs/literate/sir_statechart_with_networks.jl @@ -0,0 +1,113 @@ +# # Petri-net based ABMs +# +# First we want to load our packages with `using + +using AlgebraicABMs, Catlab, DataMigrations +using StateCharts +using Makie, CairoMakie + + +ENV["JULIA_DEBUG"] = "AlgebraicABMs"; # turn on @debug messages for this package + +# SIR Agent-based model +# This model is the same as the AnyLogic example model "SIR Agent-Based Networks.". +# This model currently support three types of networks: +# 1. Random (if the parameter "p_random_connect" = 1.0) +# 2. Ring lattice (if the parameter "p_random_connect" = 0.0) +# 3. Small world (if the parameter "p_random_connect": 0.0 < p_random_connect < 1.0 ) +# it indicates that there are "p_random_connect" percentage connections are radom connections, and +# ( 1 - "p_random_connect" ) percentage connections are connections are Ring lattice connections. +# This example model in Anylogic can be found: https://cloud.anylogic.com/model/7088e817-1dab-42b1-89fe-b19bc0a823e1?mode=SETTINGS + +# Step 1: define parameters +# parameters of state transition +total_population = 10 # Unit: persons +frac_initial_infected = 0.1 +contact_rate = 5 # Unit: contacts per day +infectivity = 0.01 +average_illness_duration = 15 # days +# parameters define the network +average_connections = 2 +p_random_connect = 0.9 + +# Step 2: create the StateChart + +# function UnitStateChartF takes in 3 input argummapents: +# states: tuple/array of all states +# transitions: tuple/array of all transitions. The syntax of each transition: +# transition_name => (source_state => target_state, tranisiton_type => values) +# alternative transitions: tuple/array of all alternative transitions. The syntax of each altanative transition: +# :source_transition_name=>((alternative_transition_name, probability_value) => target_state) +states = [:S, :I, :R] # define the states +transitions = [ :Infection => (:S => :I, :Pattern => 1), :Recovery => (:I => :R, :Rate => 1.0 / average_illness_duration)] +alternatives = [] # this model does not include alternative transitions +# 2.1 create the Statechart +SIRStatechart = UnitStateChartF(states, transitions, alternatives) +# 2.2 Visualization of the Statechart +stateColors = Dict(:S => "green",:I => "red",:R => "gray"); # argument define colors of each state in StateChart +StateCharts.Graph(SIRStatechart, stateColors = stateColors) + +# Step 3: create the model schema +### generate the model schema by composing auto state chart schema and the network schema + +# 3.1 define the persons object named :V, because we plan to compose the schema_statechart with graph schema by identifying ":V" +schema_statechart = StateChartABMSchema_MultipleObjects(SIRStatechart,:V) |> schemaACSet |> schemaPresent +# 3.2 compose the state chart schema with the network schema (symmetric reflective graph) as the ABM model schema (without equations) +schema_model_without_equations = compose(Open([:S],schemaACSet(schema_statechart),[:V]), Open([:V],schemaACSet(SchUndirectedReflectiveNetwork),[:E])) |> apex |> schemaPresent +# 3.3 add the composition equations of the model schema, since those equations disapper after composition. This would be fixed in the future using GATlab +@present schema_model <: schema_model_without_equations begin + compose(inv,inv) == id(E) + compose(inv,src) == tgt + compose(inv,tgt) == src + compose(refl, src) == id(V) + compose(refl, tgt) == id(V) + compose(refl, inv) == refl +end + +@acset_type SIRModelStaticNet(schema_model, index=[:src,:tgt]){Symbol} <: AbstractSymmetricReflexiveGraph +# show the model schema +to_graphviz(schema_model; prog="dot") + +# Step 4: Data Migration from the ACSet with schema -- schema_statechart to the ACSet with schema -- schema_model +# 4.1 define the data migration rule to automatically migrate ACSets from pure state chart schema to the composed model shema +const migrate_rule = @migration schema_model schema_statechart begin + V => V; + ID=>ID; VID => VID + S => S; SV => SV + I => I; IV => IV + R => R; RV => RV + + E => @product begin + p1::V + p2::V + end + src => p1 + tgt => p2 + inv => begin + p1 => p2 + p2 => p1 + end +end + +# 4.2 define the user_defined rewrite rules: rewrite rules for contacts of Infectives +# Note that each transitions_rule has a pair: timer=>rule +transition_rules = [ ContinuousHazard( 1.0 / (contact_rate * infectivity)) => make_infectious_rule_MultipleObjects(SIRStatechart, [[:S],[:I]],[:I],[[:I],[:I]], :V; use_DataMigration = true, acset = SIRModelStaticNet, migration_rule = migrate_rule)] + +# Step 5: define the network +network = smallworldNetWork(Int(total_population), average_connections, p_random_connect); + +# Step 6: Initialization +# In the initialization, randomly assign "total_population * frac_initial_infected" persons as Infectives (state I), and the rest persons are all Susceptibles (state S) +init = radomlyAssignInitialInfectives(SIRModelStaticNet(), Int(total_population), Int(frac_initial_infected * total_population); obn=:V, network = network) + +# Step 7: run the ABM model +abm = make_ABM(SIRStatechart, transition_rules, :V; is_schema_singObject=false, use_DataMigration=true, acset = SIRModelStaticNet, migration_rule = migrate_rule) +res = run!(abm, init; maxtime=50.0) + +# Step 8: Visualization of results +# 8.1 plot out the time series of each state +Makie.plot(res; Dict(o=>X->nparts(X,o) for o in states)...) +# 8.2 show the networks at time t +t = 5 +StateCharts.Graph(res,t,stateColors) + diff --git a/ext/StateChartsInterface.jl b/ext/StateChartsInterface.jl new file mode 100644 index 0000000..9895e3d --- /dev/null +++ b/ext/StateChartsInterface.jl @@ -0,0 +1,256 @@ +module StateChartsInterface + +#export StateChartABMSchema_SingleObject, StateChartABMSchema_MultipleObjects, StateChartCset_SingleObject, StateChartCset_MultipleObjects +#, make_rule_SingleObject, make_infectious_rule_SingleObject, make_rule_MultipleObjects, make_infectious_rule_MultipleObjects, get_rule, +#get_timer, make_ABM + +using StateCharts +using AlgebraicABMs, AlgebraicRewriting, Catlab +import AlgebraicABMs: StateChartABMSchema_MultipleObjects, make_infectious_rule_MultipleObjects, radomlyAssignInitialInfectives, make_ABM + +const DST = :TimeOut # the transtition type for discrete hazard +const CTT = :Rate # the transtition type for continueous hazard +const RLT = :Pattern # the transtition type for user defined rewrite rule + +""" automatically generate the ABM model schema from the state chart. This schema contains one object, and the name of this object is given by the input +argument: obn::Symbol. And each unit state chart (with one entrance transition) contributes one attributes of the object "obn" in the schema + +""" + +homs_name(cod::Symbol,dom::Symbol)=Symbol(String(cod)*String(dom)) + +vectorify(n::Vector) = collect(n) +vectorify(n::Tuple) = length(n) == 1 ? [n] : collect(n) +vectorify(n::SubArray) = collect(n) +vectorify(n) = [n] + + +# function return the ABM model Schema based on the input state charts +# the schema's structure is: +# 1. one object: usually is Person +# 2. each unit state chart contributes to one attribute + +function StateChartABMSchema_SingleObject(ss::AbstractUnitStateChart,obn::Symbol=:P) + attrsType = dnames(ss) # the attributes' type is the name of each unit state chart + obns = repeat([obn], length(attrsType)) + attrs = map((x,y)->homs_name(x,y),obns,attrsType) # define the name of the map from the object "obn" to each attribute type + attrsTuple = map((x,y,z)->(x,y,z),attrs,obns,attrsType) + schema = BasicSchema([obn], [], attrsType, attrsTuple) + schema +end + +# function return ABM model schema based on the input state charts +# the schema's structure is: +# 1. Objects: total number of states in state charts + 1 (the coproduct of the objects of states), the name of the coproduct +# object is given by input argument "obn" +# 2. Morphisms: each object of state has a morphism to the "1" -- the coproduct of states + +function StateChartABMSchema_MultipleObjects(ss::AbstractUnitStateChart,obn::Symbol=:P) + # generate the list of objects + obs = vcat(snames(ss),[obn]) + # generate the list of morphisms + morphisms_names = map((x,y)-> homs_name(x,y), snames(ss), repeat([obn],length(snames(ss)))) + morphisms = map((x,y,z)->(x,y,z),morphisms_names,snames(ss),repeat([obn],length(snames(ss)))) + # add the attribute of ID to track the individuals + attrsType = [:ID] + attrsTuple = [(homs_name(obn,:ID), obn, :ID)] + # create the schema + schema = BasicSchema(obs, morphisms, attrsType, attrsTuple) + # return the cset generated + schema +end + +function StateChartCset_SingleObject(ss::AbstractUnitStateChart,obn::Symbol=:P) + attrsType = dnames(ss) + schema = StateChartABMSchema_SingleObject(ss,obn) + tp_asgn = Dict(zip(attrsType,repeat([Symbol], length(attrsType)))) + acset = AnonACSet(schema;type_assignment=tp_asgn) + acset +end + +StateChartCset_MultipleObjects(ss::AbstractUnitStateChart,obn::Symbol=:P) = AnonACSet(StateChartABMSchema_MultipleObjects(ss,obn);type_assignment=Dict(:ID=>Symbol)) + +##### TODO: need to test the functions of "singObject" and combination of DataMigration + +## create a representable of single object schema +## this is an instance with one person, and the attributes values given by a vecoter of pairs e.g.,[:Inf=>:R,:Healthy=>:H] (Attr=>state). Then the rest attributes +## are not included in this vector are all given variables. +function representable_SingleObject(ss::AbstractUnitStateChart,as=[],obn::Symbol=:P) + rep = StateChartCset_SingleObject(ss,obn) + as = vectorify(as) + # add the attributes have vaules + [add_part!(rep,obn;x[1]=>x[2]) for x in as] + # add the attributes without values + asall = attrs(acset_schema(rep);just_names=true) + asvar = asall[asall .∉ Ref(map(first,as))] # ∉ is the unicode of \notin + # based on the schema design logic, the Attr name is the cocatenate of the codomain (obn) and domain (AttrType), so we can get the AttrType simple chop the obn from the Attr + ats = map(x->Symbol(chopprefix(String(x),String(obn))),asvar) + map((attrn,attrt)->add_part!(rep,obn;attrn=>AttrVar(add_part!(rep,attrt))),asvar,ats) + return rep +end + +# this function defines rewrite rules generate by single transition +# here i should be []. and it generate only one single person without state +function singleTransitionLIR_SingleObject(ss::AbstractUnitStateChart, l, i, r, obn::Symbol=:P) + L,R = [representable_SingleObject(ss,x,obn) for x in [l,r]] + I = representable_SingleObject(ss,i,obn) + return [L,I,R] +end + +#inf_net = Rule(only(homomorphisms(inf_m_net,inf_l_net;monic=[:V])),only(homomorphisms(inf_m_net,inf_r_net;monic=[:V],initial=(I=[2],)))) + +function mk_rule(LIR; obn::Symbol=:P, use_DataMigration::Bool = false, acset=nothing, migration_rule=nothing) + L,I,R = LIR + if !use_DataMigration + return Rule(homomorphism(I,L;monic=true),homomorphism(I,R;monic=true)) + else + Lm,Im,Rm = [migrate(acset, lir, migration_rule) for lir in LIR] + return Rule(homomorphism(Im,Lm;monic=true),homomorphism(Im,Rm;monic=true)) + end +end + +# this function defines rewrite rules generate by single transition +# here i should be []. and it generate only one single person without state +function make_rule_SingleObject(ss::AbstractUnitStateChart, l, i, r, obn::Symbol=:P; use_DataMigration::Bool = false, acset=nothing, migration_rule=nothing) + LIR = singleTransitionLIR_SingleObject(ss,l,i,r,obn) + mk_rule(LIR;use_DataMigration = use_DataMigration, obn=obn, acset = acset, migration_rule=migration_rule) +end + +function infectiousruleLIR_SingleObject(ss::AbstractUnitStateChart, l, i, r, obn::Symbol=:P) + L = coproduct([representable_SingleObject(ss,p,obn) for p in l]) |> apex + I = coproduct([representable_SingleObject(ss,[],obn),representable_SingleObject(ss,i,obn)]) |> apex + R = coproduct([representable_SingleObject(ss,p,obn) for p in r]) |> apex + return [L,I,R] +end + +function make_infectious_rule_SingleObject(ss::AbstractUnitStateChart, l, i, r, obn::Symbol=:P; use_DataMigration::Bool = false, acset=nothing, migration_rule=nothing) + LIR = infectiousruleLIR_SingleObject(ss, l, i, r, obn) + mk_rule(LIR;use_DataMigration = use_DataMigration, acset = acset, migration_rule=migration_rule) +end +############################### end of the functions of "singleObject" ################################################################ + +## create a representable of multiple object schema +## this is an instance with one person, and each state (may accross multiple state charts) is an object point to the person +## as a set of states. e.g., as =[:S,:H], which indicates this person is Susceptible and Healthy. ID attributes is a variable +## NOTE: if represent a single person, without any states, it can write as []. But this indicates a person with a variable attributes ID +function representable_MultipleObjects(ss::AbstractUnitStateChart,as,obn::Symbol=:P) + println(as) + rep = StateChartCset_MultipleObjects(ss,obn) + as = vectorify(as) + + # add the object of person + add_part!(rep,obn;homs_name(obn,:ID)=>AttrVar(add_part!(rep,:ID))) + # add the (state) objects to L, I and R + [add_part!(rep, s; homs_name(s,obn)=>i) for (i,s) in enumerate(as)] + + return rep +end + +# this function defines rewrite rules generate by single transition +# here i should be []. and it generate only one single person without state +function make_rule_MultipleObjects(ss::AbstractUnitStateChart, l, i, r, obn::Symbol=:P; use_DataMigration::Bool = false, acset=nothing, migration_rule=nothing) + LIR = [representable_MultipleObjects(ss,x,obn) for x in [l,i,r]] + mk_rule(LIR;obn=obn, use_DataMigration = use_DataMigration, acset = acset, migration_rule=migration_rule) +end + +function infectiousruleLIR_MultipleObjects(ss::AbstractUnitStateChart, l, i, r, obn::Symbol=:P) + L = coproduct([representable_MultipleObjects(ss,p,obn) for p in l]) |> apex + I = coproduct([representable_MultipleObjects(ss,[],obn),representable_MultipleObjects(ss,i,obn)]) |> apex + R = coproduct([representable_MultipleObjects(ss,p,obn) for p in r]) |> apex + return [L,I,R] +end + +function make_infectious_rule_MultipleObjects(ss::AbstractUnitStateChart, l, i, r, obn::Symbol=:P; use_DataMigration::Bool = false, acset=nothing, migration_rule=nothing) + LIR = infectiousruleLIR_MultipleObjects(ss,l,i,r,obn) + L,I,R = use_DataMigration ? [migrate(acset, lir, migration_rule) for lir in LIR] : LIR + return Rule(homomorphism(I,L;monic=[obn]),homomorphism(I,R;monic=[obn])) +end + +# generate the rewrite rule +# t is the index of a transition +# transitions_rules is an array with pair of (timer, rule) as its elements +function get_rule(ss::AbstractUnitStateChart,t,transitions_rules,obn::Symbol=:P;is_schema_singObject::Bool=true, use_DataMigration::Bool = false, acset=nothing, migration_rule=nothing) + tt = ttype(ss,t) + if tt == RLT + tidx=Int(texpr(ss,t)) + return transitions_rules[tidx][2] + elseif tt in [DST, CTT] + # we assume all the states of other state charts are all unkown (set as variables) + l = sname(ss,tsource(ss,t)) + r = sname(ss,ttarget(ss,t)) +# statechart = dname(ss,spartition(ss,tsource(ss,t))) + statechart = dname(ss,t) + attrtype = homs_name(obn,statechart) + if is_schema_singObject + if !use_DataMigration + return make_rule_SingleObject(ss,[attrtype=>l],[],[attrtype=>r],obn) + else + return make_rule_SingleObject(ss,[attrtype=>l],[],[attrtype=>r],obn; use_DataMigration = true, acset = acset, migration_rule=migration_rule) + end + else # multiple object schema + if !use_DataMigration + return make_rule_MultipleObjects(ss,[l],[],[r],obn) + else + return make_rule_MultipleObjects(ss,[l],[],[r],obn; use_DataMigration = true, acset = acset, migration_rule=migration_rule) + end + end + else + tn = tname(ss,t) + throw("transition:$tn does not have valid transition type!") + end +end + +function get_timer(ss::AbstractUnitStateChart, transitions_rules,t) + tt=ttype(ss,t) + if tt == DST + timeout = texpr(ss,t) + timer = DiscreteHazard(timeout) + elseif tt == CTT + r = texpr(ss,t) + timer = ContinuousHazard(1/r) + elseif tt == RLT + tidx=Int(texpr(ss,t)) + timer = transitions_rules[tidx][1] + else + tn = tname(ss,t) + throw("transition:$tn does not have valid transition type!") + end + timer +end + +# randomly assign "ninf" persons as Infectives, and the rest are all Susceptibles +function radomlyAssignInitialInfectives(acset, totalPopulation::Int64, ninf::Int64; I::Symbol = :I, S::Symbol = :S, obn::Symbol = :P, network = nothing) + init = acset # create an instance of the model based on model schema + + # get the radom samples of infectives + infs = rand(1:totalPopulation, ninf) + # the rest are all Susceptibles + spts = findall(!in(infs),1:totalPopulation) + + if isnothing(network) + add_parts!(init,obn,totalPopulation;homs_name(obn,:ID)=>map(id->:V*Symbol(id),1:Int(totalPopulation))) + else + copy_parts!(init, network) # set the vertices' ID + set_subpart!(init, 1:Int(totalPopulation), homs_name(obn,:ID), map(id->:V*Symbol(id),1:Int(totalPopulation))) + end + + add_parts!(init,S,Int(totalPopulation-ninf);homs_name(S,obn)=>spts) + add_parts!(init,I,ninf;homs_name(I,obn)=>infs) + + return init +end + + +make_ABMRule(ss::AbstractUnitStateChart,t,transitions_rules, obn::Symbol=:P; is_schema_singObject::Bool=true, use_DataMigration::Bool = false, acset=nothing, migration_rule=nothing) = begin + tn = tname(ss, t) + rule = get_rule(ss, t,transitions_rules,obn; is_schema_singObject=is_schema_singObject, use_DataMigration=use_DataMigration, acset = acset, migration_rule=migration_rule) + timer = get_timer(ss,transitions_rules,t) + ABMRule(tn, rule, timer) +end + +make_ABM(ss::AbstractUnitStateChart, transitions_rules, obn=:P;is_schema_singObject::Bool=true, use_DataMigration::Bool = false, acset=nothing, migration_rule=nothing) = begin + abmrules = [make_ABMRule(ss,t,transitions_rules, obn; is_schema_singObject=is_schema_singObject, use_DataMigration=use_DataMigration, acset=acset, migration_rule=migration_rule) for t in 1:nt(ss)] + ABM(abmrules) +end + +end diff --git a/src/AlgebraicABMs.jl b/src/AlgebraicABMs.jl index 986be80..ffd31e0 100644 --- a/src/AlgebraicABMs.jl +++ b/src/AlgebraicABMs.jl @@ -2,7 +2,7 @@ """ module AlgebraicABMs -export PetriNetCSet +export PetriNetCSet, StateChartABMSchema_MultipleObjects, make_infectious_rule_MultipleObjects, radomlyAssignInitialInfectives, make_ABM using Reexport @@ -19,4 +19,12 @@ include("Visualization.jl") # Methods to be implemented by extensions function PetriNetCSet end +# statecharts +function StateChartABMSchema_MultipleObjects end +function make_infectious_rule_MultipleObjects end +function radomlyAssignInitialInfectives end +function make_ABM end + + + end # module diff --git a/test/Project.toml b/test/Project.toml index 1b51266..c1e6250 100644 --- a/test/Project.toml +++ b/test/Project.toml @@ -11,4 +11,5 @@ DifferentialEquations = "0c46a032-eb83-5123-abaf-570d42b7fbaa" Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f" Makie = "ee78f7c6-11fb-53f2-987a-cfe4a2b5a57a" ProfileView = "c46f51b8-102a-5cf2-8d2c-8597cb0e0da7" +StateCharts = "bfc335c0-d564-425b-98aa-d5c699c8b7ab" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" diff --git a/test/StateChartsInterface.jl b/test/StateChartsInterface.jl new file mode 100644 index 0000000..9224571 --- /dev/null +++ b/test/StateChartsInterface.jl @@ -0,0 +1,59 @@ +module TestStateChartsInterface + +using Test +using Makie, CairoMakie +using StateCharts +using AlgebraicABMs +using AlgebraicABMs.ABMs: RepresentableP, RegularP +using Catlab + + +totalPopulation=5.0 +β = 0.1 +p_becomingInfective = 5.0 + +# create the pertussis state chart +pertussis = UnitStateChartF((:S,:E,:I,:R₄,:V₄), # states +(:newExposure=>(:S=>:E,:Pattern=>1),:becomingInfective=>(:E=>:I,:TimeOut=>p_becomingInfective),:recovery=>(:I=>:R₄,:Rate=>1.0/21.0),:vaccinated=>(:S=>:V₄,:Rate=>0.01)), #non-start transitions +() # alternatives for non-start transitions +) + +StateCharts.Graph(pertussis) +schema = StateChartABMSchema_MultipleObjects(pertussis) |> schemaACSet |> schemaPresent +to_graphviz(schema; prog="dot") + + +########### use ABM model schema of one object P ##################### +########### the states of state charts are attributes ################ + +# define the new_infectious rule +transition_rules = [ ContinuousHazard( 1.0 / β) => make_infectious_rule_MultipleObjects(pertussis, [[:S],[:I]],[:I],[[:E],[:I]])] + +# Initial state for schema of multiple object +@acset_type PertussisNet(schema, index=[:src,:tgt]){Symbol} +init = radomlyAssignInitialInfectives(PertussisNet(), Int(totalPopulation), Int(1)) + +## test ABM rules +abm = make_ABM(pertussis, transition_rules; is_schema_singObject=false) + +@test abm.rules[1].pattern_type isa RepresentableP +@test abm.rules[2].pattern_type == RegularP() +@test abm.rules[3].pattern_type isa RepresentableP +@test abm.rules[4].pattern_type isa RepresentableP +# transition "newExposure" +@test nparts.(Ref(codom(Catlab.left(abm.rules[1].rule))), [:S,:E,:I,:R₄,:V₄,:P]) == [1,0,1,0,0,2] # S->p1, I->p2 +@test nparts.(Ref(dom(Catlab.left(abm.rules[1].rule))), [:S,:E,:I,:R₄,:V₄,:P]) == [0,0,1,0,0,2] # p1, I->p2 +@test nparts.(Ref(codom(Catlab.right(abm.rules[1].rule))), [:S,:E,:I,:R₄,:V₄,:P]) == [0,1,1,0,0,2] # E->p1, I->p2 +# transition "recovery" +@test nparts.(Ref(codom(Catlab.left(abm.rules[3].rule))), [:S,:E,:I,:R₄,:V₄,:P]) == [0,0,1,0,0,1] # I->p +@test nparts.(Ref(dom(Catlab.left(abm.rules[3].rule))), [:S,:E,:I,:R₄,:V₄,:P]) == [0,0,0,0,0,1] # p +@test nparts.(Ref(codom(Catlab.right(abm.rules[3].rule))), [:S,:E,:I,:R₄,:V₄,:P]) == [0,0,0,1,0,1] # R4->p +# transition "vaccinated" +@test nparts.(Ref(codom(Catlab.left(abm.rules[4].rule))), [:S,:E,:I,:R₄,:V₄,:P]) == [1,0,0,0,0,1] # S->p +@test nparts.(Ref(dom(Catlab.left(abm.rules[4].rule))), [:S,:E,:I,:R₄,:V₄,:P]) == [0,0,0,0,0,1] # p +@test nparts.(Ref(codom(Catlab.right(abm.rules[4].rule))), [:S,:E,:I,:R₄,:V₄,:P]) == [0,0,0,0,1,1] # V4->p + +#res = run!(abm, init; maxtime=50.0) +#Makie.plot(res; Dict(o=>X->nparts(X,o) for o in [:S,:E,:I,:R₄,:V₄])...) + +end \ No newline at end of file diff --git a/test/runtests.jl b/test/runtests.jl index 51b6afd..f958400 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -11,3 +11,7 @@ end @testset "PetriInterface" begin include("PetriInterface.jl") end + +@testset "StateChartsInterface" begin + include("TestStateChartsInterface") +end From 1b203829b0a71ca027ef37b1c05d8632602cd60a Mon Sep 17 00:00:00 2001 From: Xiaoyan Li Date: Mon, 11 Nov 2024 11:54:39 -0600 Subject: [PATCH 2/3] fix a bug --- ext/StateChartsInterface.jl | 13 ++++++++++--- test/StateChartsInterface.jl | 2 +- 2 files changed, 11 insertions(+), 4 deletions(-) diff --git a/ext/StateChartsInterface.jl b/ext/StateChartsInterface.jl index 9895e3d..5bb4d81 100644 --- a/ext/StateChartsInterface.jl +++ b/ext/StateChartsInterface.jl @@ -134,7 +134,6 @@ end ## as a set of states. e.g., as =[:S,:H], which indicates this person is Susceptible and Healthy. ID attributes is a variable ## NOTE: if represent a single person, without any states, it can write as []. But this indicates a person with a variable attributes ID function representable_MultipleObjects(ss::AbstractUnitStateChart,as,obn::Symbol=:P) - println(as) rep = StateChartCset_MultipleObjects(ss,obn) as = vectorify(as) @@ -160,10 +159,18 @@ function infectiousruleLIR_MultipleObjects(ss::AbstractUnitStateChart, l, i, r, return [L,I,R] end -function make_infectious_rule_MultipleObjects(ss::AbstractUnitStateChart, l, i, r, obn::Symbol=:P; use_DataMigration::Bool = false, acset=nothing, migration_rule=nothing) +# RdoubleI = false, indicates the R (right object) include one Infective person. E.g., the Susceptible transite to Exposed (E) if infected +# RdoubleI = true, indicates the R (right object) include two Infective persons. E.g., the Susceptible transite to Infective (I) if infected +# the basic idea is whether +function make_infectious_rule_MultipleObjects(ss::AbstractUnitStateChart, l, i, r, obn::Symbol=:P; RdoubleI::Bool = true, use_DataMigration::Bool = false, acset=nothing, migration_rule=nothing) LIR = infectiousruleLIR_MultipleObjects(ss,l,i,r,obn) L,I,R = use_DataMigration ? [migrate(acset, lir, migration_rule) for lir in LIR] : LIR - return Rule(homomorphism(I,L;monic=[obn]),homomorphism(I,R;monic=[obn])) + if RdoubleI + rule = Rule(homomorphism(I,L;monic=[obn]),homomorphism(I,R;monic=[obn],initial=(I=[2],))) + else + rule = Rule(homomorphism(I,L;monic=[obn]),homomorphism(I,R;monic=[obn])) + end + return rule end # generate the rewrite rule diff --git a/test/StateChartsInterface.jl b/test/StateChartsInterface.jl index 9224571..0a293d7 100644 --- a/test/StateChartsInterface.jl +++ b/test/StateChartsInterface.jl @@ -27,7 +27,7 @@ to_graphviz(schema; prog="dot") ########### the states of state charts are attributes ################ # define the new_infectious rule -transition_rules = [ ContinuousHazard( 1.0 / β) => make_infectious_rule_MultipleObjects(pertussis, [[:S],[:I]],[:I],[[:E],[:I]])] +transition_rules = [ ContinuousHazard( 1.0 / β) => make_infectious_rule_MultipleObjects(pertussis, [[:S],[:I]],[:I],[[:E],[:I]], RdoubleI = false)] # Initial state for schema of multiple object @acset_type PertussisNet(schema, index=[:src,:tgt]){Symbol} From d55234b18ba0de4208d9e2ee2d3e5c2aef857d33 Mon Sep 17 00:00:00 2001 From: Xiaoyan Li Date: Wed, 13 Nov 2024 15:26:55 -0600 Subject: [PATCH 3/3] StateCharts as weak deps tested work well --- docs/literate/sir_statechart_with_networks.jl | 25 +++++++++++++++---- test/StateChartsInterface.jl | 2 +- test/runtests.jl | 2 +- 3 files changed, 22 insertions(+), 7 deletions(-) diff --git a/docs/literate/sir_statechart_with_networks.jl b/docs/literate/sir_statechart_with_networks.jl index 2ff3a3f..64a4bdf 100644 --- a/docs/literate/sir_statechart_with_networks.jl +++ b/docs/literate/sir_statechart_with_networks.jl @@ -21,8 +21,8 @@ ENV["JULIA_DEBUG"] = "AlgebraicABMs"; # turn on @debug messages for this packag # Step 1: define parameters # parameters of state transition -total_population = 10 # Unit: persons -frac_initial_infected = 0.1 +total_population = 50 # Unit: persons +frac_initial_infected = 0.02 contact_rate = 5 # Unit: contacts per day infectivity = 0.01 average_illness_duration = 15 # days @@ -101,13 +101,28 @@ network = smallworldNetWork(Int(total_population), average_connections, p_random init = radomlyAssignInitialInfectives(SIRModelStaticNet(), Int(total_population), Int(frac_initial_infected * total_population); obn=:V, network = network) # Step 7: run the ABM model +ts = 50.0 # the total running time abm = make_ABM(SIRStatechart, transition_rules, :V; is_schema_singObject=false, use_DataMigration=true, acset = SIRModelStaticNet, migration_rule = migrate_rule) -res = run!(abm, init; maxtime=50.0) +res = run!(abm, init; maxtime=ts); # Step 8: Visualization of results # 8.1 plot out the time series of each state Makie.plot(res; Dict(o=>X->nparts(X,o) for o in states)...) # 8.2 show the networks at time t -t = 5 -StateCharts.Graph(res,t,stateColors) +t0 = 0 # initial state +StateCharts.Graph(res,t0,stateColors) +t10 = 10 # initial state +StateCharts.Graph(res,t10,stateColors) + +t20 = 20 +StateCharts.Graph(res,t20,stateColors) + +t30 = 30 +StateCharts.Graph(res,t30,stateColors) + +t40 = 40 +StateCharts.Graph(res,t40,stateColors) + +t50 = Int(ts) # final state +StateCharts.Graph(res,t50,stateColors) \ No newline at end of file diff --git a/test/StateChartsInterface.jl b/test/StateChartsInterface.jl index 0a293d7..8098182 100644 --- a/test/StateChartsInterface.jl +++ b/test/StateChartsInterface.jl @@ -53,7 +53,7 @@ abm = make_ABM(pertussis, transition_rules; is_schema_singObject=false) @test nparts.(Ref(dom(Catlab.left(abm.rules[4].rule))), [:S,:E,:I,:R₄,:V₄,:P]) == [0,0,0,0,0,1] # p @test nparts.(Ref(codom(Catlab.right(abm.rules[4].rule))), [:S,:E,:I,:R₄,:V₄,:P]) == [0,0,0,0,1,1] # V4->p -#res = run!(abm, init; maxtime=50.0) +res = run!(abm, init; maxtime=10.0) #Makie.plot(res; Dict(o=>X->nparts(X,o) for o in [:S,:E,:I,:R₄,:V₄])...) end \ No newline at end of file diff --git a/test/runtests.jl b/test/runtests.jl index f958400..b7d3a3c 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -13,5 +13,5 @@ end end @testset "StateChartsInterface" begin - include("TestStateChartsInterface") + include("StateChartsInterface.jl") end