updates test

Author: TorkelE

src/expression_utils.jl

@@ -94,18 +94,6 @@ function is_escaped_expr(expr)
 end
 
 
-### Generic Expression Handling ###
-
-# Convert an expression that is a vector with symbols that have values assigned using `=` 
-# (e.g. :([a=1.0, b=2.0])) to a vector where the assignment instead uses symbols and pairs 
-# (e.g. :([a=>1.0, b=>2.0])). Used to e.g. convert default reaction metadata to a form that can be 
-# evaluated as actual code.
-function expr_equal_vector_to_pairs(expr_vec)
-    pair_vector = :([])
-    foreach(arg -> push!(pair_vector.args, arg.args[1] => arg.args[2]), expr_vec.args) 
-    return pair_vector
-end
-
 ### Old Stuff ###
 
 #This will be called whenever a function stored in funcdict is called.

src/reaction_network.jl

@@ -360,7 +360,6 @@ function make_reaction_system(ex::Expr; name = :(gensym(:ReactionSystem)))
     check_default_noise_scaling!(default_reaction_metadata, options)
     compound_expr, compound_species = read_compound_options(options)
     check_default_noise_scaling!(default_reaction_metadata, options)
-    default_reaction_metadata = expr_equal_vector_to_pairs(default_reaction_metadata)
 
     # Parses reactions, species, and parameters.
     reactions = get_reactions(reaction_lines)
@@ -806,7 +805,7 @@ function check_default_noise_scaling!(default_reaction_metadata, options)
         if (length(options[:default_noise_scaling].args) != 3) # Becasue of how expressions are, 3 is used.
             error("@default_noise_scaling should only have a single input, this appears not to be the case: \"$(options[:default_noise_scaling])\"")
         end
-        push!(default_reaction_metadata.args, :(noise_scaling=$(options[:default_noise_scaling].args[3])))
+        push!(default_reaction_metadata.args, :(:noise_scaling => $(options[:default_noise_scaling].args[3])))
     end
 end
 

src/reactionsystem.jl

@@ -1538,7 +1538,7 @@ function Base.convert(::Type{<:NonlinearSystem}, rs::ReactionSystem; name = name
                       kwargs...)
     iscomplete(rs) || error(COMPLETENESS_ERROR)
     spatial_convert_err(rs::ReactionSystem, NonlinearSystem)
