From 9ca9e43cbb0f6eded4f597200d1c169732f90553 Mon Sep 17 00:00:00 2001 From: Fredrik Bagge Carlson Date: Wed, 10 Dec 2025 16:24:35 +0100 Subject: [PATCH 1/7] add LinearMPC extension --- Project.toml | 6 ++ ext/RobustAndOptimalControlLinearMPCExt.jl | 95 ++++++++++++++++++++++ src/ExtendedStateSpace.jl | 2 +- test/runtests.jl | 5 ++ test/test_linearmpc_ext.jl | 93 +++++++++++++++++++++ 5 files changed, 200 insertions(+), 1 deletion(-) create mode 100644 ext/RobustAndOptimalControlLinearMPCExt.jl create mode 100644 test/test_linearmpc_ext.jl diff --git a/Project.toml b/Project.toml index 3a656890..2df64c4a 100644 --- a/Project.toml +++ b/Project.toml @@ -23,6 +23,12 @@ Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" UUIDs = "cf7118a7-6976-5b1a-9a39-7adc72f591a4" UnPack = "3a884ed6-31ef-47d7-9d2a-63182c4928ed" +[weakdeps] +LinearMPC = "82e1c212-e1a2-49d2-b26a-a31d6968e3bd" + +[extensions] +RobustAndOptimalControlLinearMPCExt = "LinearMPC" + [compat] ChainRulesCore = "1" ComponentArrays = "0.9, 0.10, 0.11, 0.12, 0.13, 0.14, 0.15" diff --git a/ext/RobustAndOptimalControlLinearMPCExt.jl b/ext/RobustAndOptimalControlLinearMPCExt.jl new file mode 100644 index 00000000..63d641ad --- /dev/null +++ b/ext/RobustAndOptimalControlLinearMPCExt.jl @@ -0,0 +1,95 @@ +""" + RobustAndOptimalControlLinearMPCExt + +This extension allows you to convert an `LQGProblem` from RobustAndOptimalControl.jl to a `LinearMPC.MPC` object from LinearMPC.jl. +""" +module RobustAndOptimalControlLinearMPCExt + +using LinearMPC +using LinearAlgebra +using RobustAndOptimalControl: LQGProblem + + +""" + LinearMPC.MPC(prob::LQGProblem; N, Nc=N, Q3=nothing, Qf=nothing, + umin=nothing, umax=nothing, ymin=nothing, ymax=nothing, kwargs...) + +Convert an `LQGProblem` from RobustAndOptimalControl.jl to a `LinearMPC.MPC` object. + +# Arguments +- `prob::LQGProblem`: The LQG problem to convert +- `N::Int`: Prediction horizon (required) +- `Nc::Int`: Control horizon (default: `N`) +- `Q3`: Input rate cost matrix (default: zeros) +- `Qf`: Terminal state cost matrix (default: none) +- `umin`, `umax`: Input bounds (default: none) +- `ymin`, `ymax`: Output bounds (default: none) +- `kwargs...`: Additional arguments passed to LinearMPC.MPC + +# Notes +- Only discrete-time systems are supported +- The Kalman filter/observer from LQGProblem is not converted +- Uses C1 (performance output) as the controlled output matrix +- Cost matrices Q1, Q2 from LQGProblem map to Q, R in LinearMPC + +# Example +```julia +using RobustAndOptimalControl, LinearMPC, LinearAlgebra + +sys = ss([1 0.1; 0 1], [0; 0.1], [1 0], 0, 0.1) +lqg = LQGProblem(sys, I(2), I(1), I(2), 0.01*I(1)) +mpc = LinearMPC.MPC(lqg; N=20, umin=[-0.3], umax=[0.3]) + +sim = LinearMPC.Simulation(mpc, N=100, r=[1.0, 0]) + +using Plots +plot(sim) +``` +""" +function LinearMPC.MPC(prob::LQGProblem; + N::Int, + Nc::Int = N, + Q3 = zeros(0,0), + Qf = zeros(0,0), + umin = zeros(0), + umax = zeros(0), + ymin = zeros(0), + ymax = zeros(0), + kwargs... +) + # Validate discrete-time system + if !ControlSystemsBase.isdiscrete(prob) + error("Only discrete-time systems are supported. Got a continuous-time system.") + end + + # Extract system matrices + F = Matrix(prob.A) # State transition matrix + G = Matrix(prob.B2) # Control input matrix + C = Matrix(prob.C1) # Performance output matrix + + # Get sampling time + Ts = prob.Ts + + # Get dimensions + nu = size(G, 2) + ny = size(C, 1) + + # Create the LinearMPC.MPC object + mpc = LinearMPC.MPC(F, G; Ts, C, Np=N, Nc, kwargs...) + + # Set objective + Q = Matrix(prob.Q1) + R = Matrix(prob.Q2) + Rr = Matrix(Q3) + + LinearMPC.set_objective!(mpc; Q, R, Rr, Qf) + LinearMPC.set_bounds!(mpc; umin, umax, ymin, ymax) + + # Set labels from system names + # sys = prob.sys + # set_labels!(mpc; x=Symbol.(state_names(sys)), u=Symbol.(input_names(sys))) + + return mpc +end + +end # module \ No newline at end of file diff --git a/src/ExtendedStateSpace.jl b/src/ExtendedStateSpace.jl index 6c74588a..084953e6 100644 --- a/src/ExtendedStateSpace.jl +++ b/src/ExtendedStateSpace.jl @@ -599,7 +599,7 @@ The conversion from a regular statespace object to an `ExtendedStateSpace` creat i.e., the system and performance mappings are identical, `system_mapping(se) == performance_mapping(se) == s`. However, all matrices `B1, B2, C1, C2; D11, D12, D21, D22` are overridable by a corresponding keyword argument. In this case, the controlled outputs are the same as measured outputs. -Related: `se = convert(ExtendedStateSpace, s::StateSpace)` produces an `ExtendedStateSpace` with empty `performance_mapping` from w->z such that `ss(se) == s`. +Related: `se = convert(ExtendedStateSpace, s::StateSpace)` produces an `ExtendedStateSpace` with empty `performance_mapping` from w->z such that `ss(se) == s`. `ExtendedStateSpace(sys, B1=I, C1=I)` leads to a system where all state variables are affected by noise, and all are considered performance outputs, this corresponds to the traditional use of the functions `lqr` and `kalman`. """ function ExtendedStateSpace(s::AbstractStateSpace; A = s.A, diff --git a/test/runtests.jl b/test/runtests.jl index ebad5ece..ec429020 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -131,4 +131,9 @@ using Test include("test_mcm_nugap.jl") end + @testset "LinearMPC extension" begin + @info "Testing LinearMPC extension" + include("test_linearmpc_ext.jl") + end + end diff --git a/test/test_linearmpc_ext.jl b/test/test_linearmpc_ext.jl new file mode 100644 index 00000000..4926fea5 --- /dev/null +++ b/test/test_linearmpc_ext.jl @@ -0,0 +1,93 @@ +using Test +using LinearAlgebra +using LinearMPC +using RobustAndOptimalControl + +# Test basic LQGProblem conversion (from docstring example) +Ts = 0.1 +sys = ss([1 Ts; 0 1], [0; Ts], [1 0], 0, Ts) + +# Create LQGProblem: Q1=state cost, Q2=input cost, R1=process noise, R2=measurement noise +lqg = LQGProblem(sys, I(2), I(1), I(2), 0.01*I(1)) + +# Convert to LinearMPC with constraints +mpc = LinearMPC.MPC(lqg; N=20, umin=[-0.3], umax=[0.3]) + +# Check basic properties +@test mpc.model.Ts == Ts +@test mpc.model.F == sys.A +@test mpc.model.G == sys.B +@test mpc.settings.Np == 20 + +# Simulate and verify it works +sim = LinearMPC.Simulation(mpc; N=100, r=[1.0, 0]) + +# Check that simulation ran +@test size(sim.xs, 2) == 101 # 100 steps + initial state +@test size(sim.us, 2) == 100 + +# Check constraint satisfaction +@test all(sim.us .>= -0.3 - 1e-6) +@test all(sim.us .<= 0.3 + 1e-6) + +# Check that output converges toward reference +final_output = (sys.C * sim.xs[:, end])[1] +@test abs(final_output - 1.0) < 0.1 # Should be close to reference + +@testset "LQGProblem with Q3 (input rate penalty)" begin + # Test with input rate cost + Q3 = 0.5 * I(1) + mpc_q3 = LinearMPC.MPC(lqg; N=15, Q3, umin=[-0.5], umax=[0.5]) + + sim_q3 = LinearMPC.Simulation(mpc_q3; N=50, r=[0.5, 0]) + + # With Q3, control should be smoother (smaller rate of change) + du = diff(sim_q3.us, dims=2) + max_rate = maximum(abs.(du)) + @test max_rate < 0.3 # Rate should be limited due to Q3 penalty +end + +@testset "LQGProblem with terminal cost Qf" begin + # Test with terminal cost + Qf = 10.0 * I(2) # Higher terminal cost + mpc_qf = LinearMPC.MPC(lqg; N=10, Qf, umin=[-1.0], umax=[1.0]) + + sim_qf = LinearMPC.Simulation(mpc_qf; N=30, r=[0.3, 0]) + + # Should still work and converge + final_output = (sys.C * sim_qf.xs[:, end])[1] + @test abs(final_output - 0.3) < 0.1 +end + +@testset "LQGProblem MIMO system" begin + # Test with 2-input 2-output system + A = [0.9 0.1; 0.05 0.95] + B = [1.0 0.0; 0.0 1.0] + C = [1.0 0.0; 0.0 1.0] + D = zeros(2, 2) + sys_mimo = ss(A, B, C, D, Ts) + + lqg_mimo = LQGProblem(sys_mimo, I(2), I(2), I(2), 0.01*I(2)) + mpc_mimo = LinearMPC.MPC(lqg_mimo; N=15, + umin=[-0.5, -0.5], umax=[0.5, 0.5]) + + sim_mimo = LinearMPC.Simulation(mpc_mimo; N=50, r=[0.3, -0.2]) + + # Check dimensions + @test size(sim_mimo.us, 1) == 2 + @test size(sim_mimo.xs, 1) == 2 + + # Check constraints on both inputs + @test all(sim_mimo.us[1, :] .>= -0.5 - 1e-6) + @test all(sim_mimo.us[1, :] .<= 0.5 + 1e-6) + @test all(sim_mimo.us[2, :] .>= -0.5 - 1e-6) + @test all(sim_mimo.us[2, :] .<= 0.5 + 1e-6) +end + +@testset "LQGProblem continuous-time error" begin + # Test that continuous-time systems throw an error + sys_cont = ss([0 1; -1 -1], [0; 1], [1 0], 0) # Continuous-time + lqg_cont = LQGProblem(sys_cont, I(2), I(1), I(2), 0.01*I(1)) + + @test_throws ErrorException LinearMPC.MPC(lqg_cont; N=10) +end From 3905dfc5b03aa23b1abd46bc47f904697015a393 Mon Sep 17 00:00:00 2001 From: Fredrik Bagge Carlson Date: Wed, 10 Dec 2025 17:05:37 +0100 Subject: [PATCH 2/7] add test dep --- Project.toml | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index 2df64c4a..f4e844fb 100644 --- a/Project.toml +++ b/Project.toml @@ -36,6 +36,7 @@ ControlSystemsBase = "1.17" DescriptorSystems = "1.2" Distributions = "0.25" GenericSchur = "0.5.2" +LinearMPC = "0.7" MatrixEquations = "2" MatrixPencils = "1" MonteCarloMeasurements = "1.0" @@ -50,9 +51,10 @@ julia = "1.7" [extras] FiniteDiff = "6a86dc24-6348-571c-b903-95158fe2bd41" GenericLinearAlgebra = "14197337-ba66-59df-a3e3-ca00e7dcff7a" +LinearMPC = "82e1c212-e1a2-49d2-b26a-a31d6968e3bd" Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" Zygote = "e88e6eb3-aa80-5325-afca-941959d7151f" [targets] -test = ["Test", "Plots", "Zygote", "FiniteDiff", "GenericLinearAlgebra"] +test = ["Test", "Plots", "LinearMPC", "Zygote", "FiniteDiff", "GenericLinearAlgebra"] From 7edbc141982e089fc3acf7454a63d4048f1c6d7c Mon Sep 17 00:00:00 2001 From: Fredrik Bagge Carlson Date: Wed, 10 Dec 2025 17:08:25 +0100 Subject: [PATCH 3/7] add mention in docs --- docs/src/index.md | 18 ++++++++++++++++++ 1 file changed, 18 insertions(+) diff --git a/docs/src/index.md b/docs/src/index.md index 4e8b8b05..6832e3a1 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -297,3 +297,21 @@ A video tutorial on how to perform closed-loop analysis is available [here](http - [`hess_form`](@ref) - [`frequency_separation`](@ref) - [`RobustAndOptimalControl.blockdiagonalize`](@ref) + +## MPC +This package includes an extension for [LinearMPC.jl](https://darnstrom.github.io/LinearMPC.jl/stable/) that allows you to convert an [`LQGProblem`](@ref) into a linear MPC controller with added constraints. Example: +```@example lqg_mpc +using RobustAndOptimalControl, LinearMPC, LinearAlgebra + +sys = ss([1 0.1; 0 1], [0; 0.1], [1 0], 0, 0.1) +lqg = LQGProblem(sys, I(2), I(1), I(2), 0.01*I(1)) + +mpc = LinearMPC.MPC(lqg; N=20, umin=[-0.3], umax=[0.3]) + +sim = LinearMPC.Simulation(mpc, N=100, r=[1.0, 0]) + +using Plots +plot(sim) +``` + +See the [docs for LinearMPC.jl](https://darnstrom.github.io/LinearMPC.jl/stable/) for more details on how to use the MPC controller and modify its settings. \ No newline at end of file From abb7a5d78dc1ced766a22d86c171a49a954f8d31 Mon Sep 17 00:00:00 2001 From: Fredrik Bagge Carlson Date: Thu, 11 Dec 2025 06:37:42 +0100 Subject: [PATCH 4/7] add LinMPC to docs env --- docs/Project.toml | 1 + ext/RobustAndOptimalControlLinearMPCExt.jl | 3 ++- 2 files changed, 3 insertions(+), 1 deletion(-) diff --git a/docs/Project.toml b/docs/Project.toml index db590b6f..e6b72b55 100644 --- a/docs/Project.toml +++ b/docs/Project.toml @@ -3,6 +3,7 @@ ControlSystemsBase = "aaaaaaaa-a6ca-5380-bf3e-84a91bcd477e" DisplayAs = "0b91fe84-8a4c-11e9-3e1d-67c38462b6d6" Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4" ForwardDiff = "f6369f11-7733-5829-9624-2563aa707210" +LinearMPC = "82e1c212-e1a2-49d2-b26a-a31d6968e3bd" MonteCarloMeasurements = "0987c9cc-fe09-11e8-30f0-b96dd679fdca" Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80" RobustAndOptimalControl = "21fd56a4-db03-40ee-82ee-a87907bee541" diff --git a/ext/RobustAndOptimalControlLinearMPCExt.jl b/ext/RobustAndOptimalControlLinearMPCExt.jl index 63d641ad..9cf912b7 100644 --- a/ext/RobustAndOptimalControlLinearMPCExt.jl +++ b/ext/RobustAndOptimalControlLinearMPCExt.jl @@ -8,6 +8,7 @@ module RobustAndOptimalControlLinearMPCExt using LinearMPC using LinearAlgebra using RobustAndOptimalControl: LQGProblem +using ControlSystemsBase: isdiscrete """ @@ -58,7 +59,7 @@ function LinearMPC.MPC(prob::LQGProblem; kwargs... ) # Validate discrete-time system - if !ControlSystemsBase.isdiscrete(prob) + if !isdiscrete(prob) error("Only discrete-time systems are supported. Got a continuous-time system.") end From d006909d1c23ba580d75f7d386304d850ab20313 Mon Sep 17 00:00:00 2001 From: Fredrik Bagge Carlson Date: Thu, 11 Dec 2025 07:08:02 +0100 Subject: [PATCH 5/7] fix --- test/test_linearmpc_ext.jl | 6 +----- 1 file changed, 1 insertion(+), 5 deletions(-) diff --git a/test/test_linearmpc_ext.jl b/test/test_linearmpc_ext.jl index 4926fea5..67b39426 100644 --- a/test/test_linearmpc_ext.jl +++ b/test/test_linearmpc_ext.jl @@ -17,14 +17,10 @@ mpc = LinearMPC.MPC(lqg; N=20, umin=[-0.3], umax=[0.3]) @test mpc.model.Ts == Ts @test mpc.model.F == sys.A @test mpc.model.G == sys.B -@test mpc.settings.Np == 20 # Simulate and verify it works sim = LinearMPC.Simulation(mpc; N=100, r=[1.0, 0]) -# Check that simulation ran -@test size(sim.xs, 2) == 101 # 100 steps + initial state -@test size(sim.us, 2) == 100 # Check constraint satisfaction @test all(sim.us .>= -0.3 - 1e-6) @@ -44,7 +40,7 @@ final_output = (sys.C * sim.xs[:, end])[1] # With Q3, control should be smoother (smaller rate of change) du = diff(sim_q3.us, dims=2) max_rate = maximum(abs.(du)) - @test max_rate < 0.3 # Rate should be limited due to Q3 penalty + @test max_rate < 0.1 # Rate should be limited due to Q3 penalty end @testset "LQGProblem with terminal cost Qf" begin From e09543adb4e3a9b916ac8fc1a0fc4f35af21a026 Mon Sep 17 00:00:00 2001 From: Fredrik Bagge Carlson Date: Thu, 11 Dec 2025 08:16:53 +0100 Subject: [PATCH 6/7] add example file --- docs/src/index.md | 2 + examples/lqg_mpc_disturbance.jl | 106 ++++++++++++++++++++++++++++++++ src/lqg.jl | 12 ++++ test/runtests.jl | 5 ++ 4 files changed, 125 insertions(+) create mode 100644 examples/lqg_mpc_disturbance.jl diff --git a/docs/src/index.md b/docs/src/index.md index 6832e3a1..277e3362 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -314,4 +314,6 @@ using Plots plot(sim) ``` +See [Example: `lqg_mpc_disturbance.jl`](https://github.com/JuliaControl/RobustAndOptimalControl.jl/blob/master/examples/lqg_mpc_disturbance.jl) for a more detailed example. + See the [docs for LinearMPC.jl](https://darnstrom.github.io/LinearMPC.jl/stable/) for more details on how to use the MPC controller and modify its settings. \ No newline at end of file diff --git a/examples/lqg_mpc_disturbance.jl b/examples/lqg_mpc_disturbance.jl new file mode 100644 index 00000000..96c01314 --- /dev/null +++ b/examples/lqg_mpc_disturbance.jl @@ -0,0 +1,106 @@ +#= +This file demonstrates +- Adding integral action to an LQG controller by means of state augmentation +- Conversion of the LQG controller into an MPC contorller with added constraints +=# +using RobustAndOptimalControl, ControlSystemsBase, Plots, LinearAlgebra +Ts = 1 # Sample time +G = c2d(ss(tf(1, [10, 1])), Ts) # Process model + +nx = G.nx +nu = G.nu +ny = G.ny +x0 = zeros(G.nx) # Initial condition + +Q1 = 100diagm(ones(G.nx)) # state cost matrix +Q2 = 0.01diagm(ones(nu)) # control cost matrix + +R1 = 0.001I(nx) # State noise covariance +R2 = I(ny) # measurement noise covariance +prob = LQGProblem(G, Q1, Q2, R1, R2) + +disturbance = (x, t) -> t * Ts ≥ 10 # This is our load disturbance, a step at ``t = 10`` +Gcl = G_PS(prob) # This forms the transfer function from load disturbance to output +res = lsim(Gcl, disturbance, 100) +plot(res) + +Gd1 = add_low_frequency_disturbance(G, ϵ = 1e-6, measurement=false) +Gd2 = add_low_frequency_disturbance(G, ϵ = 1e-6, measurement=true) # The ϵ moves the integrator pole slightly into the stable region +plots = map([Gd1, Gd2]) do Gd + nx = Gd.nx + + C = Gd.C + Q1 = 100C'diagm(ones(G.nx)) * C # state cost matrix + x0 = zeros(nx) + + R1 = diagm([0.001, 1]) + R2 = I(ny) + prob = LQGProblem(Gd, Q1, Q2, R1, R2) + Gcl = [G_PS(prob); -comp_sensitivity(prob)] # -comp_sensitivity(prob) is the same as the transfer function from load disturbance to control signal + res = lsim(Gcl, disturbance, 100) + f1 = plot(res, ylabel=["y" "u"]); ylims!((-0.05, 0.3), sp = 1) + + w = exp10.(LinRange(-3, log10(pi / Ts), 200)) + f2 = gangoffourplot(prob, w, lab = "", legend = :bottomright) + + R1 = diagm([0.001, 0.2]) # Reduce the noise on the integrator state from 1 to 0.2 + R2 = I(ny) + prob = LQGProblem(Gd, Q1, Q2, R1, R2) + + Gcl = [G_PS(prob); -comp_sensitivity(prob)] + res = lsim(Gcl, disturbance, 100) + f3 = plot(res, ylabel=["y" "u"]); ylims!((-0.05, 0.3), sp = 1) + f4 = gangoffourplot(prob, w, lab = "", legend = :bottomright) + + plot(f1, f2, f3, f4, titlefontsize=10) +end + +plot(plots..., size=(1200,1000)) + + + +# ============================================================================== +## LinearMPC +# The example below illustrate how we can convert the LQGProblem into an MPC contorller by loading LinearMPC.jl +# We then perform a rather low-level simulation with a manual loop, where we form the observer `obs` and step the plant +# ============================================================================== + +using LinearMPC + +Gd = add_low_frequency_disturbance(G, ϵ = 1e-6, measurement=false) + +C = Gd.C +Q1 = 100diagm([1.0]) # output cost matrix +R1 = diagm([0.001, 1]) # Dynamics noise covariance +R2 = I(ny) # Measurement noise covariance +Gde = ExtendedStateSpace(Gd, B1=I) # Since B1=I, R1 has size determined by state dimension, but C1=C, so Q1 has size determined by the output dimension +prob = LQGProblem(Gde, Q1, Q2, R1, R2) +obs = observer_predictor(prob, direct=false) # If a predictor is used, the observer update should be carried out in the end of the loop, if we use the filter below, we should instead perform the observer update in the beginning of the loop directly after obtaining the new measurement but before computing a new control signal. +obs = observer_filter(prob) + +mpc = LinearMPC.MPC(prob; N=20, umin=[-3], umax=[3]) +x = zeros(G.nx) # True plant state +xh = zeros(Gd.nx) # Observer state +X = [x[]] # For storage +U = Float64[] +u_mpc = [0.0] +for i = 1:50 + y = G.C*x # Compute the true measurement output + xh = obs.A * xh + obs.B*[u_mpc; y] # Predict one step with the observer, u here is the control signal from the previous iteration, if using the predictor, use u from the current iteration and perform the observer update in the end of the loop instead + u_disturbance = i * Ts ≥ 10 ? 1.0 : 0.0 + r = [1.0] # output reference + u_mpc = LinearMPC.compute_control(mpc, xh; r) # Call MPC optimizer with estimated state + u_tot = u_mpc .+ u_disturbance # Total input is control signal + disturbance + x = G.A*x + G.B*u_tot # Advance the true plant state + push!(X, x[]) # Store data for plotting + push!(U, u_mpc[]) +end + + +plot(X*G.C', layout=2, sp=1, label="\$y\$") +plot!(U, sp=2, label="\$u\$") +hline!([1.0 3.0], linestyle=:dash, color=:black, label=["Reference" "\$u_{max}\$"], sp=[1 2]) + + +using Test +@test (G.C*X[end])[] ≈ 1 rtol=1e-4 \ No newline at end of file diff --git a/src/lqg.jl b/src/lqg.jl index d7e8b50b..f049fcfd 100644 --- a/src/lqg.jl +++ b/src/lqg.jl @@ -354,6 +354,18 @@ function ControlSystemsBase.observer_controller(l::LQGProblem, L::AbstractMatrix ss(Ac, Bc, Cc, Dc, l.timeevol) end +function ControlSystemsBase.observer_predictor(l::LQGProblem, K::Union{AbstractMatrix, Nothing} = nothing; direct = false, kwargs...) + P = system_mapping(l, identity) + if K === nothing + K = kalman(l; direct) + end + observer_predictor(P, K; kwargs...) +end + +function ControlSystemsBase.observer_filter(l::LQGProblem, K = kalman(l); kwargs...) + P = system_mapping(l, identity) + observer_filter(P, K; kwargs...) +end """ ff_controller(sys, L, K; comp_dc = true) diff --git a/test/runtests.jl b/test/runtests.jl index ec429020..1b631ccf 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -126,6 +126,11 @@ using Test include("../examples/uncertain.jl") end + @testset "lqg_mpc_disturbance" begin + @info "Testing LQG MPC example" + include("../examples/lqg_mpc_disturbance.jl") + end + @testset "mcm_nugap" begin @info "Testing mcm_nugap" include("test_mcm_nugap.jl") From 4dc08bfd17f731c2f854d02cb4a7be1a294cca27 Mon Sep 17 00:00:00 2001 From: Fredrik Bagge Carlson Date: Thu, 11 Dec 2025 09:02:22 +0100 Subject: [PATCH 7/7] global --- examples/lqg_mpc_disturbance.jl | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/examples/lqg_mpc_disturbance.jl b/examples/lqg_mpc_disturbance.jl index 96c01314..38398180 100644 --- a/examples/lqg_mpc_disturbance.jl +++ b/examples/lqg_mpc_disturbance.jl @@ -85,6 +85,7 @@ X = [x[]] # For storage U = Float64[] u_mpc = [0.0] for i = 1:50 + global x, xh, u_mpc y = G.C*x # Compute the true measurement output xh = obs.A * xh + obs.B*[u_mpc; y] # Predict one step with the observer, u here is the control signal from the previous iteration, if using the predictor, use u from the current iteration and perform the observer update in the end of the loop instead u_disturbance = i * Ts ≥ 10 ? 1.0 : 0.0 @@ -96,11 +97,10 @@ for i = 1:50 push!(U, u_mpc[]) end +using Test +@test (G.C*X[end])[] ≈ 1 rtol=1e-4 plot(X*G.C', layout=2, sp=1, label="\$y\$") plot!(U, sp=2, label="\$u\$") hline!([1.0 3.0], linestyle=:dash, color=:black, label=["Reference" "\$u_{max}\$"], sp=[1 2]) - -using Test -@test (G.C*X[end])[] ≈ 1 rtol=1e-4 \ No newline at end of file