added: reduce PredictiveController allocs with new buffer object

Author: franckgaga

docs/src/internals/predictive_control.md

@@ -34,4 +34,5 @@ ModelPredictiveControl.linconstraint!(::PredictiveController, ::LinModel)
 
 ```@docs
 ModelPredictiveControl.optim_objective!(::PredictiveController)
+ModelPredictiveControl.getinput
 ```

src/controller/execute.jl

@@ -52,9 +52,9 @@ julia> u = moveinput!(mpc, ry); round.(u, digits=3)
 function moveinput!(
     mpc::PredictiveController, 
     ry::Vector = mpc.estim.model.yop, 
-    d ::Vector = mpc.estim.buffer.empty;
-    Dhat ::Vector = repeat(d,  mpc.Hp),
-    Rhaty::Vector = repeat(ry, mpc.Hp),
+    d ::Vector = mpc.buffer.empty;
+    Dhat ::Vector = repeat!(mpc.buffer.D̂,  d,  mpc.Hp),
+    Rhaty::Vector = repeat!(mpc.buffer.R̂y, ry, mpc.Hp),
     Rhatu::Vector = mpc.Uop,
     D̂  = Dhat,
     R̂y = Rhaty,
@@ -67,13 +67,9 @@ function moveinput!(
     initpred!(mpc, mpc.estim.model, d, D̂, R̂y, R̂u)
     linconstraint!(mpc, mpc.estim.model)
     ΔŨ = optim_objective!(mpc)
-    Δu = ΔŨ[1:mpc.estim.model.nu] # receding horizon principle: only Δu(k) is used (1st one)
-    u = mpc.estim.lastu0 + mpc.estim.model.uop + Δu
-    return u
+    return getinput(mpc, ΔŨ)
 end
 
-
-
 @doc raw"""
     getinfo(mpc::PredictiveController) -> info
 
@@ -506,6 +502,24 @@ function preparestate!(mpc::PredictiveController, ym, d=mpc.estim.buffer.empty)
     return preparestate!(mpc.estim, ym, d)
 end
 
+@doc raw"""
+    getinput(mpc::PredictiveController, ΔŨ) -> u
+
+Get current manipulated input `u` from a [`PredictiveController`](@ref) solution `ΔŨ`.
+
+The first manipulated input ``\mathbf{u}(k)`` is extracted from the input increments vector
+``\mathbf{ΔŨ}`` and applied on the plant (from the receding horizon principle).
+"""
+function getinput(mpc, ΔŨ)
+    Δu  = mpc.buffer.u
+    for i in 1:mpc.estim.model.nu
+        Δu[i] = ΔŨ[i]
+    end
+    u   = Δu
+    u .+= mpc.estim.lastu0 .+ mpc.estim.model.uop
+    return u
+end
+
 """
     updatestate!(mpc::PredictiveController, u, ym, d=[]) -> x̂next
 

src/controller/explicitmpc.jl

@@ -34,6 +34,7 @@ struct ExplicitMPC{NT<:Real, SE<:StateEstimator} <: PredictiveController{NT}
     Uop::Vector{NT}
     Yop::Vector{NT}
     Dop::Vector{NT}
+    buffer::PredictiveControllerBuffer{NT}
     function ExplicitMPC{NT, SE}(
         estim::SE, Hp, Hc, M_Hp, N_Hc, L_Hp
     ) where {NT<:Real, SE<:StateEstimator}
@@ -65,6 +66,7 @@ struct ExplicitMPC{NT<:Real, SE<:StateEstimator} <: PredictiveController{NT}
         Uop, Yop, Dop = repeat(model.uop, Hp), repeat(model.yop, Hp), repeat(model.dop, Hp)
         nΔŨ = size(Ẽ, 2)
         ΔŨ = zeros(NT, nΔŨ)
+        buffer = PredictiveControllerBuffer{NT}(nu, ny, nd, Hp)
         mpc = new{NT, SE}(
             estim,
             ΔŨ, ŷ,
@@ -78,6 +80,7 @@ struct ExplicitMPC{NT<:Real, SE<:StateEstimator} <: PredictiveController{NT}
             Ks, Ps,
             d0, D̂0, D̂E,
             Uop, Yop, Dop,
+            buffer
         )
         return mpc
     end

src/controller/linmpc.jl

@@ -43,6 +43,7 @@ struct LinMPC{
     Uop::Vector{NT}
     Yop::Vector{NT}
     Dop::Vector{NT}
+    buffer::PredictiveControllerBuffer{NT}
     function LinMPC{NT, SE, JM}(
         estim::SE, Hp, Hc, M_Hp, N_Hc, L_Hp, Cwt, optim::JM
     ) where {NT<:Real, SE<:StateEstimator, JM<:JuMP.GenericModel}
@@ -74,6 +75,7 @@ struct LinMPC{
         Uop, Yop, Dop = repeat(model.uop, Hp), repeat(model.yop, Hp), repeat(model.dop, Hp)
         nΔŨ = size(Ẽ, 2)
         ΔŨ = zeros(NT, nΔŨ)
+        buffer = PredictiveControllerBuffer{NT}(nu, ny, nd, Hp)
         mpc = new{NT, SE, JM}(
             estim, optim, con,
             ΔŨ, ŷ,
@@ -86,6 +88,7 @@ struct LinMPC{
             Ks, Ps,
             d0, D̂0, D̂E,
             Uop, Yop, Dop,
+            buffer
         )
         init_optimization!(mpc, model, optim)
         return mpc

src/controller/nonlinmpc.jl

@@ -45,6 +45,7 @@ struct NonLinMPC{
     Uop::Vector{NT}
     Yop::Vector{NT}
     Dop::Vector{NT}
+    buffer::PredictiveControllerBuffer{NT}
     function NonLinMPC{NT, SE, JM, JEFunc}(
         estim::SE, Hp, Hc, M_Hp, N_Hc, L_Hp, Cwt, Ewt, JE::JEFunc, optim::JM
     ) where {NT<:Real, SE<:StateEstimator, JM<:JuMP.GenericModel, JEFunc<:Function}
@@ -75,6 +76,7 @@ struct NonLinMPC{
         Uop, Yop, Dop = repeat(model.uop, Hp), repeat(model.yop, Hp), repeat(model.dop, Hp)
         nΔŨ = size(Ẽ, 2)
         ΔŨ = zeros(NT, nΔŨ)
+        buffer = PredictiveControllerBuffer{NT}(nu, ny, nd, Hp)
         mpc = new{NT, SE, JM, JEFunc}(
             estim, optim, con,
             ΔŨ, ŷ,
@@ -87,6 +89,7 @@ struct NonLinMPC{
             Ks, Ps,
             d0, D̂0, D̂E,
             Uop, Yop, Dop,
+            buffer
         )
         init_optimization!(mpc, model, optim)
         return mpc

src/general.jl

@@ -44,6 +44,15 @@ end
 "Generate a block diagonal matrix repeating `n` times the matrix `A`."
 repeatdiag(A, n::Int) = kron(I(n), A)
 
+"In-place version of `repeat` but for vectors only."
+function repeat!(Y::Vector, a::Vector, n::Int)
+    na = length(a)
+    for i=0:n-1
+        Y[(1+na*i):(na*(i+1))] = a
+    end
+    return Y
+end
+
 "Convert 1-element vectors and normal matrices to Hermitians."
 to_hermitian(A::AbstractVector) = Hermitian(reshape(A, 1, 1), :L)
 to_hermitian(A::AbstractMatrix) = Hermitian(A, :L)

src/predictive_control.jl

@@ -19,6 +19,28 @@ julia> u = mpc([5]); round.(u, digits=3)
 """
 abstract type PredictiveController{NT<:Real} end
 
+struct PredictiveControllerBuffer{NT<:Real}
+    u ::Vector{NT}
+    R̂y::Vector{NT}
+    D̂ ::Vector{NT}
+    empty::Vector{NT}
+end
+
+@doc raw"""
+    PredictiveControllerBuffer{NT}(nu::Int, ny::Int, nd::Int)
+
+Create a buffer for `PredictiveController` objects.
+
+The buffer is used to store intermediate results during computation without allocating.
+"""
+function PredictiveControllerBuffer{NT}(nu::Int, ny::Int, nd::Int, Hp::Int) where NT <: Real
+    u  = Vector{NT}(undef, nu)
+    R̂y = Vector{NT}(undef, ny*Hp)
+    D̂  = Vector{NT}(undef, nd*Hp)
+    empty = Vector{NT}(undef, 0)
+    return PredictiveControllerBuffer{NT}(u, R̂y, D̂, empty)
+end
+
 include("controller/construct.jl")
 include("controller/execute.jl")
 include("controller/explicitmpc.jl")

src/sim_model.jl

@@ -29,7 +29,7 @@ struct SimModelBuffer{NT<:Real}
 end
 
 @doc raw"""
-    SimModelBuffer(nu::Int, nx::Int, ny::Int, nd::Int) -> SimModelBuffer{NT}
+    SimModelBuffer{NT}(nu::Int, nx::Int, ny::Int, nd::Int)
 
 Create a buffer for `SimModel` objects for inputs, states, outputs, and disturbances.
 

src/state_estim.jl

@@ -33,11 +33,11 @@ struct StateEstimatorBuffer{NT<:Real}
 end
 
 @doc raw"""
-    StateEstimatorBuffer(nx̂::Int, nym::Int) -> StateEstimatorBuffer{NT}
+    StateEstimatorBuffer{NT}(nu::Int, nx̂::Int, nym::Int, ny::Int, nd::Int)
 
 Create a buffer for `StateEstimator` objects for estimated states and measured outputs.
 
-The buffer is used to store intermediate results during simulation without allocating.
+The buffer is used to store intermediate results during estimation without allocating.
 """
 function StateEstimatorBuffer{NT}(
     nu::Int, nx̂::Int, nym::Int, ny::Int, nd::Int