Skip to content
Closed
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
29 changes: 22 additions & 7 deletions test/saveat_regression.jl
Original file line number Diff line number Diff line change
@@ -1,6 +1,5 @@
using DiffEqBase, JumpProcesses, Test
using StableRNGs
rng = StableRNG(12345)

rate_consts = [10.0]
reactant_stoich = [[1 => 1, 2 => 1]]
Expand All @@ -10,12 +9,26 @@ maj = MassActionJump(rate_consts, reactant_stoich, net_stoich)
n0 = [1, 1, 0]
tspan = (0, 0.2)
dprob = DiscreteProblem(n0, tspan)
jprob = JumpProblem(dprob, Direct(), maj, save_positions = (false, false), rng = rng)
ts = collect(0:0.002:tspan[2])
NA = zeros(length(ts))
Nsims = 10_000
sol = JumpProcesses.solve(EnsembleProblem(jprob), SSAStepper(), saveat = ts,
trajectories = Nsims)

# Give every trajectory its own independent, seeded RNG. EnsembleProblems default to
# EnsembleThreads(), so sharing one mutable RNG across trajectories is a data race that
# corrupts the RNG state. A fresh StableRNG per trajectory keeps the solve race-free and
# reproducible.
function jprob_func(save_positions)
function (prob, ctx)
rng = StableRNG(12345 + ctx.sim_id)
JumpProblem(dprob, Direct(), maj; save_positions, rng)
end
end

jprob = JumpProblem(dprob, Direct(), maj, save_positions = (false, false),
rng = StableRNG(12345))
NA = zeros(length(ts))
sol = JumpProcesses.solve(
EnsembleProblem(jprob; prob_func = jprob_func((false, false))), SSAStepper(),
saveat = ts, trajectories = Nsims)

for i in 1:length(sol.u)
NA .+= sol.u[i][1, :]
Expand All @@ -26,10 +39,12 @@ for i in 1:length(ts)
end

NA = zeros(length(ts))
jprob = JumpProblem(dprob, Direct(), maj; rng = rng)
jprob = JumpProblem(dprob, Direct(), maj; rng = StableRNG(12345))
sol = nothing;
GC.gc();
sol = JumpProcesses.solve(EnsembleProblem(jprob), SSAStepper(), trajectories = Nsims)
sol = JumpProcesses.solve(
EnsembleProblem(jprob; prob_func = jprob_func((false, true))), SSAStepper(),
trajectories = Nsims)

for i in 1:Nsims
for n in 1:length(ts)
Expand Down
Loading