Interoperating direct and indirect optimal control solvers Olivier - PowerPoint PPT Presentation

Interoperating direct and indirect optimal control solvers Olivier Cots, Joint work with J.-B. Caillau, P. Martinon and the invaluable help of Inria Sophia SED 9th International Congress on Industrial and Applied Mathematics (ICIAM 2019), July 15-19, Valencia, Spain.

Contents • Context • Direct approach • Indirect approach • Interoperating direct and indirect solvers

Scientific context - Optimal control 5 0 r 3 −5 40 40 20 20 0 0 −20 −20 −40 −40 r 2 40 r 1 2 20 1 0 r 2 r 3 0 −1 −20 −2 −40 −50 0 50 −40 −20 0 20 40 r 1 r 2 1 Norm of thrust 0.5 0 0 20 40 60 80 100 120 140 time Minimum time control of the Kepler equation (CNES / TAS / Inria / CNRS collaboration): � t f min 0 � u � d t ≡ max m ( t f ) q = − µ � q � 3 + u ¨ m , m = − β � u � , ˙ � u � ≤ T max , t ∈ [0, t f ]. Fixed initial and final Keplerian orbits, free final time t f .

Scientific context - Optimal Control Problem � t f  t 0 f 0 ( t , x ( t ), u ( t )) d t min g ( x ( t 0 ), x ( t f )) +                x ( t ) = f ( t , x ( t ), u ( t )), ˙ u ( t ) ∈ U , t ∈ [ t 0 , t f ] a.e. ,    (OCP)    h ( t , x ( t ), u ( t )) ≤ 0, t ∈ [ t 0 , t f ]             c ( x ( t 0 ), x ( t f )) = 0,    • g : boundary cost • f 0 : running cost • f : dynamics • U : control domain • h : control/state constraints • c : boundary constraints

Scientific context - Numerical methods and softwares Numerical methods: • Direct and Indirect: simple shooting, multiple shooting, collocation. – indirect: OCP − → PMP − → NLE: Newton solvers. – direct: OCP − → NLP − → KKT: Interior points solvers, SQP solvers. • HJB (value function, cf Dyn. Prog.) and LMI (to compute a lower bound). • Homotopy techniques to solve a one-parameter family of ocp’s. • Conjugate points computation. Softwares: • Bocop : direct method (collocation) and “HJB”, • HamPath : indirect method (simple and multiple shooting), homotopy and conjugate points, • . . . See appendix A of the following ref for a survey on numerical tools: S. M. Rogers, Optimal Control of Nonholonomic Mechanical Systems , PHD Thesis, 2017.

Technological context - Bocop software The Bocop project (bocop.org) began in 2010 in the framework of the Inria-Saclay initiative for an open source optimal control toolbox, and is supported by team COMMANDS. • Direct method: OCP − → NLP − → KKT • interior point solver IPOPT • automatic differentiation Adol-C/CppAD

Technological context - HamPath software The HamPath package (hampath.org) is a an open-source software developed since 2009, jointly by CNRS (Toulouse, Dijon) and Inria (Sophia). • Indirect method: OCP − → PMP − → NLE • Newton solver: hybrj from MINPACK • automatic differentiation: TAPENADE • integration: Runge-Kutta schemes with variable step-size as dopri5 , dop853 , radau5 (9, 13) Subroutines (Matlab/Python) • hfun , hvfun , dhvfun , Subroutines (Fortran90) Inputs Outputs HamPath • hfun • exphvfun , expdhvfun , • sfun • sfun , sjac , • ssolve , hampath . . . Possibilities: indirect simple and multiple shooting, homotopy and conjugate points computation

Context - Objectives Inria project (started June 1st 2019): ct (Control Toolbox) - Towards a collaborative set of reference tools to solve ocp’s Objectives. • Interoperate Bocop and HamPath (complementary approaches), • Provide a high-end interface. People involved. • McTAO (Inria Sophia): J.-B. Caillau, J.-B. Pomet • APO (CNRS): O. Cots, J. Gergaud • CAGE (Inria Paris): P. Martinon • COMMANDS (Inria Saclay): F. Bonnans • Sophia SED Interoperating of Bocop and HamPath : GUI demo, Python notebook demo.

Direct approach - Bocop software • C / C++, linux/win/mac, GUI, EPL License • OCP, model identification, delay systems • General Runge-Kutta schemes, IPOPT solver, automatic differentiation (AdolC/CppAD)

Direct approach - Time discretization, Runge Kutta schemes � t f  t 0 f 0 ( t , x ( t ), u ( t )) d t min g ( x ( t 0 ), x ( t f )) +            x ( t ) = f ( t , x ( t ), u ( t )), ˙ u ( t ) ∈ U , t ∈ [ t 0 , t f ] a.e. ,      (OCP)  h ( t , x ( t ), u ( t )) ≤ 0, t ∈ [ t 0 , t f ],            c ( x ( t 0 ), x ( t f )) = 0,     Time discretization. ( t i ) i =0,..., N , usually uniform with step h . Decision variable. new state and control variables: X := ( x i , u i ) . NLP Transcription. The OCP is reformulated in terms of the unknown X . • objective: boundary cost g ( x 0 , x N ) , running cost as sum of terms f 0 ( t i , x i , u i ) . • boundary conditions: c ( x 0 , x N ) . • path constraints: h ( t i , x i , u i ) , i = 0, · · · , N . • dynamics: general Runge Kutta formulas

