change to get_variables!

src/reaction.jl

@@ -344,7 +344,7 @@ MT.is_diff_equation(rx::Reaction) = false
 MT.is_alg_equation(rx::Reaction) = false
 
 """
-    get_variables(rx::Reaction)
+    get_symbolics(set, rx::Reaction)
 
 Returns all symbolic variables that are part of a reaction. This includes all variables 
 encountered in:
@@ -353,19 +353,31 @@ encountered in:
     - Among stoichiometries.
     - Among potential noise scaling metadata.
 """
-function ModelingToolkit.get_variables(rx::Reaction)
-    sym_vars = get_variables(rx.rate)
-    sym_vars = unique!([sym_vars; rx.substrates; rx.products])
-    for stoich in rx.substoich
-        sym_vars = unique!([sym_vars; get_variables(stoich)])
-    end
-    for stoich in rx.prodstoich
-        sym_vars = unique!([sym_vars; get_variables(stoich)])
+function get_symbolics(rx::Reaction)
+    return ModelingToolkit.get_variables!([], rx)
+end
+
+"""
+    get_variables!(set, rx::Reaction)
+
+Adds all symbolic variables that are part of a reaction to set. This includes all variables 
+encountered in:
+    - Rates.
+    - Among substrates and products.
+    - Among stoichiometries.
+    - Among potential noise scaling metadata.
+"""
+function ModelingToolkit.get_variables!(set, rx::Reaction)
+    get_variables!(set, rx.rate)
+    append!(set, rx.substrates)
+    append!(set, rx.products)
+    for stoichs in (rx.substoich, rx.prodstoich), stoich in stoichs
+        get_variables!(set, stoich)
     end
     if has_noise_scaling(rx)
-        sym_vars = unique!([sym_vars; get_variables(get_noise_scaling(rx))])
+        get_variables!(set, get_noise_scaling(rx))
     end
-    return sym_vars
+    return unique!(set)
 end
 
 ### Dependency-related Functions ###

test/reactionsystem_core/reaction.jl

@@ -2,6 +2,8 @@
 
 # Fetch packages.
 using Catalyst, Test
+using Catalyst: get_symbolics
+using ModelingToolkit: value, get_variables!
 
 # Sets the default `t` to use.
 t = default_t()
@@ -11,7 +13,7 @@ t = default_t()
 # Tests the `get_variables` function.
 let
     # Declare symbolic variables.
-    @parameters k1 k2 n1 n2 η1 η2
+    @parameters k1 k2 n1 n2 η1 η2 p
     @species X(t) Y(t) Z(t)
     @variables A(t)
 
@@ -21,11 +23,17 @@ let
     rx3 = Reaction(k1 + k2 + A, [X], [X, Y, Z], [1], [n1 + n2, 2, 1])
     rx4 = Reaction(X + t, [], [Y]; metadata = [:noise_scaling => η1 + η2])
 
-    # Test `get_variables`.
-    @test issetequal(get_variables(rx1), [k1, X])
-    @test issetequal(get_variables(rx2), [k1, k2, X, Y, n1, η1])
-    @test issetequal(get_variables(rx3), [k1, k2, A, X, Y, Z, n1, n2])
-    @test issetequal(get_variables(rx4), [X, t, Y, η1, η2])
+    # Test `get_variables!`.
+    @test issetequal(get_variables!([value(p)], rx1), [k1, X, p])
+    @test issetequal(get_variables!([value(p)], rx2), [k1, k2, X, Y, n1, η1, p])
+    @test issetequal(get_variables!([value(p)], rx3), [k1, k2, A, X, Y, Z, n1, n2, p])
+    @test issetequal(get_variables!([value(p)], rx4), [X, t, Y, η1, η2, p])
+
+    # Test `get_symbolics`.
+    @test issetequal(get_symbolics(rx1), [k1, X])
+    @test issetequal(get_symbolics(rx2), [k1, k2, X, Y, n1, η1])
+    @test issetequal(get_symbolics(rx3), [k1, k2, A, X, Y, Z, n1, n2])
+    @test issetequal(get_symbolics(rx4), [X, t, Y, η1, η2])
 end
 
 ### Tests Metadata ###