bug fix: missing V=0 upcrossing

in some cases, the maximum V is below 0mV, so we now detect upcrossings at V=-20mV

Author: mjowen

mcmcSetup.jl

@@ -26,7 +26,7 @@ Get the period of the limit cycle.
 """
 function get_period(lc::Vector{Float64}, prob::ODEProblem)::Number
     # Use callbacks to find state at start and end of a period (using upcrossings of V=0mV)
-    condition(u, _, _) = u[1]
+    condition(u, _, _) = u[1]+20
     NUM_TIMES_EFFECT_HIT::Int = 0
     function affect!(integrator)
         NUM_TIMES_EFFECT_HIT += 1 
@@ -103,7 +103,7 @@ function optimiseParameters()
     function model_simulator(p)
         prob = remake(prob, p=Tools.param_map(p))::ODEProblem
         # Converge
-        condition(u, _, _) = u[1]
+        condition(u, _, _) = u[1]+20
         STATE::Vector{Float64} = zeros(size(Model.ic))
         function affect!(integrator)
             error = STATE .- integrator.u

mhmcmc.jl

@@ -278,7 +278,7 @@ Solve the ODE until convergence starting from the default initial conditions.
 - `lc`: The converged limit cycle.
 """
 function odeSolverStandard(x, prob::ODEProblem, lc, xlc, paramMap::Function, verbose=1::Integer)
-    condition(u, _, _) = u[1]
+    condition(u, _, _) = u[1]+20
 	STATE::Vector{Float64} = zeros(size(Model.ic))
 	function affect!(integrator)
 		error = STATE .- integrator.u
@@ -310,7 +310,7 @@ Solve the ODE until convergence but starting from the previous limit cycle.
 - `lc`: The converged limit cycle.
 """
 function odeSolverTracking(x, prob::ODEProblem, lc, _, paramMap::Function, verbose::Integer)
-    condition(u, _, _) = u[1]
+    condition(u, _, _) = u[1]+20
 	STATE::Vector{Float64} = zeros(size(Model.ic))
 	function affect!(integrator)
 		error = STATE .- integrator.u

simulationTimings.jl

@@ -31,7 +31,7 @@ prob = ODEProblem(Model.ode!, Model.ic, (0.0, 10000.0), tmp, abstol=1e-10, relto
 # ODE Convergence - Standard
 params = [pSmall, pLarge]
 for i in eachindex(params)
-	condition(u, _, _) = u[1]
+	condition(u, _, _) = u[1]+20
 	STATE::Vector{Float64} = zeros(size(Model.ic))
 	function affect!(integrator)
 		error = STATE .- integrator.u
@@ -52,7 +52,7 @@ end
 
 # ODE Convergence - Tracking
 for i in eachindex(params)
-	condition(u, _, _) = u[1]
+	condition(u, _, _) = u[1]+20
 	STATE::Vector{Float64} = zeros(size(Model.ic))
 	function affect!(integrator)
 		error = STATE .- integrator.u

tools.jl

@@ -7,7 +7,7 @@ using .Model
 # Check all variables are converged automatically
 function auto_converge_check(prob::ODEProblem, ic, p::NamedTuple)::Bool
     # Use callbacks to find state at start and end of period (using upcrossings of V=0mV)
-    condition(u, _, _) = u[1]
+    condition(u, _, _) = u[1]+20
     NUM_TIMES_EFFECT_HIT::Int = 0
     function affect!(integrator)
         NUM_TIMES_EFFECT_HIT += 1 
@@ -19,7 +19,14 @@ function auto_converge_check(prob::ODEProblem, ic, p::NamedTuple)::Bool
     save_positions = (true, false))
     sol = solve(prob, Tsit5(), u0=ic, p=p, tspan=(0.0, 10.0), maxiters=1e9, 
     save_everystep=false, save_start=false, save_end=false, callback=cb)
-    error = sol[end] - sol[1]
+    if length(sol) < 2
+        # Run simulation and report max and min voltage
+        sol = solve(prob, Tsit5(), u0=ic, p=p, tspan=(0.0, 10.0), maxiters=1e9)
+        println("Max V: ", maximum(sol[1,:]))
+        println("Min V: ", minimum(sol[1,:]))
+    else
+        error = sol[end] - sol[1]
+    end
     return sum(abs.(error))<1e-6
 end
 
@@ -39,7 +46,7 @@ Align the limit cycle in the solution to start at the max of V and fixes the tim
 """
 function aligned_sol(lc, prob::ODEProblem, period::Number; save_only_V::Bool = true)
     # Use callbacks to find state at start of period (using upcrossings of V=0mV)
-    condition(u, _, _) = u[1]
+    condition(u, _, _) = u[1]+20
     function affect!(integrator)
         terminate!(integrator)
     end