Introduction

Discretise optimal control problems

An optimal control problem (OCP) with fixed initial and final times can be described as minimising the cost functional

\[g(x(t_0), x(t_f)) + \int_{t_0}^{t_f} f^{0}(t, x(t), u(t))~\mathrm{d}t\]

where the state $x$ and the control $u$ are functions of time $t$, subject for $t \in [t_0, t_f]$ to the differential constraint

\[ \dot{x}(t) = f(t, x(t), u(t))\]

and other constraints such as

\[\begin{array}{llcll} x_{\mathrm{lower}} & \le & x(t) & \le & x_{\mathrm{upper}}, \\ u_{\mathrm{lower}} & \le & u(t) & \le & u_{\mathrm{upper}}, \\ c_{\mathrm{lower}} & \le & c(t, x(t), u(t)) & \le & c_{\mathrm{upper}}, \\ b_{\mathrm{lower}} & \le & b(x(t_0), x(t_f)) & \le & b_{\mathrm{upper}}. \end{array}\]

The so-called direct approach transforms the infinite-dimensional optimal control problem (OCP) into a finite-dimensional optimisation problem (NLP). This is achieved by discretising time, typically with Runge–Kutta methods applied to the state, control variables, and dynamics equation. These methods are usually less precise than indirect methods based on Pontryagin’s Maximum Principle, but they are more robust with respect to initialisation. They are also easier to apply, which explains their widespread use in industrial applications.

In the OptimalControlProblems package, each OCP is discretised using the trapezoidal rule on a uniform grid:

\[\begin{array}{lclr} t \in [t_0,t_f] & \to & \{t_0, \ldots, t_N=t_f\} & \\[0.5em] x(\cdot),\, u(\cdot) & \to & X=\{x_0, \ldots, x_N, u_0, \ldots, u_N\} & \\[1em] \hline \\ \text{step} & \to & \displaystyle h = \frac{t_f-t_0}{N} & \\[0.5em] \text{criterion} & \to & \displaystyle g(x_0, x_N) + \frac{h}{2} \sum_{i=1}^{N} \left( f^0(t_i, x_i, u_i) + f^0(t_{i-1}, x_{i-1}, u_{i-1}) \right) & \\[1em] \text{dynamics} & \to & \displaystyle x_{i} = x_{i-1} + \frac{h}{2} \left( f(t_i, x_i, u_i) + f(t_{i-1}, x_{i-1}, u_{i-1}) \right), & i = 1:N \\[1em] \text{state constraints} & \to & x_{\mathrm{lower}} \le x_i \le x_{\mathrm{upper}}, & i = 0:N \\[1em] \text{control constraints} & \to & u_{\mathrm{lower}} \le u_i \le u_{\mathrm{upper}}, & i = 0:N \\[1em] \text{path constraints} & \to & c_{\mathrm{lower}} \le c(t_i, x_i, u_i) \le c_{\mathrm{upper}}, & i = 0:N \\[1em] \text{boundary constraints} & \to & b_{\mathrm{lower}} \le b(x_0, x_N) \le b_{\mathrm{upper}} & \end{array}\]

We therefore obtain a nonlinear programming problem (NLP) on the discretised state and control variables of the general form:

