OptimalControl news

[ocots]

OptimalControl v0.3.0

To have a look at the new version of OptimalControl, see

• basic functional
• goddard functional

Basic functional

By default, the Model is autonomous and the usage is scalar, that is the user has to treat one dimensional state, control, adjoint and others as scalar and not Vector of dimension 1.

ocp = Model() # equivalent to ocp = Model{:autonomous, :scalar}()
constraint!(ocp, :dynamics, (x, u) -> A*x + B*u)  # u and not u[1] since the control 
                                                  # is one dimensional and the usage is scalar
objective!(ocp, :lagrange, (x, u) -> 0.5u^2)      # idem for u

Remark.

• In the scalar usage, it works to replace u by u[1] in the example above.
• But you can't do things like in the following example if u is one dimensional:

ocp = Model{:autonomous, :vectorial}()
constraint!(ocp, :dynamics, (x, u) -> A*x + B*u) # u[1] is required since the usage is vectorial

If you take a look at the output of the direct solver, you can see that I have removed the Ipopt header:

sol = solve(ocp)

produces

Method = (:direct, :ADNLProblem, :ipopt)
This is Ipopt version 3.14.4, running with linear solver MUMPS 5.5.1.

Number of nonzeros in equality constraint Jacobian...:     1205
Number of nonzeros in inequality constraint Jacobian.:        0
Number of nonzeros in Lagrangian Hessian.............:    81810

See also the plot of a solution. A solution from a direct or a direct shooting method is of same Type which is OptimalControlSolution. I have written a plot of an OptimalControlSolution that you can find in CTBase.

plot(sol)

Remark. The height of the control plot depends on the number of controls and states.

[ocots]

Goddard functional

Again, since the control is one dimensional and the usage is scalar, cf. ocp = Model(), then we write the dynamics as follows:

f(x, u) = F0(x) + u*F1(x)

You can compare the height of the control plot with the basic case:

plot(direct_sol)

Remark. The direct method is very slow!

In the indirect method, we have a scalar usage too:

u0(x, p) = 0.  # instead of u0(x, p) = [0.]
u1(x, p) = 1.
us(x, p) = -H001(x, p) / H101(x, p)
ub(x, _) = -Ad(F0, g)(x) / Ad(F1, g)(x)

Remark. For the moment, when we have a saturated state constraint, we have to pass to Flow a "mixed" constraint which is a function of (x,u).

fb = Flow(ocp, ub, (x, _) -> g(x), μb)

In this new version of OptimalControl, the state, control and adjoint solutions are provided as functions of time. It implies some changes:

t = direct_sol.times
x = direct_sol.state   # it a function of time
u = direct_sol.control # it a function of time
p = direct_sol.adjoint # it a function of time
H1(t) = H1(x(t), p(t)) # the switching function

u_plot = plot(t, t -> u(t)[1], xlabel = "t", ylabel = "u", legend = false)
H1_plot = plot(t, H1, xlabel = "t", ylabel = "H1", legend = false)
g_plot = plot(t, g∘x, xlabel = "t", ylabel = "g", legend = false)
display(plot(u_plot, H1_plot, g_plot, layout=(3,1), size=(400,300)))

t13 = t[ abs.(H1.(t)) .≤ η ]
t23 = t[ 0 .≤ (g∘x).(t) .≤ η ]
p0 = p(t0)

Finally, for the moment, there is no automatic plot of a solution of a flow:

f1sb0 = f1 * (t1, fs) * (t2, fb) * (t3, f0) # concatenation of the Hamiltonian flows
flow_sol = f1sb0((t0, tf), x0, p0)
r_plot = plot(flow_sol, idxs=(0, 1), xlabel="t", label="r")
v_plot = plot(flow_sol, idxs=(0, 2), xlabel="t", label="v")
m_plot = plot(flow_sol, idxs=(0, 3), xlabel="t", label="m")
pr_plot = plot(flow_sol, idxs=(0, 4), xlabel="t", label="p_r")
pv_plot = plot(flow_sol, idxs=(0, 5), xlabel="t", label="p_v")
pm_plot = plot(flow_sol, idxs=(0, 6), xlabel="t", label="p_m")
plot(r_plot, pr_plot, v_plot, pv_plot, m_plot, pv_plot, layout=(3, 2), size=(600, 300))

[PierreMartinon]

@jbcaillau @ocots @joseph-gergaud Est-ce qu'on pourrait faire un notebook julia qui illustre la facilite d'utilisation des differentes methodes, avec

• definition fonctionnelle (unique) du probleme toto. [Plus tard definition abstraite suivant avancement parser]
• solve(toto, direct) avec plot
• solve(toto, tir direct) avec plot
• solve(toto, indirect) avec plot

Je pense que voir les resolutions avec un appel solve par methode sur le meme probleme defini une seule fois sera plutot convaincant.

[ocots]

For the moment, the direct shooting method solves only problems without path constraints, but with fixed initial condition and with a given final constraint. Besides, there is no automatic indirect shooting method. But it could be nice for the new version of OptimalControl to have the notebook you mentioned. We need to had Ipopt for the direct shooting method and we need to create a CTIndirect.jl package.

[joseph-gergaud]

The only case where it should be possible is the case of a ''smooth'' optimal control , so without singular arc, active path constraints, ...
For the moment we have the case of the Goddard problem

[jbcaillau]

@PierreMartinon Can be done on basic.jl

On this:

• regarding doc, maybe a markdown file containing julia code 😱 will do better... as it is the format accepted by Documenter.jl
• still nice to have real code (.jl, notebooks) somewhere that can be run and reproduced

[PierreMartinon]

Also the basic integrator

…

On Mar 19, 2023 at 14:24, joseph-gergaud ***@***.***> wrote: The only case where it should be possible is the case of a ''smooth'' optimal control , so without singular arc, active path constraints, ... For the moment we have the case of the Goddard problem — Reply to this email directly, view it on GitHub, or unsubscribe. You are receiving this because you commented.Message ID: ***@***.***>