From a2a063169176a817b4f68790e49f115bf9799bd0 Mon Sep 17 00:00:00 2001 From: franckgaga Date: Wed, 4 Mar 2026 12:52:54 -0500 Subject: [PATCH 1/9] doc: minor addition --- src/controller/transcription.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/controller/transcription.jl b/src/controller/transcription.jl index 5e3754bca..e273d7e1f 100644 --- a/src/controller/transcription.jl +++ b/src/controller/transcription.jl @@ -123,7 +123,7 @@ end @doc raw""" OrthogonalCollocation( - h::Int=0, no=3; f_threads=false, h_threads=false, roots=:gaussradau + h::Int=0, no::Int=3; f_threads=false, h_threads=false, roots=:gaussradau ) Construct an orthogonal collocation on finite elements [`TranscriptionMethod`](@ref). From 78e5f5f367426b36dda0d5a4032bb7416e1b8e87 Mon Sep 17 00:00:00 2001 From: franckgaga Date: Wed, 4 Mar 2026 13:06:52 -0500 Subject: [PATCH 2/9] added: precompile `OrthogonalCollocation` I expect it to be more popular than `TrapezoidalCollocation`, so let's replace this precompile workload. --- src/precompile.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/precompile.jl b/src/precompile.jl index 2408a2017..ff8e64ae2 100644 --- a/src/precompile.jl +++ b/src/precompile.jl @@ -101,7 +101,7 @@ R̂y = repeat([55; 30], 10) preparestate!(nmpc_ekf, [55, 30]) nmpc_ekf([55, 30]) - transcription = TrapezoidalCollocation() + transcription = OrthogonalCollocation() nmpc_mhe = setconstraint!(NonLinMPC( MovingHorizonEstimator(nlmodel, He=2); transcription, Hp=10, Cwt=Inf), ymin=[45, -Inf] ) From 827e04a44bebf9f9302d335fa299ed64c08af00e Mon Sep 17 00:00:00 2001 From: franckgaga Date: Wed, 4 Mar 2026 14:21:23 -0500 Subject: [PATCH 3/9] bench: `OrthogonalCollocation` unit tests --- benchmark/3_bench_predictive_control.jl | 23 +++++++++++++++++++++++ 1 file changed, 23 insertions(+) diff --git a/benchmark/3_bench_predictive_control.jl b/benchmark/3_bench_predictive_control.jl index 19e8c386d..635fbb01e 100644 --- a/benchmark/3_bench_predictive_control.jl +++ b/benchmark/3_bench_predictive_control.jl @@ -68,6 +68,17 @@ nmpc_nonlin_tc = NonLinMPC( nonlinmodel_c, transcription=TrapezoidalCollocation(), Mwt=[1], Nwt=[0.1], Lwt=[0.1], Hp=10 ) +JuMP.set_attribute(nmpc_nonlin_tc.optim, "tol", 1e-7) +nmpc_nonlin_oc = NonLinMPC( + nonlinmodel_c, transcription=OrthogonalCollocation(), + Mwt=[1], Nwt=[0.1], Lwt=[0.1], Hp=10 +) +JuMP.set_attribute(nmpc_nonlin_oc.optim, "tol", 1e-7) +nmpc_nonlin_oc_hess = NonLinMPC( + nonlinmodel_c, transcription=OrthogonalCollocation(), hessian=true, + Mwt=[1], Nwt=[0.1], Lwt=[0.1], Hp=10 +) +JuMP.set_attribute(nmpc_nonlin_oc_hess.optim, "tol", 1e-7) samples, evals, seconds = 10000, 1, 60 UNIT_MPC["NonLinMPC"]["moveinput!"]["LinModel"]["SingleShooting"] = @@ -118,6 +129,18 @@ UNIT_MPC["NonLinMPC"]["moveinput!"]["NonLinModel"]["TrapezoidalCollocation"] = setup=preparestate!($nmpc_nonlin_tc, $y_c, $d_c), samples=samples, evals=evals, seconds=seconds ) +UNIT_MPC["NonLinMPC"]["moveinput!"]["NonLinModel"]["OrthogonalCollocation"] = + @benchmarkable( + moveinput!($nmpc_nonlin_oc, $y_c, $d_c), + setup=preparestate!($nmpc_nonlin_oc, $y_c, $d_c), + samples=samples, evals=evals, seconds=seconds + ) +UNIT_MPC["NonLinMPC"]["moveinput!"]["NonLinModel"]["OrthogonalCollocationHessian"] = + @benchmarkable( + moveinput!($nmpc_nonlin_oc_hess, $y_c, $d_c), + setup=preparestate!($nmpc_nonlin_oc_hess, $y_c, $d_c), + samples=samples, evals=evals, seconds=seconds + ) UNIT_MPC["NonLinMPC"]["getinfo!"]["NonLinModel"] = @benchmarkable( getinfo($nmpc_nonlin_ss), From b11d6245d42fc351bb3862703e3594b0d65a0510 Mon Sep 17 00:00:00 2001 From: franckgaga Date: Wed, 4 Mar 2026 14:22:24 -0500 Subject: [PATCH 4/9] bench: increase tolerance in MPC unit tests --- benchmark/3_bench_predictive_control.jl | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/benchmark/3_bench_predictive_control.jl b/benchmark/3_bench_predictive_control.jl index 635fbb01e..4027a14e0 100644 --- a/benchmark/3_bench_predictive_control.jl +++ b/benchmark/3_bench_predictive_control.jl @@ -44,26 +44,32 @@ nmpc_lin_ss = NonLinMPC( linmodel, transcription=SingleShooting(), Mwt=[1, 1], Nwt=[0.1, 0.1], Lwt=[0.1, 0.1], Hp=10 ) +JuMP.set_attribute(nmpc_lin_ss.optim, "tol", 1e-7) nmpc_lin_ms = NonLinMPC( linmodel, transcription=MultipleShooting(), Mwt=[1, 1], Nwt=[0.1, 0.1], Lwt=[0.1, 0.1], Hp=10 ) +JuMP.set_attribute(nmpc_lin_ms.optim, "tol", 1e-7) nmpc_nonlin_ss = NonLinMPC( nonlinmodel, transcription=SingleShooting(), Mwt=[1, 1], Nwt=[0.1, 0.1], Lwt=[0.1, 0.1], Hp=10 ) +JuMP.set_attribute(nmpc_nonlin_ss.optim, "tol", 1e-7) nmpc_nonlin_ss_hess = NonLinMPC( nonlinmodel, transcription=SingleShooting(), hessian=true, Mwt=[1, 1], Nwt=[0.1, 0.1], Lwt=[0.1, 0.1], Hp=10 ) +JuMP.set_attribute(nmpc_nonlin_ss_hess.optim, "tol", 1e-7) nmpc_nonlin_ms = NonLinMPC( nonlinmodel, transcription=MultipleShooting(), Mwt=[1, 1], Nwt=[0.1, 0.1], Lwt=[0.1, 0.1], Hp=10 ) +JuMP.set_attribute(nmpc_nonlin_ms.optim, "tol", 1e-7) nmpc_nonlin_ms_hess = NonLinMPC( nonlinmodel, transcription=MultipleShooting(), hessian=true, Mwt=[1, 1], Nwt=[0.1, 0.1], Lwt=[0.1, 0.1], Hp=10 ) +JuMP.set_attribute(nmpc_nonlin_ms_hess.optim, "tol", 1e-7) nmpc_nonlin_tc = NonLinMPC( nonlinmodel_c, transcription=TrapezoidalCollocation(), Mwt=[1], Nwt=[0.1], Lwt=[0.1], Hp=10 From 68825c5dabd042dd0c19bf6768db7ce12ad136dd Mon Sep 17 00:00:00 2001 From: franckgaga Date: Wed, 4 Mar 2026 14:54:02 -0500 Subject: [PATCH 5/9] bench: MadNLP with `MultipleShooting` and Hessian --- benchmark/3_bench_predictive_control.jl | 25 ++++++++++++------------- 1 file changed, 12 insertions(+), 13 deletions(-) diff --git a/benchmark/3_bench_predictive_control.jl b/benchmark/3_bench_predictive_control.jl index 4027a14e0..3d89cd6b0 100644 --- a/benchmark/3_bench_predictive_control.jl +++ b/benchmark/3_bench_predictive_control.jl @@ -379,19 +379,13 @@ nmpc_uno_ms_hess = NonLinMPC(estim; Hp, Hc, Mwt, Nwt, Cwt, optim, transcription, nmpc_uno_ms_hess = setconstraint!(nmpc_uno_ms_hess; umin, umax) JuMP.unset_time_limit_sec(nmpc_uno_ms_hess.optim) -# TODO: does not work well with MadNLP and MultipleShooting or TrapezoidalCollocation, -# figure out why. Current theory: -# MadNLP LBFGS approximation is less robust than Ipopt version. Re-test when exact Hessians -# will be supported in ModelPredictiveControl.jl. The following attributes kinda work with -# the MadNLP LBFGS approximation but super slow (~1000 times slower than Ipopt): -# optim = JuMP.Model(MadNLP.Optimizer) -# transcription = MultipleShooting() -# nmpc_madnlp_ms = NonLinMPC(estim; Hp, Hc, Mwt, Nwt, Cwt, optim, transcription) -# nmpc_madnlp_ms = setconstraint!(nmpc_madnlp_ms; umin, umax) -# JuMP.unset_time_limit_sec(nmpc_madnlp_ms.optim) -# JuMP.set_attribute(nmpc_madnlp_ms.optim, "hessian_approximation", MadNLP.CompactLBFGS) -# MadNLP_QNopt = MadNLP.QuasiNewtonOptions(; max_history=42) -# JuMP.set_attribute(nmpc_madnlp_ms.optim, "quasi_newton_options", MadNLP_QNopt) +# skip MadNLP.jl with MultipleShooting and hessian=false, their LBFGS does not work well + +optim = JuMP.Model(MadNLP.Optimizer) +transcription, hessian = MultipleShooting(), true +nmpc_madnlp_ms_hess = NonLinMPC(estim; Hp, Hc, Mwt, Nwt, Cwt, optim, transcription, hessian) +nmpc_madnlp_ms_hess = setconstraint!(nmpc_madnlp_ms_hess; umin, umax) +JuMP.unset_time_limit_sec(nmpc_madnlp_ms_hess.optim) samples, evals, seconds = 100, 1, 15*60 CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["Ipopt"]["SingleShooting"] = @@ -439,6 +433,11 @@ CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["MadNLP"]["SingleShooting"] = sim!($nmpc_madnlp_ss, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), samples=samples, evals=evals, seconds=seconds ) +CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["MadNLP"]["MultipleShooting (Hessian)"] = + @benchmarkable( + sim!($nmpc_madnlp_ms_hess, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), + samples=samples, evals=evals, seconds=seconds + ) # CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["Uno"]["MultipleShooting (Hessian)"] = # @benchmarkable( # sim!($nmpc_uno_ms_hess, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), From 8eefb3cf60c4c94384dad93f61bb932d9fd7b35e Mon Sep 17 00:00:00 2001 From: franckgaga Date: Wed, 4 Mar 2026 15:01:47 -0500 Subject: [PATCH 6/9] bench: re-enabling `UnoSolver.jl` --- benchmark/3_bench_predictive_control.jl | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/benchmark/3_bench_predictive_control.jl b/benchmark/3_bench_predictive_control.jl index 3d89cd6b0..c3efeb808 100644 --- a/benchmark/3_bench_predictive_control.jl +++ b/benchmark/3_bench_predictive_control.jl @@ -438,11 +438,11 @@ CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["MadNLP"]["MultipleShooting (He sim!($nmpc_madnlp_ms_hess, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), samples=samples, evals=evals, seconds=seconds ) -# CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["Uno"]["MultipleShooting (Hessian)"] = -# @benchmarkable( -# sim!($nmpc_uno_ms_hess, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), -# samples=samples, evals=evals, seconds=seconds -# ) +CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["Uno"]["MultipleShooting (Hessian)"] = + @benchmarkable( + sim!($nmpc_uno_ms_hess, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), + samples=samples, evals=evals, seconds=seconds + ) # ----------------- Case study: Pendulum economic -------------------------------- model2, p = pendulum_model2, pendulum_p2 From 28139ef104042bb6be2f95e5711a164c51db013d Mon Sep 17 00:00:00 2001 From: franckgaga Date: Wed, 4 Mar 2026 16:22:59 -0500 Subject: [PATCH 7/9] bench: `OrthogonalCollocation` on pendulum Also called `GC.gc()` in `setup` function to try to avoid OOM. --- benchmark/3_bench_predictive_control.jl | 138 ++++++++++++++++-------- 1 file changed, 94 insertions(+), 44 deletions(-) diff --git a/benchmark/3_bench_predictive_control.jl b/benchmark/3_bench_predictive_control.jl index c3efeb808..493a4a3d8 100644 --- a/benchmark/3_bench_predictive_control.jl +++ b/benchmark/3_bench_predictive_control.jl @@ -268,12 +268,8 @@ optim = JuMP.Model(DAQP.Optimizer, add_bridges=false) transcription = SingleShooting() mpc_d_daqp_ss = setconstraint!(LinMPC(model_d; optim, transcription), ymin=[45, -Inf]) -# # Skip DAQP with MultipleShooting, it is not designed for sparse Hessians. Kind of works -# # with "eps_prox" configured to 1e-6, but not worth it. -# optim = JuMP.Model(DAQP.Optimizer, add_bridges=false) -# transcription = MultipleShooting() -# mpc_d_daqp_ms = setconstraint!(LinMPC(model_d; optim, transcription), ymin=[45, -Inf]) -# JuMP.set_attribute(mpc_d_daqp_ms.optim, "eps_prox", 1e-6) +# Skip DAQP with MultipleShooting, it is not designed for sparse Hessians. Kind of works +# with "eps_prox" configured to 1e-6, but not worth it. optim = JuMP.Model(optimizer_with_attributes(Ipopt.Optimizer,"sb"=>"yes"), add_bridges=false) transcription = SingleShooting() @@ -367,11 +363,25 @@ nmpc_ipopt_tct = NonLinMPC(estim; Hp, Hc, Mwt, Nwt, Cwt, optim, transcription) nmpc_ipopt_tct = setconstraint!(nmpc_ipopt_tct; umin, umax) JuMP.unset_time_limit_sec(nmpc_ipopt_tct.optim) +optim = JuMP.Model(optimizer_with_attributes(Ipopt.Optimizer,"sb"=>"yes"), add_bridges=false) +transcription = OrthogonalCollocation() +nmpc_ipopt_oc = NonLinMPC(estim; Hp, Hc, Mwt, Nwt, Cwt, optim, transcription) +nmpc_ipopt_oc = setconstraint!(nmpc_ipopt_oc; umin, umax) +JuMP.unset_time_limit_sec(nmpc_ipopt_oc.optim) + +optim = JuMP.Model(optimizer_with_attributes(Ipopt.Optimizer,"sb"=>"yes"), add_bridges=false) +transcription, hessian = OrthogonalCollocation(), true +nmpc_ipopt_oc_hess = NonLinMPC(estim; Hp, Hc, Mwt, Nwt, Cwt, optim, transcription, hessian) +nmpc_ipopt_oc_hess = setconstraint!(nmpc_ipopt_oc_hess; umin, umax) +JuMP.unset_time_limit_sec(nmpc_ipopt_oc_hess.optim) + +# skip MadNLP.jl with hessian=false, their LBFGS does not work well + optim = JuMP.Model(MadNLP.Optimizer, add_bridges=false) -transcription = SingleShooting() -nmpc_madnlp_ss = NonLinMPC(estim; Hp, Hc, Mwt, Nwt, Cwt, optim, transcription) -nmpc_madnlp_ss = setconstraint!(nmpc_madnlp_ss; umin, umax) -JuMP.unset_time_limit_sec(nmpc_madnlp_ss.optim) +transcription, hessian = SingleShooting(), true +nmpc_madnlp_ss_hess = NonLinMPC(estim; Hp, Hc, Mwt, Nwt, Cwt, optim, transcription, hessian) +nmpc_madnlp_ss_hess = setconstraint!(nmpc_madnlp_ss_hess; umin, umax) +JuMP.unset_time_limit_sec(nmpc_madnlp_ss_hess.optim) optim = JuMP.Model(()->UnoSolver.Optimizer(preset="filtersqp"), add_bridges=false) transcription, hessian = MultipleShooting(), true @@ -379,8 +389,6 @@ nmpc_uno_ms_hess = NonLinMPC(estim; Hp, Hc, Mwt, Nwt, Cwt, optim, transcription, nmpc_uno_ms_hess = setconstraint!(nmpc_uno_ms_hess; umin, umax) JuMP.unset_time_limit_sec(nmpc_uno_ms_hess.optim) -# skip MadNLP.jl with MultipleShooting and hessian=false, their LBFGS does not work well - optim = JuMP.Model(MadNLP.Optimizer) transcription, hessian = MultipleShooting(), true nmpc_madnlp_ms_hess = NonLinMPC(estim; Hp, Hc, Mwt, Nwt, Cwt, optim, transcription, hessian) @@ -391,58 +399,69 @@ samples, evals, seconds = 100, 1, 15*60 CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["Ipopt"]["SingleShooting"] = @benchmarkable( sim!($nmpc_ipopt_ss, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["Ipopt"]["SingleShooting (Hessian)"] = @benchmarkable( sim!($nmpc_ipopt_ss_hess, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["Ipopt"]["MultipleShooting"] = @benchmarkable( sim!($nmpc_ipopt_ms, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["Ipopt"]["MultipleShooting (Hessian)"] = @benchmarkable( sim!($nmpc_ipopt_ms_hess, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["Ipopt"]["MultipleShooting (threaded)"] = @benchmarkable( sim!($nmpc_ipopt_mst, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["Ipopt"]["TrapezoidalCollocation"] = @benchmarkable( sim!($nmpc_ipopt_tc, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["Ipopt"]["TrapezoidalCollocation (Hessian)"] = @benchmarkable( sim!($nmpc_ipopt_tc_hess, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["Ipopt"]["TrapezoidalCollocation (threaded)"] = @benchmarkable( sim!($nmpc_ipopt_tct, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) -CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["MadNLP"]["SingleShooting"] = +CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["Ipopt"]["OrthogonalCollocation"] = + @benchmarkable( + sim!($nmpc_ipopt_oc, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() + ) +CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["Ipopt"]["OrthogonalCollocation (Hessian)"] = + @benchmarkable( + sim!($nmpc_ipopt_oc_hess, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() + ) +CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["MadNLP"]["SingleShooting (Hessian)"] = @benchmarkable( sim!($nmpc_madnlp_ss, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["MadNLP"]["MultipleShooting (Hessian)"] = @benchmarkable( sim!($nmpc_madnlp_ms_hess, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["Uno"]["MultipleShooting (Hessian)"] = @benchmarkable( sim!($nmpc_uno_ms_hess, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds - ) + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() + ) + # ----------------- Case study: Pendulum economic -------------------------------- model2, p = pendulum_model2, pendulum_p2 @@ -493,49 +512,80 @@ empc_ipopt_tc_hess = NonLinMPC(estim2; Hp, Hc, Nwt, Mwt=Mwt2, Cwt, JE, Ewt, opti empc_ipopt_tc_hess = setconstraint!(empc_ipopt_tc_hess; umin, umax) JuMP.unset_time_limit_sec(empc_ipopt_tc_hess.optim) +optim = JuMP.Model(optimizer_with_attributes(Ipopt.Optimizer,"sb"=>"yes"), add_bridges=false) +transcription = OrthogonalCollocation() +empc_ipopt_oc = NonLinMPC(estim2; Hp, Hc, Nwt, Mwt=Mwt2, Cwt, JE, Ewt, optim, transcription, p) +empc_ipopt_oc = setconstraint!(empc_ipopt_oc; umin, umax) +JuMP.unset_time_limit_sec(empc_ipopt_oc.optim) + +optim = JuMP.Model(optimizer_with_attributes(Ipopt.Optimizer,"sb"=>"yes"), add_bridges=false) +transcription, hessian = OrthogonalCollocation(), true +empc_ipopt_oc_hess = NonLinMPC(estim2; Hp, Hc, Nwt, Mwt=Mwt2, Cwt, JE, Ewt, optim, transcription, p) +empc_ipopt_oc_hess = setconstraint!(empc_ipopt_oc_hess; umin, umax) +JuMP.unset_time_limit_sec(empc_ipopt_oc_hess.optim) + optim = JuMP.Model(MadNLP.Optimizer, add_bridges=false) -transcription = SingleShooting() -empc_madnlp_ss = NonLinMPC(estim2; Hp, Hc, Nwt, Mwt=Mwt2, Cwt, JE, Ewt, optim, transcription, p) -empc_madnlp_ss = setconstraint!(empc_madnlp_ss; umin, umax) -JuMP.unset_time_limit_sec(empc_madnlp_ss.optim) +transcription, hessian = SingleShooting(), true +empc_madnlp_ss_hess = NonLinMPC(estim2; Hp, Hc, Nwt, Mwt=Mwt2, Cwt, JE, Ewt, optim, transcription, hessian, p) +empc_madnlp_ss_hess = setconstraint!(empc_madnlp_ss_hess; umin, umax) +JuMP.unset_time_limit_sec(empc_madnlp_ss_hess.optim) -# TODO: test EMPC with MadNLP and MultipleShooting and TrapezoidalCollocation, see comment above. +optim = JuMP.Model(()->UnoSolver.Optimizer(preset="filtersqp"), add_bridges=false) +transcription, hessian = MultipleShooting(), true +empc_uno_ss_hess = NonLinMPC(estim2; Hp, Hc, Nwt, Mwt=Mwt2, Cwt, JE, Ewt, optim, transcription, hessian, p) +empc_uno_ss_hess = setconstraint!(empc_uno_ss_hess; umin, umax) +JuMP.unset_time_limit_sec(empc_uno_ss_hess.optim) samples, evals, seconds = 100, 1, 15*60 CASE_MPC["Pendulum"]["NonLinMPC"]["Economic"]["Ipopt"]["SingleShooting"] = @benchmarkable( sim!($empc_ipopt_ss, $N, $ry; plant=$plant2, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_MPC["Pendulum"]["NonLinMPC"]["Economic"]["Ipopt"]["SingleShooting (Hessian)"] = @benchmarkable( sim!($empc_ipopt_ss_hess, $N, $ry; plant=$plant2, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_MPC["Pendulum"]["NonLinMPC"]["Economic"]["Ipopt"]["MultipleShooting"] = @benchmarkable( sim!($empc_ipopt_ms, $N, $ry; plant=$plant2, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_MPC["Pendulum"]["NonLinMPC"]["Economic"]["Ipopt"]["MultipleShooting (Hessian)"] = @benchmarkable( sim!($empc_ipopt_ms_hess, $N, $ry; plant=$plant2, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_MPC["Pendulum"]["NonLinMPC"]["Economic"]["Ipopt"]["TrapezoidalCollocation"] = @benchmarkable( sim!($empc_ipopt_tc, $N, $ry; plant=$plant2, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_MPC["Pendulum"]["NonLinMPC"]["Economic"]["Ipopt"]["TrapezoidalCollocation (Hessian)"] = @benchmarkable( sim!($empc_ipopt_tc_hess, $N, $ry; plant=$plant2, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) -CASE_MPC["Pendulum"]["NonLinMPC"]["Economic"]["MadNLP"]["SingleShooting"] = +CASE_MPC["Pendulum"]["NonLinMPC"]["Economic"]["Ipopt"]["OrthogonalCollocation"] = @benchmarkable( - sim!($empc_madnlp_ss, $N, $ry; plant=$plant2, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + sim!($empc_ipopt_oc, $N, $ry; plant=$plant2, x_0=$x_0, x̂_0=$x̂_0, progress=false), + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() + ) +CASE_MPC["Pendulum"]["NonLinMPC"]["Economic"]["Ipopt"]["OrthogonalCollocation (Hessian)"] = + @benchmarkable( + sim!($empc_ipopt_oc_hess, $N, $ry; plant=$plant2, x_0=$x_0, x̂_0=$x̂_0, progress=false), + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() + ) +CASE_MPC["Pendulum"]["NonLinMPC"]["Economic"]["MadNLP"]["SingleShooting (Hessian)"] = + @benchmarkable( + sim!($empc_madnlp_ss_hess, $N, $ry; plant=$plant2, x_0=$x_0, x̂_0=$x̂_0, progress=false), + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() + ) +CASE_MPC["Pendulum"]["NonLinMPC"]["Economic"]["Uno"]["MultipleShooting (Hessian)"] = + @benchmarkable( + sim!($empc_uno_ss_hess, $N, $ry; plant=$plant2, x_0=$x_0, x̂_0=$x̂_0, progress=false), + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) # -------------- Case study: Pendulum custom constraints -------------------------- @@ -598,27 +648,27 @@ samples, evals, seconds = 100, 1, 15*60 CASE_MPC["Pendulum"]["NonLinMPC"]["Custom constraints"]["Ipopt"]["SingleShooting"] = @benchmarkable( sim!($nmpc2_ipopt_ss, $N, $ry; plant=$plant2, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_MPC["Pendulum"]["NonLinMPC"]["Custom constraints"]["Ipopt"]["MultipleShooting"] = @benchmarkable( sim!($nmpc2_ipopt_ms, $N, $ry; plant=$plant2, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_MPC["Pendulum"]["NonLinMPC"]["Custom constraints"]["Ipopt"]["MultipleShooting (Hessian)"] = @benchmarkable( sim!($nmpc2_ipopt_ms_hess, $N, $ry; plant=$plant2, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_MPC["Pendulum"]["NonLinMPC"]["Custom constraints"]["Ipopt"]["TrapezoidalCollocation"] = @benchmarkable( sim!($nmpc2_ipopt_tc, $N, $ry; plant=$plant2, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_MPC["Pendulum"]["NonLinMPC"]["Custom constraints"]["Ipopt"]["TrapezoidalCollocation (Hessian)"] = @benchmarkable( sim!($nmpc2_ipopt_tc_hess, $N, $ry; plant=$plant2, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) # ----------------- Case study: Pendulum successive linearization ------------------------- From f1a937234a2b762fae3165830f8db64cf0ff46c4 Mon Sep 17 00:00:00 2001 From: franckgaga Date: Wed, 4 Mar 2026 16:28:11 -0500 Subject: [PATCH 8/9] bench: debug last commit --- benchmark/3_bench_predictive_control.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/benchmark/3_bench_predictive_control.jl b/benchmark/3_bench_predictive_control.jl index 493a4a3d8..f2250ee88 100644 --- a/benchmark/3_bench_predictive_control.jl +++ b/benchmark/3_bench_predictive_control.jl @@ -448,7 +448,7 @@ CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["Ipopt"]["OrthogonalCollocation ) CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["MadNLP"]["SingleShooting (Hessian)"] = @benchmarkable( - sim!($nmpc_madnlp_ss, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), + sim!($nmpc_madnlp_ss_hess, $N, $ry; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_MPC["Pendulum"]["NonLinMPC"]["Noneconomic"]["MadNLP"]["MultipleShooting (Hessian)"] = From 532110e49be121b9b31a955bb9a0ce7de2224ed7 Mon Sep 17 00:00:00 2001 From: franckgaga Date: Wed, 4 Mar 2026 16:40:16 -0500 Subject: [PATCH 9/9] bench: MHE with `MadNLP` with Hessians --- benchmark/2_bench_state_estim.jl | 42 +++++++++++++++++--------------- 1 file changed, 22 insertions(+), 20 deletions(-) diff --git a/benchmark/2_bench_state_estim.jl b/benchmark/2_bench_state_estim.jl index e52e9d15b..982aeeaa2 100644 --- a/benchmark/2_bench_state_estim.jl +++ b/benchmark/2_bench_state_estim.jl @@ -343,50 +343,52 @@ JuMP.set_attribute(mhe_pendulum_ipopt_predh.optim, "tol", 1e-7) optim = JuMP.Model(MadNLP.Optimizer, add_bridges=false) direct = true -mhe_pendulum_madnlp_curr = MovingHorizonEstimator( - model; He, σQ, σR, nint_u, σQint_u, optim, direct +hessian = true +mhe_pendulum_madnlp_currh = MovingHorizonEstimator( + model; He, σQ, σR, nint_u, σQint_u, optim, direct, hessian ) -mhe_pendulum_madnlp_curr = setconstraint!(mhe_pendulum_madnlp_curr; v̂min, v̂max) -JuMP.unset_time_limit_sec(mhe_pendulum_madnlp_curr.optim) -JuMP.set_attribute(mhe_pendulum_madnlp_curr.optim, "tol", 1e-7) +mhe_pendulum_madnlp_currh = setconstraint!(mhe_pendulum_madnlp_currh; v̂min, v̂max) +JuMP.unset_time_limit_sec(mhe_pendulum_madnlp_currh.optim) +JuMP.set_attribute(mhe_pendulum_madnlp_currh.optim, "tol", 1e-7) optim = JuMP.Model(MadNLP.Optimizer, add_bridges=false) direct = false -mhe_pendulum_madnlp_pred = MovingHorizonEstimator( - model; He, σQ, σR, nint_u, σQint_u, optim, direct +hessian = true +mhe_pendulum_madnlp_predh = MovingHorizonEstimator( + model; He, σQ, σR, nint_u, σQint_u, optim, direct, hessian ) -mhe_pendulum_madnlp_pred = setconstraint!(mhe_pendulum_madnlp_pred; v̂min, v̂max) -JuMP.unset_time_limit_sec(mhe_pendulum_madnlp_pred.optim) -JuMP.set_attribute(mhe_pendulum_madnlp_pred.optim, "tol", 1e-7) +mhe_pendulum_madnlp_pred = setconstraint!(mhe_pendulum_madnlp_predh; v̂min, v̂max) +JuMP.unset_time_limit_sec(mhe_pendulum_madnlp_predh.optim) +JuMP.set_attribute(mhe_pendulum_madnlp_predh.optim, "tol", 1e-7) samples, evals, seconds = 25, 1, 15*60 CASE_ESTIM["Pendulum"]["MovingHorizonEstimator"]["Ipopt"]["Current form"] = @benchmarkable( sim!($mhe_pendulum_ipopt_curr, $N, $u; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_ESTIM["Pendulum"]["MovingHorizonEstimator"]["Ipopt"]["Current form (Hessian)"] = @benchmarkable( sim!($mhe_pendulum_ipopt_currh, $N, $u; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_ESTIM["Pendulum"]["MovingHorizonEstimator"]["Ipopt"]["Prediction form"] = @benchmarkable( sim!($mhe_pendulum_ipopt_pred, $N, $u; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) CASE_ESTIM["Pendulum"]["MovingHorizonEstimator"]["Ipopt"]["Prediction form (Hessian)"] = @benchmarkable( sim!($mhe_pendulum_ipopt_predh, $N, $u; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) -CASE_ESTIM["Pendulum"]["MovingHorizonEstimator"]["MadNLP"]["Current form"] = +CASE_ESTIM["Pendulum"]["MovingHorizonEstimator"]["MadNLP"]["Current form (Hessian)"] = @benchmarkable( - sim!($mhe_pendulum_madnlp_curr, $N, $u; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + sim!($mhe_pendulum_madnlp_currh, $N, $u; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) -CASE_ESTIM["Pendulum"]["MovingHorizonEstimator"]["MadNLP"]["Prediction form"] = +CASE_ESTIM["Pendulum"]["MovingHorizonEstimator"]["MadNLP"]["Prediction form (Hessian)"] = @benchmarkable( - sim!($mhe_pendulum_madnlp_pred, $N, $u; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), - samples=samples, evals=evals, seconds=seconds + sim!($mhe_pendulum_madnlp_predh, $N, $u; plant=$plant, x_0=$x_0, x̂_0=$x̂_0, progress=false), + samples=samples, evals=evals, seconds=seconds, setup=GC.gc() ) \ No newline at end of file