Switch to OrdinaryDiffEqCore.jl (#141)

JoshuaLampert · web-flow · commit 4dc6f737cd4b · 2025-02-24T21:12:48.000+01:00
* switch to OrdinaryDiffEqCore.jl

* bump compat of SciMLBase.j to v2.68

* bump version of SymbolicIndexingInterface.jl to v0.3.36

* drop support for FastBroadcast v0.3

* do not import OrdinaryDiffEq

* add missing compat bounds

* set compat bounds of stdlibs to 1

* use OrdinaryDiffEq.jl subpackages in tests

* use OrdinaryDiffEq.jl subpackages

* fix another using OrdinaryDiffEq

* use OrdinaryDiffEq.jl subpackages in docs

* add compat for OrdinaryDiffEqRosenbrock.jl in docs

* fix typo in docs

* fix step!

* use ADTypes.jl to fix deprecation warnings

* format

* apply new formatter version

* bump compat of FastBroadcast.jl to v0.3.1 (v0.3.0 is yanked)

* fix typo

* bump FastBroadcast to v0.3.4

* lower compat bound of OrdinaryDiffEqCore.jl to v1.6

* bump compat bound of Reexport.jl to v1.2.2

* lower compat bound of OrdinaryDiffEqRosenbrock.jl v1.3

* bump compat bound of MuladdMacro.jl to v0.2.4

* bump compat for StaticArrays.jl to v1.9.7

* bump compat for StaticArrays.jl to v1.9.7

* bump compat for FastBroadcast.jl v0.3.5

* fix get_tmp_cache

* fix solve!

* bump compat of LinearSolve.jl to v2.32

* also bump LinearSolve.jl to v2.32 in docs an test

* bump compat of Plots.jl to v1.25.11

* bump compat of OrdinaryDiffEqCore.jl to v1.15.1

* bump compat of OrdinaryDiffEqCore.jl to v1.16

* bump compat of OrdinaryDiffEqRosenbrock to v1.4

* bump compat of OrdinaryDiffEqSDIRK.jl to v1.2
diff --git a/.github/workflows/Downgrade.yml b/.github/workflows/Downgrade.yml
@@ -38,7 +38,7 @@ jobs:
       - uses: julia-actions/cache@v2
       - uses: julia-actions/julia-downgrade-compat@v1
         with:
-          skip: LinearAlgebra,SparseArrays,Statistics
+          skip: LinearAlgebra,SparseArrays,Statistics,Test
           projects: ., test
       - uses: julia-actions/julia-buildpkg@v1
         env:
diff --git a/.github/workflows/FormatCheck.yml b/.github/workflows/FormatCheck.yml
@@ -28,7 +28,7 @@ jobs:
         #
         # julia  -e 'using Pkg; Pkg.add(PackageSpec(name="JuliaFormatter", version="0.13.0"))'
         run: |
-          julia  -e 'using Pkg; Pkg.add(PackageSpec(name = "JuliaFormatter", version="1.0.45"))'
+          julia  -e 'using Pkg; Pkg.add(PackageSpec(name = "JuliaFormatter", version="1.0.60"))'
           julia  -e 'using JuliaFormatter; format(".", verbose = true)'
       - name: Format check
         run: |
diff --git a/Project.toml b/Project.toml
@@ -8,7 +8,7 @@ FastBroadcast = "7034ab61-46d4-4ed7-9d0f-46aef9175898"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 LinearSolve = "7ed4a6bd-45f5-4d41-b270-4a48e9bafcae"
 MuladdMacro = "46d2c3a1-f734-5fdb-9937-b9b9aeba4221"
-OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
+OrdinaryDiffEqCore = "bbf590c4-e513-4bbe-9b18-05decba2e5d8"
 Reexport = "189a3867-3050-52da-a836-e630ba90ab69"
 SciMLBase = "0bca4576-84f4-4d90-8ffe-ffa030f20462"
 SimpleUnPack = "ce78b400-467f-4804-87d8-8f486da07d0a"
@@ -17,15 +17,15 @@ StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
 SymbolicIndexingInterface = "2efcf032-c050-4f8e-a9bb-153293bab1f5"
 
 [compat]
-FastBroadcast = "0.2.5, 0.3"
-LinearAlgebra = "1.7"
-LinearSolve = "2.21"
-MuladdMacro = "0.2.1"
-OrdinaryDiffEq = "6.59"
-Reexport = "1"
-SciMLBase = "2"
+FastBroadcast = "0.3.5"
+LinearAlgebra = "1"
+LinearSolve = "2.32"
+MuladdMacro = "0.2.4"
+OrdinaryDiffEqCore = "1.16"
+Reexport = "1.2.2"
+SciMLBase = "2.68"
 SimpleUnPack = "1"
-SparseArrays = "1.7"
-StaticArrays = "1.5"
-SymbolicIndexingInterface = "0.2, 0.3"
+SparseArrays = "1"
+StaticArrays = "1.9.7"
+SymbolicIndexingInterface = "0.3.36"
 julia = "1.10"
diff --git a/docs/Project.toml b/docs/Project.toml
@@ -1,23 +1,27 @@
 [deps]
 BenchmarkTools = "6e4b80f9-dd63-53aa-95a3-0cdb28fa8baf"
+DiffEqBase = "2b5f629d-d688-5b77-993f-72d75c75574e"
+DiffEqCallbacks = "459566f4-90b8-5000-8ac3-15dfb0a30def"
 Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
 InteractiveUtils = "b77e0a4c-d291-57a0-90e8-8db25a27a240"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
 LinearSolve = "7ed4a6bd-45f5-4d41-b270-4a48e9bafcae"
-OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
+OrdinaryDiffEqRosenbrock = "43230ef6-c299-4910-a778-202eb28ce4ce"
 Pkg = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f"
 Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
 SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
 StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
 
 [compat]
 BenchmarkTools = "1"
+DiffEqBase = "6.160"
+DiffEqCallbacks = "4"
 Documenter = "1"
 InteractiveUtils = "1"
-LinearAlgebra = "1.7"
-LinearSolve = "2.21"
-OrdinaryDiffEq = "6.59"
+LinearAlgebra = "1"
+LinearSolve = "2.32"
+OrdinaryDiffEqRosenbrock = "1.4"
 Pkg = "1"
 Plots = "1"
-SparseArrays = "1.7"
-StaticArrays = "1.5"
+SparseArrays = "1"
+StaticArrays = "1.9.7"
diff --git a/docs/make.jl b/docs/make.jl
@@ -79,13 +79,13 @@ makedocs(modules = [PositiveIntegrators],
                  "Stratospheric reaction problem" => "stratospheric_reaction.md",
                  "Linear Advection" => "linear_advection.md",
                  "Heat Equation, Neumann BCs" => "heat_equation_neumann.md",
-                 "Heat Equation, Dirichlet BCs" => "heat_equation_dirichlet.md",
+                 "Heat Equation, Dirichlet BCs" => "heat_equation_dirichlet.md"
              ],
              "Troubleshooting, FAQ" => "faq.md",
              "API reference" => "api_reference.md",
              "Contributing" => "contributing.md",
              "Code of conduct" => "code_of_conduct.md",
-             "License" => "license.md",
+             "License" => "license.md"
          ])
 
 deploydocs(repo = "github.com/SKopecz/PositiveIntegrators.jl",
diff --git a/docs/src/heat_equation_dirichlet.md b/docs/src/heat_equation_dirichlet.md
@@ -235,7 +235,7 @@ versioninfo()
 println()
 
 using Pkg
-Pkg.status(["PositiveIntegrators", "SparseArrays", "KLU", "LinearSolve", "OrdinaryDiffEq"],
+Pkg.status(["PositiveIntegrators", "SparseArrays", "KLU", "LinearSolve"],
            mode=PKGMODE_MANIFEST)
 nothing # hide
 ```
diff --git a/docs/src/heat_equation_neumann.md b/docs/src/heat_equation_neumann.md
@@ -213,7 +213,7 @@ versioninfo()
 println()
 
 using Pkg
-Pkg.status(["PositiveIntegrators", "SparseArrays", "KLU", "LinearSolve", "OrdinaryDiffEq"],
+Pkg.status(["PositiveIntegrators", "SparseArrays", "KLU", "LinearSolve"],
            mode=PKGMODE_MANIFEST)
 nothing # hide
 ```
diff --git a/docs/src/index.md b/docs/src/index.md
@@ -72,7 +72,7 @@ Consequently the production matrix ``\mathbf P`` and destruction vector ``\mathb
 ```
 
 ```@setup LotkaVolterra
-import Pkg; Pkg.add("OrdinaryDiffEq");  Pkg.add("Plots")
+import Pkg; Pkg.add("Plots")
 ```
 To solve this PDS together with initial values ``u_1(0)=u_2(0)=2`` on the time domain ``(0,10)``, we first need to create a `PDSProblem`.
 ```@example LotkaVolterra
@@ -134,10 +134,6 @@ The corresponding production matrix ``\mathbf P`` is
 
 The following example shows how to implement the above SIR model with ``\beta=0.4, \gamma=0.04``, initial conditions ``S(0)=997, I(0)=3, R(0)=0`` and time domain ``(0, 100)`` using `ConservativePDSProblem` from [PositiveIntegrators.jl](https://github.com/SKopecz/PositiveIntegrators.jl).
 
-```@setup SIR
-import Pkg; Pkg.add("OrdinaryDiffEq");
-```
-
 ```@example SIR
 using PositiveIntegrators
 
diff --git a/docs/src/linear_advection.md b/docs/src/linear_advection.md
@@ -164,7 +164,7 @@ versioninfo()
 println()
 
 using Pkg
-Pkg.status(["PositiveIntegrators", "SparseArrays", "KLU", "LinearSolve", "OrdinaryDiffEq"],
+Pkg.status(["PositiveIntegrators", "SparseArrays", "KLU", "LinearSolve"],
            mode=PKGMODE_MANIFEST)
 nothing # hide
 ```
diff --git a/docs/src/npzd_model.md b/docs/src/npzd_model.md
@@ -14,7 +14,7 @@ Z' &= 0.5(1 - e^{-1.21P^2})Z - 0.01Z - 0.02Z,\\
 D' &= 0.05P + 0.02Z - 0.003D,
 \end{aligned}
 ```
-and we consider the initial conditions ``N=8``, ``P=2``, ``Z=1`` and ``D=4``. The time domain of interest is ``t\in[0,10]``. 
+and we consider the initial conditions ``N=8``, ``P=2``, ``Z=1`` and ``D=4``. The time domain of interest is ``t\in[0,10]``.
 
 The model can be represented as a conservative PDS with production terms
 ```math
@@ -28,7 +28,7 @@ whereby production terms not listed have the value zero. Since the PDS is conser
 
 ## Solution of the production-destruction system
 
-Now we are ready to define a [`ConservativePDSProblem`](@ref) and to solve this problem with a method of [PositiveIntegrators.jl](https://github.com/SKopecz/PositiveIntegrators.jl) or [OrdinaryDiffEq.jl](https://docs.sciml.ai/OrdinaryDiffEq/stable/). 
+Now we are ready to define a [`ConservativePDSProblem`](@ref) and to solve this problem with a method of [PositiveIntegrators.jl](https://github.com/SKopecz/PositiveIntegrators.jl) or [OrdinaryDiffEq.jl](https://docs.sciml.ai/OrdinaryDiffEq/stable/).
 
 As mentioned above, we will try different approaches to solve this PDS and compare their efficiency. These are
 1. an out-of-place implementation with standard (dynamic) matrices and vectors,
@@ -70,7 +70,7 @@ sol_oop = solve(prob_oop, MPRK43I(1.0, 0.5))
 
 nothing #hide
 ```
-Plotting the solution shows that the components ``N`` and ``P`` are in danger of becoming negative. 
+Plotting the solution shows that the components ``N`` and ``P`` are in danger of becoming negative.
 ```@example NPZD
 using Plots
 
@@ -203,7 +203,7 @@ versioninfo()
 println()
 
 using Pkg
-Pkg.status(["PositiveIntegrators", "StaticArrays", "LinearSolve", "OrdinaryDiffEq"],
+Pkg.status(["PositiveIntegrators", "StaticArrays", "LinearSolve"],
            mode=PKGMODE_MANIFEST)
 nothing # hide
 ```
diff --git a/docs/src/robertson.md b/docs/src/robertson.md
@@ -1,6 +1,6 @@
 # [Tutorial: Solution of Robertson problem](@id tutorial-robertson)
 
-In this tutorial we show that MPRK schemes can be used to integrate stiff problems. 
+In this tutorial we show that MPRK schemes can be used to integrate stiff problems.
 We also show how callbacks can be used to change the time step in non-adaptive schemes.
 
 ## Definition of the production-destruction system
@@ -13,7 +13,7 @@ u_2' &=  0.04u_1-10^4 u_2u_3-3⋅10^7 u_2^2, & u_2(0)&=0 \\
 u_3' &= 3⋅10^7 u_2^2, & u_3(0)&=0.
 \end{aligned}
 ```
-The time domain of interest is ``t\in[0,10^{11}]``, because of which some kind of adaptive time stepping is required. 
+The time domain of interest is ``t\in[0,10^{11}]``, because of which some kind of adaptive time stepping is required.
 
 The model can be represented as a conservative PDS with production terms
 ```math
@@ -35,8 +35,8 @@ Since this PDS consists of only three differential equations we provide an out-o
 using PositiveIntegrators, StaticArrays
 
 function prod(u, p, t)
-    @SMatrix [0.0 1e4*u[2]*u[3] 0.0; 
-              4e-2*u[1] 0.0 0.0; 
+    @SMatrix [0.0 1e4*u[2]*u[3] 0.0;
+              4e-2*u[1] 0.0 0.0;
               0.0 3e7*u[2]^2 0.0]
 end
 u0 = @SVector [1.0, 0.0, 0.0]
@@ -60,11 +60,12 @@ isnonnegative(sol)
 
 ### Using callbacks to solve the Robertson problem with non-adatpive schemes
 
-The `SSPMPRK43()` scheme is only available with fixed time stepping. With a scheme like this, it would take a huge amount of time to solve the Robertson problem, since the time step must be chosen very small to accurately solve the problem in its initial phase. However, the use of a `callback` allows us to modify the time step size after each step, which makes a solution with a fixed step method possible. 
+The `SSPMPRK43()` scheme is only available with fixed time stepping. With a scheme like this, it would take a huge amount of time to solve the Robertson problem, since the time step must be chosen very small to accurately solve the problem in its initial phase. However, the use of a `callback` allows us to modify the time step size after each step, which makes a solution with a fixed step method possible.
 
 In the following example the `callback` increases the time step size by a factor of 1.5 after each time step.
 ```@example robertson
-using OrdinaryDiffEq
+using DiffEqCallbacks
+using DiffEqBase
 
 stepsize_callback = DiscreteCallback(
     Returns(true), # adapt the step size after every time step
@@ -93,7 +94,7 @@ versioninfo()
 println()
 
 using Pkg
-Pkg.status(["PositiveIntegrators", "StaticArrays", "LinearSolve", "OrdinaryDiffEq"],
+Pkg.status(["PositiveIntegrators", "StaticArrays", "LinearSolve", "DiffEqCallbacks", "DiffEqBase"],
            mode=PKGMODE_MANIFEST)
 nothing # hide
 ```
diff --git a/docs/src/stratospheric_reaction.md b/docs/src/stratospheric_reaction.md
@@ -1,6 +1,6 @@
 # [Tutorial: Solution of a stratospheric reaction problem](@id tutorial-stratos)
 
-This tutorial is about the efficient solution of a stiff non-autonomous and non-conservative production-destruction systems (PDS) with a small number of differential equations. 
+This tutorial is about the efficient solution of a stiff non-autonomous and non-conservative production-destruction systems (PDS) with a small number of differential equations.
 We will compare the use of standard arrays and static arrays from [StaticArrays.jl](https://juliaarrays.github.io/StaticArrays.jl/stable/) and assess their efficiency.
 
 ## Definition of the production-destruction system
@@ -54,7 +54,7 @@ In addition, all production and destruction terms not listed have the value zero
 
 ## Solution of the production-destruction system
 
-Now we are ready to define a [`PDSProblem`](@ref) and to solve this problem with a method of [PositiveIntegrators.jl](https://github.com/SKopecz/PositiveIntegrators.jl) or [OrdinaryDiffEq.jl](https://docs.sciml.ai/OrdinaryDiffEq/stable/). 
+Now we are ready to define a [`PDSProblem`](@ref) and to solve this problem with a method of [PositiveIntegrators.jl](https://github.com/SKopecz/PositiveIntegrators.jl) or [OrdinaryDiffEq.jl](https://docs.sciml.ai/OrdinaryDiffEq/stable/).
 
 As mentioned above, we will try different approaches to solve this PDS and compare their efficiency. These are
 1. an out-of-place implementation with standard (dynamic) matrices and vectors,
@@ -138,7 +138,7 @@ sol_oop = solve(prob_oop, MPRK43I(1.0, 0.5))
 
 nothing #hide
 ```
-Plotting the solution shows that the components ``O¹ᴰ``, ``O`` and ``NO`` are in danger of becoming negative. 
+Plotting the solution shows that the components ``O¹ᴰ``, ``O`` and ``NO`` are in danger of becoming negative.
 ```@example stratreac
 using Plots
 
@@ -149,7 +149,7 @@ plot(sol_oop,
     xticks = (range(first(tspan), last(tspan), 4), range(12.0, 84.0, 4)),
     yguide=["O¹ᴰ" "O" "O₃" "O₂" "NO" "NO₂"],
     tickfontsize = 7,
-    legend = :none, 
+    legend = :none,
     widen = true
     )
 ```
@@ -253,7 +253,7 @@ plot(sol_ip,
     xticks = (range(first(tspan), last(tspan), 4), range(12.0, 84.0, 4)),
     yguide=["O¹ᴰ" "O" "O₃" "O₂" "NO" "NO₂"],
     tickfontsize = 7,
-    legend = :none, 
+    legend = :none,
     widen = true
     )
 ```
@@ -264,7 +264,7 @@ sol_oop.t ≈ sol_ip.t && sol_oop.u ≈ sol_ip.u
 ```
 
 ### Using static arrays
-For PDS with a small number of differential equations like the stratospheric reaction model the use of static arrays will be more efficient. To create a function which computes the production matrix and returns a static matrix, we only need to add the `@SMatrix` macro. Accordingly, we use the `@SVector` macro for the destruction vector. 
+For PDS with a small number of differential equations like the stratospheric reaction model the use of static arrays will be more efficient. To create a function which computes the production matrix and returns a static matrix, we only need to add the `@SMatrix` macro. Accordingly, we use the `@SVector` macro for the destruction vector.
 
 ```@example stratreac
 using StaticArrays
@@ -353,7 +353,7 @@ plot(sol_static,
     xticks = (range(first(tspan), last(tspan), 4), range(12.0, 84.0, 4)),
     yguide=["O¹ᴰ" "O" "O₃" "O₂" "NO" "NO₂"],
     tickfontsize = 7,
-    legend = :none, 
+    legend = :none,
     widen = true
     )
 ```
@@ -376,15 +376,15 @@ a2 = [0; 0; 0; 0; 1; 1] # second linear invariant
 
 p1 = plot(sol_oop.t, relerr_lininv(a1, u0, sol_oop))
 p2 = plot(sol_oop.t, relerr_lininv(a2, u0, sol_oop))
-plot(p1, p2, 
+plot(p1, p2,
     xticks = (range(first(tspan), last(tspan), 4), range(12.0, 84.0, 4)),
     legend = :none)
 ```
 
 In contrast to MPRK schemes, Runge-Kutta and Rosenbrock methods preserve all linear invariants, but are not guaranteed to generate nonnegative solutions.
 One way to enforce nonnegative solutions of such schemes is passing [`isnegative`](@ref) to the solver option [`isoutofdomain`](https://docs.sciml.ai/DiffEqDocs/stable/basics/common_solver_opts/). We show this using the Rosenbrock scheme `Rosenbrock23` as an example.
 ```@example stratreac
-using OrdinaryDiffEq
+using OrdinaryDiffEqRosenbrock
 
 sol_tmp = solve(prob_oop, Rosenbrock23());
 isnonnegative(sol_tmp)
@@ -430,7 +430,7 @@ versioninfo()
 println()
 
 using Pkg
-Pkg.status(["PositiveIntegrators", "StaticArrays", "LinearSolve", "OrdinaryDiffEq"],
+Pkg.status(["PositiveIntegrators", "StaticArrays", "LinearSolve", "OrdnaryDiffEqRosenbrock"],
            mode=PKGMODE_MANIFEST)
 nothing # hide
 ```
diff --git a/src/PositiveIntegrators.jl b/src/PositiveIntegrators.jl
@@ -18,25 +18,25 @@ using SciMLBase: AbstractODEFunction, NullParameters, FullSpecialize,
                  isinplace
 
 # TODO: Check imports and using statements below, reduce if possible
-using OrdinaryDiffEq: OrdinaryDiffEq, OrdinaryDiffEqAlgorithm, ODESolution
+using OrdinaryDiffEqCore: OrdinaryDiffEqCore, OrdinaryDiffEqAlgorithm, ODESolution
 
 using SymbolicIndexingInterface: SymbolicIndexingInterface
 
 using LinearSolve: LinearSolve, LinearProblem, LUFactorization, solve!
 
 import SciMLBase: interp_summary
 
-using OrdinaryDiffEq: @cache,
-                      OrdinaryDiffEqAdaptiveAlgorithm,
-                      OrdinaryDiffEqConstantCache, OrdinaryDiffEqMutableCache,
-                      False,
-                      _vec
-import OrdinaryDiffEq: alg_order, isfsal,
-                       calculate_residuals, calculate_residuals!,
-                       alg_cache, get_tmp_cache,
-                       initialize!, perform_step!,
-                       _ode_interpolant, _ode_interpolant!,
-                       get_fsalfirstlast
+using OrdinaryDiffEqCore: @cache,
+                          OrdinaryDiffEqAdaptiveAlgorithm,
+                          OrdinaryDiffEqConstantCache, OrdinaryDiffEqMutableCache,
+                          False,
+                          _vec
+import OrdinaryDiffEqCore: alg_order, isfsal,
+                           calculate_residuals, calculate_residuals!,
+                           alg_cache, get_tmp_cache,
+                           initialize!, perform_step!,
+                           _ode_interpolant, _ode_interpolant!,
+                           get_fsalfirstlast
 
 # 2. Export functionality defining the public API
 export PDSFunction, PDSProblem
diff --git a/test/Project.toml b/test/Project.toml
diff --git a/test/runtests.jl b/test/runtests.jl

Original file line number	Diff line number	Diff line change
`@@ -28,7 +28,7 @@ jobs:`
`28`	`28`	`#`
`29`	`29`	`# julia -e 'using Pkg; Pkg.add(PackageSpec(name="JuliaFormatter", version="0.13.0"))'`
`30`	`30`	`run: \|`
`31`		`- julia -e 'using Pkg; Pkg.add(PackageSpec(name = "JuliaFormatter", version="1.0.45"))'`
	`31`	`+ julia -e 'using Pkg; Pkg.add(PackageSpec(name = "JuliaFormatter", version="1.0.60"))'`
`32`	`32`	`julia -e 'using JuliaFormatter; format(".", verbose = true)'`
`33`	`33`	`- name: Format check`
`34`	`34`	`run: \|`