Support online parameter estimation with an API (a.k.a. adaptive state estimation)

Supporting adaptive state estimation is quite straightforward on any state estimator that works with NonLinModel: you just augment the state vector with the estimated parameters $p_1$, $p_2$, etc. and augment the $\mathbf{f}$ equation with a time-invariant assumption $dp_1/dt = 0$, $dp_2/dt = 0$, etc.

Because of this, I first thought that adding an explicit API that would do this automatically is a very low priority feature. But, while writing this simple example of online parameter estimation with the MHE on the inverted pendulum:

using ModelPredictiveControl, Plots, JuMP
theme(:dark); default(fontfamily="Computer Modern"); scalefontsizes(1.1)
function f!(ẋ, x, u, _ , p)
    g, L, m = p                 # [m/s²], [m], [kg]
    θ, ω, K = x[1], x[2], x[3]  # [rad], [rad/s], [kg/s]
    τ  = u[1]                   # [Nm]
    ẋ[1] = ω
    ẋ[2] = -g/L*sin(θ) - K/m*ω + τ/m/L^2
    ẋ[3] = 0                    # dK/dt = 0 (the estimated parameter)
end
h!(y, x, _ , _ ) = (y[1] = 180/π*x[1]; nothing) # [°]
p = [9.8, 0.4, 0.3]
nu, nx, ny, Ts = 1, 3, 1, 0.1
vu, vx, vy = ["\$τ\$ (Nm)"], ["\$θ\$ (rad)", "\$ω\$ (rad/s)", "\$K\$ (kg/s)"], ["\$θ\$ (°)"]
model = setname!(NonLinModel(f!, h!, Ts, nu, nx, ny; p=p); u=vu, x=vx, y=vy)
σQ=[0.1, 1.0, 1.0]; σR=[5.0]; nint_ym=0; nint_u=0; 
transcription = OrthogonalCollocation(); hessian = true;
mhe = MovingHorizonEstimator(model; σQ, σR, nint_ym, nint_u, He=10, transcription, hessian)
setconstraint!(mhe, x̂min=[-Inf, -Inf, 1.0], x̂max=[+Inf, +Inf, 2.0], v̂min=[-2.5], v̂max=[2.5])
mhe |> display
N = 30; u = [1.10];
unset_time_limit_sec(mhe.optim)
res = sim!(mhe, N, u, x_0=[0, 0, 1.5], x̂_0=[0, 0, 1.2], y_noise=[0.5])
T = @elapsed sim!(mhe, N, u, x_0=[0, 0, 1.5], x̂_0=[0, 0, 1.2], y_noise=[0.5])
println("Sampling time: $Ts s, Average compute time: $(T/N) s")
plot(res, plotu=false, plotxwithx̂=true) |> display

I noticed that the states associated to the estimated parameter are in fact LTI models (the 0 in the RHS, that is, an integrator). It means that we could specialized the equality constraint as linear one on these states, for non-SingleShooting transcription. This is very similar to what is already done with the generic integrating states from nint_ym and nint_u options. The goal is to reduces the dimensions of the Jacobian thus more efficient AD. So this feature is more interesting to implement than I initially thought.

Start by tracing how NonLinModel states and MovingHorizonEstimator build equality constraints for OrthogonalCollocation, then compare the existing nint_ym and nint_u handling. Check how SingleShooting is excluded and identify where an API for parameter states would fit. Done means parameter estimation can be configured explicitly and non-SingleShooting Jacobian dimensions are reduced without changing the illustrated behavior.