Direct approach - Time discretization, Runge Kutta schemes i − 1 s General RK formula. s stages, coefficients a ij , b i , c i such that i =1 b i = 1 and c i = j =1 a ij . � � Butcher form: The formula for one step is: c 1 a 11 · · · a 1 s s . . . k i = f ( t 0 + c i h , x 0 + h j =1 a ij k j ), for i ∈ � 1 , s � , ... � . . . . . . c s a s 1 · · · a ss s x 1 = x 0 + h i =1 b i k i , � b 1 · · · b s Setting the unknown X := ( x i , u i , k i ) , we finally obtain the discretized problem (NLP) min F ( X ), C LB ≤ C ( X ) ≤ C UB . → A direct method solves (NLP) as an approximation of (OCP). Remark. KKT for NLP tends to PMP for OCP when h → 0 .

Direct approach - Bocop - Goddard problem 1D rocket ascent with maximal final mass: • State: q = ( h , v , m ) • Objective: max m ( t f ) • Dynamics: ˙  = v h     = 1    ˙ m ( c u − D ( h , v )) − g ( h ) v     m = − b u ˙    • Constraints: t f free     u ∈ [0, 1]     h (0) = 1, v (0) = 0, m (0) = 1     h ( t f ) = 1.01   

Direct approach - Bocop - Goddard problem 1D rocket ascent with maximal final mass: solution of the form B + SB 0 . ++ no a priori on the structure, robust to the initial guess, – accuracy can be limited by discretization/structure. New Bocop : micro-swimmer example.

Indirect approach - HamPath software Fortran routines S p y,λ q H λ p z q Fortran core AD AD AD AD Ý Ñ B S B y p y,λ q B S B λ p y, λ q H λ p z q QR Newton d Ý Ñ T p r p s qq RK H λ p z q RK + Newton RK ssolve hampath exphvfun expdhvfun Fortran Fortran Fortran Fortran Interface routines ssolve hampath exphvfun expdhvfun Interface Interface Interface Interface

Indirect approach - Simple shooting - Introduction Let consider: J ( u ( · )) := 1 � t f  0 u ( t ) 2 d t − → min     2     (P 1 ) x ( t ) = − x ( t ) + u ( t ), ˙ u ( t ) ∈ R , t ∈ [0 , t f ] p.p. , x (0) = x 0 ,      x ( t f ) = x f ,    with t f := 1 , x 0 := − 1 , x f := 0 and ∀ t ∈ [0 , t f ] , x ( t ) ∈ R .

Indirect approach - Simple shooting - Introduction Let consider: J ( u ( · )) := 1 � t f  0 u ( t ) 2 d t − → min     2     (P 1 ) x ( t ) = − x ( t ) + u ( t ), ˙ u ( t ) ∈ R , t ∈ [0 , t f ] p.p. , x (0) = x 0 ,      x ( t f ) = x f ,    with t f := 1 , x 0 := − 1 , x f := 0 and ∀ t ∈ [0 , t f ] , x ( t ) ∈ R . p 0 = − 1 (cas normal) • Pseudo-Hamiltonian: H ( x , p , p 0 , u ) := p ( − x + u ) + p 0 1 2 u 2 , • Maximization condition: u s ( x , p ) := − p / p 0 = p , ∂ u 2 = p 0 = − 1 < 0 Å ∂ 2 H ã , • We get the boundary value problem (BVP) :  x ( t ) = ˙ ∂ p H [ t ] = − x ( t ) + u s ( x ( t ), p ( t )) = − x ( t ) + p ( t ),        p ( t ) = − ∂ x H [ t ] = ˙ p ( t ), (BVP 1 )    x (0) = x 0 , x ( t f ) = x f ,     where [ t ] := ( x ( t ), p ( t ), p 0 , u s ( x ( t ), p ( t ))) .

Indirect approach - Simple shooting - Introduction We want to solve (BVP 1 ):  x ( t ) = ˙ ∂ p H [ t ] = − x ( t ) + p ( t ),        p ( t ) = − ∂ x H [ t ] = ˙ p ( t ), (BVP 1 )     x (0) = x 0 , x ( t f ) = x f .    • We introduce : Ñ ∂ H ∂ p ( z , p 0 , u ), − ∂ H é # — ∂ x ( z , p 0 , u ) H ( z , u ) := , z := ( x , p ). z ( t ) = # — • We denote by z ( · , x 0 , p 0 ) the solution of ˙ H ( z ( t ), u s ( z ( t ))) , z (0) = ( x 0 , p 0 ) .

Indirect approach - Simple shooting - Introduction We want to solve (BVP 1 ):  x ( t ) = ˙ ∂ p H [ t ] = − x ( t ) + p ( t ),        p ( t ) = − ∂ x H [ t ] = ˙ p ( t ), (BVP 1 )     x (0) = x 0 , x ( t f ) = x f .    • We introduce : Ñ ∂ H ∂ p ( z , p 0 , u ), − ∂ H é # — ∂ x ( z , p 0 , u ) H ( z , u ) := , z := ( x , p ). z ( t ) = # — • We denote by z ( · , x 0 , p 0 ) the solution of ˙ H ( z ( t ), u s ( z ( t ))) , z (0) = ( x 0 , p 0 ) . • We define the shooting function by: S : R − → R → S ( p 0 ) := Π x ( z ( t f , x 0 , p 0 )) − x f , p 0 �− où Π x ( x , p ) = x . • Solve (BVP 1 ) amount to solve S ( p 0 ) = 0. This is the simple shooting method.