[Dev] Dual variables

[ocots]

If we have this:

ocp = @def begin

    t ∈ [0, 1], time
    x ∈ R², state
    u ∈ R, control

    x₂(t) ≤ 1.2
    x₂(t) ≤ 2.0
    x₂(t) ≤ 3.0

    x(0) == [-1, 0]
    x(1) == [0, 0]

    ẋ(t) == [x₂(t), u(t)]

    ∫( 0.5u(t)^2 ) → min

end;

when we solve the problem and call build_solution, do we pass in state_constraints_ub_dual a function that returns a vector of dimension 3 because there are 3 constraints on x₂ or a vector of dimension 2 since x is of dimension 2.

@PierreMartinon?

[jbcaillau]

@ocots @PierreMartinon not sure to understand: the dimension of the multiplier depends on the number of constraints, not on the state dimension; I am missing sth?

[ocots]

Usually yes, but for the solver I guess that we cannot give more than one bound and so I don't know this situation is handled. And if there is no box constraints I think that CTDirect returns anyway a multiplier but it is 0.

[ocots]

I need to know the dimensions of

    path_constraints_dual::TPCD=__constraints(),
    boundary_constraints_dual::Union{Vector{Float64},Nothing}=__constraints(),
    state_constraints_lb_dual::Union{Matrix{Float64},Nothing}=__constraints(),
    state_constraints_ub_dual::Union{Matrix{Float64},Nothing}=__constraints(),
    control_constraints_lb_dual::Union{Matrix{Float64},Nothing}=__constraints(),
    control_constraints_ub_dual::Union{Matrix{Float64},Nothing}=__constraints(),
    variable_constraints_lb_dual::Union{Vector{Float64},Nothing}=__constraints(),
    variable_constraints_ub_dual::Union{Vector{Float64},Nothing}=__constraints(),

when CTDirect calls

CTModels.jl/src/solution.jl

Lines 34 to 54 in 7b5021d

	function build_solution(
	ocp::Model,
	T::Vector{Float64},
	X::TX,
	U::TU,
	v::Vector{Float64},
	P::TP;
	objective::Float64,
	iterations::Int,
	constraints_violation::Float64,
	message::String,
	stopping::Symbol,
	success::Bool,
	path_constraints_dual::TPCD=__constraints(),
	boundary_constraints_dual::Union{Vector{Float64},Nothing}=__constraints(),
	state_constraints_lb_dual::Union{Matrix{Float64},Nothing}=__constraints(),
	state_constraints_ub_dual::Union{Matrix{Float64},Nothing}=__constraints(),
	control_constraints_lb_dual::Union{Matrix{Float64},Nothing}=__constraints(),
	control_constraints_ub_dual::Union{Matrix{Float64},Nothing}=__constraints(),
	variable_constraints_lb_dual::Union{Vector{Float64},Nothing}=__constraints(),
	variable_constraints_ub_dual::Union{Vector{Float64},Nothing}=__constraints(),

[jbcaillau]

@ocots @PierreMartinon OK. in the implementation for ExaModels in CTParser v0.3, box constraints are passed to the modeller when declaring the (discretised) variable. several bounds can be set, the last one replacing the previous (no check), see e.g. the case of boxes for state (aka state range): https://github.com/control-toolbox/CTParser.jl/blob/314399f5d0bd806c00c5e98d2d2a8c2f01926c2a/src/onepass.jl#L635

[ocots]

But when calling build_solution, I think that it is always the length of the state (for boxes for state).

[jbcaillau]

@PierreMartinon ?

[ocots]

@PierreMartinon any idea?

[PierreMartinon]

Ok, let me see. Keep in mind that everything in the DOCP is discretized, be it the variables or constraints.

• what I call in the code the 'box multipliers', i.e. the dual variables for the lower and upper bounds on the state / control / optimization variables, have the same dimension as the corresponding variables. This is why they are handled in parse_docp_primal, using the same getters as the variables. More precisely
  • state_constraints_lb_dual and state_constraints_ub_dual have the same dimension and layout as the discrete state variables [X_1(t0) .. X_n(t0), X_1(t1) .. X_n(t1), ..., X_1(tf) .. X_n(tf)]
  • control_constraints_lb_dual and control_constraints_ub_dual follow the discrete controls [U_1(t0) .. U_m(t0), U_1(t1) .. U_m(t1), ..., U_1(tf) .. U_m(tf)] for Trapeze, with Euler and Midpoint stopping at the penultimate time step instead.
  • variable_constraints_lb_dual and variables_constraints_ub_dual are just like v.
    Multipliers for missing bounds should be zero.

The dual variables for the path and boundary constraints are handled in parse_docp_dual.

• the multipliers for the path constraints follow the discretized path constraints, so path_constraints_lb_dual and path_constraints_ub_dual have the same dimension and layout as `[g_1(t0) .. g_p(t0), ..., g_1(tf) .. g_p(tf)]. By the way this is the same for all discretizations since we chose to evaluate the path constraints at the time steps (and not time stages).
• finally for the boundary conditions we simply have boundary_constraints_lb_dual and boundary_constraints_ub_dual organized as [bc(x0,xf)].

Where do you think would be the best place to describe the dimensions ? In the description of the output arguments of the two parsing functions, and repeated in the call to build_solution ?

@ocots Edit: the arrays are 'deinterlaced', i.e.

• the states X is a (N+1) x n array, same with the state bounds multipliers
• the controls U is (N+1) x m, same with bound multipliers. Note that component N+1 is a copy of component N if the discretization has no final control
• path constraints multipliers are passed as (N+1) x dim_path_constraints

In summary, everything that corresponds to a function of time is passed as a (N+1) x (dim) array, same with the corresponding discrete multipliers for constraints, including variable bounds.

I realize that the transpose (dim) x (N+1) may be more useful since we could isolate components with the first index, such as X[1,:] for the first state variable etc. This is trivial to change on my side, so tell me which orientation you prefer.

[jbcaillau]

@PierreMartinon @ocots we should settle control-toolbox/CTParser.jl#94 (comment) before going / updating parse_docp_dual

[ocots]

@PierreMartinon We can indicate this:

• at CTDirect level: https://control-toolbox.org/CTDirect.jl/stable/solution.html#Documentation
• at CTModels level: build_ocp_solution
• in the documentation (explain that if several bounds are set, the last one replaces the previous): https://control-toolbox.org/OptimalControl.jl/stable/tutorial-abstract.html#Constraints