Use subpackages for non-ODEs

Preparing for the DifferentialEquations.jl split SciML/DifferentialEquations.jl#1086

Author: ChrisRackauckas

docs/src/examples/diffusion_implicit_heat_equation.md

@@ -16,8 +16,7 @@ First, we'll use / import some packages:
 import Plots
 import LinearAlgebra
 import DiffEqBase
-import OrdinaryDiffEq as ODE
-import OrdinaryDiffEq: SplitODEProblem, solve, IMEXEuler
+import OrdinaryDiffEqBDF as BDF
 import SciMLBase, SciMLOperators
 ```
 
@@ -179,14 +178,14 @@ params = (; n)
 
 tspan = (FT(0), N_t * FT(Δt))
 
-prob = SplitODEProblem(SciMLOperators.MatrixOperator(D),
+prob = SciMLBase.SplitODEProblem(SciMLOperators.MatrixOperator(D),
     rhs!,
     T,
     tspan,
     params)
-alg = IMEXEuler()
+alg = BDF.IMEXEuler()
 println("Solving...")
-sol = solve(prob,
+sol = SciMLBase.solve(prob,
     alg,
     dt = Δt,
     saveat = range(FT(0), N_t * FT(Δt), length = 5),

docs/src/examples/kepler_problem.md

@@ -93,7 +93,7 @@ Both Runge-Kutta-Nyström and Runge-Kutta integrator do not conserve energy nor
 In this example, we know that energy and angular momentum should be conserved. We can achieve this through manifold projection. As the name implies, it is a procedure to project the ODE solution to a manifold. Let's start with a base case, where manifold projection isn't being used.
 
 ```@example kepler
-import DiffEqCallbacks
+import DiffEqCallbacks as CB
 
 function plot_orbit2(sol)
     Plots.plot(sol, vars = (1, 2), lab = "Orbit", title = "Kepler Problem Solution")
@@ -127,7 +127,7 @@ function first_integrals_manifold(residual, u, p, t)
     residual[3:4] .= initial_first_integrals[2] - L(u[1:2], u[3:4])
 end
 
