Skip to content

Commit

Permalink
Merge pull request #157 from control-toolbox/doc
Browse files Browse the repository at this point in the history
Doc
  • Loading branch information
ocots authored Oct 5, 2023
2 parents 6abbd46 + 7eb662d commit 825c580
Show file tree
Hide file tree
Showing 4 changed files with 200 additions and 5 deletions.
2 changes: 1 addition & 1 deletion docs/make.jl
Original file line number Diff line number Diff line change
Expand Up @@ -23,12 +23,12 @@ makedocs(
"tutorial-basic-example-f.md",
"tutorial-double-integrator.md",
"tutorial-goddard.md",
"tutorial-iss.md",
"tutorial-init.md",
"tutorial-lqr-basic.md",
"tutorial-plot.md",
#"tutorial-model.md",
#"tutorial-solvers.md",
#"tutorial-iss.md",
#"tutorial-flows.md",
#"tutorial-ct_repl.md",
],
Expand Down
64 changes: 64 additions & 0 deletions docs/src/juliaopt2023.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,64 @@
```@raw html
<img width="800" alt="juliacon" src="./assets/juliacon.jpg">
```

# Solving optimal control problems with Julia
### [Jean-Baptiste Caillau](http://caillau.perso.math.cnrs.fr), [Olivier Cots](https://ocots.github.io), [Joseph Gergaud](https://scholar.google.com/citations?user=pkH4An4AAAAJ&hl=fr), [Pierre Martinon](https://www.linkedin.com/in/pierre-martinon-b4603a17), [Sophia Sed](https://iww.inria.fr/sed-sophia)

```@raw html
<img width="800" alt="affiliations" src="./assets/affil.jpg">
```

# What it's about
- Nonlinear optimal control of ODEs:
```math
g(x(t_0),x(t_f)) + \int_{t_0}^{t_f} f^0(x(t), u(t)) \to \min
```
subject to
```math
\dot{x}(t) = f(x(t), u(t)),\quad t \in [0, t_f]
```
plus boundary, control and state constraints
- Our core interests: numerical & geometrical methods in control, applications

