up tuto lqr

Author: ocots

docs/make.jl

@@ -28,7 +28,7 @@ makedocs(;
             "NLP manipulations" => "tutorial-nlp.md",
             "Indirect simple shooting" => "tutorial-iss.md",
             "Goddard: direct, indirect" => "tutorial-goddard.md",
-            "Linear–quadratic regulator" => "tutorial-lqr-basic.md",
+            "Linear–quadratic regulator" => "tutorial-lqr.md",
             "Minimal action" => "tutorial-mam.md",
             "Model Predictive Control" => "tutorial-mpc.md",
         ],

docs/src/tutorial-continuation.md

@@ -10,7 +10,7 @@ We illustrate this using a simple double integrator problem, where the fixed fin
 
 First we load the required packages:
 
-```@example main
+```@example main-cont
 using DataFrames
 using OptimalControl
 using NLPModelsIpopt
@@ -20,7 +20,7 @@ using Plots
 
 and write a function that returns the OCP for a given final time:
 
-```@example main
+```@example main-cont
 function problem(T)
 
     ocp = @def begin
@@ -49,7 +49,7 @@ nothing # hide
 
 Then we perform the continuation with a simple *for* loop, using each solution to initialize the next problem.
 
-```@example main
+```@example main-cont
 init = ()
 data = DataFrame(T=Float64[], Objective=Float64[], Iterations=Int[])
 for T ∈ range(1, 2, length=5)
@@ -67,7 +67,7 @@ As a second example, we show how to avoid redefining a new optimal control probl
 
 Let us first define the Goddard problem. Note that the formulation below illustrates all types of constraints, and the problem could be written more compactly.
 
-```@example main
+```@example main-cont
 # Parameters
 r0 = 1
 v0 = 0
@@ -127,7 +127,7 @@ sol0 = solve(goddard; display=false)
 
 Then, we perform the continuation on the maximal thrust.
 
-```@example main
+```@example main-cont
 sol = sol0 # Initialize the solution with the reference solution
 data = DataFrame(Tmax=Float64[], Objective=Float64[], Iterations=Int[])
 for Tmax_local ∈ range(Tmax_0, Tmax_f, length=6)
@@ -140,7 +140,7 @@ println(data)
 
 We plot now the objective with respect to the maximal thrust, as well as both solutions for `Tmax=3.5` and `Tmax=1`.
 
