Move GlobalDiffEq.jl into the monorepo as lib/GlobalDiffEq (repo transfer)

docs/Project.toml

@@ -2,6 +2,7 @@
 DiffEqBase = "2b5f629d-d688-5b77-993f-72d75c75574e"
 DiffEqDevTools = "f3b72e0c-5b89-59e1-b016-84e28bfd966d"
 Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
+GlobalDiffEq = "1d72d19b-84cc-4cb7-a099-7cbdb9ccc67c"
 ImplicitDiscreteSolve = "3263718b-31ed-49cf-8a0f-35a466e8af96"
 OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
 OrdinaryDiffEqAMF = "08082164-f6a1-4363-b3df-3aa6fcf571ad"
@@ -39,6 +40,7 @@ OrdinaryDiffEqVerner = "79d7bb75-1356-48c1-b8c0-6832512096c2"
 [sources]
 DiffEqBase = {path = "../lib/DiffEqBase"}
 DiffEqDevTools = {path = "../lib/DiffEqDevTools"}
+GlobalDiffEq = {path = "../lib/GlobalDiffEq"}
 ImplicitDiscreteSolve = {path = "../lib/ImplicitDiscreteSolve"}
 OrdinaryDiffEqAMF = {path = "../lib/OrdinaryDiffEqAMF"}
 OrdinaryDiffEqAdamsBashforthMoulton = {path = "../lib/OrdinaryDiffEqAdamsBashforthMoulton"}
@@ -76,6 +78,7 @@ OrdinaryDiffEqVerner = {path = "../lib/OrdinaryDiffEqVerner"}
 DiffEqBase = "7.6.1"
 DiffEqDevTools = "3, 4.0"
 Documenter = "0.27, 1"
+GlobalDiffEq = "1.2.1"
 ImplicitDiscreteSolve = "2"
 OrdinaryDiffEq = "7"
 OrdinaryDiffEqAMF = "2"

docs/make.jl

@@ -36,6 +36,7 @@ using OrdinaryDiffEqStabilizedRK
 using OrdinaryDiffEqSymplecticRK
 using OrdinaryDiffEqTsit5
 using OrdinaryDiffEqVerner
+using GlobalDiffEq
 
 cp(joinpath(@__DIR__, "Manifest.toml"), joinpath(@__DIR__, "src", "assets", "Manifest.toml"), force = true)
 cp(joinpath(@__DIR__, "Project.toml"), joinpath(@__DIR__, "src", "assets", "Project.toml"), force = true)
@@ -83,6 +84,7 @@ makedocs(
         OrdinaryDiffEqAMF,
         ImplicitDiscreteSolve,
         DiffEqDevTools,
+        GlobalDiffEq,
     ],
     linkcheck_ignore = [r"https://github.com/JuliaDiff/ForwardDiff.jl"],
     warnonly = [:docs_block, :missing_docs, :eval_block],

docs/pages.jl

