From aa48eb1315e4864f435e3ab6025aa3217fb7e92e Mon Sep 17 00:00:00 2001 From: franckgaga Date: Fri, 7 Nov 2025 16:21:38 -0500 Subject: [PATCH 1/6] changed: better `jacobian` and `hessian` defaults for the sparse case --- src/controller/nonlinmpc.jl | 7 ++++--- src/estimator/mhe/construct.jl | 6 ++++-- src/general.jl | 9 +++++++++ 3 files changed, 17 insertions(+), 5 deletions(-) diff --git a/src/controller/nonlinmpc.jl b/src/controller/nonlinmpc.jl index e622e9965..2703e489b 100644 --- a/src/controller/nonlinmpc.jl +++ b/src/controller/nonlinmpc.jl @@ -5,7 +5,7 @@ const DEFAULT_NONLINMPC_JACDENSE = AutoForwardDiff() const DEFAULT_NONLINMPC_JACSPARSE = AutoSparse( AutoForwardDiff(); sparsity_detector=TracerSparsityDetector(), - coloring_algorithm=GreedyColoringAlgorithm(), + coloring_algorithm=GreedyColoringAlgorithm(ALL_COLORING_ORDERS), ) const DEFAULT_NONLINMPC_HESSIAN = DEFAULT_NONLINMPC_JACSPARSE @@ -291,10 +291,11 @@ NonLinMPC controller with a sample time Ts = 10.0 s: AutoSparse( AutoForwardDiff(); sparsity_detector = TracerSparsityDetector(), - coloring_algorithm = GreedyColoringAlgorithm() + coloring_algorithm = GreedyColoringAlgorithm(ALL_COLORING_ORDERS) ) ``` - This is also the sparse backend selected for the Hessian of the Lagrangian function if + that is, it will test all coloring orders at preparation and keep the best. This is + also the sparse backend selected for the Hessian of the Lagrangian function if `oracle=true` and `hessian=true`, which is the second exception. Second order derivatives are only supported with `oracle=true` option. diff --git a/src/estimator/mhe/construct.jl b/src/estimator/mhe/construct.jl index 6069b6ea4..b17225426 100644 --- a/src/estimator/mhe/construct.jl +++ b/src/estimator/mhe/construct.jl @@ -5,7 +5,7 @@ const DEFAULT_NONLINMHE_JACOBIAN = AutoForwardDiff() const DEFAULT_NONLINMHE_HESSIAN = AutoSparse( AutoForwardDiff(); sparsity_detector=TracerSparsityDetector(), - coloring_algorithm=GreedyColoringAlgorithm(), + coloring_algorithm=GreedyColoringAlgorithm(ALL_COLORING_ORDERS), ) @doc raw""" @@ -382,9 +382,11 @@ MovingHorizonEstimator estimator with a sample time Ts = 10.0 s: AutoSparse( AutoForwardDiff(); sparsity_detector = TracerSparsityDetector(), - coloring_algorithm = GreedyColoringAlgorithm() + coloring_algorithm = GreedyColoringAlgorithm(ALL_COLORING_ORDERS) ) ``` + that is, it will test all coloring orders at preparation and keep the best. + The slack variable ``ε`` relaxes the constraints if enabled, see [`setconstraint!`](@ref). It is disabled by default for the MHE (from `Cwt=Inf`) but it should be activated for problems with two or more types of bounds, to ensure feasibility (e.g. on the estimated diff --git a/src/general.jl b/src/general.jl index 3972d6c32..4ce717392 100644 --- a/src/general.jl +++ b/src/general.jl @@ -6,6 +6,15 @@ const DEFAULT_LWT = 0.0 const DEFAULT_CWT = 1e5 const DEFAULT_EWT = 0.0 +"All deterministic algorithms for matrix coloring order in `SparseMatrixColoring.jl`." +const ALL_COLORING_ORDERS = ( + NaturalOrder(), + LargestFirst(), + SmallestLast(), + IncidenceDegree(), + DynamicLargestFirst(), +) + "Termination status that means 'no solution available'." const ERROR_STATUSES = ( JuMP.INFEASIBLE, JuMP.DUAL_INFEASIBLE, JuMP.LOCALLY_INFEASIBLE, From 0bbb0b6803ee525438d1aae124ebad024c618cea Mon Sep 17 00:00:00 2001 From: franckgaga Date: Fri, 7 Nov 2025 16:28:50 -0500 Subject: [PATCH 2/6] debug: import symbols from `SparseMatrixColoring` --- src/ModelPredictiveControl.jl | 2 ++ 1 file changed, 2 insertions(+) diff --git a/src/ModelPredictiveControl.jl b/src/ModelPredictiveControl.jl index 423475c5f..8bb73005f 100644 --- a/src/ModelPredictiveControl.jl +++ b/src/ModelPredictiveControl.jl @@ -14,6 +14,8 @@ using DifferentiationInterface: hessian!, value_gradient_and_hessian!, prepare_h using DifferentiationInterface: Constant, Cache using SparseConnectivityTracer: TracerSparsityDetector using SparseMatrixColorings: GreedyColoringAlgorithm, sparsity_pattern +using SparseMatrixColorings: NaturalOrder, LargestFirst, SmallestLast +using SparseMatrixColorings: IncidenceDegree, DynamicLargestFirst import ProgressLogging From 083dedec9ad155883534097cbe3838aaa6f3df1d Mon Sep 17 00:00:00 2001 From: franckgaga Date: Fri, 7 Nov 2025 16:51:22 -0500 Subject: [PATCH 3/6] doc: add estimation horizon in `README.md` --- README.md | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/README.md b/README.md index b35894ded..b80c4bce5 100644 --- a/README.md +++ b/README.md @@ -51,14 +51,14 @@ Our goal is controlling the first output $y_1$, but the second one $y_2$ should 35: ```julia -mhe = MovingHorizonEstimator(model) +mhe = MovingHorizonEstimator(model, He=10) mpc = LinMPC(mhe, Mwt=[1, 0], Nwt=[0.1]) mpc = setconstraint!(mpc, ymax=[Inf, 35]) ``` -The keyword arguments `Mwt` and `Nwt` are the output setpoint tracking and move suppression -weights, respectively. A setpoint step change of five tests `mpc` controller in closed-loop. -The result is displayed with [`Plots.jl`](https://github.com/JuliaPlots/Plots.jl): +The keyword arguments `He`, `Mwt` and `Nwt` are the estimation horizon, the output setpoint +tracking and move suppression weights, respectively. A setpoint step change of five tests +`mpc` controller in closed-loop. The result is displayed with [`Plots.jl`](https://github.com/JuliaPlots/Plots.jl): ```julia using Plots From 364ac2b0722372c53130a6436c5b124038dc24c9 Mon Sep 17 00:00:00 2001 From: franckgaga Date: Fri, 7 Nov 2025 17:07:24 -0500 Subject: [PATCH 4/6] revert: going back to old README file --- README.md | 9 ++++----- 1 file changed, 4 insertions(+), 5 deletions(-) diff --git a/README.md b/README.md index b80c4bce5..d160d35ea 100644 --- a/README.md +++ b/README.md @@ -51,14 +51,13 @@ Our goal is controlling the first output $y_1$, but the second one $y_2$ should 35: ```julia -mhe = MovingHorizonEstimator(model, He=10) -mpc = LinMPC(mhe, Mwt=[1, 0], Nwt=[0.1]) +mpc = LinMPC(model, Mwt=[1, 0], Nwt=[0.1]) mpc = setconstraint!(mpc, ymax=[Inf, 35]) ``` -The keyword arguments `He`, `Mwt` and `Nwt` are the estimation horizon, the output setpoint -tracking and move suppression weights, respectively. A setpoint step change of five tests -`mpc` controller in closed-loop. The result is displayed with [`Plots.jl`](https://github.com/JuliaPlots/Plots.jl): +The keyword arguments `Mwt` and `Nwt` are the output setpoint tracking and move suppression +weights, respectively. A setpoint step change of five tests `mpc` controller in closed-loop. +The result is displayed with [`Plots.jl`](https://github.com/JuliaPlots/Plots.jl): ```julia using Plots From 7e3ab467ca53ebb1077ee22bfa4b41649fc58fff Mon Sep 17 00:00:00 2001 From: franckgaga Date: Sat, 8 Nov 2025 10:43:58 -0500 Subject: [PATCH 5/6] doc: explicit color orders --- src/controller/nonlinmpc.jl | 10 ++++++++-- src/estimator/mhe/construct.jl | 8 +++++++- 2 files changed, 15 insertions(+), 3 deletions(-) diff --git a/src/controller/nonlinmpc.jl b/src/controller/nonlinmpc.jl index 2703e489b..7ac5ddc75 100644 --- a/src/controller/nonlinmpc.jl +++ b/src/controller/nonlinmpc.jl @@ -291,10 +291,16 @@ NonLinMPC controller with a sample time Ts = 10.0 s: AutoSparse( AutoForwardDiff(); sparsity_detector = TracerSparsityDetector(), - coloring_algorithm = GreedyColoringAlgorithm(ALL_COLORING_ORDERS) + coloring_algorithm = GreedyColoringAlgorithm(( + NaturalOrder(), + LargestFirst(), + SmallestLast(), + IncidenceDegree(), + DynamicLargestFirst() + )) ) ``` - that is, it will test all coloring orders at preparation and keep the best. This is + that is, it will test many coloring orders at preparation and keep the best. This is also the sparse backend selected for the Hessian of the Lagrangian function if `oracle=true` and `hessian=true`, which is the second exception. Second order derivatives are only supported with `oracle=true` option. diff --git a/src/estimator/mhe/construct.jl b/src/estimator/mhe/construct.jl index b17225426..08e2a8d85 100644 --- a/src/estimator/mhe/construct.jl +++ b/src/estimator/mhe/construct.jl @@ -382,7 +382,13 @@ MovingHorizonEstimator estimator with a sample time Ts = 10.0 s: AutoSparse( AutoForwardDiff(); sparsity_detector = TracerSparsityDetector(), - coloring_algorithm = GreedyColoringAlgorithm(ALL_COLORING_ORDERS) + coloring_algorithm = GreedyColoringAlgorithm(( + NaturalOrder(), + LargestFirst(), + SmallestLast(), + IncidenceDegree(), + DynamicLargestFirst() + )) ) ``` that is, it will test all coloring orders at preparation and keep the best. From 63369c19ec5ca94682214b25e08855754098e27c Mon Sep 17 00:00:00 2001 From: franckgaga Date: Sat, 8 Nov 2025 10:44:08 -0500 Subject: [PATCH 6/6] doc: idem --- src/estimator/mhe/construct.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/estimator/mhe/construct.jl b/src/estimator/mhe/construct.jl index 08e2a8d85..5d7c9120a 100644 --- a/src/estimator/mhe/construct.jl +++ b/src/estimator/mhe/construct.jl @@ -391,7 +391,7 @@ MovingHorizonEstimator estimator with a sample time Ts = 10.0 s: )) ) ``` - that is, it will test all coloring orders at preparation and keep the best. + that is, it will test many coloring orders at preparation and keep the best. The slack variable ``ε`` relaxes the constraints if enabled, see [`setconstraint!`](@ref). It is disabled by default for the MHE (from `Cwt=Inf`) but it should be activated for