Added benchmark in Synapse

Author: sivasathyaseeelan

benchmarks/HybridJumps/Synapse.jmd

@@ -1894,21 +1894,22 @@ First, we define the problem for the `PDMP` package.
 
 ```julia
 function SynapseProblem(
-        xc,
-        xd,
-        t1,
-        t2,
-        events_bap,
-        bap_by_epsp,
-        glu,
-        p_synapse,
-        nu,
-        algo::T,
-        agg = nothing;
-        saveat = [],
-        save_everystep = isempty(saveat),
-        kwargs...
-) where {T <: CHV}
+    xc,
+    xd,
+    t1,
+    t2,
+    events_bap,
+    bap_by_epsp,
+    glu,
+    p_synapse,
+    nu,
+    algo::T,
+    agg = nothing
+    vr_agg = VR_FRM();
+    saveat = [],
+    save_everystep = isempty(saveat),
+    kwargs...,
+) where {T<:CHV}
     problem = PDMP.PDMPProblem(
         (dxc, xc, xd, p, t) -> F_synapse(dxc, xc, xd, p, t, events_bap, bap_by_epsp),
         (rate, xc, xd, p, t, sum_rate) -> R_synapse(rate, xc, xd, p, t, sum_rate, glu),
@@ -2066,22 +2067,23 @@ function buildRxDependencyGraph(nu)
 end
 
 function SynapseProblem(
-        xc,
-        xd,
-        t1,
-        t2,
-        events_bap,
-        bap_by_epsp,
-        glu,
-        p_synapse,
-        nu,
-        algo,
-        agg;
-        jumps = nothing,
-        save_positions = (false, true),
-        saveat = [],
-        save_everystep = isempty(saveat),
-        kwargs...
+    xc,
+    xd,
+    t1,
+    t2,
+    events_bap,
+    bap_by_epsp,
+    glu,
+    p_synapse,
+    nu,
+    algo,
+    agg
+    vr_agg;
+    jumps = nothing,
+    save_positions = (false, true),
+    saveat = [],
+    save_everystep = isempty(saveat),
+    kwargs...,
 )
     p = (
         xd0 = copy(xd),
@@ -2107,6 +2109,7 @@ function SynapseProblem(
         oprob,
         agg,
         jumps;
+        vr_aggregator = vr_agg,
         dep_graph,
         save_positions,
         saveat,
@@ -2125,22 +2128,23 @@ We define functions to run the evolution of the whole synapse. The evolution inc
 ```julia
 # adapted from SynapseElife/src/SynapseModel.jl
 function evolveSynapse(
-        xc0::Vector{T},
-        xd0,
-        p_synapse::SynapseParams,
-        events_sorted_times,
-        is_pre_or_post_event,
-        bap_by_epsp,
-        is_glu_released,
-        nu,
-        algos,
-        agg = nothing;
-        progress = false,
-        abstol = 1e-8,
-        reltol = 1e-7,
-        save_positions = (false, true),
-        saveat = [],
-        kwargs...
+    xc0::Vector{T},
+    xd0,
+    p_synapse::SynapseParams,
+    events_sorted_times,
+    is_pre_or_post_event,
+    bap_by_epsp,
+    is_glu_released,
+    nu,
+    algos,
+    agg = nothing
+    vr_agg = VR_FRM();
+    progress = false,
+    abstol = 1e-8,
+    reltol = 1e-7,
+    save_positions = (false, true),
+    saveat = [],
+    kwargs...,
 ) where {T}
     tt, XC,
     XD = evolveSynapse_noformat(
@@ -2153,7 +2157,8 @@ function evolveSynapse(
         is_glu_released,
         nu,
         algos,
-        agg;
+        agg
+        vr_agg;
         progress,
         abstol,
         reltol,
@@ -2167,22 +2172,23 @@ end
 
 # adapted from SynapseElife/src/SynapseModel.jl
 function evolveSynapse_noformat(
-        xc0::Vector{T},
-        xd0,
-        p_synapse::SynapseParams,
-        events_sorted_times,
-        is_pre_or_post_event,
-        bap_by_epsp,
-        is_glu_released,
-        nu,
-        algos,
-        agg = nothing;
-        progress = false,
-        abstol = 1e-8,
-        reltol = 1e-7,
-        save_positions = (false, true),
-        saveat = [],
-        kwargs...
+    xc0::Vector{T},
+    xd0,
+    p_synapse::SynapseParams,
+    events_sorted_times,
+    is_pre_or_post_event,
+    bap_by_epsp,
+    is_glu_released,
+    nu,
+    algos,
+    agg = nothing
+    vr_agg = VR_FRM();
+    progress = false,
+    abstol = 1e-8,
+    reltol = 1e-7,
+    save_positions = (false, true),
+    saveat = [],
+    kwargs...,
 ) where {T}
     if save_positions isa Tuple{Bool, Bool}
         save_positionsON = save_positions
@@ -2209,53 +2215,50 @@ function evolveSynapse_noformat(
     events_bap = events_sorted_times[is_pre_or_post_event .== false]
 
     # function to simulate the synapse when Glutamate is ON
-    SimGluON = (xc,
-        xd,
-        t1,
-        t2,
-        glu) -> SynapseProblem(
-        xc,
-        xd,
-        t1,
-        t2,
-        events_bap,
-        bap_by_epsp,
-        glu,
-        p_synapse,
-        nu,
-        algos[1],
-        agg;
-        jumps,
-        reltol,
-        abstol,
-        saveat,
-        save_positions = save_positionsON,
-        kwargs...
-    )
+    SimGluON =
+        (xc, xd, t1, t2, glu) -> SynapseProblem(
+            xc,
+            xd,
+            t1,
+            t2,
+            events_bap,
+            bap_by_epsp,
+            glu,
+            p_synapse,
+            nu,
+            algos[1],
+            agg
+            vr_agg;
+            jumps,
+            reltol,
+            abstol,
+            saveat,
+            save_positions = save_positionsON,
+            kwargs...,
+        )
 
     # function to simulate the synapse when Glutamate is OFF
-    SimGluOFF = (xc,
-        xd,
-        t1,
-        t2) -> SynapseProblem(
-        xc,
-        xd,
-        t1,
-        t2,
-        events_bap,
-        bap_by_epsp,
-        zero(T),
-        p_synapse,
-        nu,
-        algos[2],
-        agg;
-        jumps,
-        reltol,
-        abstol,
-        saveat,
-        save_positions = save_positionsOFF,
-        kwargs...
-    )
+    SimGluOFF =
+        (xc, xd, t1, t2) -> SynapseProblem(
+            xc,
+            xd,
+            t1,
+            t2,
+            events_bap,
+            bap_by_epsp,
+            zero(T),
+            p_synapse,
+            nu,
+            algos[2],
+            agg
+            vr_agg;
+            jumps,
+            reltol,
+            abstol,
+            saveat,
+            save_positions = save_positionsOFF,
+            kwargs...,
+        )
 
     # random variable for Glutamate concentration
     gluDist = Gamma(1 / p_synapse.glu_cv^2, p_synapse.glu_cv^2)
@@ -2368,15 +2371,38 @@ const algorithms = (
     (
         label = "PDMP",
         agg = nothing,
+        vr_agg = VR_FRM(),
         solver = (CHV(solver), CHV(solver)),
-        saveat = []
+        saveat = [],
     ),
     (
         label = "Coevolve",
         agg = Coevolve(),
+        vr_agg = VR_FRM(),
         solver = (solver, solver),
-        saveat = 1 / p_synapse.sampling_rate
-    )
+        saveat = 1 / p_synapse.sampling_rate,
+    ),
+    (
+        label = "VR_Direct",
+        agg = Coevolve(),
+        vr_agg = VR_Direct(),
+        solver = (solver, solver),
+        saveat = 1 / p_synapse.sampling_rate,
+    ),
+    (
+        label = "VR_DirectFW",
+        agg = Coevolve(),
+        vr_agg = VR_DirectFW(),
+        solver = (solver, solver),
+        saveat = 1 / p_synapse.sampling_rate,
+    ),
+    (
+        label = "VR_FRM",
+        agg = Coevolve(),
+        vr_agg = VR_FRM(),
+        solver = (solver, solver),
+        saveat = 1 / p_synapse.sampling_rate,
+    ),
 );
 ```
 
@@ -2398,7 +2424,8 @@ for algo in algorithms
             [true],
             nu,
             algo.solver,
-            algo.agg;
+            algo.agg
+            alg.vr_agg;
             save_positions = (false, true),
             saveat = algo.saveat,
             save_everystep = false
@@ -2443,7 +2470,8 @@ for algo in algorithms
                 [true],
                 nu,
                 $(algo).solver,
-                $(algo).agg;
+                $(algo).agg
+                $(algo).vr_agg;
                 save_positions = (false, true),
                 saveat = $(algo).saveat,
                 save_everystep = false