-```@example main
+```@example main-cont
 using Plots.PlotMeasures # for leftmargin
 
 plt_obj = plot(data.Tmax, data.Objective;

docs/src/tutorial-discretisation.md

@@ -5,7 +5,7 @@ These methods are used to convert a continuous-time optimal control problem (OCP
 
 Let us import the necessary packages and define the optimal control problem ([Goddard problem](@ref tutorial-goddard)) we will use as an example throughout this tutorial.
 
-```@example main
+```@example main-disc
 using BenchmarkTools  # for benchmarking
 using DataFrames      # to store the results
 using OptimalControl  # to define the optimal control problem and more
@@ -64,14 +64,14 @@ When calling `solve`, the option `disc_method=...` can be used to specify the di
 
 Let us first solve the problem with the default `:trapeze` method and display the solution.
 
-```@example main
+```@example main-disc
 sol = solve(ocp; disc_method=:trapeze, display=false)
 plot(sol; size=(800, 800))
 ```
 
 Let us now compare different discretization schemes to evaluate their accuracy and performance.
 
-```@example main
+```@example main-disc
 # Solve the problem with different discretization methods
 solutions = []
 data = DataFrame(
@@ -102,7 +102,7 @@ end
 println(data)
 ```
 
-```@example main
+```@example main-disc
 # Plot the results
 x_style = (legend=:none,)
 p_style = (legend=:none,)
@@ -120,14 +120,14 @@ plot(plt; size=(800, 800))
 
 For some large problems, you may notice that the solver takes a long time before the iterations actually begin. This is due to the computation of sparse derivatives — specifically, the Jacobian of the constraints and the Hessian of the Lagrangian — which can be quite costly. One possible alternative is to set the option `adnlp_backend=:manual`, which uses simpler sparsity patterns. The resulting matrices are faster to compute but are also less sparse, so this represents a trade-off between automatic differentiation (AD) preparation time and the efficiency of the optimization itself.
 
-```@example main
+```@example main-disc
 solve(ocp; disc_method=:gauss_legendre_3, grid_size=1000, adnlp_backend=:manual)
 nothing # hide
 ```
 
 Let us now compare the performance of the two backends. We will use the `@btimed` macro from the `BenchmarkTools` package to measure the time taken for both the preparation of the NLP problem and the execution of the solver. Besides, we will also collect the number of non-zero elements in the Jacobian and Hessian matrices, which can be useful to understand the sparsity of the problem, thanks to the functions `get_nnzo`, `get_nnzj`, and `get_nnzj` from the `NLPModels` package. The problem is first discretized with the `direct_transcription` method and then solved with the `ipopt` solver, see the [tutorial on direct transcription](@ref tutorial-nlp) for more details.
 
-```@example main
+```@example main-disc
 # DataFrame to store the results
 data = DataFrame(
     Backend=Symbol[],
@@ -191,7 +191,7 @@ println(data)
 
 The option `time_grid=...` allows you to provide the full time grid vector `t0, t1, ..., tf`, which is especially useful if a non-uniform grid is desired. In the case of a free initial and/or final time, you should provide a normalized grid ranging from 0 to 1. Note that `time_grid` overrides `grid_size` if both options are specified.
 
-```@example main
+```@example main-disc
 sol = solve(ocp; time_grid=[0, 0.1, 0.5, 0.9, 1], display=false)
 println(time_grid(sol))
 ```

docs/src/tutorial-goddard.md

@@ -39,7 +39,7 @@ We import the [Plots.jl](https://docs.juliaplots.org) package to plot the soluti
 The [OrdinaryDiffEq.jl](https://docs.sciml.ai/OrdinaryDiffEq) package is used to 
 define the shooting function for the indirect method and the [MINPACK.jl](https://github.com/sglyon/MINPACK.jl) package permits to solve the shooting equation.
 
-```@example main
+```@example main-goddard
 using OptimalControl  # to define the optimal control problem and more
 using NLPModelsIpopt  # to solve the problem via a direct method
 using OrdinaryDiffEq  # to get the Flow function from OptimalControl
@@ -51,7 +51,7 @@ using Plots           # to plot the solution
 
 We define the problem
 
-```@example main
+```@example main-goddard
 const t0 = 0      # initial time
 const r0 = 1      # initial altitude
 const v0 = 0      # initial speed
@@ -101,14 +101,14 @@ nothing # hide
 
 We then solve it
 
-```@example main
+```@example main-goddard
 direct_sol = solve(ocp; grid_size=100)
 nothing # hide
 ```
 
 and plot the solution
 
-```@example main
+```@example main-goddard
 plt = plot(direct_sol, label="direct", size=(800, 800))
 ```
 
@@ -119,7 +119,7 @@ bang arc with maximal control, followed by a singular arc, then by a boundary ar
 arc is with zero control. Note that the switching function vanishes along the singular and
 boundary arcs.
 
-```@example main
+```@example main-goddard
 t = time_grid(direct_sol)   # the time grid as a vector
 x = state(direct_sol)       # the state as a function of time
 u = control(direct_sol)     # the control as a function of time
@@ -191,7 +191,7 @@ as well as the associated multiplier for the *order one* state constraint on the
 
 With the help of differential geometry primitives, these expressions are straightforwardly translated into Julia code:
 
-```@example main
+```@example main-goddard
 # Controls
 const u0 = 0                            # off control
 const u1 = 1                            # bang control
@@ -218,7 +218,7 @@ nothing # hide
 Then, we define the shooting function according to the optimal structure we have determined,
 that is a concatenation of four arcs.
 
-```@example main
+```@example main-goddard
 x0 = [r0, v0, m0] # initial state
 
 function shoot!(s, p0, t1, t2, t3, tf)
@@ -245,7 +245,7 @@ To solve the problem by an indirect shooting method, we then need a good initial
 that is a good approximation of the initial costate, the three switching times and the
 final time.
 
-```@example main
+```@example main-goddard
 η = 1e-3
 t13 = t[ abs.(φ.(t)) .≤ η ]
 t23 = t[ 0 .≤ (g ∘ x).(t) .≤ η ]
@@ -272,7 +272,7 @@ println("\nNorm of the shooting function: ‖s‖ = ", norm(s), "\n")
 
 We aggregate the data to define the initial guess vector.
 
-```@example main
+```@example main-goddard
 ξ = [p0..., t1, t2, t3, tf] # initial guess
 ```
 
@@ -304,7 +304,7 @@ function fsolve(f, j, x; kwargs...)
 end
 ```
 
-```@example main
+```@example main-goddard
 using DifferentiationInterface
 import ForwardDiff
 backend = AutoForwardDiff()
@@ -313,7 +313,7 @@ nothing # hide
 
 Let us define the problem to solve.
 
-```@example main
+```@example main-goddard
 # auxiliary function with aggregated inputs
 nle!  = ( s, ξ) -> shoot!(s, ξ[1:3], ξ[4], ξ[5], ξ[6], ξ[7])
 
@@ -326,7 +326,7 @@ We are now in position to solve the problem with the `hybrj` solver from MINPACK
 function, providing the Jacobian.
 Let us solve the problem and retrieve the initial costate solution.
 
-```@example main
+```@example main-goddard
 # resolution of S(ξ) = 0
 indirect_sol = fsolve(nle!, jnle!, ξ, show_trace=true)
 
@@ -355,7 +355,7 @@ println("\nNorm of the shooting function: ‖s‖ = ", norm(s), "\n")
 We plot the solution of the indirect solution (in red) over the solution of the direct method 
 (in blue).
 
-```@example main
+```@example main-goddard
 f = f1 * (t1, fs) * (t2, fb) * (t3, f0) # concatenation of the flows
 flow_sol = f((t0, tf), x0, p0)          # compute the solution: state, costate, control...