Skip to content

Commit

Permalink
Merge pull request #823 from JuliaReach/schillic/ode
Browse files Browse the repository at this point in the history
Replace DifferentialEquations by OrdinaryDiffEq
  • Loading branch information
schillic authored Apr 10, 2024
2 parents ff7dc77 + 1881eaf commit 4f8737a
Show file tree
Hide file tree
Showing 11 changed files with 33 additions and 30 deletions.
4 changes: 2 additions & 2 deletions docs/Project.toml
Original file line number Diff line number Diff line change
@@ -1,5 +1,4 @@
[deps]
DifferentialEquations = "0c46a032-eb83-5123-abaf-570d42b7fbaa"
DisplayAs = "0b91fe84-8a4c-11e9-3e1d-67c38462b6d6"
Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"
ExponentialUtilities = "d4d017d3-3776-5f7e-afef-a10c40355c18"
Expand All @@ -9,12 +8,12 @@ LaTeXStrings = "b964fa9f-0449-5b57-a5c2-d3ea65f4040f"
LazySets = "b4f0291d-fe17-52bc-9479-3d1a343d9043"
Literate = "98b081ad-f1c9-55d3-8b20-4c87d4299306"
MathematicalSystems = "d14a8603-c872-5ed3-9ece-53e0e82e39da"
OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
ReachabilityBase = "379f33d0-9447-4353-bd03-d664070e549f"
Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7"

[compat]
DifferentialEquations = "7"
DisplayAs = "0.1"
Documenter = "1"
ExponentialUtilities = "1"
Expand All @@ -24,6 +23,7 @@ LaTeXStrings = "1"
LazySets = "2"
Literate = "2"
MathematicalSystems = "0.11 - 0.13"
OrdinaryDiffEq = "6"
Plots = "1"
ReachabilityBase = "0.2.3"
Symbolics = "5"
6 changes: 3 additions & 3 deletions docs/src/man/faq.md
Original file line number Diff line number Diff line change
Expand Up @@ -123,20 +123,20 @@ DisplayAs.Text(DisplayAs.PNG(fig)) # hide

