Skip to content

Commit

Permalink
Merge branch 'master' into relaxation_method
Browse files Browse the repository at this point in the history
  • Loading branch information
Theozeud authored May 27, 2024
2 parents e35d9fa + 1e8716a commit 7fada13
Show file tree
Hide file tree
Showing 62 changed files with 3,334 additions and 995 deletions.
8 changes: 6 additions & 2 deletions Project.toml
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
name = "OrdinaryDiffEq"
uuid = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
authors = ["Chris Rackauckas <[email protected]>", "Yingbo Ma <[email protected]>"]
version = "6.74.1"
version = "6.80.1"

[deps]
ADTypes = "47edcb42-4c32-4615-8424-f2b9edc5f35b"
Expand All @@ -10,6 +10,7 @@ ArrayInterface = "4fba245c-0d91-5ea0-9b3e-6abc04ee57a9"
DataStructures = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8"
DiffEqBase = "2b5f629d-d688-5b77-993f-72d75c75574e"
DocStringExtensions = "ffbed154-4ef7-542d-bbb7-c09d3a79fcae"
EnumX = "4e289a0a-7415-4d19-859d-a7e5c4648b56"
ExponentialUtilities = "d4d017d3-3776-5f7e-afef-a10c40355c18"
FastBroadcast = "7034ab61-46d4-4ed7-9d0f-46aef9175898"
FastClosures = "9aa1b823-49e4-5ca5-8b0f-3971ec8bab6a"
Expand All @@ -34,6 +35,7 @@ RecursiveArrayTools = "731186ca-8d62-57ce-b412-fbd966d074cd"
Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
SciMLBase = "0bca4576-84f4-4d90-8ffe-ffa030f20462"
SciMLOperators = "c0aeaf25-5076-4817-a8d5-81caf7dfa961"
SciMLStructures = "53ae85a6-f571-4167-b2af-e1d143709226"
SimpleNonlinearSolve = "727e6d20-b764-4bd8-a329-72de5adea6c7"
SimpleUnPack = "ce78b400-467f-4804-87d8-8f486da07d0a"
SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
Expand All @@ -49,8 +51,9 @@ ArrayInterface = "7"
DataStructures = "0.18"
DiffEqBase = "6.147"
DocStringExtensions = "0.9"
EnumX = "1"
ExponentialUtilities = "1.22"
FastBroadcast = "0.2"
FastBroadcast = "0.2, 0.3"
FastClosures = "0.3"
FillArrays = "1.9"
FiniteDiff = "2"
Expand All @@ -74,6 +77,7 @@ RecursiveArrayTools = "2.36, 3"
Reexport = "1.0"
SciMLBase = "2.27.1"
SciMLOperators = "0.3"
SciMLStructures = "1"
SimpleNonlinearSolve = "1"
SimpleUnPack = "1"
SparseArrays = "1.9"
Expand Down
3 changes: 2 additions & 1 deletion docs/Project.toml
Original file line number Diff line number Diff line change
Expand Up @@ -3,4 +3,5 @@ Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"

[compat]
Documenter = "0.27"
Documenter = "0.27, 1"
OrdinaryDiffEq = "6"
15 changes: 7 additions & 8 deletions docs/src/stiff/rosenbrock.md
Original file line number Diff line number Diff line change
Expand Up @@ -3,13 +3,6 @@
## Standard Rosenbrock Methods

```@docs
ROS2
ROS3
ROS2PR
ROS3PR
Scholz47
ROS3PRL
ROS3PRL2
ROS3P
Rodas3
Rodas3P
Expand All @@ -25,7 +18,6 @@ Rodas4P
Rodas4P2
Rodas5
Rodas5P
GeneralRosenbrock
```

