Merge pull request #270 from SciML/myb/maxiters

Handle maxiters properly

Author: YingboMa

Project.toml

@@ -1,7 +1,7 @@
 name = "Sundials"
 uuid = "c3572dad-4567-51f8-b174-8c6c989267f4"
 authors = ["Chris Rackauckas <[email protected]>"]
-version = "4.2.3"
+version = "4.2.4"
 
 [deps]
 CEnum = "fa961155-64e5-5f13-b03f-caf6b980ea82"

src/common_interface/integrator_types.jl

@@ -16,6 +16,7 @@ mutable struct DEOptions{SType,TstopType,CType,reltolType,abstolType,F5}
     progress::Bool
     progress_name::String
     progress_message::F5
+    maxiters::Int
 end
 
 abstract type AbstractSundialsIntegrator{algType} <: DiffEqBase.AbstractODEIntegrator{algType,true,Vector{Float64},Float64} end

src/common_interface/solve.jl

@@ -315,7 +315,7 @@ function DiffEqBase.__init(
     opts = DEOptions(saveat_internal,tstops_internal,save_everystep,dense,
                      timeseries_errors,dense_errors,save_on,save_end,
                      callbacks_internal,abstol,reltol,verbose,advance_to_tstop,stop_at_next_tstop,
-                     progress,progress_name,progress_message)
+                     progress,progress_name,progress_message,maxiters)
     integrator = CVODEIntegrator(u0,prob.p,t0,t0,mem,_LS,_A,sol,alg,f!,userfun,jac,opts,
                        tout,tdir,sizeu,false,tmp,uprev,Cint(flag),false,0,1,callback_cache,0.)
 
@@ -744,7 +744,7 @@ function DiffEqBase.__init(
     opts = DEOptions(saveat_internal,tstops_internal,save_everystep,dense,
                      timeseries_errors,dense_errors,save_on,save_end,
                      callbacks_internal,abstol,reltol,verbose,advance_to_tstop,stop_at_next_tstop,
-                     progress,progress_name,progress_message)
+                     progress,progress_name,progress_message,maxiters)
     integrator = ARKODEIntegrator(utmp,prob.p,t0,t0,mem,_LS,_A,_MLS,_M,sol,alg,f!,userfun,jac,opts,
                        tout,tdir,sizeu,false,tmp,uprev,Cint(flag),false,0,1,callback_cache,0.)
 
@@ -1094,7 +1094,7 @@ function DiffEqBase.__init(
     opts = DEOptions(saveat_internal,tstops_internal,save_everystep,dense,
                     timeseries_errors,dense_errors,save_on,save_end,
                     callbacks_internal,abstol,reltol,verbose,advance_to_tstop,stop_at_next_tstop,
-                    progress,progress_name,progress_message)
+                    progress,progress_name,progress_message,maxiters)
 
     integrator = IDAIntegrator(utmp,dutmp,prob.p,t0,t0,mem,_LS,_A,sol,alg,f!,userfun,jac,opts,
                    tout,tdir,sizeu,sizedu,false,tmp,uprev,Cint(flag),false,0,1,callback_cache,0.)
@@ -1155,8 +1155,19 @@ function set_stop_time(integrator::IDAIntegrator,tstop)
     IDASetStopTime(integrator.mem,tstop)
 end
 
+function get_iters!(integrator::CVODEIntegrator,iters)
+    CVodeGetNumSteps(integrator.mem,iters)
+end
+function get_iters!(integrator::ARKODEIntegrator,iters)
+    ARKStepGetNumSteps(integrator.mem,iters)
+end
+function get_iters!(integrator::IDAIntegrator,iters)
+    IDAGetNumSteps(integrator.mem,iters)
+end
+
 function DiffEqBase.solve!(integrator::AbstractSundialsIntegrator)
     uType = eltype(integrator.sol.u)
+    iters = Ref(-1)
     while !isempty(integrator.opts.tstops)
         # Sundials can have floating point issues approaching a tstop if
         # there is a modifying event each
@@ -1176,8 +1187,13 @@ function DiffEqBase.solve!(integrator::AbstractSundialsIntegrator)
             if isempty(integrator.opts.tstops)
               break
             end
+            get_iters!(integrator, iters)
+            if iters[] + 1 > integrator.opts.maxiters
+                integrator.flag = -1 # MaxIters: -1
+                break
+            end
         end
-        (integrator.flag < 0) && break
+        integrator.flag < 0 && break
         handle_tstop!(integrator)
     end
 
@@ -1209,11 +1225,11 @@ function DiffEqBase.solve!(integrator::AbstractSundialsIntegrator)
         dense_errors=integrator.opts.dense_errors)
     end
 
-    if integrator.sol.retcode != :Default
-      return integrator.sol
+    if integrator.sol.retcode === :Default
+        integrator.sol = DiffEqBase.solution_new_retcode(integrator.sol,interpret_sundials_retcode(integrator.flag))
     end
-    integrator.sol = DiffEqBase.solution_new_retcode(integrator.sol,interpret_sundials_retcode(integrator.flag))
-    nothing
+
+    return integrator.sol
 end
 
 function handle_tstop!(integrator::AbstractSundialsIntegrator)

test/common_interface/cvode.jl

@@ -107,3 +107,6 @@ u0 = 0.0;
 tspan = (0.0, 10)
 prob = ODEProblem(number_test,u0,tspan,(2.0,0.01))
 sol = solve(prob,CVODE_BDF())
+
+sol = solve(prob2,CVODE_BDF(), maxiters=1)
+@test length(sol.t) == 2