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.

Beware!

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.OptimalControlModule

OptimalControl module.

List of all the exported names:

source

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
source
*(
    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 to G at t_switch.

Example

julia> F * (1.0, G)
source
*(
    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 from F to G.

Example

julia> F * (1.0, η, G)
source
CTFlows.FlowFunction
Flow(
    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)
source
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)
source
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)
source
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 using CTModels.
  • u::CTFlows.ControlLaw: A feedback control law generated by ControlLaw(...) 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 from t0 to tf.
  • f((t0, tf), x0, p0) returns the full trajectory over the interval.

Example

julia> u = (x, p) -> p
julia> f = Flow(ocp, ControlLaw(u))
source
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)
  • 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)
source
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)
Warning

All input functions must match the autonomous/non-autonomous nature of the problem.

source
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.

source
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 of autonomous and variable.
  • autonomous: Whether the dynamics are time-independent (false by default).
  • variable: Whether the dynamics depend on a control or parameter v.
  • alg, abstol, reltol, saveat, internalnorm: Solver settings passed to OrdinaryDiffEq.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])
source
CTFlows.HamiltonianType
struct 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 or NonAutonomous
  • VD: Fixed or NonFixed

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])
source
CTFlows.HamiltonianLiftType
struct 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])
source
CTFlows.HamiltonianVectorFieldType
struct 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])
source
CTFlows.@LieMacro

Compute 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 variable v.

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 variable v (default: false).

Bracket type detection

  • Square brackets [...] denote Lie brackets between VectorField objects.
  • Curly brackets {...} denote Poisson brackets between Hamiltonian objects or between raw functions.
  • The macro automatically dispatches to Lie or Poisson 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
source
CTFlows.LieFunction

Lie 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
source

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
source

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]
source

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]
source
CTFlows.LiftFunction
Lift(
    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 of X.

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
source
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) if autonomous=true and variable=false
  • (x, p, v) if autonomous=true and variable=true
  • (t, x, p) if autonomous=false and variable=false
  • (t, x, p, v) if autonomous=false and variable=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
source
CTModels.ModelType
struct 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}

source
CTBase.ParsingErrorType
struct 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
source
CTFlows.PoissonFunction
Poisson(
    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
source
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
source
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
source
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
source
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
source
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
source
CTModels.SolutionType
struct 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

source
CTFlows.VectorFieldType
struct 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])
source
OptimalControl.available_methodsFunction
available_methods(

) -> Tuple{Vararg{Tuple{Symbol, Symbol, Symbol}}}

Return the list of available methods that can be used to solve optimal control problems.

source
CTDirect.build_OCP_solutionFunction
build_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)

source
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)

source
ExaModels.constraintMethod
constraint(
    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.
source
CTModels.constraintsFunction
constraints(
    ocp::Model{<:CTModels.TimeDependence, <:CTModels.AbstractTimesModel, <:CTModels.AbstractStateModel, <:CTModels.AbstractControlModel, <:CTModels.AbstractVariableModel, <:Function, <:CTModels.AbstractObjectiveModel, C<:CTModels.AbstractConstraintsModel}
) -> CTModels.AbstractConstraintsModel

Return the constraints struct.

source
CTModels.controlFunction
control(
    ocp::Model{<:CTModels.TimeDependence, <:CTModels.TimesModel, <:CTModels.AbstractStateModel, T<:CTModels.AbstractControlModel}
) -> CTModels.AbstractControlModel

Return the control struct.

source
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
source
CTModels.control_componentsFunction
control_components(ocp::Model) -> Vector{String}

Return the names of the components of the control.

source
control_components(sol::Solution) -> Vector{String}

Return the names of the components of the control.

source
CTModels.control_dimensionFunction
control_dimension(ocp::Model) -> Int64

Return the control dimension.

source
control_dimension(sol::Solution) -> Int64

Return the dimension of the control.

source
CTModels.control_nameFunction
control_name(ocp::Model) -> String

Return the name of the control.

source
control_name(sol::Solution) -> String

Return the name of the control.

source
CTModels.costateFunction
costate(
    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
source
CTModels.criterionFunction
criterion(model::CTModels.MayerObjectiveModel) -> Symbol

Return the criterion (:min or :max).

source
criterion(model::CTModels.LagrangeObjectiveModel) -> Symbol

Return the criterion (:min or :max).

source
criterion(model::CTModels.BolzaObjectiveModel) -> Symbol

Return the criterion (:min or :max).

source
criterion(ocp::Model) -> Symbol

Return the type of criterion (:min or :max).

source
CTParser.@defMacro

Define 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
source
CTModels.definitionFunction
definition(ocp::Model) -> Expr

Return the model definition of the optimal control problem.

source
definition(ocp::CTModels.PreModel) -> Union{Nothing, Expr}

Return the model definition of the optimal control problem or nothing.

source
CTDirect.direct_transcriptionFunction
direct_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] or exa) 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

