Merge pull request #13 from Shreyas-Ekanathan/upstream

Upstream

Project.toml

@@ -1,7 +1,7 @@
 name = "OrdinaryDiffEq"
 uuid = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
 authors = ["Chris Rackauckas <[email protected]>", "Yingbo Ma <[email protected]>"]
-version = "6.87.0"
+version = "6.88.1"
 
 [deps]
 ADTypes = "47edcb42-4c32-4615-8424-f2b9edc5f35b"
@@ -98,19 +98,33 @@ MacroTools = "0.5"
 MuladdMacro = "0.2.1"
 NonlinearSolve = "3"
 OrdinaryDiffEqAdamsBashforthMoulton = "1"
+OrdinaryDiffEqBDF = "1"
 OrdinaryDiffEqCore = "1"
+OrdinaryDiffEqDefault = "1"
 OrdinaryDiffEqDifferentiation = "1"
+OrdinaryDiffEqExplicitRK = "1"
+OrdinaryDiffEqExponentialRK = "1"
 OrdinaryDiffEqExtrapolation = "1"
+OrdinaryDiffEqFIRK = "1"
 OrdinaryDiffEqFeagin = "1"
 OrdinaryDiffEqFunctionMap = "1"
+OrdinaryDiffEqHighOrderRK = "1"
+OrdinaryDiffEqIMEXMultistep = "1"
+OrdinaryDiffEqLinear = "1"
 OrdinaryDiffEqLowOrderRK = "1"
 OrdinaryDiffEqLowStorageRK = "1"
-OrdinaryDiffEqHighOrderRK = "1"
 OrdinaryDiffEqNonlinearSolve = "1"
+OrdinaryDiffEqNordsieck = "1"
+OrdinaryDiffEqPDIRK = "1"
 OrdinaryDiffEqPRK = "1"
 OrdinaryDiffEqQPRK = "1"
 OrdinaryDiffEqRKN = "1"
+OrdinaryDiffEqRosenbrock = "1"
+OrdinaryDiffEqSDIRK = "1"
+OrdinaryDiffEqStabilizedIRK = "1"
+OrdinaryDiffEqSSPRK = "1"
 OrdinaryDiffEqStabilizedRK = "1"
+OrdinaryDiffEqSymplecticRK = "1"
 OrdinaryDiffEqTsit5 = "1"
 OrdinaryDiffEqVerner = "1"
 Polyester = "0.7"

lib/OrdinaryDiffEqAdamsBashforthMoulton/src/adams_bashforth_moulton_caches.jl

@@ -1,7 +1,9 @@
 abstract type ABMMutableCache <: OrdinaryDiffEqMutableCache end
 abstract type ABMVariableCoefficientMutableCache <: OrdinaryDiffEqMutableCache end
-get_fsalfirstlast(cache::ABMMutableCache,u) = (cache.fsalfirst, cache.k)
-get_fsalfirstlast(cache::ABMVariableCoefficientMutableCache,u) = (cache.fsalfirst, cache.k4)
+get_fsalfirstlast(cache::ABMMutableCache, u) = (cache.fsalfirst, cache.k)
+function get_fsalfirstlast(cache::ABMVariableCoefficientMutableCache, u)
+    (cache.fsalfirst, cache.k4)
+end
 @cache mutable struct AB3Cache{uType, rateType} <: ABMMutableCache
     u::uType
     uprev::uType

lib/OrdinaryDiffEqBDF/src/algorithms.jl

