control-toolbox / CTDirect.jl

Direct transcription of an optimal control problem and resolution
http://control-toolbox.org/CTDirect.jl/
MIT License
9 stars 5 forks source link

Adjoint (dynamics multipliers) parsing #184

Closed jbcaillau closed 1 month ago

jbcaillau commented 1 month ago

@PierreMartinon As discussed: probable shift in parse_DOCP_solution_dual when retrieving the multipliers associated with the dynamics (= adjoints). Bug seem to be activated for problems with nonlinear constraints (= not mere boxes on slices of state / control). See for instance this variant of Goddard that, contrary to this one, has a such a nonlinear constraint: https://github.com/control-toolbox/BVP-DAE/issues/1#issuecomment-2246348545

Similar behaviour for these other projects that also have nonlinear (control or mixed) constraints: @alesiagr https://github.com/mctao-inria/sail-therm @PhysicDev https://github.com/PhysicDev/McTao_lunar_solar_sail

jbcaillau commented 1 month ago

@PierreMartinon in the current code below, compared to the previous tags, I don't see where you would skip the multipliers associated with constraints other than the dynamics (e.g. state constraints) in order to retrieve only the dynamics:

"""
$(TYPEDSIGNATURES)

Recover OCP costate from DOCP multipliers
"""
# +++todo: parse multipliers for the path and boundary constraints (use tuples)
function parse_DOCP_solution_dual(docp, multipliers)

    # constraints, costate and constraints multipliers
    N = docp.dim_NLP_steps
    P = zeros(N, docp.dim_NLP_x)

    if !isnothing(multipliers)
        index = 1
        for i in 1:N
            # state equation multiplier for costate
            P[i,:] = multipliers[index:index+docp.dim_NLP_x-1]
            index = index + docp.dim_NLP_x
        end
    end

    return P
end
PierreMartinon commented 1 month ago

@jbcaillau
Should be fixed in 0.9.6, can you try updating ?

alesiagr commented 1 month ago

J'arrive pas à update la version 0.9.6, j'ai tjrs la 0.9.4. Sinon je pourrais tester de suite

jbcaillau commented 1 month ago

👍🏽 @PierreMartinon vu pour la ligne ci-dessous, je teste https://github.com/control-toolbox/CTDirect.jl/blob/d2321e293ae43d0187e53c901432f5916b902247/src/solution.jl#L46

jbcaillau commented 1 month ago

@PierreMartinon CTDirect 0.9.6 now available. Latency? Anyway, seems much better 🙂 @ocots @alesiagr IMG_3362

alesiagr commented 1 month ago

BINGO! @jbcaillau @PierreMartinon @ocots Merci 🙏🏼 Qu'ils sont beaux ces adjoints quand ils oscillent pas 🥹

PierreMartinon commented 1 month ago

Il faudra que je complete les tests pour le parsing de solution, verifier qu'on recupere bien ce qu'il faut. Actuellement on ne teste que l'objectif.

ocots commented 1 month ago

Top si on peut corriger ce bug. Il faudra relancer les simulations avec Paul.

ocots commented 1 month ago

Je ne sais pas si ça peut aider mais il y des choses sur l'adjoint dans le cas de contraintes d'états dans ce papier : https://inria.hal.science/hal-01367918

Il y a cette proposition page 15

Capture d’écran 2024-07-25 à 21 36 40

et cette figure juste après

Capture d’écran 2024-07-25 à 21 36 21
alesiagr commented 1 month ago

Je sais pas si ça s'adressait à moi ton message mais merci @ocots :D

ocots commented 1 month ago

Plutôt à @PierreMartinon et @jbcaillau mais tu peux regarder 😀

jbcaillau commented 1 month ago

Top si on peut corriger ce bug. Il faudra relancer les simulations avec Paul.

@PierreMartinon has corrected. Leaving the issue open so that we remember to add the relevant unitary test in CTDirect.

Check https://github.com/control-toolbox/BVP-DAE/issues/1#issuecomment-2251483987

jbcaillau commented 1 month ago

@PierreMartinon still missing unit tests before closing, right?

PierreMartinon commented 1 month ago

@PierreMartinon still missing unit tests before closing, right?

Yes. Ideally we would test the closeness of the interpolated trajectory and costate to a known functional solution (evaluate both on a given time grid then use isapprox) . Does anyone have a basic problem with simple expressions for x,u,p (preferably non-constant) ?

jbcaillau commented 1 month ago

@PierreMartinon is our basic example OK? (One of the adjoint is not constant.)

ocots commented 1 month ago

@PierreMartinon is our basic example OK? (One of the adjoint is not constant.)

There is a problem here : https://control-toolbox.org/OptimalControl.jl/stable/tutorial-double-integrator.html

PierreMartinon commented 1 month ago

@PierreMartinon is our basic example OK? (One of the adjoint is not constant.)

There is a problem here : control-toolbox.org/OptimalControl.jl/stable/tutorial-double-integrator.html

Ok, I'll take the double integrator ! @ocots, the second one does not converge on the page :D

ocots commented 1 month ago

@PierreMartinon is our basic example OK? (One of the adjoint is not constant.)

There is a problem here : control-toolbox.org/OptimalControl.jl/stable/tutorial-double-integrator.html

Ok, I'll take the double integrator ! @ocots, the second one does not converge on the page :D

Do you understand why?

PierreMartinon commented 1 month ago

I'll have a look as soon as I can get CTBase 0.12 (there is usually some latency)

ocots commented 1 month ago

I'll have a look as soon as I can get CTBase 0.12 (there is usually some latency)

It is ready. You can get it whenever you want.

PierreMartinon commented 1 month ago

When retesting locally the problem converges. sol

Let's try to rebuild the doc and check ?

PierreMartinon commented 1 month ago

On the other hand, tests for state, control and costate added and OK for the double integrator

ocots commented 1 month ago

It converged with 100 points bu not with 200. Here how many do you have?

PierreMartinon commented 1 month ago

New default in CTDirect is 250, so 250 probably

jbcaillau commented 1 month ago

New default in CTDirect is 250, so 250 probably

@PierreMartinon 👍🏽

NB. probably safer to use double integrator with $L^2$ cost (no issue of convergence, whatever the grid size)