diff --git a/docs/src/examples/solvers/Unicycle robot.jl b/docs/src/examples/solvers/Unicycle robot.jl index c1aabb706..73d0efadb 100644 --- a/docs/src/examples/solvers/Unicycle robot.jl +++ b/docs/src/examples/solvers/Unicycle robot.jl @@ -48,9 +48,10 @@ x_upp = -x_low @variable(model, -1 <= u[1:2] <= 1) # Define the dynamics -@constraint(model, ∂(x[1]) == x[1] + u[1] * cos(x[3])) -@constraint(model, ∂(x[2]) == x[2] + u[1] * sin(x[3])) -@constraint(model, ∂(x[3]) == x[3] + u[2]) +@constraint(model, Δ(x[1]) == x[1] + u[1] * cos(x[3])) +@constraint(model, Δ(x[2]) == x[2] + u[1] * sin(x[3])) +#@constraint(model, Δ(x[3]) == rem(x[3] + u[2], 2 * π)) +@constraint(model, Δ(x[3]) == x[3] + u[2]) # Define the initial and target sets x_initial = [1.0, -1.7, 0.0] @@ -129,23 +130,6 @@ end # ### Definition of the abstraction -# First we need to set the mode of the abstraction to `DIONYSOS.Optim.Abstraction.UniformGridAbstraction.DSICRETE_TIME`: -set_attribute( - model, - "approx_mode", - Dionysos.Optim.Abstraction.UniformGridAbstraction.DISCRETE_TIME, -) - -# We define the system map $f$: -function sys_map(x, u, _) - return StaticArrays.SVector{3}( - x[1] + u[1] * cos(x[3]), - x[2] + u[1] * sin(x[3]), - (x[3] + u[2]) % (2 * π), - ) -end -set_attribute(model, "system_map", sys_map) - # We define the growth bound function of $f$: function growth_bound(r, u, _) β = u[1] * r[3] diff --git a/src/MOI_wrapper.jl b/src/MOI_wrapper.jl index 31ec06bbd..f1881cf49 100644 --- a/src/MOI_wrapper.jl +++ b/src/MOI_wrapper.jl @@ -7,10 +7,13 @@ import Symbolics @enum(VariableType, INPUT, STATE, MODE) +@enum(TimeType, UNKNOWN, CONTINUOUS, DISCRETE) + mutable struct Optimizer <: MOI.AbstractOptimizer inner::Any nlp_model::Any evaluator::Union{Nothing, MOI.Nonlinear.Evaluator} + time_type::TimeType variable_values::Vector{Float64} derivative_values::Vector{Float64} variable_types::Vector{VariableType} @@ -29,6 +32,7 @@ mutable struct Optimizer <: MOI.AbstractOptimizer MOI.instantiate(Dionysos.Optim.Abstraction.UniformGridAbstraction.Optimizer), nothing, nothing, + UNKNOWN, Float64[], Float64[], VariableType[], @@ -49,6 +53,7 @@ end MOI.is_empty(model::Optimizer) = isempty(model.variable_types) function MOI.empty!(model::Optimizer) + model.time_type = UNKNOWN empty!(model.variable_values) empty!(model.derivative_values) empty!(model.variable_types) @@ -168,7 +173,21 @@ function MOI.add_constraint( if lhs isa MOI.ScalarNonlinearFunction && length(lhs.args) == 1 var = lhs.args[] if var isa MOI.VariableIndex - if lhs.head == :diff + if lhs.head == :∂ + if model.time_type == DISCRETE + error("Cannot add constraint with `∂` since you already added a constraint with `Δ`.") + end + model.time_type = CONTINUOUS + model.dynamic[var.value] = rhs + model.nonlinear_index += 1 + return MOI.ConstraintIndex{typeof(func), typeof(set)}( + model.nonlinear_index, + ) + elseif lhs.head == :Δ + if model.time_type == CONTINUOUS + error("Cannot add constraint with `Δ` since you already added a constraint with `∂`.") + end + model.time_type = DISCRETE model.dynamic[var.value] = rhs model.nonlinear_index += 1 return MOI.ConstraintIndex{typeof(func), typeof(set)}( @@ -303,13 +322,23 @@ function system( ) _U_ = Dionysos.Utils.HyperRectangle(_svec(model.lower, u_idx), _svec(model.upper, u_idx)) - return MathematicalSystems.ConstrainedBlackBoxControlContinuousSystem( - dynamic(model, x_idx, u_idx), - Dionysos.Utils.get_dims(_X_), - Dionysos.Utils.get_dims(_U_), - _X_, - _U_, - ) + if model.time_type == CONTINUOUS + return MathematicalSystems.ConstrainedBlackBoxControlContinuousSystem( + dynamic(model, x_idx, u_idx), + Dionysos.Utils.get_dims(_X_), + Dionysos.Utils.get_dims(_U_), + _X_, + _U_, + ) + else + return MathematicalSystems.ConstrainedBlackBoxControlDiscreteSystem( + dynamic(model, x_idx, u_idx), + Dionysos.Utils.get_dims(_X_), + Dionysos.Utils.get_dims(_U_), + _X_, + _U_, + ) + end end function problem(model::Optimizer) @@ -400,10 +429,13 @@ function MOI.set(model::Optimizer, attr::MOI.RawOptimizerAttribute, value) return MOI.set(model.inner, attr, value) end -export ∂, final, start +export ∂, Δ, final, start function _diff end -∂ = JuMP.NonlinearOperator(_diff, :diff) +∂ = JuMP.NonlinearOperator(_diff, :∂) + +function _delta end +Δ = JuMP.NonlinearOperator(_delta, :Δ) function _final end final = JuMP.NonlinearOperator(_final, :final) @@ -411,6 +443,9 @@ final = JuMP.NonlinearOperator(_final, :final) function _start end start = JuMP.NonlinearOperator(_start, :start) +function rem end +rem_op = JuMP.NonlinearOperator(rem, :rem) + # Type piracy function JuMP.parse_constraint_call( error_fn::Function, @@ -429,3 +464,7 @@ function JuMP.parse_constraint_call( build_call = :(build_constraint($error_fn, $f, $(OuterSet)($set))) return parse_code, build_call end + +function Base.rem(x::JuMP.AbstractJuMPScalar, d) + return rem_op(x, d) +end diff --git a/src/optim/abstraction/uniform_grid_abstraction.jl b/src/optim/abstraction/uniform_grid_abstraction.jl index f0021a335..37e53231c 100644 --- a/src/optim/abstraction/uniform_grid_abstraction.jl +++ b/src/optim/abstraction/uniform_grid_abstraction.jl @@ -4,6 +4,8 @@ module UniformGridAbstraction import StaticArrays +import MathematicalSystems + import Dionysos const DI = Dionysos const UT = DI.Utils @@ -12,7 +14,7 @@ const ST = DI.System const SY = DI.Symbolic const PR = DI.Problem -@enum ApproxMode GROWTH LINEARIZED DELTA_GAS DISCRETE_TIME +@enum ApproxMode GROWTH LINEARIZED DELTA_GAS using JuMP @@ -39,8 +41,7 @@ mutable struct Optimizer{T} <: MOI.AbstractOptimizer approx_mode::ApproxMode ## for the discrete time (no need for the jacobian_bound) - #but for the system_map, growthbound_map and sys_inv_map - system_map::Union{Nothing, Function} + # but for the `growthbound_map` and `sys_inv_map` growthbound_map::Union{Nothing, Function} sys_inv_map::Union{Nothing, Function} @@ -63,7 +64,6 @@ mutable struct Optimizer{T} <: MOI.AbstractOptimizer GROWTH, nothing, nothing, - nothing, ) end end @@ -97,23 +97,7 @@ function build_abstraction(concrete_system, model) error("Please set the `time_step`.") end - model.discretized_system = if model.approx_mode == GROWTH - if isnothing(model.jacobian_bound) - error("Please set the `jacobian_bound`.") - end - Dionysos.System.NewControlSystemGrowthRK4( - model.time_step, - concrete_system.f, - model.jacobian_bound, - noise, - noise, - model.num_sub_steps_system_map, - model.num_sub_steps_growth_bound, - ) - elseif model.approx_mode == DISCRETE_TIME - if isnothing(model.system_map) - error("Please set the `system_map`.") - end + model.discretized_system = if concrete_system isa MathematicalSystems.ConstrainedBlackBoxControlDiscreteSystem if isnothing(model.growthbound_map) error("Please set the `growthbound_map`.") end @@ -124,10 +108,23 @@ function build_abstraction(concrete_system, model) model.time_step, noise, noise, - model.system_map, + concrete_system.f, model.growthbound_map, model.sys_inv_map, ) + elseif model.approx_mode == GROWTH + if isnothing(model.jacobian_bound) + error("Please set the `jacobian_bound`.") + end + Dionysos.System.NewControlSystemGrowthRK4( + model.time_step, + concrete_system.f, + model.jacobian_bound, + noise, + noise, + model.num_sub_steps_system_map, + model.num_sub_steps_growth_bound, + ) elseif model.approx_mode == LINEARIZED Dionysos.System.NewControlSystemLinearizedRK4( model.time_step, diff --git a/src/system/controlsystem.jl b/src/system/controlsystem.jl index cd5acca73..ba23fdc58 100644 --- a/src/system/controlsystem.jl +++ b/src/system/controlsystem.jl @@ -276,7 +276,7 @@ end ############################################ ############################################ -struct SymbolicSystem{} +struct SymbolicSystem fsymbolicT::Any fsymbolic::Any Ts::Any