@@ -150,7 +150,7 @@ end
 function QNDF1(; chunk_size = Val{0}(), autodiff = Val{true}(), standardtag = Val{true}(),
         concrete_jac = nothing, diff_type = Val{:forward},
         linsolve = nothing, precs = DEFAULT_PRECS, nlsolve = NLNewton(),
-        extrapolant = :linear, kappa = -37//200,
+        extrapolant = :linear, kappa = -37 // 200,
         controller = :Standard, step_limiter! = trivial_limiter!)
     QNDF1{
         _unwrap_val(chunk_size), _unwrap_val(autodiff), typeof(linsolve), typeof(nlsolve),
@@ -233,7 +233,8 @@ function QNDF(; max_order::Val{MO} = Val{5}(), chunk_size = Val{0}(),
         diff_type = Val{:forward},
         linsolve = nothing, precs = DEFAULT_PRECS, nlsolve = NLNewton(), κ = nothing,
         tol = nothing,
-        extrapolant = :linear, kappa = (-37//200, -1//9, -823//10000, -83//2000, 0//1),
+        extrapolant = :linear, kappa = (
+            -37 // 200, -1 // 9, -823 // 10000, -83 // 2000, 0 // 1),
         controller = :Standard, step_limiter! = trivial_limiter!) where {MO}
     QNDF{MO, _unwrap_val(chunk_size), _unwrap_val(autodiff), typeof(linsolve),
         typeof(nlsolve), typeof(precs), diff_type, _unwrap_val(standardtag),

lib/OrdinaryDiffEqBDF/src/bdf_caches.jl

@@ -1,5 +1,7 @@
 abstract type BDFMutableCache <: OrdinaryDiffEqMutableCache end
-get_fsalfirstlast(cache::BDFMutableCache,u) = (cache.fsalfirst, du_alias_or_new(cache.nlsolver, cache.fsalfirst))
+function get_fsalfirstlast(cache::BDFMutableCache, u)
+    (cache.fsalfirst, du_alias_or_new(cache.nlsolver, cache.fsalfirst))
+end
 
 @cache mutable struct ABDF2ConstantCache{N, dtType, rate_prototype} <:
                       OrdinaryDiffEqConstantCache

lib/OrdinaryDiffEqBDF/src/controllers.jl

@@ -206,7 +206,7 @@ function choose_order!(alg::FBDF, integrator,
                 Val(max_order))
             local terk_tmp
             if u isa Number
-                terk_tmp = fd_weights[k - 2, 1]*u
+                terk_tmp = fd_weights[k - 2, 1] * u
                 for i in 2:(k - 2)
                     terk_tmp += fd_weights[i, k - 2] * u_history[i - 1]
                 end
@@ -215,7 +215,7 @@ function choose_order!(alg::FBDF, integrator,
                 # we need terk_tmp to be mutable.
                 # so it can be updated
                 terk_tmp = similar(u)
-                @.. terk_tmp = fd_weights[k - 2, 1]*_vec(u)
+                @.. terk_tmp = fd_weights[k - 2, 1] * _vec(u)
                 for i in 2:(k - 2)
                     @.. @views terk_tmp += fd_weights[i, k - 2] * u_history[:, i - 1]
                 end

lib/OrdinaryDiffEqBDF/src/dae_caches.jl

@@ -1,5 +1,7 @@
 abstract type DAEBDFMutableCache <: OrdinaryDiffEqMutableCache end
-get_fsalfirstlast(cache::DAEBDFMutableCache,u) = (cache.fsalfirst, du_alias_or_new(cache.nlsolver, cache.fsalfirst))
+function get_fsalfirstlast(cache::DAEBDFMutableCache, u)
+    (cache.fsalfirst, du_alias_or_new(cache.nlsolver, cache.fsalfirst))
+end
 
 @cache mutable struct DImplicitEulerCache{uType, rateType, uNoUnitsType, N} <:
                       DAEBDFMutableCache
@@ -13,7 +15,7 @@ get_fsalfirstlast(cache::DAEBDFMutableCache,u) = (cache.fsalfirst, du_alias_or_n
 end
 
 # Not FSAL
-get_fsalfirstlast(cache::DImplicitEulerCache,u) = (u,u)
+get_fsalfirstlast(cache::DImplicitEulerCache, u) = (u, u)
 
 mutable struct DImplicitEulerConstantCache{N} <: OrdinaryDiffEqConstantCache
     nlsolver::N

lib/OrdinaryDiffEqBDF/src/dae_perform_step.jl

@@ -160,8 +160,6 @@ end
 end
 
 function initialize!(integrator, cache::DABDF2Cache)
-
-
     integrator.kshortsize = 2
     @unpack k₁, k₂ = cache.eulercache
     resize!(integrator.k, integrator.kshortsize)

lib/OrdinaryDiffEqCore/Project.toml

@@ -1,7 +1,7 @@
 name = "OrdinaryDiffEqCore"
 uuid = "bbf590c4-e513-4bbe-9b18-05decba2e5d8"
 authors = ["ParamThakkar123 <[email protected]>"]
-version = "1.2.0"
+version = "1.3.0"
 
 [deps]
 ADTypes = "47edcb42-4c32-4615-8424-f2b9edc5f35b"

lib/OrdinaryDiffEqCore/ext/OrdinaryDiffEqCoreEnzymeCoreExt.jl

@@ -1,12 +1,20 @@
 module OrdinaryDiffEqCoreEnzymeCoreExt
 import OrdinaryDiffEqCore, EnzymeCore
 
-Enzyme.EnzymeCore.EnzymeRules.inactive(::typeof(OrdinaryDiffEqCore.increment_nf!), args...) = true
-Enzyme.EnzymeCore.EnzymeRules.inactive(::typeof(OrdinaryDiffEqCore.fixed_t_for_floatingpoint_error!), args...) = true
-Enzyme.EnzymeCore.EnzymeRules.inactive(::typeof(OrdinaryDiffEqCore.increment_accept!), args...) = true
-Enzyme.EnzymeCore.EnzymeRules.inactive(::typeof(OrdinaryDiffEqCore.increment_reject!), args...) = true
-Enzyme.EnzymeCore.EnzymeRules.inactive(::typeof(OrdinaryDiffEqCore.increment_nf_perform_step!), args...) = true
-Enzyme.EnzymeCore.EnzymeRules.inactive(::typeof(OrdinaryDiffEqCore.check_error!), args...) = true
-Enzyme.EnzymeCore.EnzymeRules.inactive(::typeof(OrdinaryDiffEqCore.log_step!), args...) = true
+EnzymeCore.EnzymeRules.inactive(::typeof(OrdinaryDiffEqCore.increment_nf!), args...) = true
+function EnzymeCore.EnzymeRules.inactive(
+        ::typeof(OrdinaryDiffEqCore.fixed_t_for_floatingpoint_error!), args...)
+    true
+end
+function EnzymeCore.EnzymeRules.inactive(
+        ::typeof(OrdinaryDiffEqCore.increment_accept!), args...)
+    true
+end
+function EnzymeCore.EnzymeRules.inactive(
+        ::typeof(OrdinaryDiffEqCore.increment_reject!), args...)
+    true
+end
+EnzymeCore.EnzymeRules.inactive(::typeof(OrdinaryDiffEqCore.check_error!), args...) = true
+EnzymeCore.EnzymeRules.inactive(::typeof(OrdinaryDiffEqCore.log_step!), args...) = true
 
-end
+end

lib/OrdinaryDiffEqCore/src/OrdinaryDiffEqCore.jl

@@ -41,7 +41,7 @@ using SimpleUnPack, RecursiveArrayTools, DataStructures, ArrayInterface
 import TruncatedStacktraces
 
 import StaticArraysCore: SArray, MVector, SVector, StaticArray, MMatrix,
-                     StaticMatrix
+                         StaticMatrix
 
 # Integrator Interface
 import DiffEqBase: resize!, deleteat!, addat!, full_cache, user_cache, u_cache, du_cache,
@@ -54,7 +54,6 @@ import DiffEqBase: resize!, deleteat!, addat!, full_cache, user_cache, u_cache,
                    isautodifferentiable,
                    get_tstops, get_tstops_array, get_tstops_max
 
-
 using DiffEqBase: check_error!, @def, _vec, _reshape
 
 using FastBroadcast: @.., True, False

lib/OrdinaryDiffEqCore/src/alg_utils.jl

@@ -345,6 +345,18 @@ qsteady_max_default(alg::OrdinaryDiffEqImplicitAlgorithm) = isadaptive(alg) ? 1
 #DiffEqBase.nlsolve_default(::QNDF, ::Val{κ}) = 1//2
 
 # SSP coefficients
+
+"""
+    ssp_coefficient(alg)
+Return the SSP coefficient of the ODE algorithm `alg`. If one time step of size
+`dt` with `alg` can be written as a convex combination of explicit Euler steps
+with step sizes `cᵢ * dt`, the SSP coefficient is the minimal value of `1/cᵢ`.
+# Examples
+```julia-repl
+julia> ssp_coefficient(SSPRK104())
+6
+```
+"""
 ssp_coefficient(alg) = error("$alg is not a strong stability preserving method.")
 
 # We shouldn't do this probably.

lib/OrdinaryDiffEqCore/src/caches/basic_caches.jl

@@ -5,15 +5,15 @@ struct ODEEmptyCache <: OrdinaryDiffEqConstantCache end
 struct ODEChunkCache{CS} <: OrdinaryDiffEqConstantCache end
 
 # Don't worry about the potential alloc on a constant cache
-get_fsalfirstlast(cache::OrdinaryDiffEqConstantCache,u) = (zero(u), zero(u))
+get_fsalfirstlast(cache::OrdinaryDiffEqConstantCache, u) = (zero(u), zero(u))
 
 mutable struct CompositeCache{T, F} <: OrdinaryDiffEqCache
     caches::T
     choice_function::F
     current::Int
 end
 
-get_fsalfirstlast(cache::CompositeCache,u) = get_fsalfirstlast(cache.caches[1],u)
+get_fsalfirstlast(cache::CompositeCache, u) = get_fsalfirstlast(cache.caches[1], u)
 
 mutable struct DefaultCache{T1, T2, T3, T4, T5, T6, A, F, uType} <: OrdinaryDiffEqCache
     args::A
@@ -28,12 +28,13 @@ mutable struct DefaultCache{T1, T2, T3, T4, T5, T6, A, F, uType} <: OrdinaryDiff
     cache6::T6
     function DefaultCache{T1, T2, T3, T4, T5, T6, F, uType}(
             args, choice_function, current, u) where {T1, T2, T3, T4, T5, T6, F, uType}
-        new{T1, T2, T3, T4, T5, T6, typeof(args), F, uType}(args, choice_function, current, u)
+        new{T1, T2, T3, T4, T5, T6, typeof(args), F, uType}(
+            args, choice_function, current, u)
     end
 end
 
-function get_fsalfirstlast(cache::DefaultCache,u) 
-    (cache.u,cache.u) # will be overwritten by the cache choice
+function get_fsalfirstlast(cache::DefaultCache, u)
+    (cache.u, cache.u) # will be overwritten by the cache choice
 end
 
 function alg_cache(alg::CompositeAlgorithm, u, rate_prototype, ::Type{uEltypeNoUnits},

lib/OrdinaryDiffEqCore/src/integrators/integrator_utils.jl

@@ -255,8 +255,8 @@ function _loopfooter!(integrator)
     end
     if integrator.opts.progress && integrator.iter % integrator.opts.progress_steps == 0
         log_step!(integrator.opts.progress_name, integrator.opts.progress_id,
-                  integrator.opts.progress_message, integrator.dt, integrator.u,
-                  integrator.p, integrator.t, integrator.sol.prob.tspan)
+            integrator.opts.progress_message, integrator.dt, integrator.u,
+            integrator.p, integrator.t, integrator.sol.prob.tspan)
     end
 
     # Take value because if t is dual then maxeig can be dual
@@ -280,18 +280,18 @@ end
 
 function log_step!(progress_name, progress_id, progress_message, dt, u, p, t, tspan)
     t1, t2 = tspan
-    @logmsg(LogLevel(-1),progress_name,
-    _id=progress_id,
-    message=progress_message(dt, u, p, t),
-    progress=(t - t1) / (t2 - t1))
+    @logmsg(LogLevel(-1), progress_name,
+        _id=progress_id,
+        message=progress_message(dt, u, p, t),
+        progress=(t - t1) / (t2 - t1))
 end
 
 function fixed_t_for_floatingpoint_error!(integrator, ttmp)
     if has_tstop(integrator)
         tstop = integrator.tdir * first_tstop(integrator)
         if abs(ttmp - tstop) <
-        100eps(float(max(integrator.t, tstop) / oneunit(integrator.t))) *
-        oneunit(integrator.t)
+           100eps(float(max(integrator.t, tstop) / oneunit(integrator.t))) *
+           oneunit(integrator.t)
             tstop
         else
             ttmp

lib/OrdinaryDiffEqCore/src/interp_func.jl

@@ -62,7 +62,6 @@ end
 
 # strip interpolation of function information
 function SciMLBase.strip_interpolation(id::InterpolationData)
-
     cache = strip_cache(id.cache)
 
     InterpolationData(nothing, id.timeseries,
@@ -78,7 +77,7 @@ end
 function strip_cache(cache)
     if hasfield(typeof(cache), :jac_config) || hasfield(typeof(cache), :grad_config)
         fieldnums = length(fieldnames(typeof(cache)))
-        noth_list = fill(nothing,fieldnums)
+        noth_list = fill(nothing, fieldnums)
         cache_type_name = Base.typename(typeof(cache)).wrapper
         cache_type_name(noth_list...)
     else