# Where it comes from
- [BOCOP: the optimal control solver](https://www.bocop.org)
- [HamPath: indirect and Hamiltonian pathfollowing](http://www.hampath.org)
- [Coupling direct and indirect solvers, examples](https://ct.gitlabpages.inria.fr/gallery//notebooks.html)

# OptimalControl.jl
- [Basic example: double integrator](https://control-toolbox.org/docs/optimalcontrol/dev/tutorial-basic-example-f.html)
- [Basic example: double integrator (cont'ed)](https://control-toolbox.org/docs/optimalcontrol/dev/tutorial-basic-example.html)
- [Advanced example: Goddard problem](https://control-toolbox.org/docs/optimalcontrol/dev/tutorial-goddard.html)

# Wrap up
- [X] High level modelling of optimal control problems
- [X] Efficient numerical resolution coupling direct and indirect methods
- [X] Collection of examples

# Future
- [ct_repl](./assets/repl.mp4)
- Additional solvers: direct shooting, collocation for BVP, Hamiltonian pathfollowing...
- ... and open to contributions!
- [CTProblems.jl](https://control-toolbox.org/docs/ctproblems/stable/problems-list.html)

# control-toolbox.org
- Open toolbox
- Collection of Julia Packages rooted at [OptimalControl.jl](https://control-toolbox.org/docs/optimalcontrol)

```@raw html
<a href="https://control-toolbox.org"><img width="800" alt="control-toolbox.org" src="./assets/control-toolbox.jpg"></a>
```

# Credits (not exhaustive!)
- [DifferentialEquations.jl](https://github.com/SciML/DifferentialEquations.jl)
- [JuMP](https://jump.dev/JuMP.jl),
[InfiniteOpt.jl](https://docs.juliahub.com/InfiniteOpt/p3GvY/0.4.1),
[ADNLPModels.jl](https://jso.dev/ADNLPModels.jl)
- [Ipopt](https://github.com/coin-or/ipopt)
- [JuliaDiff](https://juliadiff.org)
([FowardDiff.jl](https://juliadiff.org/ForwardDiff.jl),
[Zygote.jl](https://fluxml.ai/Zygote.jl))
- [MLStyle.jl](https://thautwarm.github.io/MLStyle.jl)
- [REPLMaker.jl](https://docs.juliahub.com/ReplMaker)

8 changes: 5 additions & 3 deletions docs/src/tutorial-goddard.md
Original file line number Diff line number Diff line change
Expand Up @@ -277,9 +277,11 @@ using MINPACK # NLE solver
nle = (s, ξ) -> shoot!(s, ξ[1:3], ξ[4], ξ[5], ξ[6], ξ[7]) # auxiliary function
# with aggregated inputs
ξ = [ p0 ; t1 ; t2 ; t3 ; tf ] # initial guess
indirect_sol = @suppress_err begin # hide
fsolve(nle, ξ)
end # hide
global indirect_sol = # hide
@suppress_err begin # hide
fsolve(nle, ξ) # hide
indirect_sol = fsolve(nle, ξ) # resolution of S(ξ) = 0
end # hide
# we retrieve the costate solution together with the times
p0 = indirect_sol.x[1:3]
Expand Down
131 changes: 130 additions & 1 deletion docs/src/tutorial-iss.md
Original file line number Diff line number Diff line change
@@ -1,3 +1,132 @@
# Indirect simple shooting

In construction.
In this tutorial we present the indirect simple shooting method on a simple example.

```@setup main
using Suppressor # to suppress warnings
```

Let us start by importing the necessary packages.

```@example main
using OptimalControl
using MINPACK # NLE solver
```

Let us consider the following optimal control problem.

```@example main
t0 = 0
tf = 1
x0 = -1
xf = 0
α = 1.5
@def ocp begin
t ∈ [ t0, tf ], time
x ∈ R, state
u ∈ R, control
x(t0) == x0
x(tf) == xf
ẋ(t) == -x(t) + α * x(t)^2 + u(t)
∫( 0.5u(t)^2 ) → min
end;
nothing # hide
```

From the Pontryagin Maximum Principle, we have that the maximising control is given by

```math
u(x, p) = p
```

and we can define the flow of the associated Hamiltonian vector field

```@example main
u(x, p) = p
f = Flow(ocp, u)
nothing # hide
```

We define an auxiliary exponential map.

```@example main
exp(p0; saveat=[]) = f((t0, tf), x0, p0, saveat=saveat).ode_sol
nothing # hide
```

We are now in position to define the shooting function and solve the shooting equations.

```@example main
S(p0) = exp(p0)(tf)[1] - xf; # shooting function
nle = (s, ξ) -> s[1] = S(ξ[1]) # auxiliary function
ξ = [ 0.0 ] # initial guess
global indirect_sol = # hide
@suppress_err begin # hide
fsolve(nle, ξ) # hide
indirect_sol = fsolve(nle, ξ) # resolution of S(p0) = 0
end # hide
p0_sol = indirect_sol.x[1] # costate solution
println("costate: p0 = ", p0_sol)
@suppress_err begin # hide
println("shoot: |S(p0)| = ", abs(S(p0_sol)), "\n")
end # hide
nothing # hide
```

We get:

```@example main
# parameters
times = range(t0, tf, length=2)
p0min = -0.5
p0max = 2
# plot of the flow
plt_flow = plot()
p0s = range(p0min, p0max, length=20)
for i ∈ 1:length(p0s) # plot some trajectories
sol = exp(p0s[i])
x = [sol(t)[1] for t ∈ sol.t]
p = [sol(t)[2] for t ∈ sol.t]
label = i==1 ? "extremals" : false
plot!(plt_flow, x, p, color=:blue, label=label)
end
p0s = range(p0min, p0max, length=200)
xs = zeros(length(p0s), length(times))
ps = zeros(length(p0s), length(times))
for i ∈ 1:length(p0s)
sol = exp(p0s[i], saveat=times)
xs[i, :] = [z[1] for z ∈ sol.(times)]
ps[i, :] = [z[2] for z ∈ sol.(times)]
end
for j ∈ 1:length(times) # plot intermediate points
label = j==1 ? "flow at times" : false
plot!(plt_flow, xs[:, j], ps[:, j], color=:green, linewidth=2, label=label)
end
plot!(plt_flow, xlims = (-1.1, 1), ylims = (p0min, p0max))
plot!(plt_flow, [0, 0], [p0min, p0max], color=:black, xlabel="x", ylabel="p", label="x=xf")
# solution
sol = exp(p0_sol)
x = [sol(t)[1] for t ∈ sol.t]
p = [sol(t)[2] for t ∈ sol.t]
plot!(plt_flow, x, p, color=:red, linewidth=2, label="extremal solution")
plot!(plt_flow, [x[end]], [p[end]], seriestype=:scatter, color=:green, label=false) # solution
# plot of the shooting function
plt_shoot = plot(xlims=(p0min, p0max), ylims=(-2, 4), xlabel="p₀", ylabel="y")
plot!(plt_shoot, p0s, S, linewidth=2, label="S(p₀)", color=:green)
plot!(plt_shoot, [p0min, p0max], [0, 0], color=:black, label="y=0")
plot!(plt_shoot, [p0_sol, p0_sol], [-2, 0], color=:black, label="p₀ solution", linestyle=:dash)
plot!(plt_shoot, [p0_sol], [0], seriestype=:scatter, color=:green, label=false) # solution
# plots
plot(plt_flow, plt_shoot, layout=(1,2), size=(800, 450))
```

0 comments on commit 825c580

Please sign in to comment.