It is often necessary to plot a "bunch" of trajectories starting from a set of initial conditions.
The [parallel ensemble simulations](https://diffeq.sciml.ai/stable/features/ensemble/) capabilities
from the Julia [DifferentialEquations.jl](https://github.com/SciML/DifferentialEquations.jl) suite
from the Julia [OrdinaryDiffEq.jl](https://github.com/SciML/OrdinaryDiffEq.jl) suite
can be used to numerically simulate a given number of trajectories and it counts with state-of-the-art
algorithms for stiff and non-stiff ODEs as well as many other advanced features,
such as distributed computing, multi-threading and GPU support.
See [EnsembleAlgorithms](https://diffeq.sciml.ai/stable/features/ensemble/#EnsembleAlgorithms) for details.

As a simple example consider the scalar ODE ``x'(t) = 1.01x(t)`` with initial condition on the
interval ``x(0) ∈ [0, 0.5]``. To solve it using ensemble simulations, pass the `ensemble=true`
keyword argument to the solve function (if the library `DifferentialEquations` was not loaded in
keyword argument to the solve function (if the library `OrdinaryDiffEq` was not loaded in
your current session, an error is triggered). The number of trajectories can be specified
with the `trajectories` keyword argument.

```@example
using ReachabilityAnalysis, DifferentialEquations
using ReachabilityAnalysis, OrdinaryDiffEq
# formulate initial-value problem
prob = @ivp(x' = 1.01x, x(0) ∈ 0 .. 0.5)
Expand Down
2 changes: 1 addition & 1 deletion docs/src/tutorials/linear_methods/introduction.md
Original file line number Diff line number Diff line change
Expand Up @@ -72,7 +72,7 @@ plot!(trange, 2.0 * exp.(-trange), vars=(0, 1), c=:magenta, lab="")
```

The user-facing interface is designed to be intuitive and interactive, and it is
inspired by of other Julia packages such as DifferentialEquations.jl.
inspired by of other Julia packages such as [OrdinaryDiffEq.jl](https://github.com/SciML/OrdinaryDiffEq.jl).
Moreover, the library internals are written in a modular and composable way,
such that advanced users are able to modify and changed easily, or to compose with
other algorithms, different steps of the solution process.
Expand Down
4 changes: 2 additions & 2 deletions docs/src/tutorials/taylor_methods/introduction.md
Original file line number Diff line number Diff line change
Expand Up @@ -18,10 +18,10 @@ x'(t) = -x(t) ~ \sin(t),\qquad t ≥ 0.
```

Standard integration schemes fail to produce helpful solutions if the initial state is an interval. We illustrate this point
by solving the given differential equation with the `Tsit5` algorithm from [DifferentialEquations.jl](https://diffeq.sciml.ai/stable/) suite.
by solving the given differential equation with the `Tsit5` algorithm from [OrdinaryDiffEq.jl](https://github.com/SciML/OrdinaryDiffEq.jl) suite.

```@example nonlinear_univariate
using DifferentialEquations, IntervalArithmetic
using OrdinaryDiffEq, IntervalArithmetic
# initial condition
x₀ = [-1 .. 1]
Expand Down
2 changes: 1 addition & 1 deletion src/Continuous/solve.jl
Original file line number Diff line number Diff line change
Expand Up @@ -67,7 +67,7 @@ function solve(ivp::IVP{<:AbstractContinuousSystem_}, args...; kwargs...)
got_ensemble = get(kwargs, :ensemble, false)
dict = Dict{Symbol,Any}(:ensemble => nothing)
if got_ensemble
@requires DifferentialEquations
@requires OrdinaryDiffEq
ensemble_sol = _solve_ensemble(ivp, args...; kwargs...)
dict[:ensemble] = ensemble_sol
end
Expand Down
2 changes: 1 addition & 1 deletion src/Hybrid/solve.jl
Original file line number Diff line number Diff line change
Expand Up @@ -170,7 +170,7 @@ function solve(ivp::IVP{<:AbstractHybridSystem}, args...;

got_ensemble = get(kwargs, :ensemble, false)
if got_ensemble
@requires DifferentialEquations
@requires OrdinaryDiffEq
ensemble_sol = _solve_ensemble(ivp, args...; kwargs...)
dict = Dict{Symbol,Any}(:ensemble => ensemble_sol)
sol = ReachSolution(HFout, cpost, dict)
Expand Down
2 changes: 1 addition & 1 deletion src/Initialization/init.jl
Original file line number Diff line number Diff line change
Expand Up @@ -113,7 +113,7 @@ using Requires

function __init__()
# numerical differential equations suite
@require DifferentialEquations = "0c46a032-eb83-5123-abaf-570d42b7fbaa" include("init_DifferentialEquations.jl")
@require OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed" include("init_OrdinaryDiffEq.jl")

# sparse dynamic representation of multivariate polynomials
@require DynamicPolynomials = "7c1d4256-1411-5781-91ec-d7bc3513ac07" include("init_DynamicPolynomials.jl")
Expand Down
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
import .DifferentialEquations as DE
import .OrdinaryDiffEq as ODE
import Random

const DEFAULT_TRAJECTORIES = 10
Expand All @@ -10,8 +10,8 @@ const DEFAULT_TRAJECTORIES = 10
# =====================================

function _solve_ensemble(ivp::InitialValueProblem, args...;
trajectories_alg=DE.Tsit5(),
ensemble_alg=DE.EnsembleThreads(),
trajectories_alg=ODE.Tsit5(),
ensemble_alg=ODE.EnsembleThreads(),
inplace=true,
initial_states=nothing,
kwargs...)
Expand Down Expand Up @@ -39,8 +39,8 @@ function _solve_ensemble(ivp::InitialValueProblem, args...;
end

# formulate ensemble ODE problem
ensemble_prob = DE.ODEProblem(field, first(initial_states), tspan)
_prob_func(prob, i, repeat) = DE.remake(prob; u0=initial_states[i])
ensemble_prob = ODE.ODEProblem(field, first(initial_states), tspan)
_prob_func(prob, i, repeat) = ODE.remake(prob; u0=initial_states[i])

# choose tolerances
reltol = get(kwargs, :reltol, 1e-3)
Expand All @@ -49,14 +49,14 @@ function _solve_ensemble(ivp::InitialValueProblem, args...;
callback = get(kwargs, :callback, nothing)
dtmax = get(kwargs, :dtmax, Inf)

ensemble_prob = DE.EnsembleProblem(ensemble_prob; prob_func=_prob_func)
ensemble_prob = ODE.EnsembleProblem(ensemble_prob; prob_func=_prob_func)

if isnothing(callback)
result = DE.solve(ensemble_prob, trajectories_alg, ensemble_alg;
result = ODE.solve(ensemble_prob, trajectories_alg, ensemble_alg;
trajectories=trajectories, reltol=reltol,
abstol=abstol, dtmax=dtmax)
else
result = DE.solve(ensemble_prob, trajectories_alg, ensemble_alg;
result = ODE.solve(ensemble_prob, trajectories_alg, ensemble_alg;
trajectories=trajectories, reltol=reltol,
abstol=abstol, dtmax=dtmax, callback=callback)
end
Expand Down Expand Up @@ -90,7 +90,7 @@ function _solve_ensemble(ivp::InitialValueProblem{<:AbstractHybridSystem},
t0_g = tstart(tspan)
T = tend(tspan)

termination_action = (integrator) -> DE.terminate!(integrator)
termination_action = (integrator) -> ODE.terminate!(integrator)
use_discrete_callback = get(kwargs, :use_discrete_callback, false)
dtmax = get(kwargs, :dtmax, Inf)
if use_discrete_callback && isinf(dtmax)
Expand Down Expand Up @@ -118,11 +118,11 @@ function _solve_ensemble(ivp::InitialValueProblem{<:AbstractHybridSystem},
if use_discrete_callback
# use a discrete callback function based on membership
condition = (u, t, integrator) -> u I⁻
callback = DE.DiscreteCallback(condition, termination_action)
callback = ODE.DiscreteCallback(condition, termination_action)
else
# use a continuous callback function based on membership
condition = (u, t, integrator) -> u I⁻ ? -1 : 1
callback = DE.ContinuousCallback(condition, termination_action)
callback = ODE.ContinuousCallback(condition, termination_action)
end
trajectory = _solve_ensemble(ivp_loc, args...; initial_states=[x0],
callback=callback, tspan=t0 .. T,
Expand Down
6 changes: 4 additions & 2 deletions test/Project.toml
Original file line number Diff line number Diff line change
@@ -1,12 +1,13 @@
[deps]
Aqua = "4c88cf16-eb10-579e-8560-4a9242c79595"
CDDLib = "3391f64e-dcde-5f30-b752-e11513730f60"
DifferentialEquations = "0c46a032-eb83-5123-abaf-570d42b7fbaa"
Expokit = "a1e7a1ef-7a5d-5822-a38c-be74e1bb89f4"
Flowstar = "a8054ddd-9dca-4d20-8ffe-ae96ec1541f1"
IntervalArithmetic = "d1acc4aa-44c8-5952-acd4-ba5d80a2a253"
LazySets = "b4f0291d-fe17-52bc-9479-3d1a343d9043"
LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
Optim = "429524aa-4258-5aef-a3af-852621145aeb"
OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
Polyhedra = "67491407-f73d-577b-9b50-8179a7c68029"
SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf"
StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
Expand All @@ -18,11 +19,12 @@ Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
[compat]
Aqua = "0.8"
CDDLib = "0.5 - 0.9"
DifferentialEquations = "6.17, 7"
Expokit = "0.2"
Flowstar = "0.2.4"
IntervalArithmetic = "0.16 - 0.20, =0.20.9" # new versions require updates and are incompatible with dependencies
LazySets = "2.3"
Optim = "0.15 - 0.22, 1"
OrdinaryDiffEq = "6"
Polyhedra = "0.5 - 0.7"
StaticArrays = "0.12, 1"
Suppressor = "0.2"
Expand Down
12 changes: 6 additions & 6 deletions test/continuous/traces.jl
Original file line number Diff line number Diff line change
@@ -1,18 +1,18 @@
import DifferentialEquations
import OrdinaryDiffEq
using ReachabilityAnalysis: _solve_ensemble, ReachSolution

@testset "Simulation traces using DifferentialEquations.jl solvers" begin
@testset "Simulation traces using OrdinaryDiffEq.jl solvers" begin
prob = @ivp(x' = 1.0x, x(0) [1 / 2])

# `solve` extends CommonSolve.jl
sol = DifferentialEquations.solve(prob; tspan=(0.0, 1.0))
sol = OrdinaryDiffEq.solve(prob; tspan=(0.0, 1.0))
@test sol isa ReachSolution

# ensemble problems API
# ensemble problems API
sol = _solve_ensemble(prob; tspan=(0.0, 1.0))
@test sol isa DifferentialEquations.EnsembleSolution
@test sol isa OrdinaryDiffEq.EnsembleSolution
sol = _solve_ensemble(prob; T=1.0)
@test sol isa DifferentialEquations.EnsembleSolution
@test sol isa OrdinaryDiffEq.EnsembleSolution

# inplace (default)
sol = _solve_ensemble(prob; T=1.0, inplace=true)
Expand Down
1 change: 1 addition & 0 deletions test/hybrid/hybrid.jl
Original file line number Diff line number Diff line change
@@ -1,5 +1,6 @@
using ReachabilityAnalysis: _distribute, WaitingList, StateInLocation, state, location,
TimeInterval, DeterministicSwitching, NonDeterministicSwitching
import Optim

@testset "Hybrid utility functions" begin
prob, _ = bouncing_ball()
Expand Down

0 comments on commit 4f8737a

Please sign in to comment.