[Bug] Error matrix formulation

[ocots]

@jbcaillau This code produces an error:

julia> x0 = [ 0
              1 ]
2-element Vector{Int64}:
 0
 1


julia> A = [  0 1
             -1 0 ]
2×2 Matrix{Int64}:
  0  1
 -1  0


julia> B = [ 0
             1 ]
2-element Vector{Int64}:
 0
 1


julia> Q = [ 1 0
             0 1 ]
2×2 Matrix{Int64}:
 1  0
 0  1


julia> R = 1
1


julia> tf = 3
3


julia> ocp = @def begin
           t ∈ [0, tf], time
           x ∈ R², state
           u ∈ R, control
           x(0) == x0
           ẋ(t) == A * x(t) + B * u(t)
           0.5∫( x(t)' * Q * x(t) + u(t)' * R * u(t) ) → min
       end
    t ∈ [0, tf], time
    x ∈ R², state
    u ∈ R, control
    x(0) == x0
    ẋ(t) == A * x(t) + B * u(t)
    0.5 * ∫((x(t))' * Q * x(t) + (u(t))' * R * u(t)) → min

The optimal control problem is of the form:

    minimize  J(x, u) = ∫ f⁰(t, x(t), u(t)) dt, over [0, 3]

    subject to

        ẋ(t) = f(t, x(t), u(t)), t in [0, 3] a.e.,

        ϕ₋ ≤ ϕ(x(0), x(3)) ≤ ϕ₊,

    where x(t) ∈ R² and u(t) ∈ R.

Declarations (* required):
╭────────┬────────┬──────────┬──────────┬───────────┬────────────┬─────────────╮
│ times* │ state* │ control* │ variable │ dynamics* │ objective* │ constraints │
├────────┼────────┼──────────┼──────────┼───────────┼────────────┼─────────────┤
│   V    │   V    │    V     │    X     │     V     │     V      │      V      │
╰────────┴────────┴──────────┴──────────┴───────────┴────────────┴─────────────╯


julia> solve(ocp)
ERROR: Cannot determine ordering of Dual tags ForwardDiff.Tag{ReverseDiff.var"#131#134"{typeof(+), ForwardDiff.Dual{ForwardDiff.Tag{ADNLPModels.var"#ψ#72"{CTDirect.var"#14#16"{CTDirect.DOCP{CTDirect.Trapeze, Model{CTModels.TimesModel{CTModels.FixedTimeModel{Int64}, CTModels.FixedTimeModel{Int64}}, CTModels.StateModel, CTModels.ControlModel, CTModels.EmptyVariableModel, Main.var"Main".var"##484#5", CTModels.LagrangeObjectiveModel{Main.var"Main".var"##488#6"}, CTModels.ConstraintsModel{Tuple{Vector{Real}, CTModels.var"#path_cons_nl!#106"{Int64, Vector{Int64}, Tuple{}}, Vector{Real}, Vector{Symbol}}, Tuple{Vector{Real}, Main.var"Main".var"##478#4", Vector{Real}, Vector{Symbol}}, Tuple{Vector{Real}, Vector{Int64}, Vector{Real}, Vector{Symbol}}, Tuple{Vector{Real}, Vector{Int64}, Vector{Real}, Vector{Symbol}}, Tuple{Vector{Real}, Vector{Int64}, Vector{Real}, Vector{Symbol}}}}}}}, Float64}, Float64, 1}}, ForwardDiff.Dual{ForwardDiff.Tag{ADNLPModels.var"#ψ#72"{CTDirect.var"#14#16"{CTDirect.DOCP{CTDirect.Trapeze, Model{CTModels.TimesModel{CTModels.FixedTimeModel{Int64}, CTModels.FixedTimeModel{Int64}}, CTModels.StateModel, CTModels.ControlModel, CTModels.EmptyVariableModel, Main.var"Main".var"##484#5", CTModels.LagrangeObjectiveModel{Main.var"Main".var"##488#6"}, CTModels.ConstraintsModel{Tuple{Vector{Real}, CTModels.var"#path_cons_nl!#106"{Int64, Vector{Int64}, Tuple{}}, Vector{Real}, Vector{Symbol}}, Tuple{Vector{Real}, Main.var"Main".var"##478#4", Vector{Real}, Vector{Symbol}}, Tuple{Vector{Real}, Vector{Int64}, Vector{Real}, Vector{Symbol}}, Tuple{Vector{Real}, Vector{Int64}, Vector{Real}, Vector{Symbol}}, Tuple{Vector{Real}, Vector{Int64}, Vector{Real}, Vector{Symbol}}}}}}}, Float64}, Float64, 1}} and ForwardDiff.Tag{ADNLPModels.var"#ψ#72"{CTDirect.var"#14#16"{CTDirect.DOCP{CTDirect.Trapeze, Model{CTModels.TimesModel{CTModels.FixedTimeModel{Int64}, CTModels.FixedTimeModel{Int64}}, CTModels.StateModel, CTModels.ControlModel, CTModels.EmptyVariableModel, Main.var"Main".var"##484#5", CTModels.LagrangeObjectiveModel{Main.var"Main".var"##488#6"}, CTModels.ConstraintsModel{Tuple{Vector{Real}, CTModels.var"#path_cons_nl!#106"{Int64, Vector{Int64}, Tuple{}}, Vector{Real}, Vector{Symbol}}, Tuple{Vector{Real}, Main.var"Main".var"##478#4", Vector{Real}, Vector{Symbol}}, Tuple{Vector{Real}, Vector{Int64}, Vector{Real}, Vector{Symbol}}, Tuple{Vector{Real}, Vector{Int64}, Vector{Real}, Vector{Symbol}}, Tuple{Vector{Real}, Vector{Int64}, Vector{Real}, Vector{Symbol}}}}}}}, Float64}

[jbcaillau]

@ocots yet another issue with reverse diff. using forward over forward (can be changed, as pointed by @amontoison and checked by @PierreMartinon) should slow down the computation but solve the issue. pénible

[amontoison]

@jbcaillau @ocots I finally understand the issue and can (probably) fix it in ADNLPModels.jl. I need to force the same AD context for the two components of the Lagrangian (f(x) and lambda^T c(x)).

-> https://github.com/JuliaDiff/ForwardDiff.jl/blob/382363a196a34edc6f0916e782b30a403f5a8fa9/src/dual.jl#L37-L38

[ocots]

Writing a time minimum problem like this does not work neither:

 ∫( 1 ) → min

[jbcaillau]

Writing a time minimum problem like this does not work neither:

∫( 1 ) → min

same error? (dual tags...) what about ∫( 1 + 0u(t))?

[jbcaillau]

#481