source
CTModels.dualFunction
dual(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.

source
CTModels.dynamicsFunction
dynamics(
    ocp::Model{<:CTModels.TimeDependence, <:CTModels.AbstractTimesModel, <:CTModels.AbstractStateModel, <:CTModels.AbstractControlModel, <:CTModels.AbstractVariableModel, D<:Function}
) -> Function

Return the dynamics.

source
CTModels.export_ocp_solutionFunction
export_ocp_solution(
    sol::CTModels.AbstractSolution;
    format,
    filename
)

Export a solution in JLD or JSON formats. Redirect to one of the methods:

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
source
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"
source
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"
source
CTModels.final_timeFunction
final_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.

source
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.

source
final_time(ocp::CTModels.AbstractModel) -> Real
source
final_time(
    ocp::CTModels.AbstractModel,
    variable::AbstractVector
) -> Any
source
final_time(
    ocp::Model{<:CTModels.TimeDependence, <:CTModels.TimesModel{<:CTModels.AbstractTimeModel, CTModels.FixedTimeModel{T<:Real}}}
) -> Real

Return the final time, for a fixed final time.

source
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.

source
final_time(
    ocp::Model{<:CTModels.TimeDependence, <:CTModels.TimesModel{<:CTModels.AbstractTimeModel, CTModels.FreeTimeModel}},
    variable::Real
) -> Real

Return the final time, for a free final time.

source
CTModels.final_time_nameFunction
final_time_name(model::CTModels.TimesModel) -> String

Get the name of the final time from the times model.

source
final_time_name(ocp::Model) -> String

Return the name of the final time.

source
final_time_name(sol::Solution) -> String

Return the name of the final time.

source
CTModels.get_build_examodelFunction
get_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.

source
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.

source
CTModels.has_fixed_final_timeFunction
has_fixed_final_time(
    times::CTModels.TimesModel{<:CTModels.AbstractTimeModel, <:CTModels.FixedTimeModel{T<:Real}}
) -> Bool

Check if the final time is fixed. Return true.

source
has_fixed_final_time(
    times::CTModels.TimesModel{<:CTModels.AbstractTimeModel, CTModels.FreeTimeModel}
) -> Bool

Check if the final time is free. Return false.

source
has_fixed_final_time(ocp::Model) -> Bool

Check if the final time is fixed.

source
CTModels.has_fixed_initial_timeFunction
has_fixed_initial_time(
    times::CTModels.TimesModel{<:CTModels.FixedTimeModel{T<:Real}}
) -> Bool

Check if the initial time is fixed. Return true.

source
has_fixed_initial_time(
    times::CTModels.TimesModel{CTModels.FreeTimeModel}
) -> Bool

Check if the initial time is free. Return false.

source
has_fixed_initial_time(ocp::Model) -> Bool

Check if the initial time is fixed.

source
CTModels.has_free_final_timeFunction
has_free_final_time(times::CTModels.TimesModel) -> Bool

Check if the final time is free.

source
has_free_final_time(ocp::Model) -> Bool

Check if the final time is free.

source
CTModels.has_free_initial_timeFunction
has_free_initial_time(times::CTModels.TimesModel) -> Bool

Check if the final time is free.

source
has_free_initial_time(ocp::Model) -> Bool

Check if the initial time is free.

source
CTModels.has_lagrange_costFunction
has_lagrange_cost(_::CTModels.MayerObjectiveModel) -> Bool

Return false.

source
has_lagrange_cost(
    _::CTModels.LagrangeObjectiveModel
) -> Bool

Return true.

source
has_lagrange_cost(_::CTModels.BolzaObjectiveModel) -> Bool

Return true.

source
has_lagrange_cost(ocp::Model) -> Bool

Check if the model has a Lagrange cost.

source
CTModels.has_mayer_costFunction
has_mayer_cost(_::CTModels.MayerObjectiveModel) -> Bool

Return true.

source
has_mayer_cost(_::CTModels.LagrangeObjectiveModel) -> Bool

Return false.

source
has_mayer_cost(_::CTModels.BolzaObjectiveModel) -> Bool

Return true.

source
has_mayer_cost(ocp::Model) -> Bool

Check if the model has a Mayer cost.

source
CTModels.import_ocp_solutionFunction
import_ocp_solution(
    ocp::CTModels.AbstractModel;
    format,
    filename
) -> Any

Import a solution from a JLD or JSON file. Redirect to one of the methods:

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
source
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")
source
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")
source
CTModels.infosFunction
infos(sol::Solution) -> Dict{Symbol, Any}

Return a dictionary of additional infos depending on the solver or nothing.

source
CTModels.initial_timeFunction
initial_time(
    model::CTModels.TimesModel{<:CTModels.FixedTimeModel{T<:Real}}
) -> Real

Get the initial time from the times model, from a fixed initial time model.

source
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.

source
initial_time(
    ocp::Model{<:CTModels.TimeDependence, <:CTModels.TimesModel{CTModels.FixedTimeModel{T<:Real}}}
) -> Real

Return the initial time, for a fixed initial time.

source
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.

source
initial_time(
    ocp::Model{<:CTModels.TimeDependence, <:CTModels.TimesModel{CTModels.FreeTimeModel}},
    variable::Real
) -> Real

Return the initial time, for a free initial time.

source
CTModels.initial_time_nameFunction
initial_time_name(model::CTModels.TimesModel) -> String

Get the name of the initial time from the times model.

source
initial_time_name(ocp::Model) -> String

Return the name of the initial time.

source
initial_time_name(sol::Solution) -> String

Return the name of the initial time.

source
CTModels.iterationsFunction
iterations(sol::Solution) -> Int64

Return the number of iterations (if solved by an iterative method).

source
CTModels.lagrangeFunction
lagrange(
    model::CTModels.LagrangeObjectiveModel{L<:Function}
) -> Function

Return the Lagrange function.

source
lagrange(
    model::CTModels.BolzaObjectiveModel{<:Function, L<:Function}
) -> Function

Return the Lagrange function.

source
lagrange(
    ocp::Model{<:CTModels.TimeDependence, <:CTModels.AbstractTimesModel, <:CTModels.AbstractStateModel, <:CTModels.AbstractControlModel, <:CTModels.AbstractVariableModel, <:Function, CTModels.LagrangeObjectiveModel{L<:Function}}
) -> Function

Return the Lagrange cost.

source
lagrange(
    ocp::Model{<:CTModels.TimeDependence, <:CTModels.AbstractTimesModel, <:CTModels.AbstractStateModel, <:CTModels.AbstractControlModel, <:CTModels.AbstractVariableModel, <:Function, <:CTModels.BolzaObjectiveModel{<:Function, L<:Function}}
) -> Any

Return the Lagrange cost.

source
CTModels.mayerFunction
mayer(
    model::CTModels.MayerObjectiveModel{M<:Function}
) -> Function

Return the Mayer function.

source
mayer(
    model::CTModels.BolzaObjectiveModel{M<:Function}
) -> Function

Return the Mayer function.

source
mayer(
    ocp::Model{<:CTModels.TimeDependence, <:CTModels.AbstractTimesModel, <:CTModels.AbstractStateModel, <:CTModels.AbstractControlModel, <:CTModels.AbstractVariableModel, <:Function, <:CTModels.MayerObjectiveModel{M<:Function}}
) -> Any

Return the Mayer cost.

source
mayer(
    ocp::Model{<:CTModels.TimeDependence, <:CTModels.AbstractTimesModel, <:CTModels.AbstractStateModel, <:CTModels.AbstractControlModel, <:CTModels.AbstractVariableModel, <:Function, <:CTModels.BolzaObjectiveModel{M<:Function}}
) -> Any

Return the Mayer cost.

source
CTModels.messageFunction
message(sol::Solution) -> String

Return the message associated to the status criterion.

source
CTDirect.modelFunction
model(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.

source
RecipesBase.plotMethod
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,
    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)       
source
RecipesBase.plot!Method
plot!(
    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.

source
CommonSolve.solveMethod
solve(
    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)
Note

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)
Note

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))
Note

