OptimalControl.jl
OptimalControl.jl is the core package of the control-toolbox ecosystem. Below, we group together the documentation of all the functions and types exported by OptimalControl.
Even if the following functions are prefixed by another package, such as CTFlows.Lift
, they can all be used with OptimalControl. In fact, all functions prefixed with another package are simply reexported. For example, Lift
is defined in CTFlows but accessible from OptimalControl.
julia> using OptimalControl
julia> F(x) = 2x
julia> H = Lift(F)
julia> x = 1
julia> p = 2
julia> H(x, p)
4
Exported functions and types
OptimalControl.OptimalControl
— ModuleOptimalControl module.
List of all the exported names:
*
Flow
Hamiltonian
HamiltonianLift
HamiltonianVectorField
@Lie
Lie
Lift
Model
ParsingError
Poisson
Solution
VectorField
available_methods
build_OCP_solution
constraint
constraints
constraints_violation
control
control_components
control_dimension
control_name
costate
criterion
@def
definition
direct_transcription
dual
dynamics
export_ocp_solution
final_time
final_time_name
get_build_examodel
has_fixed_final_time
has_fixed_initial_time
has_free_final_time
has_free_initial_time
has_lagrange_cost
has_mayer_cost
import_ocp_solution
infos
initial_time
initial_time_name
is_autonomous
iterations
lagrange
mayer
message
model
objective
set_initial_guess
solve
state
state_components
state_dimension
state_name
status
successful
time_grid
time_name
times
variable
variable_components
variable_dimension
variable_name
⋅
Documentation
Base.:*
— Method*(x, y...)
Multiplication operator.
Infix x*y*z*...
calls this function with all arguments, i.e. *(x, y, z, ...)
, which by default then calls (x*y) * z * ...
starting from the left.
Juxtaposition such as 2pi
also calls *(2, pi)
. Note that this operation has higher precedence than a literal *
. Note also that juxtaposition "0x..." (integer zero times a variable whose name starts with x
) is forbidden as it clashes with unsigned integer literals: 0x01 isa UInt8
.
Note that overflow is possible for most integer types, including the default Int
, when multiplying large numbers.
Examples
julia> 2 * 7 * 8
112
julia> *(2, 7, 8)
112
julia> [2 0; 0 3] * [1, 10] # matrix * vector
2-element Vector{Int64}:
2
30
julia> 1/2pi, 1/2*pi # juxtaposition has higher precedence
(0.15915494309189535, 1.5707963267948966)
julia> x = [1, 2]; x'x # adjoint vector * vector
5
*(
F::CTFlowsODE.AbstractFlow,
g::Tuple{Real, TF<:CTFlowsODE.AbstractFlow}
) -> Any
Shorthand for concatenate(F, g)
when g
is a tuple (t_switch, G)
.
Arguments
F::AbstractFlow
: The first flow.g::Tuple{ctNumber, AbstractFlow}
: Tuple containing the switching time and second flow.
Returns
- A new flow that switches from
F
toG
att_switch
.
Example
julia> F * (1.0, G)
*(
F::CTFlowsODE.AbstractFlow,
g::Tuple{Real, Any, TF<:CTFlowsODE.AbstractFlow}
) -> Any
Shorthand for concatenate(F, g)
when g
is a tuple (t_switch, η_switch, G)
including a jump.
Arguments
F::AbstractFlow
: The first flow.g::Tuple{ctNumber, Any, AbstractFlow}
: Tuple with switching time, jump value, and second flow.
Returns
- A flow with a jump at
t_switch
and a switch fromF
toG
.
Example
julia> F * (1.0, η, G)
CTFlows.Flow
— FunctionFlow(
vf::VectorField;
alg,
abstol,
reltol,
saveat,
internalnorm,
kwargs_Flow...
) -> CTFlowsODE.VectorFieldFlow
Constructs a flow object for a classical (non-Hamiltonian) vector field.
This creates a VectorFieldFlow
that integrates the ODE system dx/dt = vf(t, x, v)
using DifferentialEquations.jl. It handles both fixed and parametric dynamics, as well as jump discontinuities and event stopping.
Keyword Arguments
alg
,abstol
,reltol
,saveat
,internalnorm
: Solver options.kwargs_Flow...
: Additional arguments passed to the solver configuration.
Example
julia> vf(t, x, v) = -v * x
julia> flow = CTFlows.Flow(CTFlows.VectorField(vf))
julia> x1 = flow(0.0, 1.0, 1.0)
Flow(
h::CTFlows.AbstractHamiltonian;
alg,
abstol,
reltol,
saveat,
internalnorm,
kwargs_Flow...
) -> CTFlowsODE.HamiltonianFlow
Constructs a Hamiltonian flow from a scalar Hamiltonian.
This method builds a numerical integrator that simulates the evolution of a Hamiltonian system given a Hamiltonian function h(t, x, p, l)
or h(x, p)
.
Internally, it computes the right-hand side of Hamilton’s equations via automatic differentiation and returns a HamiltonianFlow
object.
Keyword Arguments
alg
,abstol
,reltol
,saveat
,internalnorm
: solver options.kwargs_Flow...
: forwarded to the solver.
Example
julia> H(x, p) = dot(p, p) + dot(x, x)
julia> flow = CTFlows.Flow(CTFlows.Hamiltonian(H))
julia> xf, pf = flow(0.0, x0, p0, 1.0)
Flow(
hv::HamiltonianVectorField;
alg,
abstol,
reltol,
saveat,
internalnorm,
kwargs_Flow...
) -> CTFlowsODE.HamiltonianFlow
Constructs a Hamiltonian flow from a precomputed Hamiltonian vector field.
This method assumes you already provide the Hamiltonian vector field (dx/dt, dp/dt)
instead of deriving it from a scalar Hamiltonian.
Returns a HamiltonianFlow
object that integrates the given system.
Keyword Arguments
alg
,abstol
,reltol
,saveat
,internalnorm
: solver options.kwargs_Flow...
: forwarded to the solver.
Example
julia> hv(t, x, p, l) = (∇ₚH, -∇ₓH)
julia> flow = CTFlows.Flow(CTFlows.HamiltonianVectorField(hv))
julia> xf, pf = flow(0.0, x0, p0, 1.0, l)
Flow(
ocp::Model,
u::CTFlows.ControlLaw;
alg,
abstol,
reltol,
saveat,
internalnorm,
kwargs_Flow...
) -> CTFlowsODE.OptimalControlFlow
Construct a flow for an optimal control problem using a given control law.
This method builds the Hamiltonian system associated with the optimal control problem (ocp
) and integrates the corresponding state–costate dynamics using the specified control law u
.
Arguments
ocp::CTModels.Model
: An optimal control problem defined usingCTModels
.u::CTFlows.ControlLaw
: A feedback control law generated byControlLaw(...)
or similar.alg
: Integration algorithm (default inferred).abstol
: Absolute tolerance for the ODE solver.reltol
: Relative tolerance for the ODE solver.saveat
: Time points at which to save the solution.internalnorm
: Optional norm function used by the integrator.kwargs_Flow
: Additional keyword arguments passed to the solver.
Returns
A flow object f
such that:
f(t0, x0, p0, tf)
integrates the state and costate fromt0
totf
.f((t0, tf), x0, p0)
returns the full trajectory over the interval.
Example
julia> u = (x, p) -> p
julia> f = Flow(ocp, ControlLaw(u))
Flow(
ocp::Model,
u::Function;
autonomous,
variable,
alg,
abstol,
reltol,
saveat,
internalnorm,
kwargs_Flow...
) -> CTFlowsODE.OptimalControlFlow
Construct a flow for an optimal control problem using a control function in feedback form.
This method constructs the Hamiltonian and integrates the associated state–costate dynamics using a raw function u
. It automatically wraps u
as a control law.
Arguments
ocp::CTModels.Model
: The optimal control problem.u::Function
: A feedback control function:- If
ocp
is autonomous:u(x, p)
- If non-autonomous:
u(t, x, p)
- If
autonomous::Bool
: Whether the control law depends on time.variable::Bool
: Whether the OCP involves variable time (e.g., free final time).alg
,abstol
,reltol
,saveat
,internalnorm
: ODE solver parameters.kwargs_Flow
: Additional options.
Returns
A Flow
object compatible with function call interfaces for state propagation.
Example
julia> u = (t, x, p) -> t + p
julia> f = Flow(ocp, u)
Flow(
ocp::Model,
u::Union{CTFlows.ControlLaw{<:Function, T, V}, CTFlows.FeedbackControl{<:Function, T, V}},
g::Union{CTFlows.MixedConstraint{<:Function, T, V}, CTFlows.StateConstraint{<:Function, T, V}},
μ::CTFlows.Multiplier{<:Function, T, V};
alg,
abstol,
reltol,
saveat,
internalnorm,
kwargs_Flow...
) -> CTFlowsODE.OptimalControlFlow
Construct a flow for an optimal control problem with control and constraint multipliers in feedback form.
This variant constructs a Hamiltonian system incorporating both the control law and a multiplier law (e.g., for enforcing state or mixed constraints). All inputs must be consistent in time dependence.
Arguments
ocp::CTModels.Model
: The optimal control problem.u::ControlLaw or FeedbackControl
: Feedback control.g::StateConstraint or MixedConstraint
: Constraint function.μ::Multiplier
: Multiplier function.alg
,abstol
,reltol
,saveat
,internalnorm
: Solver settings.kwargs_Flow
: Additional options.
Returns
A Flow
object that integrates the constrained Hamiltonian dynamics.
Example
julia> f = Flow(ocp, (x, p) -> p[1], (x, u) -> x[1] - 1, (x, p) -> x[1]+p[1])
For non-autonomous cases:
julia> f = Flow(ocp, (t, x, p) -> t + p, (t, x, u) -> x - 1, (t, x, p) -> x+p)
Flow(
ocp::Model,
u::Function,
g::Function,
μ::Function;
autonomous,
variable,
alg,
abstol,
reltol,
saveat,
internalnorm,
kwargs_Flow...
) -> CTFlowsODE.OptimalControlFlow
Construct a flow from a raw feedback control, constraint, and multiplier.
This version is for defining flows directly from user functions without wrapping them into ControlLaw
, Constraint
, or Multiplier
types. Automatically wraps and adapts them based on time dependence.
Arguments
ocp::CTModels.Model
: The optimal control problem.u::Function
: Control law.g::Function
: Constraint.μ::Function
: Multiplier.autonomous::Bool
: Whether the system is autonomous.variable::Bool
: Whether time is a free variable.alg
,abstol
,reltol
,saveat
,internalnorm
: Solver parameters.kwargs_Flow
: Additional options.
Returns
A Flow
object ready for trajectory integration.
Flow(
dyn::Function;
autonomous,
variable,
alg,
abstol,
reltol,
saveat,
internalnorm,
kwargs_Flow...
) -> CTFlowsODE.ODEFlow
Constructs a Flow
from a user-defined dynamical system given as a Julia function.
This high-level interface handles:
- autonomous and non-autonomous systems,
- presence or absence of additional variables (
v
), - selection of ODE solvers and tolerances,
- and integrates with the CTFlows event system (e.g., jumps, callbacks).
Arguments
dyn
: A function defining the vector field. Its signature must match the values ofautonomous
andvariable
.autonomous
: Whether the dynamics are time-independent (false
by default).variable
: Whether the dynamics depend on a control or parameterv
.alg
,abstol
,reltol
,saveat
,internalnorm
: Solver settings passed toOrdinaryDiffEq.solve
.kwargs_Flow
: Additional keyword arguments passed to the solver.
Returns
An ODEFlow
object, wrapping both the full solver and its right-hand side (RHS).
Supported Function Signatures for dyn
Depending on the (autonomous, variable)
flags:
(false, false)
:dyn(x)
(false, true)
:dyn(x, v)
(true, false)
:dyn(t, x)
(true, true)
:dyn(t, x, v)
Example
julia> dyn(t, x, v) = [-x[1] + v[1] * sin(t)]
julia> flow = CTFlows.Flow(dyn; autonomous=true, variable=true)
julia> xT = flow((0.0, 1.0), [1.0], [0.1])
CTFlows.Hamiltonian
— Typestruct Hamiltonian{TF<:Function, TD<:CTFlows.TimeDependence, VD<:CTFlows.VariableDependence} <: CTFlows.AbstractHamiltonian{TD<:CTFlows.TimeDependence, VD<:CTFlows.VariableDependence}
Encodes the Hamiltonian function H = ⟨p, f⟩ + L
in optimal control.
Fields
f
: a callable of the form:f(x, p)
f(t, x, p)
f(x, p, v)
f(t, x, p, v)
Type Parameters
TD
:Autonomous
orNonAutonomous
VD
:Fixed
orNonFixed
Example
julia> Hf(x, p) = dot(p, [x[2], -x[1]])
julia> H = Hamiltonian{typeof(Hf), Autonomous, Fixed}(Hf)
julia> H([1.0, 0.0], [1.0, 1.0])
CTFlows.HamiltonianLift
— Typestruct HamiltonianLift{TV<:VectorField, TD<:CTFlows.TimeDependence, VD<:CTFlows.VariableDependence} <: CTFlows.AbstractHamiltonian{TD<:CTFlows.TimeDependence, VD<:CTFlows.VariableDependence}
Lifts a vector field X
into a Hamiltonian function using the canonical symplectic structure.
This is useful to convert a vector field into a Hamiltonian via the identity: H(x, p) = ⟨p, X(x)⟩
.
Constructor
Use HamiltonianLift(X::VectorField)
where X
is a VectorField{...}
.
Example
f(x) = [x[2], -x[1]]
julia> X = VectorField{typeof(f), Autonomous, Fixed}(f)
julia> H = HamiltonianLift(X)
julia> H([1.0, 0.0], [0.5, 0.5])
CTFlows.HamiltonianVectorField
— Typestruct HamiltonianVectorField{TF<:Function, TD<:CTFlows.TimeDependence, VD<:CTFlows.VariableDependence} <: CTFlows.AbstractVectorField{TD<:CTFlows.TimeDependence, VD<:CTFlows.VariableDependence}
Represents the Hamiltonian vector field associated to a Hamiltonian function, typically defined as (∂H/∂p, -∂H/∂x)
.
Fields
f
: a callable implementing the Hamiltonian vector field.
Example
julia> f(x, p) = [p[2], -p[1], -x[1], -x[2]]
julia> XH = HamiltonianVectorField{typeof(f), Autonomous, Fixed}(f)
julia> XH([1.0, 0.0], [0.5, 0.5])
CTFlows.@Lie
— MacroCompute Lie or Poisson brackets.
This macro provides a unified notation to define recursively nested Lie brackets (for vector fields) or Poisson brackets (for Hamiltonians).
Syntax
@Lie [F, G]
: computes the Lie bracket[F, G]
of two vector fields.@Lie [[F, G], H]
: supports arbitrarily nested Lie brackets.@Lie {H, K}
: computes the Poisson bracket{H, K}
of two Hamiltonians.@Lie {{H, K}, L}
: supports arbitrarily nested Poisson brackets.@Lie expr autonomous = false
: specifies a non-autonomous system.@Lie expr variable = true
: indicates presence of an auxiliary variablev
.
Keyword-like arguments can be provided to control the evaluation context for Poisson brackets with raw functions:
autonomous = Bool
: whether the system is time-independent (default:true
).variable = Bool
: whether the system depends on an extra variablev
(default:false
).
Bracket type detection
- Square brackets
[...]
denote Lie brackets betweenVectorField
objects. - Curly brackets
{...}
denote Poisson brackets betweenHamiltonian
objects or between raw functions. - The macro automatically dispatches to
Lie
orPoisson
depending on the input pattern.
Return
A callable object representing the specified Lie or Poisson bracket expression. The returned function can be evaluated like any other vector field or Hamiltonian.
Examples
■ Lie brackets with VectorField
(autonomous)
julia> F1 = VectorField(x -> [0, -x[3], x[2]])
julia> F2 = VectorField(x -> [x[3], 0, -x[1]])
julia> L = @Lie [F1, F2]
julia> L([1.0, 2.0, 3.0])
3-element Vector{Float64}:
2.0
-1.0
0.0
■ Lie brackets with VectorField
(non-autonomous, with auxiliary variable)
julia> F1 = VectorField((t, x, v) -> [0, -x[3], x[2]]; autonomous=false, variable=true)
julia> F2 = VectorField((t, x, v) -> [x[3], 0, -x[1]]; autonomous=false, variable=true)
julia> L = @Lie [F1, F2]
julia> L(0.0, [1.0, 2.0, 3.0], 1.0)
3-element Vector{Float64}:
2.0
-1.0
0.0
■ Poisson brackets with Hamiltonian
(autonomous)
julia> H1 = Hamiltonian((x, p) -> x[1]^2 + p[2]^2)
julia> H2 = Hamiltonian((x, p) -> x[2]^2 + p[1]^2)
julia> P = @Lie {H1, H2}
julia> P([1.0, 1.0], [3.0, 2.0])
-4.0
■ Poisson brackets with Hamiltonian
(non-autonomous, with variable)
julia> H1 = Hamiltonian((t, x, p, v) -> x[1]^2 + p[2]^2 + v; autonomous=false, variable=true)
julia> H2 = Hamiltonian((t, x, p, v) -> x[2]^2 + p[1]^2 + v; autonomous=false, variable=true)
julia> P = @Lie {H1, H2}
julia> P(1.0, [1.0, 3.0], [4.0, 2.0], 3.0)
8.0
■ Poisson brackets from raw functions
julia> H1 = (x, p) -> x[1]^2 + p[2]^2
julia> H2 = (x, p) -> x[2]^2 + p[1]^2
julia> P = @Lie {H1, H2}
julia> P([1.0, 1.0], [3.0, 2.0])
-4.0
■ Poisson bracket with non-autonomous raw functions
julia> H1 = (t, x, p) -> x[1]^2 + p[2]^2 + t
julia> H2 = (t, x, p) -> x[2]^2 + p[1]^2 + t
julia> P = @Lie {H1, H2} autonomous = false
julia> P(3.0, [1.0, 2.0], [4.0, 1.0])
-8.0
■ Nested brackets
julia> F = VectorField(x -> [-x[1], x[2], x[3]])
julia> G = VectorField(x -> [x[3], -x[2], 0])
julia> H = VectorField(x -> [0, 0, -x[1]])
julia> nested = @Lie [[F, G], H]
julia> nested([1.0, 2.0, 3.0])
3-element Vector{Float64}:
2.0
0.0
-6.0
julia> H1 = (x, p) -> x[2]*x[1]^2 + p[1]^2
julia> H2 = (x, p) -> x[1]*p[2]^2
julia> H3 = (x, p) -> x[1]*p[2] + x[2]*p[1]
julia> nested_poisson = @Lie {{H1, H2}, H3}
julia> nested_poisson([1.0, 2.0], [0.5, 1.0])
14.0
■ Mixed expressions with arithmetic
julia> F1 = VectorField(x -> [0, -x[3], x[2]])
julia> F2 = VectorField(x -> [x[3], 0, -x[1]])
julia> x = [1.0, 2.0, 3.0]
julia> @Lie [F1, F2](x) + 3 * [F1, F2](x)
3-element Vector{Float64}:
8.0
-4.0
0.0
julia> H1 = (x, p) -> x[1]^2
julia> H2 = (x, p) -> p[1]^2
julia> H3 = (x, p) -> x[1]*p[1]
julia> x = [1.0, 2.0, 3.0]
julia> p = [3.0, 2.0, 1.0]
julia> @Lie {H1, H2}(x, p) + 2 * {H2, H3}(x, p)
24.0
CTFlows.Lie
— FunctionLie derivative of a scalar function along a vector field.
Example:
julia> φ = x -> [x[2], -x[1]]
julia> X = VectorField(φ)
julia> f = x -> x[1]^2 + x[2]^2
julia> Lie(X,f)([1, 2])
0
julia> φ = (t, x, v) -> [t + x[2] + v[1], -x[1] + v[2]]
julia> X = VectorField(φ, NonAutonomous, NonFixed)
julia> f = (t, x, v) -> t + x[1]^2 + x[2]^2
julia> Lie(X, f)(1, [1, 2], [2, 1])
10
Lie derivative of a scalar function along a function with specified dependencies.
Example:
julia> φ = x -> [x[2], -x[1]]
julia> f = x -> x[1]^2 + x[2]^2
julia> Lie(φ,f)([1, 2])
0
julia> φ = (t, x, v) -> [t + x[2] + v[1], -x[1] + v[2]]
julia> f = (t, x, v) -> t + x[1]^2 + x[2]^2
julia> Lie(φ, f, autonomous=false, variable=true)(1, [1, 2], [2, 1])
10
Lie bracket of two vector fields in the autonomous case.
Example:
julia> f = x -> [x[2], 2x[1]]
julia> g = x -> [3x[2], -x[1]]
julia> X = VectorField(f)
julia> Y = VectorField(g)
julia> Lie(X, Y)([1, 2])
[7, -14]
Lie bracket of two vector fields in the nonautonomous case.
Example:
julia> f = (t, x, v) -> [t + x[2] + v, -2x[1] - v]
julia> g = (t, x, v) -> [t + 3x[2] + v, -x[1] - v]
julia> X = VectorField(f, NonAutonomous, NonFixed)
julia> Y = VectorField(g, NonAutonomous, NonFixed)
julia> Lie(X, Y)(1, [1, 2], 1)
[-7, 12]
CTFlows.Lift
— FunctionLift(
X::VectorField
) -> HamiltonianLift{VectorField{TF, TD, VD}} where {TF<:Function, TD<:CTFlows.TimeDependence, VD<:CTFlows.VariableDependence}
Construct the Hamiltonian lift of a VectorField
.
Arguments
X::VectorField
: The vector field to lift. Its signature determines if it is autonomous and/or variable.
Returns
- A
HamiltonianLift
callable object representing the Hamiltonian lift ofX
.
Examples
julia> HL = Lift(VectorField(x -> [x[1]^2, x[2]^2], autonomous=true, variable=false))
julia> HL([1, 0], [0, 1]) # returns 0
julia> HL2 = Lift(VectorField((t, x, v) -> [t + x[1]^2, x[2]^2 + v], autonomous=false, variable=true))
julia> HL2(1, [1, 0], [0, 1], 1) # returns 1
julia> H = Lift(x -> 2x)
julia> H(1, 1) # returns 2
julia> H2 = Lift((t, x, v) -> 2x + t - v, autonomous=false, variable=true)
julia> H2(1, 1, 1, 1) # returns 2
# Alternative syntax using symbols for autonomy and variability
julia> H3 = Lift((t, x, v) -> 2x + t - v, NonAutonomous, NonFixed)
julia> H3(1, 1, 1, 1) # returns 2
Lift(
X::Function;
autonomous,
variable
) -> CTFlows.var"#19#23"{<:Function}
Construct the Hamiltonian lift of a function.
Arguments
X::Function
: The function representing the vector field.autonomous::Bool=true
: Whether the function is autonomous (time-independent).variable::Bool=false
: Whether the function depends on an additional variable argument.
Returns
- A callable function computing the Hamiltonian lift,
(and variants depending on autonomous
and variable
).
Details
Depending on the autonomous
and variable
flags, the returned function has one of the following call signatures:
(x, p)
ifautonomous=true
andvariable=false
(x, p, v)
ifautonomous=true
andvariable=true
(t, x, p)
ifautonomous=false
andvariable=false
(t, x, p, v)
ifautonomous=false
andvariable=true
Examples
julia> H = Lift(x -> 2x)
julia> H(1, 1) # returns 2
julia> H2 = Lift((t, x, v) -> 2x + t - v, autonomous=false, variable=true)
julia> H2(1, 1, 1, 1) # returns 2
CTModels.Model
— Typestruct Model{TD<:CTModels.TimeDependence, TimesModelType<:CTModels.AbstractTimesModel, StateModelType<:CTModels.AbstractStateModel, ControlModelType<:CTModels.AbstractControlModel, VariableModelType<:CTModels.AbstractVariableModel, DynamicsModelType<:Function, ObjectiveModelType<:CTModels.AbstractObjectiveModel, ConstraintsModelType<:CTModels.AbstractConstraintsModel, BuildExaModelType<:Union{Nothing, Function}} <: CTModels.AbstractModel
Fields
times::CTModels.AbstractTimesModel
state::CTModels.AbstractStateModel
control::CTModels.AbstractControlModel
variable::CTModels.AbstractVariableModel
dynamics::Function
objective::CTModels.AbstractObjectiveModel
constraints::CTModels.AbstractConstraintsModel
definition::Expr
build_examodel::Union{Nothing, Function}
CTBase.ParsingError
— Typestruct ParsingError <: CTBase.CTException
Exception thrown during parsing when a syntax error or invalid structure is detected.
Fields
var::String
: A message describing the parsing error.
Example
julia> throw(ParsingError("unexpected token"))
ERROR: ParsingError: unexpected token
CTFlows.Poisson
— FunctionPoisson(
f::CTFlows.AbstractHamiltonian{CTFlows.Autonomous, V<:CTFlows.VariableDependence},
g::CTFlows.AbstractHamiltonian{CTFlows.Autonomous, V<:CTFlows.VariableDependence}
) -> Any
Poisson bracket of two Hamiltonian functions (subtype of AbstractHamiltonian). Autonomous case.
Returns a Hamiltonian representing the Poisson bracket {f, g}
of two autonomous Hamiltonian functions f
and g
.
Example
julia> f = (x, p) -> x[2]^2 + 2x[1]^2 + p[1]^2
julia> g = (x, p) -> 3x[2]^2 - x[1]^2 + p[2]^2 + p[1]
julia> F = Hamiltonian(f)
julia> G = Hamiltonian(g)
julia> Poisson(f, g)([1, 2], [2, 1]) # -20
julia> Poisson(f, G)([1, 2], [2, 1]) # -20
julia> Poisson(F, g)([1, 2], [2, 1]) # -20
Poisson(
f::CTFlows.AbstractHamiltonian{CTFlows.NonAutonomous, V<:CTFlows.VariableDependence},
g::CTFlows.AbstractHamiltonian{CTFlows.NonAutonomous, V<:CTFlows.VariableDependence}
) -> Any
Poisson bracket of two Hamiltonian functions. Non-autonomous case.
Returns a Hamiltonian representing {f, g}
where f
and g
are time-dependent.
Example
julia> f = (t, x, p, v) -> t*v[1]*x[2]^2 + 2x[1]^2 + p[1]^2 + v[2]
julia> g = (t, x, p, v) -> 3x[2]^2 - x[1]^2 + p[2]^2 + p[1] + t - v[2]
julia> F = Hamiltonian(f, autonomous=false, variable=true)
julia> G = Hamiltonian(g, autonomous=false, variable=true)
julia> Poisson(F, G)(2, [1, 2], [2, 1], [4, 4]) # -76
julia> Poisson(f, g, NonAutonomous, NonFixed)(2, [1, 2], [2, 1], [4, 4]) # -76
Poisson(
f::HamiltonianLift{T<:CTFlows.TimeDependence, V<:CTFlows.VariableDependence},
g::HamiltonianLift{T<:CTFlows.TimeDependence, V<:CTFlows.VariableDependence}
)
Poisson bracket of two HamiltonianLift vector fields.
Returns the HamiltonianLift corresponding to the Lie bracket of vector fields f.X
and g.X
.
Example
julia> f = x -> [x[1]^2 + x[2]^2, 2x[1]^2]
julia> g = x -> [3x[2]^2, x[2] - x[1]^2]
julia> F = Lift(f)
julia> G = Lift(g)
julia> Poisson(F, G)([1, 2], [2, 1]) # -64
julia> f = (t, x, v) -> [t*v[1]*x[2]^2, 2x[1]^2 + v[2]]
julia> g = (t, x, v) -> [3x[2]^2 - x[1]^2, t - v[2]]
julia> F = Lift(f, NonAutonomous, NonFixed)
julia> G = Lift(g, NonAutonomous, NonFixed)
julia> Poisson(F, G)(2, [1, 2], [2, 1], [4, 4]) # 100
Poisson(
f::Function,
g::Function;
autonomous,
variable
) -> Hamiltonian
Poisson bracket of two functions. The time and variable dependence are specified with keyword arguments.
Returns a Hamiltonian computed from the functions promoted as Hamiltonians.
Example
julia> f = (x, p) -> x[2]^2 + 2x[1]^2 + p[1]^2
julia> g = (x, p) -> 3x[2]^2 - x[1]^2 + p[2]^2 + p[1]
julia> Poisson(f, g)([1, 2], [2, 1]) # -20
julia> f = (t, x, p, v) -> t*v[1]*x[2]^2 + 2x[1]^2 + p[1]^2 + v[2]
julia> g = (t, x, p, v) -> 3x[2]^2 - x[1]^2 + p[2]^2 + p[1] + t - v[2]
julia> Poisson(f, g, autonomous=false, variable=true)(2, [1, 2], [2, 1], [4, 4]) # -76
Poisson(
f::Function,
g::CTFlows.AbstractHamiltonian{TD<:CTFlows.TimeDependence, VD<:CTFlows.VariableDependence}
) -> Hamiltonian
Poisson bracket of a function and a Hamiltonian.
Returns a Hamiltonian representing {f, g}
where g
is already a Hamiltonian.
Example
julia> f = (x, p) -> x[2]^2 + 2x[1]^2 + p[1]^2
julia> g = (x, p) -> 3x[2]^2 - x[1]^2 + p[2]^2 + p[1]
julia> G = Hamiltonian(g)
julia> Poisson(f, G)([1, 2], [2, 1]) # -20
julia> f = (t, x, p, v) -> t*v[1]*x[2]^2 + 2x[1]^2 + p[1]^2 + v[2]
julia> g = (t, x, p, v) -> 3x[2]^2 - x[1]^2 + p[2]^2 + p[1] + t - v[2]
julia> G = Hamiltonian(g, autonomous=false, variable=true)
julia> Poisson(f, G)(2, [1, 2], [2, 1], [4, 4]) # -76
Poisson(
f::CTFlows.AbstractHamiltonian{TD<:CTFlows.TimeDependence, VD<:CTFlows.VariableDependence},
g::Function
) -> Hamiltonian
Poisson bracket of a Hamiltonian and a function.
Returns a Hamiltonian representing {f, g}
where f
is already a Hamiltonian.
Example
julia> f = (x, p) -> x[2]^2 + 2x[1]^2 + p[1]^2
julia> g = (x, p) -> 3x[2]^2 - x[1]^2 + p[2]^2 + p[1]
julia> F = Hamiltonian(f)
julia> Poisson(F, g)([1, 2], [2, 1]) # -20
julia> f = (t, x, p, v) -> t*v[1]*x[2]^2 + 2x[1]^2 + p[1]^2 + v[2]
julia> g = (t, x, p, v) -> 3x[2]^2 - x[1]^2 + p[2]^2 + p[1] + t - v[2]
julia> F = Hamiltonian(f, autonomous=false, variable=true)
julia> Poisson(F, g)(2, [1, 2], [2, 1], [4, 4]) # -76
CTModels.Solution
— Typestruct Solution{TimeGridModelType<:CTModels.AbstractTimeGridModel, TimesModelType<:CTModels.AbstractTimesModel, StateModelType<:CTModels.AbstractStateModel, ControlModelType<:CTModels.AbstractControlModel, VariableModelType<:CTModels.AbstractVariableModel, CostateModelType<:Function, ObjectiveValueType<:Real, DualModelType<:CTModels.AbstractDualModel, SolverInfosType<:CTModels.AbstractSolverInfos, ModelType<:CTModels.AbstractModel} <: CTModels.AbstractSolution
Fields
time_grid::CTModels.AbstractTimeGridModel
times::CTModels.AbstractTimesModel
state::CTModels.AbstractStateModel
control::CTModels.AbstractControlModel
variable::CTModels.AbstractVariableModel
costate::Function
objective::Real
dual::CTModels.AbstractDualModel
solver_infos::CTModels.AbstractSolverInfos
model::CTModels.AbstractModel
CTFlows.VectorField
— Typestruct VectorField{TF<:Function, TD<:CTFlows.TimeDependence, VD<:CTFlows.VariableDependence} <: CTFlows.AbstractVectorField{TD<:CTFlows.TimeDependence, VD<:CTFlows.VariableDependence}
Represents a dynamical system dx/dt = f(...)
as a vector field.
Fields
f
: a callable of the form:f(x)
f(t, x)
f(x, v)
f(t, x, v)
Example
f(x) = [x[2], -x[1]]
vf = VectorField{typeof(f), Autonomous, Fixed}(f)
vf([1.0, 0.0])
OptimalControl.available_methods
— Functionavailable_methods(
) -> Tuple{Vararg{Tuple{Symbol, Symbol, Symbol}}}
Return the list of available methods that can be used to solve optimal control problems.
CTDirect.build_OCP_solution
— Functionbuild_OCP_solution(
docp,
docp_solution
) -> Solution{TimeGridModelType, TimesModelType, StateModelType, ControlModelType, VariableModelType, CostateModelType, Float64, DualModelType, CTModels.SolverInfos{Dict{Symbol, Any}}, ModelType} where {TimeGridModelType<:CTModels.TimeGridModel, TimesModelType<:CTModels.TimesModel, StateModelType<:Union{CTModels.StateModelSolution{TS} where TS<:CTModels.var"#114#136", CTModels.StateModelSolution{TS} where TS<:CTModels.var"#115#137"}, ControlModelType<:Union{CTModels.ControlModelSolution{TS} where TS<:CTModels.var"#116#138", CTModels.ControlModelSolution{TS} where TS<:CTModels.var"#117#139"}, VariableModelType<:Union{CTModels.VariableModelSolution{Vector{Float64}}, CTModels.VariableModelSolution{Float64}}, CostateModelType<:Union{CTModels.var"#118#140", CTModels.var"#119#141"}, DualModelType<:(CTModels.DualModel{PC_Dual, BC_Dual, SC_LB_Dual, SC_UB_Dual, CC_LB_Dual, CC_UB_Dual, VC_LB_Dual, VC_UB_Dual} where {PC_Dual<:Union{Nothing, CTModels.var"#121#143", CTModels.var"#122#144"}, BC_Dual<:Union{Nothing, Vector{Float64}}, SC_LB_Dual<:Union{Nothing, CTModels.var"#124#146", CTModels.var"#125#147"}, SC_UB_Dual<:Union{Nothing, CTModels.var"#127#149", CTModels.var"#128#150"}, CC_LB_Dual<:Union{Nothing, CTModels.var"#130#152", CTModels.var"#131#153"}, CC_UB_Dual<:Union{Nothing, CTModels.var"#133#155", CTModels.var"#134#156"}, VC_LB_Dual<:Union{Nothing, Vector{Float64}}, VC_UB_Dual<:Union{Nothing, Vector{Float64}}}), ModelType<:Model}
Build OCP functional solution from DOCP discrete solution (given as a SolverCore.GenericExecutionStats)
build_OCP_solution(docp; primal, dual, mult_LB, mult_UB)
Build OCP functional solution from DOCP discrete solution (given as array for primal variables, optionally dual variables and bounds multipliers)
ExaModels.constraint
— Methodconstraint(
ocp::Model,
label::Symbol
) -> Tuple{Symbol, Any, Any, Any}
Get a labelled constraint from the model. Returns a tuple of the form (type, f, lb, ub)
where type
is the type of the constraint, f
is the function, lb
is the lower bound and ub
is the upper bound.
The function returns an exception if the label is not found in the model.
Arguments
model
: The model from which to retrieve the constraint.label
: The label of the constraint to retrieve.
Returns
Tuple
: A tuple containing the type, function, lower bound, and upper bound of the constraint.
CTModels.constraints
— Functionconstraints(
ocp::Model{<:CTModels.TimeDependence, <:CTModels.AbstractTimesModel, <:CTModels.AbstractStateModel, <:CTModels.AbstractControlModel, <:CTModels.AbstractVariableModel, <:Function, <:CTModels.AbstractObjectiveModel, C<:CTModels.AbstractConstraintsModel}
) -> CTModels.AbstractConstraintsModel
Return the constraints struct.
CTModels.constraints_violation
— Functionconstraints_violation(sol::Solution) -> Float64
Return the constraints violation.
CTModels.control
— Functioncontrol(
ocp::Model{<:CTModels.TimeDependence, <:CTModels.TimesModel, <:CTModels.AbstractStateModel, T<:CTModels.AbstractControlModel}
) -> CTModels.AbstractControlModel
Return the control struct.
control(
sol::Solution{<:CTModels.AbstractTimeGridModel, <:CTModels.AbstractTimesModel, <:CTModels.AbstractStateModel, <:CTModels.ControlModelSolution{TS<:Function}}
) -> Function
Return the control as a function of time.
julia> u = control(sol)
julia> t0 = time_grid(sol)[1]
julia> u0 = u(t0) # control at the initial time
CTModels.control_components
— Functioncontrol_components(ocp::Model) -> Vector{String}
Return the names of the components of the control.
control_components(sol::Solution) -> Vector{String}
Return the names of the components of the control.
CTModels.control_dimension
— Functioncontrol_dimension(ocp::Model) -> Int64
Return the control dimension.
control_dimension(sol::Solution) -> Int64
Return the dimension of the control.
CTModels.control_name
— Functioncontrol_name(ocp::Model) -> String
Return the name of the control.
control_name(sol::Solution) -> String
Return the name of the control.
CTModels.costate
— Functioncostate(
sol::Solution{<:CTModels.AbstractTimeGridModel, <:CTModels.AbstractTimesModel, <:CTModels.AbstractStateModel, <:CTModels.AbstractControlModel, <:CTModels.AbstractVariableModel, Co<:Function}
) -> Function
Return the costate as a function of time.
julia> p = costate(sol)
julia> t0 = time_grid(sol)[1]
julia> p0 = p(t0) # costate at the initial time
CTModels.criterion
— Functioncriterion(model::CTModels.MayerObjectiveModel) -> Symbol
Return the criterion (:min or :max).
criterion(model::CTModels.LagrangeObjectiveModel) -> Symbol
Return the criterion (:min or :max).
criterion(model::CTModels.BolzaObjectiveModel) -> Symbol
Return the criterion (:min or :max).
criterion(ocp::Model) -> Symbol
Return the type of criterion (:min or :max).
CTParser.@def
— MacroDefine an optimal control problem. One pass parsing of the definition. Can be used writing either ocp = @def begin ... end
or @def ocp begin ... end
. In the second case, setting log
to true
will display the parsing steps.
Example
ocp = @def begin
tf ∈ R, variable
t ∈ [ 0, tf ], time
x ∈ R², state
u ∈ R, control
tf ≥ 0
-1 ≤ u(t) ≤ 1
q = x₁
v = x₂
q(0) == 1
v(0) == 2
q(tf) == 0
v(tf) == 0
0 ≤ q(t) ≤ 5, (1)
-2 ≤ v(t) ≤ 3, (2)
ẋ(t) == [ v(t), u(t) ]
tf → min
end
@def ocp begin
tf ∈ R, variable
t ∈ [ 0, tf ], time
x ∈ R², state
u ∈ R, control
tf ≥ 0
-1 ≤ u(t) ≤ 1
q = x₁
v = x₂
q(0) == 1
v(0) == 2
q(tf) == 0
v(tf) == 0
0 ≤ q(t) ≤ 5, (1)
-2 ≤ v(t) ≤ 3, (2)
ẋ(t) == [ v(t), u(t) ]
tf → min
end true # final boolean to show parsing log
CTModels.definition
— Functiondefinition(ocp::Model) -> Expr
Return the model definition of the optimal control problem.
definition(ocp::CTModels.PreModel) -> Union{Nothing, Expr}
Return the model definition of the optimal control problem or nothing
.
CTDirect.direct_transcription
— Functiondirect_transcription(
ocp::Model,
description...;
grid_size,
disc_method,
time_grid,
init,
kwargs...
) -> CTDirect.DOCP{_A, Model{TD, TimesModelType, StateModelType, ControlModelType, VariableModelType, DynamicsModelType, ObjectiveModelType, ConstraintsModelType, BuildExaModelType}} where {_A<:CTDirect.Discretization, TD<:CTModels.TimeDependence, TimesModelType<:CTModels.AbstractTimesModel, StateModelType<:CTModels.AbstractStateModel, ControlModelType<:CTModels.AbstractControlModel, VariableModelType<:CTModels.AbstractVariableModel, DynamicsModelType<:Function, ObjectiveModelType<:CTModels.AbstractObjectiveModel, ConstraintsModelType<:CTModels.AbstractConstraintsModel, BuildExaModelType<:Union{Nothing, Function}}
Discretize an optimal control problem into a nonlinear optimization problem.
Arguments
- ocp: optimal control problem as defined in
CTModels
- [description]: set the NLP model ([
:adnlp
] orexa
) and / or solver ([:ipopt
], :madnlp or :knitro)
Keyword arguments (optional)
grid_size
: number of time steps for the discretized problem ([250])disc_method
: discretization method ([:trapeze
],:euler
,:euler_implicit
,:midpoint
,gauss_legendre_2
,gauss_legendre_3
)time_grid
: explicit time grid (can be non uniform)init
: info for the starting guess (values as named tuple or existing solution)
Other kewwords arguments are passed down to the NLP modeler
CTModels.dual
— Functiondual(sol::Solution, model::Model, label::Symbol) -> Any
Return the dual variable associated with a constraint identified by its label
.
Searches through all constraint types (path, boundary, state, control, and variable constraints) defined in the model and returns the corresponding dual value from the solution.
Arguments
sol::Solution
: Solution object containing dual variables.model::Model
: Model containing constraint definitions.label::Symbol
: Symbol corresponding to a constraint label.
Returns
A function of time t
for time-dependent constraints, or a scalar/vector for time-invariant duals. If the label is not found, throws an IncorrectArgument
exception.
CTModels.dynamics
— Functiondynamics(
ocp::Model{<:CTModels.TimeDependence, <:CTModels.AbstractTimesModel, <:CTModels.AbstractStateModel, <:CTModels.AbstractControlModel, <:CTModels.AbstractVariableModel, D<:Function}
) -> Function
Return the dynamics.
CTModels.export_ocp_solution
— Functionexport_ocp_solution(
sol::CTModels.AbstractSolution;
format,
filename
)
Export a solution in JLD or JSON formats. Redirect to one of the methods:
export_ocp_solution(JLD2Tag(), sol, filename=filename)
export_ocp_solution(JSON3Tag(), sol, filename=filename)
Examples
julia> using JSON3
julia> export_ocp_solution(sol; filename="solution", format=:JSON)
julia> using JLD2
julia> export_ocp_solution(sol; filename="solution", format=:JLD) # JLD is the default
export_ocp_solution(
::CTModels.JSON3Tag,
sol::Solution;
filename
)
Export an optimal control solution to a .json
file using the JSON3 format.
This function serializes a CTModels.Solution
into a structured JSON dictionary, including all primal and dual information, which can be read by external tools.
Arguments
::CTModels.JSON3Tag
: A tag used to dispatch the export method for JSON3.sol::CTModels.Solution
: The solution to be saved.
Keyword Arguments
filename::String = "solution"
: Base filename. The.json
extension is automatically appended.
Notes
The exported JSON includes the time grid, state, control, costate, objective, solver info, and all constraint duals (if available).
Example
julia> using JSON3
julia> export_ocp_solution(JSON3Tag(), sol; filename="mysolution")
# → creates "mysolution.json"
export_ocp_solution(
::CTModels.JLD2Tag,
sol::Solution;
filename
)
Export an optimal control solution to a .jld2
file using the JLD2 format.
This function serializes and saves a CTModels.Solution
object to disk, allowing it to be reloaded later.
Arguments
::CTModels.JLD2Tag
: A tag used to dispatch the export method for JLD2.sol::CTModels.Solution
: The optimal control solution to be saved.
Keyword Arguments
filename::String = "solution"
: Base name of the file. The.jld2
extension is automatically appended.
Example
julia> using JLD2
julia> export_ocp_solution(JLD2Tag(), sol; filename="mysolution")
# → creates "mysolution.jld2"
CTModels.final_time
— Functionfinal_time(
model::CTModels.TimesModel{<:CTModels.AbstractTimeModel, <:CTModels.FixedTimeModel{T<:Real}}
) -> Real
Get the final time from the times model, from a fixed final time model.
final_time(
model::CTModels.TimesModel{<:CTModels.AbstractTimeModel, CTModels.FreeTimeModel},
variable::AbstractArray{T<:Real, 1}
) -> Any
Get the final time from the times model, from a free final time model.
final_time(ocp::CTModels.AbstractModel) -> Real
final_time(
ocp::CTModels.AbstractModel,
variable::AbstractVector
) -> Any
final_time(
ocp::Model{<:CTModels.TimeDependence, <:CTModels.TimesModel{<:CTModels.AbstractTimeModel, CTModels.FixedTimeModel{T<:Real}}}
) -> Real
Return the final time, for a fixed final time.
final_time(
ocp::Model{<:CTModels.TimeDependence, <:CTModels.TimesModel{<:CTModels.AbstractTimeModel, CTModels.FreeTimeModel}},
variable::AbstractArray{T<:Real, 1}
) -> Any
Return the final time, for a free final time.
final_time(
ocp::Model{<:CTModels.TimeDependence, <:CTModels.TimesModel{<:CTModels.AbstractTimeModel, CTModels.FreeTimeModel}},
variable::Real
) -> Real
Return the final time, for a free final time.
CTModels.final_time_name
— Functionfinal_time_name(model::CTModels.TimesModel) -> String
Get the name of the final time from the times model.
final_time_name(ocp::Model) -> String
Return the name of the final time.
final_time_name(sol::Solution) -> String
Return the name of the final time.
CTModels.get_build_examodel
— Functionget_build_examodel(
ocp::Model{<:CTModels.TimeDependence, <:CTModels.AbstractTimesModel, <:CTModels.AbstractStateModel, <:CTModels.AbstractControlModel, <:CTModels.AbstractVariableModel, <:Function, <:CTModels.AbstractObjectiveModel, <:CTModels.AbstractConstraintsModel, BE<:Function}
) -> Function
Return the build_examodel.
get_build_examodel(
_::Model{<:CTModels.TimeDependence, <:CTModels.AbstractTimesModel, <:CTModels.AbstractStateModel, <:CTModels.AbstractControlModel, <:CTModels.AbstractVariableModel, <:Function, <:CTModels.AbstractObjectiveModel, <:CTModels.AbstractConstraintsModel, <:Nothing}
)
Return an error (UnauthorizedCall) since the model is not built with the :exa backend.
CTModels.has_fixed_final_time
— Functionhas_fixed_final_time(
times::CTModels.TimesModel{<:CTModels.AbstractTimeModel, <:CTModels.FixedTimeModel{T<:Real}}
) -> Bool
Check if the final time is fixed. Return true.
has_fixed_final_time(
times::CTModels.TimesModel{<:CTModels.AbstractTimeModel, CTModels.FreeTimeModel}
) -> Bool
Check if the final time is free. Return false.
has_fixed_final_time(ocp::Model) -> Bool
Check if the final time is fixed.
CTModels.has_fixed_initial_time
— Functionhas_fixed_initial_time(
times::CTModels.TimesModel{<:CTModels.FixedTimeModel{T<:Real}}
) -> Bool
Check if the initial time is fixed. Return true.
has_fixed_initial_time(
times::CTModels.TimesModel{CTModels.FreeTimeModel}
) -> Bool
Check if the initial time is free. Return false.
has_fixed_initial_time(ocp::Model) -> Bool
Check if the initial time is fixed.
CTModels.has_free_final_time
— Functionhas_free_final_time(times::CTModels.TimesModel) -> Bool
Check if the final time is free.
has_free_final_time(ocp::Model) -> Bool
Check if the final time is free.
CTModels.has_free_initial_time
— Functionhas_free_initial_time(times::CTModels.TimesModel) -> Bool
Check if the final time is free.
has_free_initial_time(ocp::Model) -> Bool
Check if the initial time is free.
CTModels.has_lagrange_cost
— Functionhas_lagrange_cost(_::CTModels.MayerObjectiveModel) -> Bool
Return false.
has_lagrange_cost(
_::CTModels.LagrangeObjectiveModel
) -> Bool
Return true.
has_lagrange_cost(_::CTModels.BolzaObjectiveModel) -> Bool
Return true.
has_lagrange_cost(ocp::Model) -> Bool
Check if the model has a Lagrange cost.
CTModels.has_mayer_cost
— Functionhas_mayer_cost(_::CTModels.MayerObjectiveModel) -> Bool
Return true.
has_mayer_cost(_::CTModels.LagrangeObjectiveModel) -> Bool
Return false.
has_mayer_cost(_::CTModels.BolzaObjectiveModel) -> Bool
Return true.
has_mayer_cost(ocp::Model) -> Bool
Check if the model has a Mayer cost.
CTModels.import_ocp_solution
— Functionimport_ocp_solution(
ocp::CTModels.AbstractModel;
format,
filename
) -> Any
Import a solution from a JLD or JSON file. Redirect to one of the methods:
import_ocp_solution(JLD2Tag(), ocp, filename=filename)
import_ocp_solution(JSON3Tag(), ocp, filename=filename)
Examples
julia> using JSON3
julia> sol = import_ocp_solution(ocp; filename="solution", format=:JSON)
julia> using JLD2
julia> sol = import_ocp_solution(ocp; filename="solution", format=:JLD) # JLD is the default
import_ocp_solution(
::CTModels.JSON3Tag,
ocp::Model;
filename
)
Import an optimal control solution from a .json
file exported with export_ocp_solution
.
This function reads the JSON contents and reconstructs a CTModels.Solution
object, including the discretized primal and dual trajectories.
Arguments
::CTModels.JSON3Tag
: A tag used to dispatch the import method for JSON3.ocp::CTModels.Model
: The model associated with the optimal control problem. Used to rebuild the full solution.
Keyword Arguments
filename::String = "solution"
: Base filename. The.json
extension is automatically appended.
Returns
CTModels.Solution
: A reconstructed solution instance.
Notes
Handles both vector and matrix encodings of signals. If dual fields are missing or null
, the corresponding attributes are set to nothing
.
Example
julia> using JSON3
julia> sol = import_ocp_solution(JSON3Tag(), model; filename="mysolution")
import_ocp_solution(
::CTModels.JLD2Tag,
ocp::Model;
filename
)
Import an optimal control solution from a .jld2
file.
This function loads a previously saved CTModels.Solution
from disk.
Arguments
::CTModels.JLD2Tag
: A tag used to dispatch the import method for JLD2.ocp::CTModels.Model
: The associated model (used for dispatch consistency; not used internally).
Keyword Arguments
filename::String = "solution"
: Base name of the file. The.jld2
extension is automatically appended.
Returns
CTModels.Solution
: The loaded solution object.
Example
julia> using JLD2
julia> sol = import_ocp_solution(JLD2Tag(), model; filename="mysolution")
CTModels.infos
— Functioninfos(sol::Solution) -> Dict{Symbol, Any}
Return a dictionary of additional infos depending on the solver or nothing
.
CTModels.initial_time
— Functioninitial_time(
model::CTModels.TimesModel{<:CTModels.FixedTimeModel{T<:Real}}
) -> Real
Get the initial time from the times model, from a fixed initial time model.
initial_time(
model::CTModels.TimesModel{CTModels.FreeTimeModel},
variable::AbstractArray{T<:Real, 1}
) -> Any
Get the initial time from the times model, from a free initial time model.
initial_time(
ocp::Model{<:CTModels.TimeDependence, <:CTModels.TimesModel{CTModels.FixedTimeModel{T<:Real}}}
) -> Real
Return the initial time, for a fixed initial time.
initial_time(
ocp::Model{<:CTModels.TimeDependence, <:CTModels.TimesModel{CTModels.FreeTimeModel}},
variable::AbstractArray{T<:Real, 1}
) -> Any
Return the initial time, for a free initial time.
initial_time(
ocp::Model{<:CTModels.TimeDependence, <:CTModels.TimesModel{CTModels.FreeTimeModel}},
variable::Real
) -> Real
Return the initial time, for a free initial time.
CTModels.initial_time_name
— Functioninitial_time_name(model::CTModels.TimesModel) -> String
Get the name of the initial time from the times model.
initial_time_name(ocp::Model) -> String
Return the name of the initial time.
initial_time_name(sol::Solution) -> String
Return the name of the initial time.
CTModels.is_autonomous
— Methodis_autonomous(
_::Model{CTModels.Autonomous, <:CTModels.TimesModel}
) -> Bool
Return true
.
CTModels.iterations
— Functioniterations(sol::Solution) -> Int64
Return the number of iterations (if solved by an iterative method).
CTModels.lagrange
— Functionlagrange(
model::CTModels.LagrangeObjectiveModel{L<:Function}
) -> Function
Return the Lagrange function.
lagrange(
model::CTModels.BolzaObjectiveModel{<:Function, L<:Function}
) -> Function
Return the Lagrange function.
lagrange(
ocp::Model{<:CTModels.TimeDependence, <:CTModels.AbstractTimesModel, <:CTModels.AbstractStateModel, <:CTModels.AbstractControlModel, <:CTModels.AbstractVariableModel, <:Function, CTModels.LagrangeObjectiveModel{L<:Function}}
) -> Function
Return the Lagrange cost.
lagrange(
ocp::Model{<:CTModels.TimeDependence, <:CTModels.AbstractTimesModel, <:CTModels.AbstractStateModel, <:CTModels.AbstractControlModel, <:CTModels.AbstractVariableModel, <:Function, <:CTModels.BolzaObjectiveModel{<:Function, L<:Function}}
) -> Any
Return the Lagrange cost.
CTModels.mayer
— Functionmayer(
model::CTModels.MayerObjectiveModel{M<:Function}
) -> Function
Return the Mayer function.
mayer(
model::CTModels.BolzaObjectiveModel{M<:Function}
) -> Function
Return the Mayer function.
mayer(
ocp::Model{<:CTModels.TimeDependence, <:CTModels.AbstractTimesModel, <:CTModels.AbstractStateModel, <:CTModels.AbstractControlModel, <:CTModels.AbstractVariableModel, <:Function, <:CTModels.MayerObjectiveModel{M<:Function}}
) -> Any
Return the Mayer cost.
mayer(
ocp::Model{<:CTModels.TimeDependence, <:CTModels.AbstractTimesModel, <:CTModels.AbstractStateModel, <:CTModels.AbstractControlModel, <:CTModels.AbstractVariableModel, <:Function, <:CTModels.BolzaObjectiveModel{M<:Function}}
) -> Any
Return the Mayer cost.
CTModels.message
— Functionmessage(sol::Solution) -> String
Return the message associated to the status criterion.
CTDirect.model
— Functionmodel(docp::CTDirect.DOCP) -> Any
Extracts the NLP model from a DOCP object.
Arguments
docp::DOCP
: The DOCP object containing the NLP model.
Returns
The NLP model stored in the DOCP object.
ExaModels.objective
— Methodobjective(ocp::Model) -> CTModels.AbstractObjectiveModel
See CTModels.objective.
ExaModels.objective
— Methodobjective(sol::Solution) -> Real
Return the objective value.
RecipesBase.plot
— Methodplot(
sol::Solution,
description::Symbol...;
layout,
control,
time,
state_style,
state_bounds_style,
control_style,
control_bounds_style,
costate_style,
time_style,
path_style,
path_bounds_style,
dual_style,
size,
kwargs...
) -> Plots.Plot
Plot the components of an optimal control solution.
This is the main user-facing function to visualise the solution of an optimal control problem solved with the control-toolbox ecosystem.
It generates a set of subplots showing the evolution of the state, control, costate, path constraints, and dual variables over time, depending on the problem and the user’s choices.
Arguments
sol::CTModels.Solution
: The optimal control solution to visualise.description::Symbol...
: A variable number of symbols indicating which components to include in the plot. Common values include::state
– plot the state.:costate
– plot the costate (adjoint).:control
– plot the control.:path
– plot the path constraints.:dual
– plot the dual variables (or Lagrange multipliers) associated with path constraints.
If no symbols are provided, a default set is used based on the problem and styles.
Keyword Arguments (Optional)
layout::Symbol = :group
: Specifies how to arrange plots.:group
: Fewer plots, grouping similar variables together (e.g., all states in one subplot).:split
: One plot per variable component, stacked in a layout.
control::Symbol = :components
: Defines how to represent control inputs.:components
: One curve per control component.:norm
: Single curve showing the Euclidean norm ‖u(t)‖.:all
: Plot both components and norm.
time::Symbol = :default
: Time normalisation for plots.:default
: Real time scale.:normalize
or:normalise
: Normalised to the interval [0, 1].
Style Options (Optional)
All style-related keyword arguments can be either a NamedTuple
of plotting attributes or the Symbol
:none
referring to not plot the associated element. These allow you to customise color, line style, markers, etc.
time_style
: Style for vertical lines at initial and final times.state_style
: Style for state components.costate_style
: Style for costate components.control_style
: Style for control components.path_style
: Style for path constraint values.dual_style
: Style for dual variables.
Bounds Decorations (Optional)
Use these options to customise bounds on the plots if applicable and defined in the model. Set to :none
to hide.
state_bounds_style
: Style for state bounds.control_bounds_style
: Style for control bounds.path_bounds_style
: Style for path constraint bounds.
Returns
- A
Plots.Plot
object, which can be displayed, saved, or further customised.
Example
# basic plot
julia> plot(sol)
# plot only the state and control
julia> plot(sol, :state, :control)
# customise layout and styles, no costate
julia> plot(sol;
layout = :group,
control = :all,
state_style = (color=:blue, linestyle=:solid),
control_style = (color=:red, linestyle=:dash),
costate_style = :none)
RecipesBase.plot!
— Methodplot!(
p::Plots.Plot,
sol::Solution,
description::Symbol...;
layout,
control,
time,
state_style,
state_bounds_style,
control_style,
control_bounds_style,
costate_style,
time_style,
path_style,
path_bounds_style,
dual_style,
kwargs...
) -> Plots.Plot
Modify Plot p
with the optimal control solution sol
.
See plot
for full behavior and keyword arguments.
CTDirect.set_initial_guess
— Functionset_initial_guess(docp::CTDirect.DOCP, init) -> Any
Set initial guess in the DOCP
CommonSolve.solve
— Methodsolve(
ocp::Model,
description::Symbol...;
display,
kwargs...
) -> Solution{TimeGridModelType, TimesModelType, StateModelType, ControlModelType, VariableModelType, CostateModelType, Float64, DualModelType, CTModels.SolverInfos{Dict{Symbol, Any}}, ModelType} where {TimeGridModelType<:CTModels.TimeGridModel, TimesModelType<:CTModels.TimesModel, StateModelType<:Union{CTModels.StateModelSolution{TS} where TS<:CTModels.var"#114#136", CTModels.StateModelSolution{TS} where TS<:CTModels.var"#115#137"}, ControlModelType<:Union{CTModels.ControlModelSolution{TS} where TS<:CTModels.var"#116#138", CTModels.ControlModelSolution{TS} where TS<:CTModels.var"#117#139"}, VariableModelType<:Union{CTModels.VariableModelSolution{Vector{Float64}}, CTModels.VariableModelSolution{Float64}}, CostateModelType<:Union{CTModels.var"#118#140", CTModels.var"#119#141"}, DualModelType<:(CTModels.DualModel{PC_Dual, Vector{Float64}, SC_LB_Dual, SC_UB_Dual, CC_LB_Dual, CC_UB_Dual, Vector{Float64}, Vector{Float64}} where {PC_Dual<:Union{CTModels.var"#121#143"{CTModels.var"#120#142"}, CTModels.var"#122#144"{CTModels.var"#120#142"}}, SC_LB_Dual<:Union{CTModels.var"#124#146"{CTModels.var"#123#145"}, CTModels.var"#125#147"{CTModels.var"#123#145"}}, SC_UB_Dual<:Union{CTModels.var"#127#149"{CTModels.var"#126#148"}, CTModels.var"#128#150"{CTModels.var"#126#148"}}, CC_LB_Dual<:Union{CTModels.var"#130#152"{CTModels.var"#129#151"}, CTModels.var"#131#153"{CTModels.var"#129#151"}}, CC_UB_Dual<:Union{CTModels.var"#133#155"{CTModels.var"#132#154"}, CTModels.var"#134#156"{CTModels.var"#132#154"}}}), ModelType<:Model}
Solve the optimal control problem ocp
by the method given by the (optional) description. The get the list of available methods:
julia> available_methods()
The higher in the list, the higher is the priority. The keyword arguments are specific to the chosen method and represent the options of the solver.
Arguments
ocp::OptimalControlModel
: the optimal control problem to solve.description::Symbol...
: the description of the method used to solve the problem.kwargs...
: the options of the method.
Examples
The simplest way to solve the optimal control problem is to call the function without any argument.
julia> sol = solve(ocp)
The method description is a list of Symbols. The default is
julia> sol = solve(ocp, :direct, :adnlp, :ipopt)
You can provide a partial description, the function will find the best match.
julia> sol = solve(ocp, :direct)
See the resolution methods section for more details about the available methods.
The keyword arguments are specific to the chosen method and correspond to the options of the different solvers. For example, the keyword max_iter
is an Ipopt option that may be used to set the maximum number of iterations.
julia> sol = solve(ocp, :direct, :ipopt, max_iter=100)
See the direct method section for more details about associated options. These options also detailed in the CTDirect.solve
documentation. This main solve
method redirects to CTDirect.solve
when the :direct
Symbol is given in the description. See also the NLP solvers section for more details about Ipopt or MadNLP options.
To help the solve converge, an initial guess can be provided within the keyword init
. You can provide the initial guess for the state, control, and variable.
julia> sol = solve(ocp, init=(state=[-0.5, 0.2], control=0.5))
See how to set an initial guess for more details.
CTModels.state
— Functionstate(
ocp::Model{<:CTModels.TimeDependence, <:CTModels.TimesModel, T<:CTModels.AbstractStateModel}
) -> CTModels.AbstractStateModel
Return the state struct.
state(
sol::Solution{<:CTModels.AbstractTimeGridModel, <:CTModels.AbstractTimesModel, <:CTModels.StateModelSolution{TS<:Function}}
) -> Function
Return the state as a function of time.
julia> x = state(sol)
julia> t0 = time_grid(sol)[1]
julia> x0 = x(t0) # state at the initial time
CTModels.state_components
— Functionstate_components(ocp::Model) -> Vector{String}
Return the names of the components of the state.
state_components(sol::Solution) -> Vector{String}
Return the names of the components of the state.
CTModels.state_dimension
— Functionstate_dimension(ocp::CTModels.PreModel) -> Int64
state_dimension(ocp::Model) -> Int64
Return the state dimension.
state_dimension(sol::Solution) -> Int64
Return the dimension of the state.
CTModels.state_name
— Functionstate_name(ocp::Model) -> String
Return the name of the state.
state_name(sol::Solution) -> String
Return the name of the state.
CTModels.status
— Functionstatus(sol::Solution) -> Symbol
Return the status criterion (a Symbol).
CTModels.successful
— Functionsuccessful(sol::Solution) -> Bool
Return the successful status.
CTModels.time_grid
— Functiontime_grid(
sol::Solution{<:CTModels.TimeGridModel{T<:Union{StepRangeLen, AbstractVector{<:Real}}}}
) -> Union{StepRangeLen, AbstractVector{<:Real}}
Return the time grid.
CTModels.time_name
— Functiontime_name(model::CTModels.TimesModel) -> String
Get the name of the time variable from the times model.
time_name(ocp::Model) -> String
Return the name of the time.
time_name(sol::Solution) -> String
Return the name of the time component.
CTModels.times
— Functiontimes(
ocp::Model{<:CTModels.TimeDependence, T<:CTModels.TimesModel}
) -> CTModels.TimesModel
Return the times struct.
ExaModels.variable
— Methodvariable(ocp::Model) -> CTModels.AbstractVariableModel
See CTModels.variable.
ExaModels.variable
— Methodvariable(
sol::Solution
) -> Union{Real, AbstractVector{<:Real}}
Return the variable or nothing
.
julia> v = variable(sol)
CTModels.variable_components
— Functionvariable_components(ocp::Model) -> Vector{String}
Return the names of the components of the variable.
variable_components(sol::Solution) -> Vector{String}
Return the names of the components of the variable.
CTModels.variable_dimension
— Functionvariable_dimension(ocp::Model) -> Int64
Return the variable dimension.
variable_dimension(sol::Solution) -> Int64
Return the dimension of the variable.
CTModels.variable_name
— Functionvariable_name(ocp::Model) -> String
Return the name of the variable.
variable_name(sol::Solution) -> String
Return the name of the variable.
CTFlows.:⋅
— FunctionLie derivative of a scalar function along a vector field in the autonomous case.
Example:
julia> φ = x -> [x[2], -x[1]]
julia> X = VectorField(φ)
julia> f = x -> x[1]^2 + x[2]^2
julia> (X⋅f)([1, 2])
0
Lie derivative of a scalar function along a vector field in the nonautonomous case.
Example:
julia> φ = (t, x, v) -> [t + x[2] + v[1], -x[1] + v[2]]
julia> X = VectorField(φ, NonAutonomous, NonFixed)
julia> f = (t, x, v) -> t + x[1]^2 + x[2]^2
julia> (X⋅f)(1, [1, 2], [2, 1])
10
Lie derivative of a scalar function along a function (considered autonomous and non-variable).
Example:
julia> φ = x -> [x[2], -x[1]]
julia> f = x -> x[1]^2 + x[2]^2
julia> (φ⋅f)([1, 2])
0
julia> φ = (t, x, v) -> [t + x[2] + v[1], -x[1] + v[2]]
julia> f = (t, x, v) -> t + x[1]^2 + x[2]^2
julia> (φ⋅f)(1, [1, 2], [2, 1])
MethodError