added: be more efficient with weights and conversion matrices for NonLinMPC

[franckgaga]

Following #198 and discussion at #193, this is a summary of the revision for structural sparsity in PredictiveController types.

The big picture here is as follows: since it is meant for real-time applications, the computational performance of matrix operations has a higher priority than the memory footprint for storing the matrices.

Conversion of decision vector Z̃ to input increment ΔŨ

Since the $\mathbf{\tilde{P}_{\Delta u}}$ matrix is required at controller construction in the case of linear plant model, it is still constructed as a normal dense matrix. But I no longer use the matrix in getΔŨ! function (which is called multiple times in the objective function of NonLinMPC). I now use indexing with a @views.

Conversion of decision vector Z̃ to input U

I made no change here. The $\mathbf{\tilde{P}_{u}}$ matrix contains around 30 to 40% 1.0 values with a SingleShooting transcription (the rest is zeros, it depends of $H_p$ and $H_c$ parameters). It is very unlikely that any special sparse type will outperform computations with a dense matrix (a matrix product). Also note that I tested in the past storing the matrix as dense array of Bool or a BitArray. The matrix product was always slower than a dense array of floats.

Moreover, this matrix will have a more weird and complex structure with the new upcoming move blocking feature. In this case, doing the conversion manually with a for loop and conditionals in getU0! would be prone to bugs, and presumably less efficient than a product with a dense matrix with a nice contiguous structure in memory.

Objective function weights

VERY frequently, the weights will be diagonal matrices. But in very rare case, they will be generic positive semi-definite hermitian matrices (e.g. tuning rules based on optimization or frequency weighting). Ideally, we need to support both cases, but without loosing the performance advantage of matrix operation with pure diagonal matrices.

Before this PR, the weight matrices were simply stored as Hermitian{NT, Matrix{NT}} in the ControllerWeights object to be as generic as possible. Doing so, we loose the performance boost of matrix products with Diagonal types, in the case of diagonal weights. This also allows the user to specify block-diagonal weights as e.g. SparseMatrixCSC{Float64, Int64} and the type will be preserved in mpc.weights structure.

As suggested by @gdalle in his PR, I introduce 3 new parameters in the parametric struct ControllerWeights to store the type of the 3 weights, in order to preserve special types like Diagonal{NT, Vector{NT}}. The performance advantage will mainly visible in NonLinMPC based on nonlinear plant models, since the objective value is computed with e.g. dot(Ȳ, mpc.weights.M_Hp, Ȳ), which is faster when mpc.weights.M_Hp isa Diagonal:

julia> using LinearAlgebra; A = Hermitian(diagm(1.0:1000.0)); x=rand(1000);

julia> using BenchmarkTools; @btime dot($x, $A, $x)
  60.053 μs (0 allocations: 0 bytes)
175696.98591124764

julia> A2 = Hermitian(Diagonal(1.0:1000.0));

julia> using BenchmarkTools; @btime dot($x, $A2, $x)
  1.439 μs (0 allocations: 0 bytes)
175696.98591124764

The default values of the weights are now also specialized e.g. Diagonal(repeat(Mwt, Hp)) instead of diagm(repeat(Mwt, Hp)).

[codecov-commenter]

Codecov Report

All modified and coverable lines are covered by tests ✅

Project coverage is 98.84%. Comparing base (5c85eeb) to head (52fe019).
Report is 9 commits behind head on main.

Additional details and impacted files

@@            Coverage Diff             @@
##             main     #202      +/-   ##
==========================================
+ Coverage   98.81%   98.84%   +0.02%     
==========================================
  Files          25       25              
  Lines        4223     4225       +2     
==========================================
+ Hits         4173     4176       +3     
+ Misses         50       49       -1     

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

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.