The solve function
In this manual, we explain the solve
function from OptimalControl.jl package.
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.
Basic usage
Let us define a basic optimal control problem.
using OptimalControl
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
We can now solve the problem:
using NLPModelsIpopt
solve(ocp)
▫ This is OptimalControl version v1.1.0 running with: direct, adnlp, ipopt.
▫ The optimal control problem is solved with CTDirect version v0.15.1.
┌─ The NLP is modelled with ADNLPModels and solved with NLPModelsIpopt.
│
├─ Number of time steps⋅: 250
└─ Discretisation scheme: trapeze
▫ This is Ipopt version 3.14.17, running with linear solver MUMPS 5.8.0.
Number of nonzeros in equality constraint Jacobian...: 3005
Number of nonzeros in inequality constraint Jacobian.: 0
Number of nonzeros in Lagrangian Hessian.............: 251
Total number of variables............................: 1004
variables with only lower bounds: 0
variables with lower and upper bounds: 0
variables with only upper bounds: 0
Total number of equality constraints.................: 755
Total number of inequality constraints...............: 0
inequality constraints with only lower bounds: 0
inequality constraints with lower and upper bounds: 0
inequality constraints with only upper bounds: 0
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
0 1.0000000e-01 1.10e+00 2.73e-14 0.0 0.00e+00 - 0.00e+00 0.00e+00 0
1 -5.0000000e-03 7.36e-02 2.66e-15 -11.0 6.08e+00 - 1.00e+00 1.00e+00h 1
2 6.0003829e+00 1.78e-15 1.78e-15 -11.0 6.01e+00 - 1.00e+00 1.00e+00h 1
Number of Iterations....: 2
(scaled) (unscaled)
Objective...............: 6.0003828724303263e+00 6.0003828724303263e+00
Dual infeasibility......: 1.7763568394002505e-15 1.7763568394002505e-15
Constraint violation....: 1.7763568394002505e-15 1.7763568394002505e-15
Variable bound violation: 0.0000000000000000e+00 0.0000000000000000e+00
Complementarity.........: 0.0000000000000000e+00 0.0000000000000000e+00
Overall NLP error.......: 1.7763568394002505e-15 1.7763568394002505e-15
Number of objective function evaluations = 3
Number of objective gradient evaluations = 3
Number of equality constraint evaluations = 3
Number of inequality constraint evaluations = 0
Number of equality constraint Jacobian evaluations = 3
Number of inequality constraint Jacobian evaluations = 0
Number of Lagrangian Hessian evaluations = 2
Total seconds in IPOPT = 1.840
EXIT: Optimal Solution Found.
Note that we must import NLPModelsIpopt.jl before calling solve
. This is because the default method uses a direct approach, which transforms the optimal control problem into a nonlinear program (NLP) of the form:
\[\text{minimize}\quad F(y), \quad\text{subject to the constraints}\quad g(y) \le 0, \quad h(y) = 0. \]
Calling solve
without loading a NLP solver package first will notify the user:
julia> solve(ocp)
ERROR: ExtensionError. Please make: julia> using NLPModelsIpopt
Resolution methods and algorithms
OptimalControl.jl offers a list of methods. To get it, simply call available_methods
.
available_methods()
(:direct, :adnlp, :ipopt)
(:direct, :adnlp, :madnlp)
(:direct, :adnlp, :knitro)
(:direct, :exa, :ipopt)
(:direct, :exa, :madnlp)
(:direct, :exa, :knitro)
Each line is a method, with priority going from top to bottom. This means that
solve(ocp)
is equivalent to
solve(ocp, :direct, :adnlp, :ipopt)
- The first symbol refers to the general class of method. The only possible value is:
:direct
: currently, only the so-called direct approach is implemented. Direct methods discretise the original optimal control problem and solve the resulting NLP. In this case, the mainsolve
method redirects toCTDirect.solve
.
- The second symbol refers to the NLP modeler. The possible values are:
:adnlp
: the NLP problem is modeled by aADNLPModels.ADNLPModel
. It provides automatic differentiation (AD)-based models that follow the NLPModels.jl API.:exa
: the NLP problem is modeled by aExaModels.ExaModel
. It provides automatic differentiation and SIMD abstraction.
- The third symbol specifies the NLP solver. Possible values are:
:ipopt
: callsNLPModelsIpopt.ipopt
to solve the NLP problem.:madnlp
: creates a MadNLP.MadNLPSolver instance from the NLP problem and solve it. MadNLP.jl is an open-source solver in Julia implementing a filter line-search interior-point algorithm like Ipopt.:knitro
: uses the Knitro solver (license required).
For instance, let us try MadNLP solver with ExaModel modeller.
using MadNLP
ocp = @def begin
t ∈ [ t0, tf ], time
x = (q, v) ∈ R², state
u ∈ R, control
x(t0) == x0
x(tf) == [0, 0]
∂(q)(t) == v(t)
∂(v)(t) == u(t)
0.5∫( u(t)^2 ) → min
end
solve(ocp, :exa, :madnlp)
▫ This is OptimalControl version v1.1.0 running with: direct, exa, madnlp.
▫ The optimal control problem is solved with CTDirect version v0.15.1.
┌─ The NLP is modelled with ExaModels and solved with MadNLP.
│
├─ Number of time steps⋅: 250
└─ Discretisation scheme: trapeze
▫ This is MadNLP version v0.8.7, running with umfpack
Number of nonzeros in constraint Jacobian............: 2004
Number of nonzeros in Lagrangian Hessian.............: 1751
Total number of variables............................: 753
variables with only lower bounds: 0
variables with lower and upper bounds: 0
variables with only upper bounds: 0
Total number of equality constraints.................: 504
Total number of inequality constraints...............: 0
inequality constraints with only lower bounds: 0
inequality constraints with lower and upper bounds: 0
inequality constraints with only upper bounds: 0
iter objective inf_pr inf_du lg(mu) ||d|| lg(rg) alpha_du alpha_pr ls
0 5.0000000e-03 1.10e+00 5.55e-17 -1.0 0.00e+00 - 0.00e+00 0.00e+00 0
1 6.0003829e+00 1.85e-15 2.15e-13 -1.0 6.08e+00 - 1.00e+00 1.00e+00h 1
Number of Iterations....: 1
(scaled) (unscaled)
Objective...............: 6.0003828724303929e+00 6.0003828724303929e+00
Dual infeasibility......: 2.1493917756743031e-13 2.1493917756743031e-13
Constraint violation....: 1.8492152253912764e-15 1.8492152253912764e-15
Complementarity.........: 0.0000000000000000e+00 0.0000000000000000e+00
Overall NLP error.......: 2.1493917756743031e-13 2.1493917756743031e-13
Number of objective function evaluations = 2
Number of objective gradient evaluations = 2
Number of constraint evaluations = 2
Number of constraint Jacobian evaluations = 2
Number of Lagrangian Hessian evaluations = 1
Total wall-clock secs in solver (w/o fun. eval./lin. alg.) = 13.144
Total wall-clock secs in linear solver = 0.311
Total wall-clock secs in NLP function evaluations = 0.000
Total wall-clock secs = 13.455
EXIT: Optimal Solution Found (tol = 1.0e-08).
Note that you can provide a partial description. If multiple full descriptions contain it, priority is given to the first one in the list. Hence, all of the following calls are equivalent:
solve(ocp)
solve(ocp, :direct )
solve(ocp, :adnlp )
solve(ocp, :ipopt)
solve(ocp, :direct, :adnlp )
solve(ocp, :direct, :ipopt)
solve(ocp, :direct, :adnlp, :ipopt)
Direct method
The main options for the direct method, with their [default] values, are:
display
([true
],false
): settingdisplay = false
disables output.init
: information for the initial guess. It can be given as numerical values, functions, or an existing solution. See how to set an initial guess.grid_size
([250
]): number of time steps in the (uniform) time discretization grid. More precisely, ifN = grid_size
and the initial and final times aret0
andtf
, then the step lengthΔt = (tf - t0) / N
.time_grid
([nothing
]): explicit time grid (can be non-uniform). Iftime_grid = nothing
, a uniform grid of lengthgrid_size
is used.disc_method
([:trapeze
],:midpoint
,:euler
,:euler_implicit
,:gauss_legendre_2
,:gauss_legendre_3
): the discretisation scheme to transform the dynamics into nonlinear equations. See the discretization method tutorial for more details.adnlp_backend
([:optimized
],:manual
,:default
): backend used for automatic differentiation to create theADNLPModels.ADNLPModel
.
For advanced usage, see:
The main solve
method from OptimalControl.jl simply redirects to CTDirect.solve
in that case.
NLP solvers specific options
In addition to these options, all remaining keyword arguments passed to solve
will be transmitted to the NLP solver used.
Please check the list of Ipopt options and the NLPModelsIpopt.jl documentation.
sol = solve(ocp; max_iter=0, tol=1e-6, display=false)
iterations(sol)
0
Similarly, please check the MadNLP.jl documentation and the list of MadNLP.jl options.
sol = solve(ocp, :madnlp; max_iter=0, tol=1e-6, display=false)
iterations(sol)
0