fix more old sde sims without seeds

Author: TorkelE

test/dsl/dsl_model_construction.jl

@@ -124,32 +124,48 @@ let
 
     for networks in identical_networks_1
         for factor in [1e-2, 1e-1, 1e0, 1e1, 1e2, 1e3]
-            u0_1 = rnd_u0(networks[1], rng; factor)
-            p_1 = rnd_ps(networks[1], rng; factor)
-            u0_2 = Pair.(unknowns(networks[2]), last.(u0_1))
-            p_2 = Pair.(parameters(networks[2]), last.(p_1))
+            u0 = rnd_u0(networks[1], rng; factor)
+            p = rnd_ps(networks[1], rng; factor)
             t = rand(rng)
 
-            @test f_eval(networks[1], u0_1, p_1, t) ≈ f_eval(networks[2], u0_2, p_2, t)
-            @test jac_eval(networks[1], u0_1, p_1, t) ≈ jac_eval(networks[2], u0_2, p_2, t)
-            @test g_eval(networks[1], u0_1, p_1, t) ≈ g_eval(networks[2], u0_2, p_2, t)
+            @test f_eval(networks[1], u0, p, t) ≈ f_eval(networks[2], u0, p, t)
+            @test jac_eval(networks[1], u0, p, t) ≈ jac_eval(networks[2], u0, p, t)
+            @test g_eval(networks[1], u0, p, t) ≈ g_eval(networks[2], u0, p, t)
         end
     end
 end
 
-# Compares networks to networks written in different ways.
+# Compares simulations for network with different species and parameter names
 let
-    identical_networks_2 = Vector{Pair}()
-
-    # Different parameter and variable names.
+    # Fetches the original network, and also declares it using alternative notation.
+    network = reaction_networks_standard[5]
     differently_written_5 = @reaction_network begin
         q, ∅ → Y1
         (l1, l2), Y1 ⟷ Y2
         (l3, l4), Y2 ⟷ Y3
         (l5, l6), Y3 ⟷ Y4
         c, Y4 → ∅
+    end    
+
+    # Checks that the networks' functions evaluates equally for various randomised inputs.
+    @unpack X1, X2, X3, X4, p, d, k1, k2, k3, k4, k5, k6 = network
+    for factor in [1e-2, 1e-1, 1e0, 1e1, 1e2, 1e3]
+        u0_1 = Dict(rnd_u0(network, rng; factor))
+        p_1 = Dict(rnd_ps(network, rng; factor))
+        u0_2 = [:Y1 => u0_1[X1], :Y2 => u0_1[X2], :Y3 => u0_1[X3], :Y4 => u0_1[X4]]
+        p_2 = [:q => p_1[p], :c => p_1[d], :l1 => p_1[k1], :l2 => p_1[k2], :l3 => p_1[k3], 
+               :l4 => p_1[k4], :l5 => p_1[k5], :l6 => p_1[k6]]
+        t = rand(rng)
+        
+        @test f_eval(network, u0_1, p_1, t) ≈ f_eval(differently_written_5, u0_2, p_2, t)
+        @test jac_eval(network, u0_1, p_1, t) ≈ jac_eval(differently_written_5, u0_2, p_2, t)
+        @test g_eval(network, u0_1, p_1, t) ≈ g_eval(differently_written_5, u0_2, p_2, t)
     end
-    push!(identical_networks_2, reaction_networks_standard[5] => differently_written_5)
+end
+
+# Compares networks to networks written in different ways.
+let
+    identical_networks_2 = Vector{Pair}()
 
     # Unfold reactions.
     differently_written_6 = @reaction_network begin
@@ -190,15 +206,13 @@ let
 
     for networks in identical_networks_2
         for factor in [1e-2, 1e-1, 1e0, 1e1, 1e2, 1e3]
-            u0_1 = rnd_u0(networks[1], rng; factor)
-            p_1 = rnd_ps(networks[1], rng; factor)
-            u0_2 = Pair.(unknowns(networks[2]), last.(u0_1))
-            p_2 = Pair.(parameters(networks[2]), last.(p_1))
+            u0 = rnd_u0(networks[1], rng; factor)
+            p = rnd_ps(networks[1], rng; factor)
             t = rand(rng)
 
-            @test f_eval(networks[1], u0_1, p_1, t) ≈ f_eval(networks[2], u0_2, p_2, t)
-            @test jac_eval(networks[1], u0_1, p_1, t) ≈ jac_eval(networks[2], u0_2, p_2, t)
-            @test g_eval(networks[1], u0_1, p_1, t) ≈ g_eval(networks[2], u0_2, p_2, t)
+            @test f_eval(networks[1], u0, p, t) ≈ f_eval(networks[2], u0, p, t)
+            @test jac_eval(networks[1], u0, p, t) ≈ jac_eval(networks[2], u0, p, t)
+            @test g_eval(networks[1], u0, p, t) ≈ g_eval(networks[2], u0, p, t)
         end
     end
 end

test/model_simulation/simulate_SDEs.jl

@@ -131,8 +131,9 @@ let
             if nameof(rn_catalyst) != :rnc9
                 sprob_1 = SDEProblem(rn_catalyst, u0_1, (0.0, 1.0), ps_1)
                 sprob_2 = SDEProblem(rn_manual.f, rn_manual.g, u0_2, (0.0, 1.0), ps_2; noise_rate_prototype = rn_manual.nrp)
-                sol1 = solve(sprob_1, ImplicitEM(); seed = rand(rng, 1:100))
-                sol2 = solve(sprob_2, ImplicitEM(); seed = rand(rng, 1:100))
+                seed =  seed = rand(rng, 1:100)
+                sol1 = solve(sprob_1, ImplicitEM(); seed)
+                sol2 = solve(sprob_2, ImplicitEM(); seed)
                 @test sol1[u0_sym] ≈ sol2.u
             end
         end

test/model_simulation/simulate_jumps.jl

@@ -111,17 +111,18 @@ let
     # Loops through all cases, checks that identical simulations are generated with/without Catalyst.
     for (rn_catalyst, rn_manual, u0_sym, ps_sym) in zip(catalyst_networks, manual_networks, u0_syms, ps_syms)
         for factor in [5, 50]
+            seed = rand(rng, 1:100)
             u0_1 = rnd_u0_Int64(rn_catalyst, rng; n = factor)
             ps_1 = rnd_ps(rn_catalyst, rng; factor = factor/100.0)
             dprob_1 = DiscreteProblem(rn_catalyst, u0_1, (0.0, 100.0), ps_1)
             jprob_1 = JumpProblem(rn_catalyst, dprob_1, Direct(); rng)
-            sol1 = solve(jprob_1, SSAStepper(); saveat = 1.0)
+            sol1 = solve(jprob_1, SSAStepper(); seed, saveat = 1.0)
 
             u0_2 = map_to_vec(u0_1, u0_sym)
             ps_2 = map_to_vec(ps_1, ps_sym)
             dprob_2 = DiscreteProblem(u0_2, (0.0, 100.0), ps_2)
             jprob_2 = JumpProblem(dprob_2, Direct(), rn_manual...; rng)
-            sol2 = solve(jprob_2, SSAStepper(); saveat = 1.0)
+            sol2 = solve(jprob_2, SSAStepper(); seed, saveat = 1.0)
 
             if nameof(rn_catalyst) == :rnh7
                 # Have spent a few hours figuring this one out. For certain seeds it actually works,