OptimalControl.jl

Appearance

Solution object ​

Everything you do after something has been computed: read the trajectories, check whether it converged, read the sensitivities. This page mirrors Inspect a problem — that page reads a model back, this one reads a solution back — and the two sections link to each other throughout.

What you get back ​

solve returns a Solution; a Flow call returns a trajectory. Both are read with the same generics — state, control, costate, objective, time_grid, plot — so everything on this page also applies to a flow's output (see Flows).

julia
using OptimalControl
using NLPModelsIpopt

t0 = 0
tf = 1
x0 = [-1, 0]

ocp = @def begin
    t ∈ [t0, tf], time
    x = (q, v) ∈ R², state
    u ∈ R, control
    x(t0) == x0
    x(tf) == [0, 0]
    ẋ(t) == [v(t), u(t)]
    0.5∫(u(t)^2) → min
end

sol = solve(ocp; display=false)

Solution and AbstractSolution are imported, not exported — they exist for dispatch, not for writing type annotations in your own code, so a signature like f(sol::Solution) = ... won't work as written; you'd need f(sol::OptimalControl.Solution) = ... or just not annotate.

Trajectories ​

state, control, variable, and costate return functions of time (except variable, which is a single vector — variables don't vary with time):

julia
x = state(sol)
u = control(sol)
p = costate(sol)
x(0.25), u(0.25), p(0.25)
([-0.8437454999279995, 1.124993999904], 3.0000480007680177, [12.00019200307206, 2.976047616761874])

These functions interpolate — they can be called anywhere in the time horizon, not just at grid points:

julia
0.25 ∈ time_grid(sol)
false
julia
x(0.25)  # still works
2-element Vector{Float64}:
 -0.8437454999279995
  1.124993999904

time_grid(sol) returns the discretization nodes; times(sol) returns the richer TimesModel struct.

1-D is a scalar here too: with a 1-D control, u(t) is a Number, not a length-1 vector — same rule as everywhere else on the site (functional-API callbacks, abstract syntax):

julia
typeof(u(0.25))
Float64

The objective ​

julia
objective(sol)
6.000096001536038

Did it work ​

julia
successful(sol), status(sol), message(sol)
(true, :first_order, "Ipopt/generic")
julia
iterations(sol), constraints_violation(sol)
(1, 2.220446049250313e-16)

infos(sol) returns a Dict of anything else the solver reported.

success is not successful

success(sol) looks like it should work but doesn't — it isn't and never was a CTModels method (it resolves to Base.success, a process-exit-status function). Calling it on a Solution now throws a migration-pointing error:

julia
julia> success(sol)
PreconditionError → success, deprecated.jl:75
│
│  `success(sol)` is deprecated
│
│  Reason  this spelling was removed in v2.1.0-beta
│
│  Hint    use successful(sol)
└─

See Migration for the full list of renamed spellings.

Dual variables ​

Dual variables (Lagrange multipliers) give sensitivity information. A richer problem to show them on:

julia
ocp = @def begin
    tf ∈ R, variable
    t ∈ [0, tf], time
    x = (q, v) ∈ R², state
    u ∈ R, control
    tf ≥ 0, (eq_tf)
    -1 ≤ u(t) ≤ 1, (eq_u)
    v(t) ≤ 0.75, (eq_v)
    x(0) == [-1, 0], (eq_x0)
    q(tf) == 0
    v(tf) == 0
    ẋ(t) == [v(t), u(t)]
    tf → min
end
sol = solve(ocp; display=false)

dual(sol, ocp, :label) returns the signed multiplier for a labeled constraint — a scalar for a variable or boundary constraint, a function of time for a path constraint:

julia
dual(sol, ocp, :eq_tf)   # variable constraint
4.8000002513961475e-12
julia
dual(sol, ocp, :eq_x0)   # boundary constraint
2-element Vector{Float64}:
 1.3298101221163192
 1.0011103794866778
julia
μ_u = dual(sol, ocp, :eq_u)   # path constraint — a function of time
μ_u(0.5)
-0.0027555370529424735

Sign convention

μ > 0 means the lower-side constraint is active, μ < 0 the upper-side, μ = 0 inactive. For box constraints the solver reports separate non-negative lower/upper multipliers internally; dual combines them as μ = μ_lb − μ_ub per component. The raw, unsigned, non-negative versions are available separately — see below.

The box-constraint duals, without going through a label — one accessor per group, no per-label lookup needed:

julia
state_constraints_lb_dual(sol), state_constraints_ub_dual(sol)
(CoercedTrajectory(identity), CoercedTrajectory(identity))
julia
control_constraints_lb_dual(sol), control_constraints_ub_dual(sol)
(CoercedTrajectory(only), CoercedTrajectory(only))
julia
variable_constraints_lb_dual(sol), variable_constraints_ub_dual(sol)
([4.8000002513961475e-12], [0.0])

with matching dimension accessors:

julia
dim_dual_state_constraints_box(sol), dim_dual_control_constraints_box(sol), dim_dual_variable_constraints_box(sol)
(2, 1, 1)

And the nonlinear (path/boundary) duals as a whole, plus their counts:

julia
path_constraints_dual(sol), boundary_constraints_dual(sol)
(CoercedTrajectory(identity), [1.3298101221163192, 1.0011103794866778, -1.3298101221163192, 1.0009448000722094])
julia
dim_path_constraints_nl(sol), dim_boundary_constraints_nl(sol)
(0, 4)

Back to the model ​

julia
model(sol) === ocp
true

model(sol) gives back the exact OCP the solution was computed from — everything on Inspect a problem works on it.

Empty solutions ​

julia
is_empty_time_grid(sol)
false

false for anything a normal solve produced — it's a defensive check for the placeholder case (an uninitialized or degenerate solution object carrying no time grid at all), not something a real solve leaves you needing to handle.

See also ​