See how to set an initial guess for more details.

source
CTModels.stateFunction
state(
    ocp::Model{<:CTModels.TimeDependence, <:CTModels.TimesModel, T<:CTModels.AbstractStateModel}
) -> CTModels.AbstractStateModel

Return the state struct.

source
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
source
CTModels.state_componentsFunction
state_components(ocp::Model) -> Vector{String}

Return the names of the components of the state.

source
state_components(sol::Solution) -> Vector{String}

Return the names of the components of the state.

source
CTModels.state_dimensionFunction
state_dimension(ocp::CTModels.PreModel) -> Int64
source
state_dimension(ocp::Model) -> Int64

Return the state dimension.

source
state_dimension(sol::Solution) -> Int64

Return the dimension of the state.

source
CTModels.state_nameFunction
state_name(ocp::Model) -> String

Return the name of the state.

source
state_name(sol::Solution) -> String

Return the name of the state.

source
CTModels.statusFunction
status(sol::Solution) -> Symbol

Return the status criterion (a Symbol).

source
CTModels.time_gridFunction
time_grid(
    sol::Solution{<:CTModels.TimeGridModel{T<:Union{StepRangeLen, AbstractVector{<:Real}}}}
) -> Union{StepRangeLen, AbstractVector{<:Real}}

Return the time grid.

