From 8ec6544eb1b41b08765ec63bb6a2d1d4854f65ad Mon Sep 17 00:00:00 2001 From: franckgaga Date: Fri, 31 Jul 2026 01:38:56 -0400 Subject: [PATCH 1/3] debug: `AutoSparse` hessian with `MovingHorizonEstimator` now works --- src/estimator/mhe/execute.jl | 17 ++++++++++++++--- 1 file changed, 14 insertions(+), 3 deletions(-) diff --git a/src/estimator/mhe/execute.jl b/src/estimator/mhe/execute.jl index 7d2fbcb80..da99b775c 100644 --- a/src/estimator/mhe/execute.jl +++ b/src/estimator/mhe/execute.jl @@ -810,9 +810,6 @@ function obj_nonlinprog(estim::MovingHorizonEstimator, ::SimModel, x̄, V̂, Ŵ nŴ, nYm = Nk*estim.nx̂, Nk*estim.nym Ŵ, V̂ = Ŵ[1:nŴ], V̂[1:nYm] end - if any(isnan, V̂) # ignore NaN values in V̂ for the objective function: - V̂ = [isnan(v) ? 0 : v for v in V̂] - end Jε = estim.nε > 0 ? estim.C*Z̃[begin]^2 : 0 return dot(x̄, invP̄, x̄) + dot(Ŵ, invQ̂_Nk, Ŵ) + dot(V̂, invR̂_Nk, V̂) + Jε end @@ -889,6 +886,13 @@ function predict_mhe!( d0next = @views estim.D0[(1 + nd*j):(nd*(j+1))] ĥ!(ŷ0next, estim, model, x̂0next, d0next) ŷ0nextm = @views ŷ0next[estim.i_ym] + if any(isnan, y0nextm) + for i in eachindex(y0nextm) + if isnan(y0nextm[i]) + y0nextm[i] = ŷ0nextm[i] + end + end + end V̂[(1 + nym*(j-1)):(nym*j)] .= y0nextm .- ŷ0nextm x̂0, d0 = x̂0next, d0next end @@ -900,6 +904,13 @@ function predict_mhe!( ŵ = @views Ŵ[(1 + nŵ*(j-1)):(nŵ*j)] ĥ!(ŷ0, estim, model, x̂0, d0) ŷ0m = @views ŷ0[estim.i_ym] + if any(isnan, y0m) + for i in eachindex(y0m) + if isnan(y0m[i]) + y0m[i] = ŷ0m[i] + end + end + end V̂[(1 + nym*(j-1)):(nym*j)] .= y0m .- ŷ0m x̂0next = @views X̂0[(1 + nx̂ *(j-1)):(nx̂ *j)] f̂!(x̂0next, û0, k, estim, model, x̂0, u0, d0) From 0ed29f8aae20a91065a5f092abba67bf33cbc98c Mon Sep 17 00:00:00 2001 From: franckgaga Date: Fri, 31 Jul 2026 01:39:36 -0400 Subject: [PATCH 2/3] bump --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 12249dea9..ab3ed796b 100644 --- a/Project.toml +++ b/Project.toml @@ -1,6 +1,6 @@ name = "ModelPredictiveControl" uuid = "61f9bdb8-6ae4-484a-811f-bbf86720c31c" -version = "2.8.2" +version = "2.8.3" authors = ["Francis Gagnon"] [deps] From 21a11f354fe9fbbcd3408436dec644c9e2a0aba1 Mon Sep 17 00:00:00 2001 From: franckgaga Date: Fri, 31 Jul 2026 02:21:57 -0400 Subject: [PATCH 3/3] test: sparse `hessian` backend for MHE --- src/estimator/mhe/execute.jl | 12 ++---------- test/2_test_state_estim.jl | 12 +++++++++--- 2 files changed, 11 insertions(+), 13 deletions(-) diff --git a/src/estimator/mhe/execute.jl b/src/estimator/mhe/execute.jl index da99b775c..dbb54ddad 100644 --- a/src/estimator/mhe/execute.jl +++ b/src/estimator/mhe/execute.jl @@ -887,11 +887,7 @@ function predict_mhe!( ĥ!(ŷ0next, estim, model, x̂0next, d0next) ŷ0nextm = @views ŷ0next[estim.i_ym] if any(isnan, y0nextm) - for i in eachindex(y0nextm) - if isnan(y0nextm[i]) - y0nextm[i] = ŷ0nextm[i] - end - end + y0nextm = [isnan(y) ? ŷ : y for (y, ŷ) in zip(y0nextm, ŷ0nextm)] end V̂[(1 + nym*(j-1)):(nym*j)] .= y0nextm .- ŷ0nextm x̂0, d0 = x̂0next, d0next @@ -905,11 +901,7 @@ function predict_mhe!( ĥ!(ŷ0, estim, model, x̂0, d0) ŷ0m = @views ŷ0[estim.i_ym] if any(isnan, y0m) - for i in eachindex(y0m) - if isnan(y0m[i]) - y0m[i] = ŷ0m[i] - end - end + y0m = [isnan(y) ? ŷ : y for (y, ŷ) in zip(y0m, ŷ0m)] end V̂[(1 + nym*(j-1)):(nym*j)] .= y0m .- ŷ0m x̂0next = @views X̂0[(1 + nx̂ *(j-1)):(nx̂ *j)] diff --git a/test/2_test_state_estim.jl b/test/2_test_state_estim.jl index 9b8912ffa..2b870a5a4 100644 --- a/test/2_test_state_estim.jl +++ b/test/2_test_state_estim.jl @@ -1113,7 +1113,8 @@ end @testitem "MHE estimation and getinfo (NonLinModel)" setup=[SetupMPCtests] begin using .SetupMPCtests, ControlSystemsBase, LinearAlgebra, ForwardDiff - using JuMP, Ipopt + using JuMP, Ipopt, DifferentiationInterface, SparseMatrixColorings, SparseConnectivityTracer + import ForwardDiff linmodel = LinModel(sys,Ts,i_u=[1,2], i_d=[3]) linmodel = setop!(linmodel, uop=[10,50], yop=[50,30], dop=[5]) f(x,u,d,model) = model.A*x + model.Bu*u + model.Bd*d @@ -1121,7 +1122,6 @@ end nonlinmodel = NonLinModel(f, h, Ts, 2, 4, 2, 1, solver=nothing, p=linmodel) nonlinmodel = setop!(nonlinmodel, uop=[10,50], yop=[50,30], dop=[5]) - mhe1 = MovingHorizonEstimator(nonlinmodel, He=2) JuMP.set_attribute(mhe1.optim, "tol", 1e-7) preparestate!(mhe1, [50, 30], [5]) @@ -1218,7 +1218,13 @@ end @test_nowarn ModelPredictiveControl.info2debugstr(info) @test_throws ErrorException setstate!(mhe1, [1,2,3,4,5,6], diagm(.1:.1:.6)) - mhe7 = MovingHorizonEstimator(nonlinmodel, He=2) + hessian = AutoSparse( + AutoForwardDiff(); + sparsity_detector=TracerSparsityDetector(), + coloring_algorithm=GreedyColoringAlgorithm(), + ) + + mhe7 = MovingHorizonEstimator(nonlinmodel; He=2, hessian) @test_logs( (:warn, "NaN values in the MHE measurements ym: ignoring them in the objective"), preparestate!(mhe7, [50, NaN], [5])