diff --git a/lib/OrdinaryDiffEqDefault/src/OrdinaryDiffEqDefault.jl b/lib/OrdinaryDiffEqDefault/src/OrdinaryDiffEqDefault.jl index 0a4f80b53a..d8cbd101d1 100644 --- a/lib/OrdinaryDiffEqDefault/src/OrdinaryDiffEqDefault.jl +++ b/lib/OrdinaryDiffEqDefault/src/OrdinaryDiffEqDefault.jl @@ -3,6 +3,7 @@ module OrdinaryDiffEqDefault using OrdinaryDiffEq: Vern7, Vern8, Vern9, Vern6, Tsit5, Rosenbrock23, Rodas5P, FBDF, alg_stability_size, beta2_default, beta1_default, AutoSwitchCache, ODEIntegrator, trivial_limiter!, + CompositeAlgorithm, OrdinaryDiffEqAlgorithm, OrdinaryDiffEqMutableCache, AutoAlgSwitch import OrdinaryDiffEq: is_mass_matrix_alg, default_autoswitch diff --git a/src/integrators/controllers.jl b/src/integrators/controllers.jl index 19a53eba77..2fee1c65a9 100644 --- a/src/integrators/controllers.jl +++ b/src/integrators/controllers.jl @@ -1,6 +1,10 @@ abstract type AbstractController end using OrdinaryDiffEq +@inline function accept_step_controller(integrator, ::AbstractController) + return integrator.EEst <= 1 +end + @inline function stepsize_controller!(integrator, alg) stepsize_controller!(integrator, integrator.opts.controller, alg) end @@ -23,6 +27,30 @@ reset_alg_dependent_opts!(controller::AbstractController, alg1, alg2) = nothing DiffEqBase.reinit!(integrator::ODEIntegrator, controller::AbstractController) = nothing +""" + NonAdaptiveController() + +This Controller exists to match the interface when one does not want to use a controller, +basically if you want to keep a fixed time step. +""" +struct NonAdaptiveController <: AbstractController +end + +@inline function stepsize_controller!(integrator, ::NonAdaptiveController, alg) + nothing +end + +@inline function accept_step_controller(integrator, ::NonAdaptiveController) + return true +end + +function step_accept_controller!(integrator, ::NonAdaptiveController, alg, q) + integrator.dt +end + +function step_reject_controller!(integrator, ::NonAdaptiveController, alg) +end + # Standard integral (I) step size controller """ IController() @@ -311,11 +339,6 @@ end return dt_factor end -# checks whether the controller should accept a step based on the error estimate -@inline function accept_step_controller(integrator, controller::AbstractController) - return integrator.EEst <= 1 -end - @inline function accept_step_controller(integrator, controller::PIDController) return integrator.qold >= controller.accept_safety end diff --git a/src/integrators/integrator_utils.jl b/src/integrators/integrator_utils.jl index a7bf28801d..a44eecdcd1 100644 --- a/src/integrators/integrator_utils.jl +++ b/src/integrators/integrator_utils.jl @@ -222,7 +222,7 @@ function _loopfooter!(integrator) end integrator.last_stepfail = true integrator.accept_step = false - elseif integrator.opts.adaptive + else q = stepsize_controller!(integrator, integrator.alg) integrator.isout = integrator.opts.isoutofdomain(integrator.u, integrator.p, ttmp) integrator.accept_step = (!integrator.isout && @@ -255,24 +255,6 @@ function _loopfooter!(integrator) else # Reject integrator.stats.nreject += 1 end - elseif !integrator.opts.adaptive #Not adaptive - integrator.stats.naccept += 1 - integrator.tprev = integrator.t - integrator.t = if has_tstop(integrator) - tstop = integrator.tdir * first_tstop(integrator) - if abs(ttmp - tstop) < - 100eps(float(integrator.t / oneunit(integrator.t))) * oneunit(integrator.t) - tstop - else - ttmp - end - else - ttmp - end - integrator.last_stepfail = false - integrator.accept_step = true - integrator.dtpropose = integrator.dt - handle_callbacks!(integrator) end if integrator.opts.progress && integrator.iter % integrator.opts.progress_steps == 0 t1, t2 = integrator.sol.prob.tspan diff --git a/src/solve.jl b/src/solve.jl index 8b1dc8b36c..71ef62a50f 100644 --- a/src/solve.jl +++ b/src/solve.jl @@ -363,7 +363,11 @@ function DiffEqBase.__init( end if controller === nothing - controller = default_controller(_alg, cache, qoldinit, beta1, beta2) + if adaptive + controller = default_controller(_alg, cache, qoldinit, beta1, beta2) + else + controller = NonAdaptiveController() + end end save_end_user = save_end