-cb = DiffEqCallbacks.ManifoldProjection(first_integrals_manifold, autodiff = NLS.AutoForwardDiff())
+cb = CB.ManifoldProjection(first_integrals_manifold, autodiff = NLS.AutoForwardDiff())
 sol5 = ODE.solve(prob2, ODE.RK4(), dt = 1 // 5, adaptive = false, callback = cb)
 analysis_plot2(sol5, H, L)
 ```
@@ -139,7 +139,7 @@ function energy_manifold(residual, u, p, t)
     residual[1:2] .= initial_first_integrals[1] - H(u[1:2], u[3:4])
     residual[3:4] .= 0
 end
-energy_cb = DiffEqCallbacks.ManifoldProjection(energy_manifold, autodiff = NLS.AutoForwardDiff())
+energy_cb = CB.ManifoldProjection(energy_manifold, autodiff = NLS.AutoForwardDiff())
 sol6 = ODE.solve(prob2, ODE.RK4(), dt = 1 // 5, adaptive = false, callback = energy_cb)
 analysis_plot2(sol6, H, L)
 ```
@@ -151,7 +151,7 @@ function angular_manifold(residual, u, p, t)
     residual[1:2] .= initial_first_integrals[2] - L(u[1:2], u[3:4])
     residual[3:4] .= 0
 end
-angular_cb = DiffEqCallbacks.ManifoldProjection(angular_manifold, autodiff = NLS.AutoForwardDiff())
+angular_cb = CB.ManifoldProjection(angular_manifold, autodiff = NLS.AutoForwardDiff())
 sol7 = ODE.solve(prob2, ODE.RK4(), dt = 1 // 5, adaptive = false, callback = angular_cb)
 analysis_plot2(sol7, H, L)
 ```

docs/src/solvers/steady_state_solve.md

@@ -41,9 +41,9 @@ large time steps as the steady state approaches.
 Example usage:
 
 ```julia
-import DifferentialEquations as DE
+import SteadyStateDiffEq as SS, DifferentialEquations as DE
 import Sundials
-sol = DE.solve(prob, DE.SSRootfind())
-sol = DE.solve(prob, DE.DynamicSS(DE.Tsit5()))
-sol = DE.solve(prob, DE.DynamicSS(Sundials.CVODE_BDF()), dt = 1.0)
+sol = SS.solve(prob, SS.SSRootfind())
+sol = SS.solve(prob, SS.DynamicSS(DE.Tsit5()))
+sol = SS.solve(prob, SS.DynamicSS(Sundials.CVODE_BDF()), dt = 1.0)
 ```

docs/src/tutorials/dde_example.md

@@ -36,7 +36,7 @@ no delays are written as in the ODE.
 Thus, the function for this model is given by:
 
 ```@example dde
-import DifferentialEquations as DE
+import DelayDiffEq as DDE, DifferentialEquations as DE
 function bc_model(du, u, h, p, t)
     p0, q0, v0, d0, p1, q1, v1, d1, d2, beta0, beta1, tau = p
     hist3 = h(p, t - tau)[3]
@@ -50,7 +50,7 @@ end
 Now we build a `DDEProblem`. The signature
 
 ```julia
-prob = DE.DDEProblem(f, u0, h, tspan, p = SciMLBase.NullParameters();
+prob = DDE.DDEProblem(f, u0, h, tspan, p = SciMLBase.NullParameters();
     constant_lags = [], dependent_lags = [], kwargs...)
 ```
 
@@ -88,28 +88,28 @@ p = (p0, q0, v0, d0, p1, q1, v1, d1, d2, beta0, beta1, tau)
 tspan = (0.0, 10.0)
 u0 = [1.0, 1.0, 1.0]
 
-prob = DE.DDEProblem(bc_model, u0, h, tspan, p; constant_lags = lags)
+prob = DDE.DDEProblem(bc_model, u0, h, tspan, p; constant_lags = lags)
 ```
 
 An efficient way to solve this problem (given the constant lags) is with the
 MethodOfSteps solver. Through the magic that is Julia, it translates an OrdinaryDiffEq.jl
 ODE solver method into a method for delay differential equations, which is highly
-efficient due to sweet compiler magic. A good choice is the order 5 method `DE.Tsit5()`:
+efficient due to sweet compiler magic. A good choice is the order 5 method `DDE.Tsit5()`:
 
 ```@example dde
-alg = DE.MethodOfSteps(DE.Tsit5())
+alg = DDE.MethodOfSteps(DE.Tsit5())
 ```
 
-For lower tolerance solving, one can use the `DE.BS3()` algorithm to good
+For lower tolerance solving, one can use the `DDE.BS3()` algorithm to good
 effect (this combination is similar to the MATLAB `dde23`, but more efficient
-tableau), and for high tolerances the `DE.Vern6()` algorithm will give a 6th order
+tableau), and for high tolerances the `DDE.Vern6()` algorithm will give a 6th order
 solution.
 
 To solve the problem with this algorithm, we do the same thing we'd do with other
 methods on the common interface:
 
 ```@example dde
-sol = DE.solve(prob, alg)
+sol = DDE.solve(prob, alg)
 ```
 
 Note that everything available to OrdinaryDiffEq.jl can be used here, including
@@ -189,20 +189,20 @@ You might have noticed DifferentialEquations.jl allows you to solve problems
 with undeclared delays, since you can interpolate `h` at any value. This is
 a feature, but use it with caution. Undeclared delays can increase the error
 in the solution. It's recommended that you use a method with a residual control,
-such as `MethodOfSteps(DE.RK4())` whenever there are undeclared delays. With this,
+such as `MethodOfSteps(DDE.RK4())` whenever there are undeclared delays. With this,
 you can use interpolated derivatives, solve functional differential equations
 by using quadrature on the interpolant, etc. However, note that residual control
 solves with a low level of accuracy, so the tolerances should be made very small,
 and the solution should not be trusted for more than 2-3 decimal places.
 
-Note: `MethodOfSteps(DE.RK4())` with undeclared delays is similar to MATLAB's
+Note: `MethodOfSteps(DDE.RK4())` with undeclared delays is similar to MATLAB's
 `ddesd`. Thus, for example, the following is similar to solving the example
 from above with residual control:
 
 ```@example dde
-prob = DE.DDEProblem(bc_model, u0, h, tspan)
-alg = DE.MethodOfSteps(DE.RK4())
-sol = DE.solve(prob, alg)
+prob = DDE.DDEProblem(bc_model, u0, h, tspan)
+alg = DDE.MethodOfSteps(DE.RK4())
+sol = DDE.solve(prob, alg)
 ```
 
 Note that this method can solve problems with state-dependent delays.
@@ -220,9 +220,9 @@ We can solve the above problem with dependent delay tracking by declaring the
 dependent lags and solving with a `MethodOfSteps` algorithm:
 
 ```@example dde
-prob = DE.DDEProblem(bc_model, u0, h, tspan; dependent_lags = ((u, p, t) -> tau,))
-alg = DE.MethodOfSteps(DE.Tsit5())
-sol = DE.solve(prob, alg)
+prob = DDE.DDEProblem(bc_model, u0, h, tspan; dependent_lags = ((u, p, t) -> tau,))
+alg = DDE.MethodOfSteps(DE.Tsit5())
+sol = DDE.solve(prob, alg)
 ```
 
 Here, we treated the single lag `t-tau` as a state-dependent delay. Of course, you

docs/src/tutorials/rode_example.md

@@ -15,15 +15,15 @@ du = f(u,p,t,W)dt
 where ``f(u,p,t,W)=2u\sin(W)`` and ``W(t)`` is a Wiener process (Gaussian process).
 
 ```@example rode
-import DifferentialEquations as DE
+import StochasticDiffEq as SDE
 import Plots
 function f3(u, p, t, W)
     2u * sin(W)
 end
 u0 = 1.00
 tspan = (0.0, 5.0)
-prob = DE.RODEProblem(f3, u0, tspan)
-sol = DE.solve(prob, DE.RandomEM(); dt = 1 / 100)
+prob = SDE.RODEProblem(f3, u0, tspan)
+sol = SDE.solve(prob, SDE.RandomEM(); dt = 1 / 100)
 Plots.plot(sol)
 ```
 
@@ -38,16 +38,16 @@ As with the other problem types, there is an in-place version which is more
 efficient for systems. The signature is `f(du,u,p,t,W)`. For example,
 
 ```@example rode2
-import DifferentialEquations as DE
+import StochasticDiffEq as SDE
 import Plots
 function f(du, u, p, t, W)
     du[1] = 2u[1] * sin(W[1] - W[2])
     du[2] = -2u[2] * cos(W[1] + W[2])
 end
 u0 = [1.00; 1.00]
 tspan = (0.0, 5.0)
-prob = DE.RODEProblem(f, u0, tspan)
-sol = DE.solve(prob, DE.RandomEM(); dt = 1 / 100)
+prob = SDE.RODEProblem(f, u0, tspan)
+sol = SDE.solve(prob, SDE.RandomEM(); dt = 1 / 100)
 Plots.plot(sol)
 ```
 
@@ -56,15 +56,15 @@ you can use the `rand_prototype` keyword to explicitly set the size of the
 random process:
 
 ```@example rode3
-import DifferentialEquations as DE
+import StochasticDiffEq as SDE
 import Plots
 function f(du, u, p, t, W)
     du[1] = -2W[3] * u[1] * sin(W[1] - W[2])
     du[2] = -2u[2] * cos(W[1] + W[2])
 end
 u0 = [1.00; 1.00]
 tspan = (0.0, 5.0)
-prob = DE.RODEProblem(f, u0, tspan; rand_prototype = zeros(3))
-sol = DE.solve(prob, DE.RandomEM(); dt = 1 / 100)
+prob = SDE.RODEProblem(f, u0, tspan; rand_prototype = zeros(3))
+sol = SDE.solve(prob, SDE.RandomEM(); dt = 1 / 100)
 Plots.plot(sol)
 ```