@@ -53,6 +53,9 @@ pages = [
     "Misc Solvers" => [
         "misc.md",
     ],
+    "Global Error Control" => [
+        "globalerrorcontrol/GlobalDiffEq.md",
+    ],
     "APIs" => [
         "api/diffeqbase.md",
         "api/ordinarydiffeqcore.md",

docs/src/globalerrorcontrol/GlobalDiffEq.md

@@ -0,0 +1,25 @@
+# GlobalDiffEq: Global Error Estimation and Control
+
+Standard adaptive ODE solvers control the *local* error of each step. The
+local tolerances only indirectly control the *global* (accumulated) error of
+the solution, which can grow arbitrarily large over long integrations or on
+unstable problems even when every step satisfies its local tolerance. The
+GlobalDiffEq sublibrary provides solvers and solver wrappers that estimate
+the global error, and in several cases control it to a requested global
+tolerance.
+
+To use these methods:
+
+```julia
+using GlobalDiffEq
+```
+
+[`GlobalRichardson`](@ref) wraps any fixed-step method in global Richardson
+extrapolation over whole solves, interpreting `abstol` and `reltol` as global
+tolerances. It is the most robust and most expensive option.
+
+## Global error controlling wrappers
+
+```@docs
+GlobalRichardson
+```

lib/GlobalDiffEq/LICENSE.md

@@ -0,0 +1,21 @@
+MIT License
+
+Copyright (c) 2020 Chris Rackauckas <accounts@chrisrackauckas.com> and contributors
+
+Permission is hereby granted, free of charge, to any person obtaining a copy
+of this software and associated documentation files (the "Software"), to deal
+in the Software without restriction, including without limitation the rights
+to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+copies of the Software, and to permit persons to whom the Software is
+furnished to do so, subject to the following conditions:
+
+The above copyright notice and this permission notice shall be included in all
+copies or substantial portions of the Software.
+
+THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+SOFTWARE.

lib/GlobalDiffEq/Project.toml

@@ -0,0 +1,41 @@
+name = "GlobalDiffEq"
+uuid = "1d72d19b-84cc-4cb7-a099-7cbdb9ccc67c"
+authors = ["Chris Rackauckas <accounts@chrisrackauckas.com> and contributors"]
+version = "1.2.1"
+
+[deps]
+DiffEqBase = "2b5f629d-d688-5b77-993f-72d75c75574e"
+OrdinaryDiffEqTsit5 = "b1df2697-797e-41e3-8120-5422d3b24e4a"
+PrecompileTools = "aea7be01-6a6a-4083-8856-8a6e6704d82a"
+Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
+Richardson = "708f8203-808e-40c0-ba2d-98a6953ed40d"
+SciMLBase = "0bca4576-84f4-4d90-8ffe-ffa030f20462"
+
+[sources]
+DiffEqBase = {path = "../DiffEqBase"}
+OrdinaryDiffEqTsit5 = {path = "../OrdinaryDiffEqTsit5"}
+
+[compat]
+DiffEqBase = "7"
+LinearAlgebra = "1"
+OrdinaryDiffEqSSPRK = "2"
+OrdinaryDiffEqTsit5 = "2"
+PrecompileTools = "1.1"
+Reexport = "0.2, 1.0"
+Richardson = "1.2"
+SafeTestsets = "0.0.1, 0.1"
+SciMLBase = "3"
+SciMLTesting = "2.1"
+Test = "1"
+julia = "1.10"
+
+[extras]
+LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
+OrdinaryDiffEqSSPRK = "669c94d9-1f4b-4b64-b377-1aa079aa2388"
+OrdinaryDiffEqTsit5 = "b1df2697-797e-41e3-8120-5422d3b24e4a"
+SafeTestsets = "1bc83da4-3b8d-516f-aca4-4fe02f6d838f"
+SciMLTesting = "09d9d899-5365-40a9-917a-5f67fddea283"
+Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
+
+[targets]
+test = ["LinearAlgebra", "OrdinaryDiffEqSSPRK", "OrdinaryDiffEqTsit5", "SafeTestsets", "SciMLTesting", "Test"]

lib/GlobalDiffEq/README.md

@@ -0,0 +1,17 @@
+# GlobalDiffEq
+
+Differential equation solvers and solver wrappers with global (accumulated)
+error estimates and global error control for the
+[DifferentialEquations.jl common solver interface](https://docs.sciml.ai/DiffEqDocs/stable/).
+Formerly the standalone [GlobalDiffEq.jl](https://github.com/SciML/GlobalDiffEq.jl)
+repository, now maintained as a sublibrary of OrdinaryDiffEq.jl.
+
+Standard adaptive ODE solvers control the *local* error of each step, which
+does not bound the error accumulated over the whole integration. This package
+provides:
+
+  - `GlobalRichardson`: global Richardson extrapolation of whole solves of any
+    fixed-step method, interpreting `abstol`/`reltol` as global tolerances.
+
+See the [OrdinaryDiffEq.jl documentation](https://docs.sciml.ai/OrdinaryDiffEq/stable/)
+for the full API reference.

lib/GlobalDiffEq/src/GlobalDiffEq.jl

@@ -0,0 +1,32 @@
+module GlobalDiffEq
+
+using Reexport: @reexport
+@reexport using DiffEqBase
+
+import OrdinaryDiffEqTsit5, Richardson, SciMLBase
+using PrecompileTools: @setup_workload, @compile_workload
+
+abstract type GlobalDiffEqAlgorithm <: SciMLBase.AbstractODEAlgorithm end
+
+include("richardson.jl")
+
+export GlobalRichardson
+
+@setup_workload begin
+    # Simple test ODE: exponential decay du/dt = -u
+    function f!(du, u, p, t)
+        du[1] = -u[1]
+    end
+    u0 = [1.0]
+    tspan = (0.0, 1.0)
+    prob = ODEProblem(f!, u0, tspan)
+
+    @compile_workload begin
+        solve(
+            prob, GlobalRichardson(OrdinaryDiffEqTsit5.Tsit5()),
+            dt = 0.1, reltol = 1.0e-3, abstol = 1.0e-6
+        )
+    end
+end
+
+end

lib/GlobalDiffEq/src/richardson.jl

@@ -0,0 +1,38 @@
+"""
+    GlobalRichardson(alg)
+
+Wrap an ODE algorithm with global Richardson extrapolation.
+"""
+struct GlobalRichardson{A <: SciMLBase.AbstractODEAlgorithm} <: GlobalDiffEqAlgorithm
+    alg::A
+end
+
+# Forward algorithm traits to the wrapped algorithm
+# This allows GlobalRichardson to inherit capabilities from the inner algorithm
+SciMLBase.allows_arbitrary_number_types(alg::GlobalRichardson) =
+    SciMLBase.allows_arbitrary_number_types(alg.alg)
+SciMLBase.allowscomplex(alg::GlobalRichardson) =
+    SciMLBase.allowscomplex(alg.alg)
+SciMLBase.isautodifferentiable(alg::GlobalRichardson) =
+    SciMLBase.isautodifferentiable(alg.alg)
+
+function SciMLBase.__solve(
+        prob::Union{SciMLBase.AbstractODEProblem, SciMLBase.AbstractDAEProblem},
+        alg::GlobalRichardson, args...;
+        dt, kwargs...
+    )
+    opt = Dict(kwargs)
+    otheropts = delete!(copy(opt), :dt)
+    tstops = get(opt, :tstops, range(prob.tspan[1], stop = prob.tspan[2], step = dt))
+    local sol
+    val,
+        err = Richardson.extrapolate(
+        dt, rtol = get(opt, :reltol, 1.0e-3),
+        atol = get(opt, :abstol, 1.0e-6), contract = 0.5
+    ) do _dt
+        sol = solve(prob, alg.alg, args...; dt = _dt, adaptive = false, otheropts...)
+        # Convert Vector{Vector{T}} to Matrix{T} for Richardson.jl compatibility
+        reduce(hcat, sol.(tstops))
+    end
+    return sol
+end

lib/GlobalDiffEq/test/algorithm_traits_tests.jl

@@ -0,0 +1,14 @@
+using GlobalDiffEq, OrdinaryDiffEqTsit5
+using OrdinaryDiffEqSSPRK
+using Test
+import DiffEqBase: SciMLBase
+
+@testset "Algorithm traits forwarding" begin
+    alg_inner = SSPRK33()
+    alg = GlobalRichardson(alg_inner)
+
+    @test SciMLBase.allows_arbitrary_number_types(alg) ==
+        SciMLBase.allows_arbitrary_number_types(alg_inner)
+    @test SciMLBase.allowscomplex(alg) == SciMLBase.allowscomplex(alg_inner)
+    @test SciMLBase.isautodifferentiable(alg) == SciMLBase.isautodifferentiable(alg_inner)
+end

lib/GlobalDiffEq/test/basic_functionality_tests.jl

@@ -0,0 +1,29 @@
+using GlobalDiffEq, OrdinaryDiffEqTsit5, LinearAlgebra
+using OrdinaryDiffEqSSPRK
+using Test
+
+@testset "Basic functionality" begin
+    l = 1.0                             # length [m]
+    m = 1.0                             # mass[m]
+    g = 9.81                            # gravitational acceleration [m/s²]
+
+    function pendulum!(du, u, p, t)
+        du[1] = u[2]                    # θ'(t) = ω(t)
+        return du[2] = -3g / (2l) * sin(u[1]) + 3 / (m * l^2) * p(t) # ω'(t) = -3g/(2l) sin θ(t) + 3/(ml^2)M(t)
+    end
+
+    θ₀ = 0.01                           # initial angular deflection [rad]
+    ω₀ = 0.0                            # initial angular velocity [rad/s]
+    u₀ = [θ₀, ω₀]                       # initial state vector
+    tspan = (0.0, 10.0)                  # time interval
+
+    M = t -> 0.1sin(t)                    # external torque [Nm]
+
+    prob = ODEProblem(pendulum!, u₀, tspan, M)
+
+    v0 = solve(prob, Tsit5(), dt = 0.1, reltol = 1.0e-12, abstol = 0).(1:10)
+    ve = solve(
+        prob, GlobalRichardson(SSPRK33()), dt = 0.2, reltol = 1.0e-12, abstol = 0
+    ).(1:10)
+    @test norm(ve - v0) / norm(v0) < 1.0e-10
+end

lib/GlobalDiffEq/test/bigfloat_tests.jl

@@ -0,0 +1,22 @@
+using GlobalDiffEq, OrdinaryDiffEqTsit5
+using OrdinaryDiffEqSSPRK
+using Test
+
+@testset "BigFloat support" begin
+    function f_bf!(du, u, p, t)
+        du[1] = -u[1]
+    end
+
+    u0_bf = BigFloat[1.0]
+    tspan_bf = (BigFloat(0.0), BigFloat(1.0))
+    prob_bf = ODEProblem(f_bf!, u0_bf, tspan_bf)
+
+    sol_bf = solve(
+        prob_bf, GlobalRichardson(SSPRK33()),
+        dt = BigFloat(0.1), reltol = BigFloat(1.0e-3), abstol = BigFloat(1.0e-6)
+    )
+
+    @test eltype(sol_bf.u[end]) == BigFloat
+    # Check solution is reasonable (e^-1 ≈ 0.368)
+    @test isapprox(sol_bf.u[end][1], exp(BigFloat(-1)), rtol = 1.0e-4)
+end

lib/GlobalDiffEq/test/qa/Project.toml

@@ -0,0 +1,21 @@
+[deps]
+Aqua = "4c88cf16-eb10-579e-8560-4a9242c79595"
+GlobalDiffEq = "1d72d19b-84cc-4cb7-a099-7cbdb9ccc67c"
+JET = "c3a54625-cd67-489e-a8e7-0a5a0ff4e31b"
+OrdinaryDiffEqSSPRK = "669c94d9-1f4b-4b64-b377-1aa079aa2388"
+OrdinaryDiffEqTsit5 = "b1df2697-797e-41e3-8120-5422d3b24e4a"
+SciMLTesting = "09d9d899-5365-40a9-917a-5f67fddea283"
+Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
+
+[sources]
+OrdinaryDiffEqSSPRK = {path = "../../../OrdinaryDiffEqSSPRK"}
+OrdinaryDiffEqTsit5 = {path = "../../../OrdinaryDiffEqTsit5"}
+
+[compat]
+Aqua = "0.8"
+GlobalDiffEq = "1"
+JET = "0.9, 0.10, 0.11"
+OrdinaryDiffEqSSPRK = "2"
+OrdinaryDiffEqTsit5 = "2"
+SciMLTesting = "2.1"
+Test = "1"

lib/GlobalDiffEq/test/qa/qa.jl

@@ -0,0 +1,45 @@
+using SciMLTesting, GlobalDiffEq, Test
+using JET
+using OrdinaryDiffEqTsit5, OrdinaryDiffEqSSPRK
+
+# `@reexport using DiffEqBase` republishes DiffEqBase's API; those names are
+# documented and rendered at their owning packages, not in the OrdinaryDiffEq
+# manual's GlobalDiffEq page, so only GlobalDiffEq's own names are held to the
+# rendered-docs check.
+reexported_names = Tuple(
+    filter(public_api_names(GlobalDiffEq)) do name
+        object = getfield(GlobalDiffEq, name)
+        !(object isa Union{Function, Type, Module}) ||
+            parentmodule(object) !== GlobalDiffEq
+    end
+)
+
+run_qa(
+    GlobalDiffEq;
+    explicit_imports = true,
+    # GlobalDiffEq's rendered documentation lives in the monorepo docs, two
+    # directories up from the sublibrary root.
+    api_docs_kwargs = (;
+        rendered = true,
+        docs_src = joinpath(dirname(dirname(pkgdir(GlobalDiffEq))), "docs", "src"),
+        rendered_ignore = reexported_names,
+    ),
+    ei_kwargs = (;
+        all_qualified_accesses_are_public = (;
+            # `SciMLBase.__solve` is SciMLBase's internal solve entry point (not
+            # part of the public API); GlobalDiffEq overloads it via its owner
+            # SciMLBase.
+            ignore = (:__solve,),
+        ),
+    ),
+    # `@reexport using DiffEqBase` deliberately reexports DiffEqBase's API, so
+    # `DiffEqBase`, `ODEProblem`, and `solve` are inherently implicit. Tracked in
+    # https://github.com/SciML/GlobalDiffEq.jl/issues/53
+    ei_broken = (:no_implicit_imports,),
+)
+
+@testset "GlobalRichardson static analysis" begin
+    @test_opt GlobalRichardson(SSPRK33())
+    @test GlobalRichardson{typeof(SSPRK33())} <: GlobalDiffEq.GlobalDiffEqAlgorithm
+    @test GlobalRichardson(Tsit5()) isa GlobalDiffEq.GlobalDiffEqAlgorithm
+end

lib/GlobalDiffEq/test/runtests.jl

@@ -0,0 +1,26 @@
+using SafeTestsets
+using Test
+using SciMLTesting
+
+run_tests(;
+    env = "ODEDIFFEQ_TEST_GROUP",
+    core = function ()
+        @safetestset "Basic functionality" begin
+            include("basic_functionality_tests.jl")
+        end
+        @safetestset "Algorithm traits forwarding" begin
+            include("algorithm_traits_tests.jl")
+        end
+        return @safetestset "BigFloat support" begin
+            include("bigfloat_tests.jl")
+        end
+    end,
+    # QA group: Aqua + ExplicitImports + JET via SciMLTesting's run_qa. Runs in its
+    # own sub-env (test/qa/Project.toml). Kept out of the curated "All" (env-bearing),
+    # matching the prior JET group's exclusion.
+    qa = (;
+        env = joinpath(@__DIR__, "qa"),
+        body = joinpath(@__DIR__, "qa", "qa.jl"),
+    ),
+    all = ["Core"],
+)

lib/GlobalDiffEq/test/test_groups.toml

@@ -0,0 +1,5 @@
+[Core]
+versions = ["lts", "1", "pre"]
+
+[QA]
+versions = ["lts", "1"]