## Rosenbrock W-Methods
Expand All @@ -39,6 +31,13 @@ ROS34PW1b
ROS34PW2
ROS34PW3
ROS34PRw
ROS2
ROS3
ROS2S
RosenbrockW6S4OS
ROS2PR
ROS3PR
Scholz4_7
ROS3PRL
ROS3PRL2
```
46 changes: 35 additions & 11 deletions src/OrdinaryDiffEq.jl
Original file line number Diff line number Diff line change
Expand Up @@ -26,7 +26,9 @@ using LinearSolve, SimpleNonlinearSolve

using LineSearches

import FillArrays: Trues
import EnumX

import FillArrays: Trues, Falses

# Interfaces
import DiffEqBase: solve!, step!, initialize!, isadaptive
Expand Down Expand Up @@ -108,6 +110,8 @@ import ADTypes: AbstractADType,
import Polyester
using MacroTools, Adapt

using SciMLStructures: canonicalize, Tunable, isscimlstructure

const CompiledFloats = Union{Float32, Float64,
ForwardDiff.Dual{
ForwardDiff.Tag{T, W},
Expand All @@ -124,6 +128,7 @@ import Preferences

DEFAULT_PRECS(W, du, u, p, t, newW, Plprev, Prprev, solverdata) = nothing, nothing

include("doc_utils.jl")
include("misc_utils.jl")

include("algorithms.jl")
Expand All @@ -139,6 +144,7 @@ include("nlsolve/functional.jl")
include("nlsolve/newton.jl")

include("generic_rosenbrock.jl")
include("composite_algs.jl")

include("caches/basic_caches.jl")
include("caches/low_order_rk_caches.jl")
Expand Down Expand Up @@ -233,7 +239,6 @@ include("constants.jl")
include("solve.jl")
include("initdt.jl")
include("interp_func.jl")
include("composite_algs.jl")

import PrecompileTools

Expand All @@ -252,9 +257,18 @@ PrecompileTools.@compile_workload begin
Tsit5(), Vern7()
]

stiff = [Rosenbrock23(), Rosenbrock23(autodiff = false),
Rodas5P(), Rodas5P(autodiff = false),
FBDF(), FBDF(autodiff = false)
stiff = [Rosenbrock23(),
Rodas5P(),
FBDF()
]

default_ode = [
DefaultODEAlgorithm(autodiff = false),
]


default_autodiff_ode = [
DefaultODEAlgorithm(),
]

autoswitch = [
Expand All @@ -275,15 +289,24 @@ PrecompileTools.@compile_workload begin
solver_list = []
solver_list_nonadaptive = []

if Preferences.@load_preference("PrecompileNonStiff", true)

if Preferences.@load_preference("PrecompileDefault", true)
append!(solver_list, default_ode)
end

if Preferences.@load_preference("PrecompileAutodiffDefault", true)
append!(solver_list, default_autodiff_ode)
end

if Preferences.@load_preference("PrecompileNonStiff", false)
append!(solver_list, nonstiff)
end

if Preferences.@load_preference("PrecompileStiff", true)
if Preferences.@load_preference("PrecompileStiff", false)
append!(solver_list, stiff)
end

if Preferences.@load_preference("PrecompileAutoSwitch", true)
if Preferences.@load_preference("PrecompileAutoSwitch", false)
append!(solver_list, autoswitch)
end

Expand Down Expand Up @@ -355,7 +378,7 @@ export FunctionMap, Euler, Heun, Ralston, Midpoint, RK4, ExplicitRK, OwrenZen3,
OwrenZen5,
BS3, BS5, RK46NL, DP5, Tsit5, DP8, Vern6, Vern7, Vern8, TanYam7, TsitPap8,
Vern9, Feagin10, Feagin12, Feagin14, CompositeAlgorithm, Anas5, RKO65, FRK65, PFRK87,
RKM, MSRK5, MSRK6, Stepanov5, SIR54, QPRK98, Tsit5_for_relaxation
RKM, MSRK5, MSRK6, Stepanov5, SIR54, QPRK98, PSRK4p7q6, PSRK3p6q5, PSRK3p5q4, Tsit5_for_relaxation

export SSPRK22, SSPRK33, KYKSSPRK42, SSPRK53, SSPRK53_2N1, SSPRK53_2N2, SSPRK53_H, SSPRK63,
SSPRK73, SSPRK83, SSPRK43, SSPRK432,
Expand Down Expand Up @@ -407,7 +430,8 @@ export SplitEuler

export Nystrom4, FineRKN4, FineRKN5, Nystrom4VelocityIndependent,
Nystrom5VelocityIndependent,
IRKN3, IRKN4, DPRKN4, DPRKN5, DPRKN6, DPRKN6FM, DPRKN8, DPRKN12, ERKN4, ERKN5, ERKN7
IRKN3, IRKN4, DPRKN4, DPRKN5, DPRKN6, DPRKN6FM, DPRKN8, DPRKN12, ERKN4, ERKN5, ERKN7,
RKN4

export ROCK2, ROCK4, RKC, IRKC, ESERK4, ESERK5, SERK2

Expand Down Expand Up @@ -441,7 +465,7 @@ export KuttaPRK2p5, PDIRK44, DImplicitEuler, DABDF2, DFBDF

export ShampineCollocationInit, BrownFullBasicInit, NoInit

export NLNewton, NLAnderson, NLFunctional
export NLNewton, NLAnderson, NLFunctional, NonlinearSolveAlg

export IController, PIController, PIDController
end # module
77 changes: 50 additions & 27 deletions src/alg_utils.jl
Original file line number Diff line number Diff line change
Expand Up @@ -78,6 +78,10 @@ isfsal(alg::SSPRK932) = false
isfsal(alg::SSPRK54) = false
isfsal(alg::SSPRK104) = false

isfsal(alg::PSRK3p5q4) = false
isfsal(alg::PSRK3p6q5) = false
isfsal(alg::PSRK4p7q6) = false

get_current_isfsal(alg, cache) = isfsal(alg)

# evaluates f(t[i])
Expand Down Expand Up @@ -168,6 +172,7 @@ isimplicit(alg::CompositeAlgorithm) = any(isimplicit.(alg.algs))

isdtchangeable(alg::Union{OrdinaryDiffEqAlgorithm, DAEAlgorithm}) = true
isdtchangeable(alg::CompositeAlgorithm) = all(isdtchangeable.(alg.algs))

function isdtchangeable(alg::Union{LawsonEuler, NorsettEuler, LieEuler, MagnusGauss4,
CayleyEuler, ETDRK2, ETDRK3, ETDRK4, HochOst4, ETD2})
false
Expand Down Expand Up @@ -201,31 +206,35 @@ qmax_default(alg::CompositeAlgorithm) = minimum(qmax_default.(alg.algs))
qmax_default(alg::DP8) = 6
qmax_default(alg::Union{RadauIIA3, RadauIIA5}) = 8

function has_chunksize(alg::OrdinaryDiffEqAlgorithm)
return alg isa Union{OrdinaryDiffEqExponentialAlgorithm,
OrdinaryDiffEqAdaptiveExponentialAlgorithm,
OrdinaryDiffEqImplicitAlgorithm,
OrdinaryDiffEqAdaptiveImplicitAlgorithm,
DAEAlgorithm,
CompositeAlgorithm}
end
function get_chunksize(alg::OrdinaryDiffEqAlgorithm)
error("This algorithm does not have a chunk size defined.")
end
get_chunksize(alg::OrdinaryDiffEqAdaptiveImplicitAlgorithm{CS, AD}) where {CS, AD} = Val(CS)
get_chunksize(alg::OrdinaryDiffEqImplicitAlgorithm{CS, AD}) where {CS, AD} = Val(CS)
get_chunksize(alg::DAEAlgorithm{CS, AD}) where {CS, AD} = Val(CS)
function get_chunksize(alg::Union{OrdinaryDiffEqExponentialAlgorithm{CS, AD},
OrdinaryDiffEqAdaptiveExponentialAlgorithm{CS, AD}}) where {
CS,
AD
}
function get_chunksize(alg::Union{OrdinaryDiffEqExponentialAlgorithm{CS},
OrdinaryDiffEqAdaptiveExponentialAlgorithm{CS},
OrdinaryDiffEqImplicitAlgorithm{CS},
OrdinaryDiffEqAdaptiveImplicitAlgorithm{CS},
DAEAlgorithm{CS},
CompositeAlgorithm{CS}}) where {CS}
Val(CS)
end

function get_chunksize_int(alg::OrdinaryDiffEqAlgorithm)
error("This algorithm does not have a chunk size defined.")
end
get_chunksize_int(alg::OrdinaryDiffEqAdaptiveImplicitAlgorithm{CS, AD}) where {CS, AD} = CS
get_chunksize_int(alg::OrdinaryDiffEqImplicitAlgorithm{CS, AD}) where {CS, AD} = CS
get_chunksize_int(alg::DAEAlgorithm{CS, AD}) where {CS, AD} = CS
function get_chunksize_int(alg::Union{OrdinaryDiffEqExponentialAlgorithm{CS, AD},
OrdinaryDiffEqAdaptiveExponentialAlgorithm{CS, AD}}) where {
CS,
AD
}
function get_chunksize_int(alg::Union{OrdinaryDiffEqExponentialAlgorithm{CS},
OrdinaryDiffEqAdaptiveExponentialAlgorithm{CS},
OrdinaryDiffEqImplicitAlgorithm{CS},
OrdinaryDiffEqAdaptiveImplicitAlgorithm{CS},
DAEAlgorithm{CS},
CompositeAlgorithm{CS}}) where {CS}
CS
end
# get_chunksize(alg::CompositeAlgorithm) = get_chunksize(alg.algs[alg.current_alg])
Expand Down Expand Up @@ -494,6 +503,7 @@ alg_order(alg::DPRKN12) = 12
alg_order(alg::ERKN4) = 4
alg_order(alg::ERKN5) = 5
alg_order(alg::ERKN7) = 7
alg_order(alg::RKN4) = 4

alg_order(alg::Midpoint) = 2

Expand Down Expand Up @@ -567,6 +577,9 @@ alg_order(alg::MSRK5) = 5
alg_order(alg::MSRK6) = 6
alg_order(alg::Stepanov5) = 5
alg_order(alg::SIR54) = 5
alg_order(alg::PSRK4p7q6) = 4
alg_order(alg::PSRK3p6q5) = 3
alg_order(alg::PSRK3p5q4) = 3

alg_order(alg::BS3) = 3
alg_order(alg::BS5) = 5
Expand Down Expand Up @@ -958,10 +971,12 @@ alg_can_repeat_jac(alg::OrdinaryDiffEqNewtonAdaptiveAlgorithm) = true
alg_can_repeat_jac(alg::IRKC) = false

function unwrap_alg(alg::SciMLBase.DEAlgorithm, is_stiff)
iscomp = alg isa CompositeAlgorithm
if !iscomp
if !(alg isa CompositeAlgorithm)
return alg
elseif alg.choice_function isa AutoSwitchCache
if length(alg.algs) > 2
return alg.algs[alg.choice_function.current]
end
if is_stiff === nothing
throwautoswitch(alg)
end
Expand All @@ -978,18 +993,21 @@ end

function unwrap_alg(integrator, is_stiff)
alg = integrator.alg
iscomp = alg isa CompositeAlgorithm
if !iscomp
if !(alg isa CompositeAlgorithm)
return alg
elseif alg.choice_function isa AutoSwitchCache
if is_stiff === nothing
throwautoswitch(alg)
end
num = is_stiff ? 2 : 1
if num == 1
return alg.algs[1]
if length(alg.algs) > 2
alg.algs[alg.choice_function.current]
else
return alg.algs[2]
if is_stiff === nothing
throwautoswitch(alg)
end
num = is_stiff ? 2 : 1
if num == 1
return alg.algs[1]
else
return alg.algs[2]
end
end
else
return _eval_index(identity, alg.algs, integrator.cache.current)
Expand Down Expand Up @@ -1064,3 +1082,8 @@ is_mass_matrix_alg(alg::Union{OrdinaryDiffEqAlgorithm, DAEAlgorithm}) = false
is_mass_matrix_alg(alg::CompositeAlgorithm) = all(is_mass_matrix_alg, alg.algs)
is_mass_matrix_alg(alg::RosenbrockAlgorithm) = true
is_mass_matrix_alg(alg::NewtonAlgorithm) = !isesdirk(alg)
# hack for the default alg
function is_mass_matrix_alg(alg::CompositeAlgorithm{
<:Any, <:Tuple{Tsit5, Vern7, Rosenbrock23, Rodas5P, FBDF, FBDF}})
true
end
Loading

0 comments on commit 7fada13

Please sign in to comment.