-    fullrs = Catalyst.flatten(rs; complete = true)
+    fullrs = Catalyst.flatten(rs)
     remove_conserved && conservationlaws(fullrs)
     ists, ispcs = get_indep_sts(fullrs, remove_conserved)
     eqs = assemble_drift(fullrs, ispcs; combinatoric_ratelaws, remove_conserved,
@@ -1582,7 +1582,7 @@ function Base.convert(::Type{<:SDESystem}, rs::ReactionSystem;
     iscomplete(rs) || error(COMPLETENESS_ERROR)
     spatial_convert_err(rs::ReactionSystem, SDESystem)
 
-    flatrs = Catalyst.flatten(rs; complete = true)
+    flatrs = Catalyst.flatten(rs)
     error_if_constraints(SDESystem, flatrs)
 
     remove_conserved && conservationlaws(flatrs)
@@ -1636,7 +1636,7 @@ function Base.convert(::Type{<:JumpSystem}, rs::ReactionSystem; name = nameof(rs
     (remove_conserved !== nothing) &&
         error("Catalyst does not support removing conserved species when converting to JumpSystems.")
 
-    flatrs = Catalyst.flatten(rs; complete = true)
+    flatrs = Catalyst.flatten(rs)
     error_if_constraints(JumpSystem, flatrs)
 
     (length(MT.continuous_events(flatrs)) > 0) &&
@@ -1821,7 +1821,7 @@ function MT.flatten(rs::ReactionSystem; name = nameof(rs))
                    balanced_bc_check = false,
                    spatial_ivs = get_sivs(rs),
                    continuous_events = MT.continuous_events(rs),
-                   discrete_events = MT.discre)
+                   discrete_events = MT.continuous_events(rs))
 end
 
 """

test/dsl/custom_functions.jl

@@ -55,7 +55,7 @@ let
         mm(X, v, K), 0 --> X2
         mm(X, v, K), 0 --> X3
     end
-    f_mm = ODEFunction(convert(ODESystem, mm_network), jac = true)
+    f_mm = ODEFunction(complete(convert(ODESystem, mm_network)), jac = true)
 
     u0 = 10 * rand(rng, length(get_unknowns(mm_network)))
     p = 10 * rand(rng, length(get_ps(mm_network)))
@@ -75,7 +75,7 @@ let
         mmr(X, v, K), 0 --> X2
         mmr(X, v, K), 0 --> X3
     end
-    f_mmr = ODEFunction(convert(ODESystem, mmr_network), jac = true)
+    f_mmr = ODEFunction(complete(convert(ODESystem, mmr_network)), jac = true)
 
     u0 = 10 * rand(rng, length(get_unknowns(mmr_network)))
     p = 10 * rand(rng, length(get_ps(mmr_network)))
@@ -94,7 +94,7 @@ let
         hill(X, v, K, 2), 0 --> X1
         hill(X, v, K, 2), 0 --> X2
     end
-    f_hill = ODEFunction(convert(ODESystem, hill_network), jac = true)
+    f_hill = ODEFunction(complete(convert(ODESystem, hill_network)), jac = true)
 
     u0 = 10 * rand(rng, length(get_unknowns(hill_network)))
     p = 10 * rand(rng, length(get_ps(hill_network)))
@@ -113,7 +113,7 @@ let
         hillr(X, v, K, 2), 0 --> X1
         hillr(X, v, K, 2), 0 --> X2
     end
-    f_hillr = ODEFunction(convert(ODESystem, hillr_network), jac = true)
+    f_hillr = ODEFunction(complete(convert(ODESystem, hillr_network)), jac = true)
 
     u0 = 10 * rand(rng, length(get_unknowns(hillr_network)))
     p = 10 * rand(rng, length(get_ps(hillr_network)))
@@ -134,7 +134,7 @@ let
         hillar(X, Y, v, K, 2), 0 --> X3
         hillar(X, Y, v, K, 2), 0 --> X4
     end
-    f_hillar = ODEFunction(convert(ODESystem, hillar_network), jac = true)
+    f_hillar = ODEFunction(complete(convert(ODESystem, hillar_network)), jac = true)
 
     u0 = 10 * rand(rng, length(get_unknowns(hillar_network)))
     p = 10 * rand(rng, length(get_ps(hillar_network)))

test/dsl/dsl_basics.jl

@@ -68,54 +68,53 @@ end
 ### Test Interpolation Within the DSL ###
 
 # Tests basic interpolation cases.
-let
-    # Declares parameters and species used across the test.
-    @parameters α k k1 k2
-    @species A(t) B(t) C(t) D(t)
 
-    let
-        AA = A
-        AAA = A^2 + B
-        rn = @reaction_network rn begin
-            @parameters k
-            @species A(t) B(t) C(t) D(t)
-            k*$AAA, C --> D
-        end
-        rn2 = ReactionSystem([Reaction(k*AAA, [C], [D])], t; name=:rn)
-        @test rn == rn2
+# Declares parameters and species used across the test.
+@parameters α k k1 k2
+@species A(t) B(t) C(t) D(t)
 
-        rn = @reaction_network rn begin
-            @parameters k
-            @species A(t) C(t) D(t)
-            k, $AA + C --> D
-        end
-        rn2 = ReactionSystem([Reaction(k, [AA,C], [D])], t; name=:rn)
-        @test rn == rn2
+let
+    AA = A
+    AAA = A^2 + B
+    rn = @reaction_network rn begin
+        @parameters k
+        @species A(t) B(t) C(t) D(t)
+        k*$AAA, C --> D
     end
+    rn2 = ReactionSystem([Reaction(k*AAA, [C], [D])], t; name=:rn)
+    @test rn == rn2
 
-    let
-        BB = B; A2 = A
-        rn = @reaction_network rn begin
-            @parameters k1 k2
-            (k1,k2), C + $A2 + $BB + $A2 <--> $BB + $BB
-        end
-        rn2 = ReactionSystem([Reaction(k1, [C, A, B], [B], [1,2,1],[2]),
-                            Reaction(k2, [B], [C, A, B], [2], [1,2,1])],
-                            t; name=:rn)
-        @test rn == rn2
+    rn = @reaction_network rn begin
+        @parameters k
+        @species A(t) C(t) D(t)
+        k, $AA + C --> D
     end
+    rn2 = ReactionSystem([Reaction(k, [AA,C], [D])], t; name=:rn)
+    @test rn == rn2
+end
 
-    let
-        AA = A
-        kk1 = k^2*A
-        kk2 = k1+k2
-        rn = @reaction_network rn begin
-            @parameters α k k1 k2
-            α+$kk1*$kk2*$AA, 2*$AA + B --> $AA
-        end
-        rn2 = ReactionSystem([Reaction(α+kk1*kk2*AA, [A, B], [A], [2, 1], [1])], t; name=:rn)
-        @test rn == rn2
+let
+    BB = B; A2 = A
+    rn = @reaction_network rn begin
+        @parameters k1 k2
+        (k1,k2), C + $A2 + $BB + $A2 <--> $BB + $BB
+    end
+    rn2 = ReactionSystem([Reaction(k1, [C, A, B], [B], [1,2,1],[2]),
+                        Reaction(k2, [B], [C, A, B], [2], [1,2,1])],
+                        t; name=:rn)
+    @test rn == rn2
+end
+
+let
+    AA = A
+    kk1 = k^2*A
+    kk2 = k1+k2
+    rn = @reaction_network rn begin
+        @parameters α k k1 k2
+        α+$kk1*$kk2*$AA, 2*$AA + B --> $AA
     end
+    rn2 = ReactionSystem([Reaction(α+kk1*kk2*AA, [A, B], [A], [2, 1], [1])], t; name=:rn)
+    @test rn == rn2
 end
 
 # Miscellaneous interpolation tests. Unsure what they do here (not related to DSL).
@@ -272,6 +271,7 @@ let
     rx2 = Reaction(2*k, [B], [D], [1], [2.5])
     rx3 = Reaction(2*k, [B], [D], [2.5], [2])
     @named mixedsys = ReactionSystem([rx1,rx2,rx3],t,[B,C,D],[k])
+    mixedsys = complete(mixedsys)
     osys = convert(ODESystem, mixedsys; combinatoric_ratelaws=false)
     rn = @reaction_network mixedsys begin
         @parameters k
@@ -313,11 +313,12 @@ let
     eq = D(B) ~ -B
     @named osys = ODESystem([eq], t)
     @named rn2 = extend(osys, rn2)
+    rn2 = complete(rn2)
     @test issetequal(unknowns(rn2), species(rn2))
     @test all(isspecies, species(rn))
     @test Catalyst.isbc(ModelingToolkit.value(B))
     @test Catalyst.isbc(ModelingToolkit.value(A)) == false
-    osys2 = convert(ODESystem, rn2)
+    osys2 = complete(convert(ODESystem, rn2))
     @test issetequal(unknowns(osys2), unknowns(rn2))
     @test length(equations(osys2)) == 2
 end

test/dsl/dsl_model_construction.jl

@@ -316,7 +316,7 @@ let
     end
 end
 
-# CHeck that various symbols can be used as species/parameter names.
+# Check that various symbols can be used as species/parameter names.
 let
     @reaction_network begin
         (a, A), n ⟷ N

test/miscellaneous_tests/api.jl

@@ -537,7 +537,7 @@ let
         b23, F2 --> F3
         b31, F3 --> F1
     end
-    osys = convert(ODESystem, rn; remove_conserved = true)
+    osys = complete(convert(ODESystem, rn; remove_conserved = true))
     @unpack A, B, C, D, E, F1, F2, F3, k1, k2, m1, m2, b12, b23, b31 = osys
     u0 = [A => 10.0, B => 10.0, C => 0.0, D => 10.0, E => 0.0, F1 => 8.0, F2 => 0.0,
         F3 => 0.0]
@@ -556,7 +556,7 @@ let
         @test isapprox(sol2(tv, idxs = s), sol2(tv, idxs = s))
     end
 
-    nsys = convert(NonlinearSystem, rn; remove_conserved = true)
+    nsys = complete(convert(NonlinearSystem, rn; remove_conserved = true))
     nprob = NonlinearProblem{true}(nsys, u0, p)
     nsol = solve(nprob, NewtonRaphson(); abstol = 1e-10)
     nprob2 = ODEProblem(rn, u0, (0.0, 100.0 * tspan[2]), p)
@@ -570,7 +570,7 @@ let
 
     u0 = [A => 100.0, B => 20.0, C => 5.0, D => 10.0, E => 3.0, F1 => 8.0, F2 => 2.0,
         F3 => 20.0]
-    ssys = convert(SDESystem, rn; remove_conserved = true)
+    ssys = complete(convert(SDESystem, rn; remove_conserved = true))
     sprob = SDEProblem(ssys, u0, tspan, p)
     sprob2 = SDEProblem(rn, u0, tspan, p)
     sprob3 = SDEProblem(rn, u0, tspan, p; remove_conserved = true)

test/miscellaneous_tests/events.jl

@@ -17,12 +17,13 @@ let
     discrete_events = [1.0 => [V ~ 1.0]]
     rxs = [(@reaction $V, 0 --> A), (@reaction 1.0, A --> 0)]
     @named rs = ReactionSystem([rxs; eqs], t; discrete_events)
+    rs = complete(rs)
     @test length(ModelingToolkit.discrete_events(rs)) == 1
     @test length(ModelingToolkit.continuous_events(rs)) == 0
 
     # Tests in simulation.
     setdefaults!(rs, [:A => 0.0])
-    osys = convert(ODESystem, rs)
+    osys = complete(convert(ODESystem, rs))
     @test length(ModelingToolkit.discrete_events(osys)) == 1
     oprob = ODEProblem(osys, [], (0.0, 20.0))
     sol = solve(oprob, Tsit5())
@@ -37,6 +38,7 @@ let
     rxs = [Reaction(α, nothing, [V]), Reaction(β, [V], nothing)]
     continuous_events = [V ~ 2.5] => [α ~ 0, β ~ 0]
     @named rs = ReactionSystem(rxs, t; continuous_events)
+    rs = complete(rs)
     @test length(ModelingToolkit.discrete_events(rs)) == 0
     @test length(ModelingToolkit.continuous_events(rs)) == 1
 

test/model_simulation/simulate_SDEs.jl

@@ -255,7 +255,6 @@ let
         (p,d), 0 <--> X2, ([description="Y"],[description="Y"])
         (p,d), 0 <--> X3
     end
-    noise_scaling_network = complete(noise_scaling_network)
 
     u0 = [:X1 => 1000.0, :X2 => 1000.0, :X3 => 1000.0]
     ps = [noise_scaling_network.p => 1000.0, noise_scaling_network.d => 1.0, noise_scaling_network.η1 => 1.0]
@@ -316,8 +315,8 @@ let
 
     # Checks that systems have the correct noise scaling terms.
     rn = set_default_noise_scaling(rn, 0.5)
-    rn1_noise_scaling = [get_noise_scaling(rx) for rx in get_rxs(rn)]
-    rn2_noise_scaling = [get_noise_scaling(rx) for rx in get_rxs(Catalyst.get_systems(rn)[1])]
+    rn1_noise_scaling = [get_noise_scaling(rx) for rx in Catalyst.get_rxs(rn)]
+    rn2_noise_scaling = [get_noise_scaling(rx) for rx in Catalyst.get_rxs(Catalyst.get_systems(rn)[1])]
     rn_noise_scaling = [get_noise_scaling(rx) for rx in reactions(rn)]
     @test issetequal(rn1_noise_scaling, [2.0, 0.5])
     @test issetequal(rn2_noise_scaling, [5.0, 0.5])

test/model_simulation/u0_n_parameter_inputs.jl

@@ -19,7 +19,7 @@ include("../test_networks.jl")
 # Tests various ways to input u0 and p for various functions.
 let
     test_network = reaction_networks_standard[7]
-    test_osys = convert(ODESystem, test_network)
+    test_osys = complete(convert(ODESystem, test_network))
     @parameters p1 p2 p3 k1 k2 k3 v1 K1 d1 d2 d3 d4 d5
     @species X1(t) X2(t) X3(t) X4(t) X5(t) X6(t) X(t)