source
CTModels.time_nameFunction
time_name(model::CTModels.TimesModel) -> String

Get the name of the time variable from the times model.

source
time_name(ocp::Model) -> String

Return the name of the time.

source
time_name(sol::Solution) -> String

Return the name of the time component.

source
CTModels.timesFunction
times(
    ocp::Model{<:CTModels.TimeDependence, T<:CTModels.TimesModel}
) -> CTModels.TimesModel

Return the times struct.

source
ExaModels.variableMethod
variable(
    sol::Solution
) -> Union{Real, AbstractVector{<:Real}}

Return the variable or nothing.

julia> v = variable(sol)
source
CTModels.variable_componentsFunction
variable_components(ocp::Model) -> Vector{String}

Return the names of the components of the variable.

source
variable_components(sol::Solution) -> Vector{String}

Return the names of the components of the variable.

source
CTModels.variable_dimensionFunction
variable_dimension(ocp::Model) -> Int64

Return the variable dimension.

source
variable_dimension(sol::Solution) -> Int64

Return the dimension of the variable.

source
CTModels.variable_nameFunction
variable_name(ocp::Model) -> String

Return the name of the variable.

source
variable_name(sol::Solution) -> String

Return the name of the variable.

source
CTFlows.:⋅Function

Lie 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
source

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
source

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
source