From bae6e938425d902055f335022e15f1dd752814c8 Mon Sep 17 00:00:00 2001
From: Mahmoud Bentriou <mahmoud.bentriou@centralesupelec.fr>
Date: Fri, 11 Dec 2020 16:51:32 +0100
Subject: [PATCH] change of model names
---
bench/pkg/catalyst.jl | 79 +++++++++++++++++++++++++++++++++++++++----
models/ER.jl | 2 +-
models/SIR.jl | 2 +-
3 files changed, 75 insertions(+), 8 deletions(-)
diff --git a/bench/pkg/catalyst.jl b/bench/pkg/catalyst.jl
index 830c58b..332c5bb 100644
--- a/bench/pkg/catalyst.jl
+++ b/bench/pkg/catalyst.jl
@@ -4,14 +4,17 @@ using MarkovProcesses
using Catalyst
using DiffEqJump
+path_latex = "./"
+
+# Bench 1: A long simulation of ER
+
load_model("ER")
set_param!(ER, :k1, 0.2)
set_param!(ER, :k2, 40.0)
ER.buffer_size = 100
ER.estim_min_states = 8000
-b_pkg = @benchmark simulate(ER)
-@show minimum(b_pkg), mean(b_pkg), maximum(b_pkg)
+bench1_pkg = @benchmark simulate(ER)
rs = @reaction_network begin
c1, S + E --> SE
@@ -21,14 +24,78 @@ end c1 c2 c3
p = (0.2,40.0,1.0) # [c1,c2,c3]
tspan = (0., 100.)
u0 = [100., 100., 0., 0.] # [S,E,SE,P]
-
# solve JumpProblem
dprob = DiscreteProblem(rs, u0, tspan, p)
jprob = JumpProblem(rs, dprob, Direct())
jsol = solve(jprob, SSAStepper())
-b_catalyst = @benchmark solve(jprob, SSAStepper())
-@show minimum(b_catalyst), mean(b_catalyst), maximum(b_catalyst)
+bench1_catalyst = @benchmark solve(jprob, SSAStepper())
+
+str_latex_bench1 = "
+\\begin{tabular}{|c|c|c|c|c|}
+ \\hline
+ Bench 1 & Mean time (ms) & Max. time (ms) &
+Min. time (ms) & \\begin{tabular}[c]{@{}c@{}}Mean\\\\Memory (KB)\\end{tabular} \\\\
+ \\hline
+ Package & $(round(time(mean(bench1_pkg))/1E6, digits=2)) & $(round(time(maximum(bench1_pkg))/1E6, digits=2)) &
+ $(round(time(minimum(bench1_pkg))/1E6, digits=2)) & $(round(memory(mean(bench1_pkg))/1024, digits=2)) \\\\
+ \\hline
+ Catalyst.jl & $(round(time(mean(bench1_catalyst))/1E6, digits=2)) & $(round(time(maximum(bench1_catalyst))/1E6, digits=2)) &
+ $(round(time(minimum(bench1_catalyst))/1E6, digits=2)) & $(round(memory(mean(bench1_catalyst))/1024, digits=2)) \\\\
+ \\hline
+\\end{tabular}"
+
+# Bench 2: Creation of model + a long simulation of ER
+
+bench2_pkg = @benchmark begin
+ ER = @network_model begin
+ R1: (E+S => ES, k1*E*S)
+ R2: (ES => E+S, k2*ES)
+ R3: (ES => E+P, k3*ES)
+ end "ER"
+ set_param!(ER, :k1, 0.2)
+ set_param!(ER, :k2, 40.0)
+ ER.buffer_size = 100
+ ER.estim_min_states = 8000
+ simulate(ER)
+end
+
+bench2_catalyst = @benchmark begin
+ rs = @reaction_network begin
+ c1, S + E --> SE
+ c2, SE --> S + E
+ c3, SE --> P + E
+ end c1 c2 c3
+ p = (0.2,40.0,1.0) # [c1,c2,c3]
+ tspan = (0., 100.)
+ u0 = [100., 100., 0., 0.] # [S,E,SE,P]
+ # solve JumpProblem
+ dprob = DiscreteProblem(rs, u0, tspan, p)
+ jprob = JumpProblem(rs, dprob, Direct())
+ jsol = solve(jprob, SSAStepper())
+end
+
+str_latex_bench2 = "
+\\begin{tabular}{|c|c|c|c|c|}
+ \\hline
+ Bench 1 & Mean time (ms) & Max. time (ms) &
+Min. time (ms) & \\begin{tabular}[c]{@{}c@{}}Mean\\\\Memory (KB)\\end{tabular} \\\\
+ \\hline
+ Package & $(round(time(mean(bench2_pkg))/1E6, digits=2)) & $(round(time(maximum(bench2_pkg))/1E6, digits=2)) &
+ $(round(time(minimum(bench2_pkg))/1E6, digits=2)) & $(round(memory(mean(bench2_pkg))/1024, digits=2)) \\\\
+ \\hline
+ Catalyst.jl & $(round(time(mean(bench2_catalyst))/1E6, digits=2)) & $(round(time(maximum(bench2_catalyst))/1E6, digits=2)) &
+ $(round(time(minimum(bench2_catalyst))/1E6, digits=2)) & $(round(memory(mean(bench2_catalyst))/1024, digits=2)) \\\\
+ \\hline
+\\end{tabular}"
+
+@show minimum(bench1_pkg), mean(bench1_pkg), maximum(bench1_pkg)
+@show minimum(bench1_catalyst), mean(bench1_catalyst), maximum(bench1_catalyst)
-#plot(jsol,lw=2,title="Gillespie: Michaelis-Menten Enzyme Kinetics")
+open(path_latex * "bench1.tex", "w+") do io
+ write(io, str_latex_bench1)
+end;
+open(path_latex * "bench2.tex", "w+") do io
+ write(io, str_latex_bench2)
+end;
diff --git a/models/ER.jl b/models/ER.jl
index b2d0d0c..350707d 100644
--- a/models/ER.jl
+++ b/models/ER.jl
@@ -52,7 +52,7 @@ isabsorbing_ER(p::Vector{Float64},xn::Vector{Int}) =
@inbounds(p[1]*xn[1]*xn[2] + (p[2]+p[3])*xn[3] === 0.0)
g_ER = [:P]
-ER = ContinuousTimeModel(d,k,dict_var_ER,dict_p_ER,l_tr_ER,p_ER,x0_ER,t0_ER,ER_f!,isabsorbing_ER; g=g_ER,name="ER pkg")
+ER = ContinuousTimeModel(d,k,dict_var_ER,dict_p_ER,l_tr_ER,p_ER,x0_ER,t0_ER,ER_f!,isabsorbing_ER; g=g_ER,name="ER SSA pkg")
function create_ER(new_p::Vector{Float64})
ER_new = deepcopy(ER)
diff --git a/models/SIR.jl b/models/SIR.jl
index 7161169..e2bfadf 100644
--- a/models/SIR.jl
+++ b/models/SIR.jl
@@ -50,7 +50,7 @@ end
isabsorbing_SIR(p::Vector{Float64}, xn::Vector{Int}) = (p[1]*xn[1]*xn[2] + p[2]*xn[2]) === 0.0
g_SIR = [:I]
-SIR = ContinuousTimeModel(d,k,dict_var,dict_p,l_tr_SIR,p_SIR,x0_SIR,t0_SIR,SIR_f!,isabsorbing_SIR; g=g_SIR, name="SIR pkg")
+SIR = ContinuousTimeModel(d,k,dict_var,dict_p,l_tr_SIR,p_SIR,x0_SIR,t0_SIR,SIR_f!,isabsorbing_SIR; g=g_SIR, name="SIR SSA pkg")
function create_SIR(new_p::Vector{Float64})
SIR_new = deepcopy(SIR)
--
GitLab