added: multiple shooting transcription for LinMPC and NonLinMPC

[franckgaga]

@baggepinnen I will probably finish the multiple shooting transcription soon. LinMPC works well, still missing some work for NonLinMPC

The documentation is ready for a quick revision, do you see any mistake or something that is not clear in the docstring of LinMPC, NonLinMPC, TranscriptionMethod, SingleShooting and MultipleShooting , here (there are only minor changes in LinMPC and NonLinMPC docstring):

https://juliacontrol.github.io/ModelPredictiveControl.jl/previews/PR155

There is also the internal init_defectmat that you can look quickly, for linear plant models.

P.S. Multiple dispatch is so powerful! Now I dispatch on the type of the plant model, the predictive controller, the transcription method and the state estimator sometimes.

[codecov-commenter]

Codecov Report

Attention: Patch coverage is 83.16498% with 100 lines in your changes missing coverage. Please review.

Project coverage is 96.29%. Comparing base (01e48a4) to head (a868fc4).

Files with missing lines	Patch %	Lines
src/controller/transcription.jl	77.72%	49 Missing ⚠️
src/controller/nonlinmpc.jl	73.28%	39 Missing ⚠️
src/controller/construct.jl	88.50%	10 Missing ⚠️
src/controller/execute.jl	96.92%	2 Missing ⚠️

Additional details and impacted files

@@            Coverage Diff             @@
##             main     #155      +/-   ##
==========================================
- Coverage   99.02%   96.29%   -2.74%     
==========================================
  Files          24       25       +1     
  Lines        3790     4017     +227     
==========================================
+ Hits         3753     3868     +115     
- Misses         37      149     +112     

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

[franckgaga]

@baggepinnen I'm curious, do you know any paper that compares a single direct shooting transcription to a multiple shooting one in the context of MPC ? It seems that everyone is saying "It scales better" but no one really prove it lol.

[baggepinnen]

Do you notice any performance and/or convergence benefit from this new transcription?

There does not appear to be any tests?

---

Edit: I noticed now that some test had been added after I initiated the review yesterday

[baggepinnen]

expressed as deviation vectors from the operating point

is this the case also for nonlinear MPC?

[franckgaga]

Yes, that's also the case for nonlinear MPC. This allows NonLinModel with nonzero operating points, for example:

using Plots, ModelPredictiveControl, ControlSystemsBase
model1 = setop!(LinModel(tf(5, [10, 1]), 0.5), uop=[50], yop=[30])
f!(xnext,x,u,_,p) = (xnext .= p.A*x + p.B*u)
h!(y,x,_,p) = (y .= p.C*x)
p = (;A=model1.A, B=model1.Bu, C=model1.C)
model2 = NonLinModel(f!, h!, model1.Ts, 1, 1, 1; p, solver=nothing)
model2 = setop!(model2, uop=[50], yop=[30])
mpc = NonLinMPC(model2, Hp=10, transcription=MultipleShooting())
sim!(mpc, 100, x_0=[0], x̂_0=[0, 0]) |> plot

and it will work exactly like a linear MPC:

mpclinear = LinMPC(model1, Hp=10, transcription=MultipleShooting())
sim!(mpclinear, 100, x_0=[0], x̂_0=[0, 0]) |> plot

It's useful for the tests, maybe not so for real life applications. 🤔

edit: maybe useful for really simple empirical but nonlinear model, e.g. models with $x^2$ or $xu$ terms, but identified around a steady-state operating point. I work a bit with them during my PhD, as simple emprical model for batch fluidized bed dryer. Linear plant model are terrible for batch processes...

[franckgaga]

Tagging you since you are probably not notified since the PR is merged @baggepinnen

[franckgaga]

Will add a short comment on that in the docstring

[baggepinnen]

do you know any paper that compares a single direct shooting transcription to a multiple shooting one in the context of MPC ? It seems that everyone is saying "It scales better" but no one really prove it lol.

I don't have a reference in mind, but it aligns very well with my own experience. The scaling is not the primary benefit, it's the convergence behavior. Getting single shooting without feedback (like in iLQG) to converge to a reasonable optimum is hopeless for many problems that multiple shooting handles effortlessly.

[franckgaga]

@baggepinnen

Thanks for the review, I'll apply them in the next PR!

Ok I see, I just did not come across problems in which convergence is really hard with a single shooting transcription.

edit: FYI, I tested it on the unstable pendulum model of the manual, and it's about 5-6 times slower than the single shooting transcription. I discuss a little bit more in details my hypothesis of why on the discourse post.

[franckgaga]

@baggepinnen

I ran some additional tests on the inverted pendulum. My example in the manual uses a control horizon of $H_c=2$, which is just easy to solve in general. I stand by my view: a $H_c$ smaller than $H_p$ is an extremely useful feature to drastically ease the optimization (and also reduce the aggressiveness).

Now, If I chose $H_c=H_p=20$, the single shooting transcription is not able to converge at all anymore at the inverted position, and the multiple shooting transcription is able (but it's too slow for real-time control). My conclusion: multiple shooting is interesting for problems in which high control horizons $H_c$ are needed. I'll write an addendum on discourse.

[baggepinnen]

Are you computing dense Jacobians with ForwardDiff? If so, that requires $O((n_x+n_u)N)$ duals making the total cost scale as $O((n_x+n_u{)}^2{N}^2)$ instead of $O((n_x+n_u)N)$ so it's not really a fair comparison. In JSControl, we use sparse ad, including sparse hessian for the multiple-shooting and collocation backends and they tend to be plenty fast

[franckgaga]

Yes, all the AD computations are dense right now in the package. That's what JuMP does when user-defined operators are defined without specifying a gradient and a hessian function.

That's why I plan to integrate DI.jl in the upcoming release.