\[\text{(NLP)} \quad \left\{ \begin{array}{lr} \min \ F(X) \\[1em] X_{\mathrm{lower}} \le X \le X_{\mathrm{upper}}\\[0.5em] C_{\mathrm{lower}} \le C(X) \le C_{\mathrm{upper}} \end{array} \right.\]

JuMP and OptimalControl models

Each optimal control problem in the OptimalControlProblems package is modelled both in JuMP and in OptimalControl. The problem definitions are stored in the OptimalControlProblems.jl/ext directory:

• JuMP models are stored in the JuMPModels directory. These codes implement the NLP problem directly.
• OptimalControl models are stored in the OptimalControlModels directory. These codes represent the OCP, and the discretisation is handled by the package. The resulting NLP is represented by an ADNLPModels.ADNLPModel, which provides automatic differentiation (AD)-based models following the NLPModels.jl API.

For more specific details about the problems, see the following pages. We provide descriptions of the optimal control problems and compare the different models.

• Beam
• Bioreactor
• Cart pendulum
• Chain
• Dielectrophoretic particle
• Double oscillator
• Ducted fan
• Electric vehicle
• Glider
• Insurance
• Jackson
• Moonlander
• Robbins
• Robot
• Rocket
• Space shuttle
• Steering
• Truck trailer
• Vanderpol

Problems metadata

For each problem, additional data is provided in the MetaData directory:

julia> metadata[:my_problem][:name]
"My Problem"

To list all metadata, use metadata. To access the metadata of a specific problem, for example chain, run:

using OptimalControlProblems
metadata[:chain]

OrderedCollections.OrderedDict{Any, Any} with 8 entries:
  :name          => "chain"
  :N             => 500
  :minimise      => true
  :state_name    => ["x1", "x2", "x3"]
  :costate_name  => ["∂x1", "∂x2", "∂x3"]
  :control_name  => ["u"]
  :variable_name => nothing
  :final_time    => (:fixed, 1)

Problems characteristics

To get the list of available problems, call the problems method.

problems()

19-element Vector{Symbol}:
 :beam
 :bioreactor
 :cart_pendulum
 :chain
 :dielectrophoretic_particle
 :double_oscillator
 :ducted_fan
 :electric_vehicle
 :glider
 :insurance
 :jackson
 :moonlander
 :robbins
 :robot
 :rocket
 :space_shuttle
 :steering
 :truck_trailer
 :vanderpol

We detail below the characteristics of the optimal control problems (OCPs) and their associated nonlinear programming problems (NLPs).

using NLPModels                 # to get the number of variables and constraints
import DataFrames: DataFrame    # to store data
using OptimalControl

data_ocp = DataFrame(           # to store data of the OCPs
    Problem=Symbol[],
    State=Int[],
    Control=Int[],
    Variable=Int[],
    Cost=Symbol[],
    FinalTime=Symbol[],
    Constraints=String[],
)

data_nlp = DataFrame(           # to store data of the NLPs
    Problem=Symbol[],
    Steps=Int[],
    Variables=Int[],
    Constraints=Int[],
)

for problem in problems()

    #
    docp = eval(problem)(OptimalControlBackend())
    nlp = nlp_model(docp)
    ocp = ocp_model(docp)

    #
    cost = if has_mayer_cost(ocp) && has_lagrange_cost(ocp)
        :Bolza
    elseif has_mayer_cost(ocp)
        :Mayer
    else
        :Lagrange
    end

    #
    final_time = has_fixed_final_time(ocp) ? :fixed : :free

    #
    using CTModels # these functions should be exported by OptimalControl
    dim_state_cons_box = CTModels.dim_state_constraints_box(ocp)
    dim_control_cons_box = CTModels.dim_control_constraints_box(ocp)
    dim_variable_cons_box = CTModels.dim_variable_constraints_box(ocp)
    dim_path_cons_nl = CTModels.dim_path_constraints_nl(ocp)
    dim_boundary_cons_nl = CTModels.dim_boundary_constraints_nl(ocp)
    constraints = ""
    constraints *= dim_state_cons_box    > 0 ? "x" : ""
    constraints *= dim_control_cons_box  > 0 ? (isempty(constraints) ? "" : ", ") * "u" : ""
    constraints *= dim_variable_cons_box > 0 ? (isempty(constraints) ? "" : ", ") * "v" : ""
    constraints *= dim_path_cons_nl      > 0 ? (isempty(constraints) ? "" : ", ") * "p" : ""
    constraints *= dim_boundary_cons_nl  > 0 ? (isempty(constraints) ? "" : ", ") * "b" : ""
    constraints = if !isempty(constraints)
        "(" * constraints * ")"
    end

    #
    push!(data_ocp,
        (
            Problem=problem,
            State=state_dimension(ocp),
            Control=control_dimension(ocp),
            Variable=variable_dimension(ocp),
            Cost=cost,
            FinalTime=final_time,
            Constraints=constraints,
        )
    )

    #
    N = metadata[problem][:N] # get default number of steps

    push!(data_nlp,
        (
            Problem=problem,
            Steps=N,
            Variables=get_nvar(nlp),
            Constraints=get_ncon(nlp),
        )
    )
end

data_ocp

data_nlp