push event progress

Author: TorkelE

src/reactionsystem.jl

@@ -214,10 +214,15 @@ function find_event_vars!(ps, us, events::Vector, ivs, vars)
 end
 # For a single event, adds quantities from its condition and affect expression(s) to `ps` and `us`.
 # Applies `findvars!` to the event's condition (`event[1])` and affec (`event[2]`).
-function find_event_vars!(ps, us, event, ivs, vars)
+# Two dispatches required, and the event can be given as a MTK structure of a Pair of symbolic expressions/equations.
+function find_event_vars!(ps, us, event::Pair, ivs, vars)
     findvars!(ps, us, event[1], ivs, vars)
     findvars!(ps, us, event[2], ivs, vars)
 end
+function find_event_vars!(ps, us, event::MT.AbstractCallback, ivs, vars)
+    findvars!(ps, us, event.conditions, ivs, vars)
+    findvars!(ps, us, event.affect, ivs, vars)
+end
 
 ### ReactionSystem Structure ###
 
@@ -479,7 +484,8 @@ function ReactionSystem(eqs, iv, unknowns, ps;
     end
 
     # Creates the continuous and discrete events.
-    continuous_events, discrete_events = MT.create_symbolic_events(continuous_events, discrete_events)
+    continuous_events = create_symbolic_events(MT.SymbolicContinuousCallback, continuous_events)
+    discrete_events = create_symbolic_events(MT.SymbolicDiscreteCallback, discrete_events)
 
     # handles system metadata.
     metadata === nothing ? Base.ImmutableDict{Symbol,Any}() : metadata
@@ -501,6 +507,43 @@ function ReactionSystem(iv; kwargs...)
     ReactionSystem(Reaction[], iv, [], []; kwargs...)
 end
 
+# handles that events can be a single event or a vector.
+create_symbolic_events(type, events::Vector) = [create_symbolic_event(type, event) for event in events]
+create_symbolic_events(type, event) = [create_symbolic_event(type, event)]
+
+# Converts an input event into a form which ModelingToolkit can handle.
+create_symbolic_event(type, event::ModelingToolkit.AbstractCallback) = event
+function create_symbolic_event(type, event)
+    discrete_parameters = find_disc_pars(event[2]) 
+    return type(event; discrete_parameters)
+end
+
+# Currently, ModelingToolkit requires all parameters which are updated in an event to be specified
+# in the event using the `discrete_parameters` keyword argument. Currently, this can only be
+# done by manually creating the event using `SymbolicContinuousCallback`/`SymbolicDiscreteCallback`.
+# This routines preserves the old behaviour until a new interface have been implemented in MTK.
+function find_disc_pars(affects)
+    discrete_parameters = Set()
+    for affect in affects
+        num_disc_pars = length(discrete_parameters)
+        infer_discrete_parameters!(discrete_parameters, affect.lhs)
+        infer_discrete_parameters!(discrete_parameters, affect.rhs)
+        (length(discrete_parameters) > num_disc_pars + 1) && error("Currently, Catalyst only supports having a single parameter outside of a `Pre` statement in event affects.")
+    end
+    return collect(discrete_parameters)
+end
+
+# Find all `expr`'s parameters that occur *outside* of a Pre(...) statement. Add these to `discrete_parameters`.
+function infer_discrete_parameters!(discrete_parameters, expr)
+    expr_pre_removed = Symbolics.replacenode(expr, precall_to_1)
+    dynamic_symvars = Symbolics.get_variables(expr_pre_removed)
+    # Change this coming line to a Symbolic append type of thing.
+    union!(discrete_parameters, filter(ModelingToolkit.isparameter, dynamic_symvars))
+end
+# Functions for replacing a Pre-call with a `1.0` (removing its content from an expression).
+is_precall(expr) = iscall(expr) ? operation(expr) isa Pre : false
+precall_to_1(expr) = (is_precall(expr) ? 1.0 : expr)
+
 # Called internally (whether DSL-based or programmatic model creation is used).
 # Creates a sorted reactions + equations vector, also ensuring reaction is first in this vector.
 # Extracts potential species, variables, and parameters from the input (if not provided as part of

test/reactionsystem_core/events.jl

@@ -142,7 +142,7 @@ let
     rs2 = ReactionSystem(rxs, t; continuous_events = [ce], discrete_events = de, name = :rs)
     rs3 = ReactionSystem(rxs, t; continuous_events = ce, discrete_events = [de], name = :rs)
     rs4 = ReactionSystem(rxs, t; continuous_events = [ce], discrete_events = [de], name = :rs)
-    @test_broken rs1 == rs2 == rs3 == rs4
+    @test_broken rs1 == rs2 == rs3 == rs4 # https://github.com/SciML/ModelingToolkit.jl/issues/3907
 end
 
 # Checks that various various erroneous forms yield errors.
@@ -236,7 +236,7 @@ end
     rn_prog = complete(rn_prog)
 
     # Tests that approaches yield identical results.
-    @test isequal(rn_dsl, rn_prog)
+    @test isequal(rn_dsl, rn_prog)  # https://github.com/SciML/ModelingToolkit.jl/issues/3907
 
     u0 = [X => 6.0, Y => 4.5, Z => 5.5]
     tspan = (0.0, 20.0)