feat: create JuMPControlProblem for optimal control

[vyudu]

will add optimal control tests here when the interface for inputs and stuff is settled.

[vyudu]

@ChrisRackauckas assuming this passes stuff I think it's ready.

[vyudu]

Requires SciML/SciMLBase.jl#996

[baggepinnen]

Do I have to define the variables to optimize over using x(..)? What if I want to solve an OC problem with components defined in MTKstdlib? Those components have x(t) variables. It would be nice to have a test for that use case as well since component-based modeling would be the fundamental reason to use MTK over any other tool

[ChrisRackauckas]

No, that's just if you also want to evaluate them at a point, since this just late binds t. We could also make an evaluate operator or something that is the same as evaluate(x(t), 0.2) == x(0.2)

[vyudu]

I will add a component test case

[baggepinnen]

there appears to be no tests that the dynamics are satisfied, testing the transcription. This could be tested by simulating the system using the found optimal control included with an interpolation

[vyudu]

Sorry yeah I've converted it to a draft for now since stuff is still being added.

https://infiniteopt.github.io/InfiniteOpt.jl/stable/develop/extensions/#Derivative-Evaluation-Methods looks like we can extend InfiniteOpt to add our own derivative methods. Should I just do this for all the ODE tableaus we have (like define a InfiniteOptTsit5, InfiniteOptImplicitEuler, InfiniteOptRadauIA3, etc in the extension)? It seems like they perform a lot better than doing the collocation using pure JuMP. Maybe as a later step if not in this PR?

This also means we wouldn't have a separate JuMPDynamicOptProblem and InfiniteOptDynamicOptProblem, it would all just be InfiniteOptDynamicOptProblem

[ChrisRackauckas]

Later step, let's do one step at a time.

It seems like they perform a lot better than doing the collocation using pure JuMP

That's interesting, I wonder why it's actually better.

This also means we wouldn't have a separate JuMPDynamicOptProblem and InfiniteOptDynamicOptProblem, it would all just be InfiniteOptDynamicOptProblem

I think for benchmarking purposes it's good to have both.

[baggepinnen]

While playing with this, I noticed the following

• The silent in solve(jprob, Ipopt.Optimizer, :RadauIIA5, silent = true) appears to have no effect, causing 1000s of lines of Ipopt printouts.
• DynamicOptSolution cannot be plotted (but sol.sol can)
• when plotting plot(isol.input_sol) for the cartpole swingup, the legend entry for the input signal is θ(t) isntead of u(t)
• Something appears to be wrong with isol.input_sol for the rocket launch problem. I expect the solution to look like this
  
  and jsol.input_sol is close (it seems the default tolerance is rather relaxed?)
  
  but isol.input_sol is completely wrong
  
  I tried running this test

    # Test solution
    @parameters (T_interp::CubicSpline)(..)
    eqs = [D(h(t)) ~ v(t),
        D(v(t)) ~ (T_interp(t) - drag(h(t), v(t))) / m(t) - gravity(h(t)),
        D(m(t)) ~ -T_interp(t) / c]
    @mtkbuild rocket_ode = ODESystem(eqs, t)
    interpmap = Dict(T_interp => ctrl_to_spline(jsol.input_sol, CubicSpline))
    oprob = ODEProblem(rocket_ode, u0map, (ts, te), merge(Dict(pmap), interpmap))
    osol = solve(oprob, RadauIIA5(); adaptive = false, dt = 0.005)
    @test ≈(jsol.sol.u, osol.u, rtol = 0.02)

for the isol as well but failed

[baggepinnen]

I tried specifying a minimum-time problem, but got an error

# Double integrator minimum time
t = M.t_nounits
D = M.D_nounits
@variables x(..) v(..)
@variables u(..) [bounds = (-1.0, 1.0), input = true]
@parameters tf
constr = [v(tf) ~ 0.0, x(tf) ~ 0]
cost = [tf] # Maximize the final distance.
@named block = ODESystem(
    [D(x(t)) ~ v(t), D(v(t)) ~ u(t)], t; costs = cost, constraints = constr)
block, input_idxs = structural_simplify(block, ([u(t)], []))

u0map = [x(t) => 1.0, v(t) => 0.0]
tspan = (0.0, tf)
parammap = [u(t) => 0.0, tf=>1.0]
jprob = JuMPDynamicOptProblem(block, u0map, tspan, parammap; steps = 51)
jsol = solve(jprob, Ipopt.Optimizer, :Verner8)

UndefVarError: `NotSymbolic` not defined in `MTKJuMPDynamicOptExt`

[vyudu]

Added your free final time problem as a test case, think it should be working now.

This rocket launch error with InfiniteOpt seems pretty elusive. I can't find any error with the model specification, but the dynamics that it comes up with are impossible given the controller...

What should the output of plotting the DynamicOptSolution? Two side-by-side plots of the inputs and states? I can add an overload for it where the plot recipes go.

[ChrisRackauckas]

DynamicOptSolution

Is there a reason to not just make the solution an ODESolution? Add a residual term to that?

[vyudu]

Just thought this would make it easier for plotting purposes. Symbolic indexing is kind of weird in this case since inputs are usually non-discrete parameters and doing sol.ps[u(t)] won't get you what you want so fetching the input solution from the solution might be a bit annoying if it's all put in one ODESolution.

[vyudu]

Also useful to return the JuMP model so they can inspect the objective value and stuff like that

[vyudu]

@baggepinnen For the rocket launch input solution being wrong, It may just be that OrthogonalCollocation for InfiniteOpt is a bit strange, I tested it with implicit Euler and got this, which still doesn't look like the original solution but gives the right dynamics when I simulate the ODEProblem using it