Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
3 changes: 2 additions & 1 deletion src/ModelPredictiveControl.jl
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,6 @@ using LinearAlgebra, SparseArrays
using Random: randn

using RecipesBase
using ProgressLogging

using DifferentiationInterface: ADTypes.AbstractADType, AutoForwardDiff, AutoSparse
using DifferentiationInterface: gradient!, jacobian!, prepare_gradient, prepare_jacobian
Expand All @@ -14,6 +13,8 @@ using DifferentiationInterface: Constant, Cache
using SparseConnectivityTracer: TracerSparsityDetector
using SparseMatrixColorings: GreedyColoringAlgorithm, sparsity_pattern

import ProgressLogging

import ForwardDiff

import ControlSystemsBase
Expand Down
13 changes: 6 additions & 7 deletions src/controller/transcription.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1169,8 +1169,8 @@ The method mutates `Ŷ0`, `x̂0end`, `X̂0`, `Û0` and `K0` arguments. The aug
[`f̂!`](@ref) and [`ĥ!`](@ref) functions is called recursively in a `for` loop:
```math
\begin{aligned}
\mathbf{x̂_0}(k+1) &= \mathbf{f̂}\Big(\mathbf{x̂_0}(k), \mathbf{u_0}(k), \mathbf{d_0}(k) \Big) \\
\mathbf{ŷ_0}(k) &= \mathbf{ĥ}\Big(\mathbf{x̂_0}(k), \mathbf{d_0}(k) \Big)
\mathbf{x̂_0}(k+1) &= \mathbf{f̂}\Big(\mathbf{x̂_0}(k), \mathbf{u_0}(k), \mathbf{d̂_0}(k) \Big) \\
\mathbf{ŷ_0}(k) &= \mathbf{ĥ}\Big(\mathbf{x̂_0}(k), \mathbf{d̂_0}(k) \Big)
\end{aligned}
```
"""
Expand Down Expand Up @@ -1211,7 +1211,7 @@ Compute vectors if `model` is a [`NonLinModel`](@ref) and other [`TranscriptionM
The method mutates `Ŷ0` and `x̂0end` arguments. The augmented output function [`ĥ!`](@ref)
is called multiple times in a `for` loop:
```math
\mathbf{ŷ_0}(k) = \mathbf{ĥ}\Big(\mathbf{x̂_0}(k), \mathbf{d_0}(k) \Big)
\mathbf{ŷ_0}(k) = \mathbf{ĥ}\Big(\mathbf{x̂_0}(k), \mathbf{d̂_0}(k) \Big)
```
in which ``\mathbf{x̂_0}`` is the augmented state extracted from the decision variable `Z̃`.
"""
Expand Down Expand Up @@ -1332,7 +1332,7 @@ Nonlinear equality constrains for [`NonLinModel`](@ref) and [`MultipleShooting`]
The method mutates the `geq`, `X̂0`, `Û0` and `K0` vectors in argument. The nonlinear
equality constraints `geq` only includes the augmented state defects, computed with:
```math
\mathbf{ŝ}(k+1) = \mathbf{f̂}\Big(\mathbf{x̂_0}(k), \mathbf{u_0}(k), \mathbf{d_0}(k)\Big)
\mathbf{ŝ}(k+1) = \mathbf{f̂}\Big(\mathbf{x̂_0}(k), \mathbf{u_0}(k), \mathbf{d̂_0}(k)\Big)
- \mathbf{x̂_0}(k+1)
```
in which the augmented state ``\mathbf{x̂_0}`` are extracted from the decision variables
Expand Down Expand Up @@ -1395,8 +1395,8 @@ extracted from the decision variables `Z̃`. The ``\mathbf{k}`` coefficients are
evaluated from the continuous-time function `model.f!` and:
```math
\begin{aligned}
\mathbf{k}_1 &= \mathbf{f}\Big(\mathbf{x_0}(k), \mathbf{û_0}(k), \mathbf{d_0}(k) \Big) \\
\mathbf{k}_2 &= \mathbf{f}\Big(\mathbf{x_0}(k+1), \mathbf{û_0}(k+h), \mathbf{d_0}(k+1)\Big)
\mathbf{k}_1 &= \mathbf{f}\Big(\mathbf{x_0}(k), \mathbf{û_0}(k), \mathbf{d̂_0}(k) \Big) \\
\mathbf{k}_2 &= \mathbf{f}\Big(\mathbf{x_0}(k+1), \mathbf{û_0}(k+h), \mathbf{d̂_0}(k+1)\Big)
\end{aligned}
```
in which ``h`` is the hold order `transcription.h` and the disturbed input is:
Expand All @@ -1419,7 +1419,6 @@ function con_nonlinprogeq!(
nk = get_nk(model, transcription)
D̂0 = mpc.D̂0
X̂0_Z̃ = @views Z̃[(nΔU+1):(nΔU+nX̂)]
#TODO: allow parallel for loop or threads?
@threadsif f_threads for j=1:Hp
if j < 2
x̂0 = @views mpc.estim.x̂0[1:nx̂]
Expand Down
11 changes: 11 additions & 0 deletions src/general.jl
Original file line number Diff line number Diff line change
Expand Up @@ -133,4 +133,15 @@ macro threadsif(flag, expr)
$expr
end
end |> esc
end

"Add `ProgressLogging.@progress` with the name `name` to a `for` loop if `flag==true`"
macro progressif(flag, name, expr)
quote
if $(flag)
ProgressLogging.@progress $name $expr
else
$expr
end
end |> esc
end
16 changes: 11 additions & 5 deletions src/plot_sim.jl
Original file line number Diff line number Diff line change
Expand Up @@ -106,13 +106,16 @@ end


@doc raw"""
sim!(plant::SimModel, N::Int, u=plant.uop.+1, d=plant.dop; x_0=plant.xop) -> res
sim!(
plant::SimModel, N::Int, u=plant.uop.+1, d=plant.dop; x_0=plant.xop, progress=true
) -> res

Open-loop simulation of `plant` for `N` time steps, default to unit bump test on all inputs.

The manipulated inputs ``\mathbf{u}`` and measured disturbances ``\mathbf{d}`` are held
constant at `u` and `d` values, respectively. The plant initial state ``\mathbf{x}(0)`` is
specified by `x_0` keyword arguments. The function returns [`SimResult`](@ref) instances
specified by `x_0` keyword arguments. If `progress` is `true`, VS Code will display a
progress percentage of the simulation. The function returns [`SimResult`](@ref) instances
that can be visualized by calling `plot` on them. Note that the method mutates `plant`
internal states.

Expand All @@ -129,15 +132,16 @@ function sim!(
N::Int,
u::Vector = plant.uop.+1,
d::Vector = plant.dop;
x_0 = plant.xop
x_0 = plant.xop,
progress = true
) where {NT<:Real}
T_data = collect(plant.Ts*(0:(N-1)))
Y_data = Matrix{NT}(undef, plant.ny, N)
U_data = Matrix{NT}(undef, plant.nu, N)
D_data = Matrix{NT}(undef, plant.nd, N)
X_data = Matrix{NT}(undef, plant.nx, N)
setstate!(plant, x_0)
@progress name="$(nameof(typeof(plant))) simulation" for i=1:N
@progressif progress name="$(nameof(typeof(plant))) simulation" for i=1:N
y = evaloutput(plant, d)
Y_data[:, i] .= y
U_data[:, i] .= u
Expand Down Expand Up @@ -187,6 +191,7 @@ vectors. The simulated sensor and process noises of `plant` are specified by `y_
- `x̂_0 = nothing` or *`xhat_0`* : initial estimate ``\mathbf{x̂}(0)``, [`initstate!`](@ref)
is used if `nothing`
- `lastu = plant.uop` : last plant input ``\mathbf{u}`` for ``\mathbf{x̂}`` initialization
- `progress = true` : display a progress percentage in VS Code if `true`

# Examples
```jldoctest
Expand Down Expand Up @@ -263,6 +268,7 @@ function sim_closedloop!(
x_0 = plant.xop,
xhat_0 = nothing,
lastu = plant.uop,
progress = true,
x̂_0 = xhat_0
) where {NT<:Real}
model = estim.model
Expand All @@ -282,7 +288,7 @@ function sim_closedloop!(
lastd, lasty = d, evaloutput(plant, d)
initstate!(est_mpc, lastu, lasty[estim.i_ym], lastd)
isnothing(x̂_0) || setstate!(est_mpc, x̂_0)
@progress name="$(nameof(typeof(est_mpc))) simulation" for i=1:N
@progressif progress name="$(nameof(typeof(est_mpc))) simulation" for i=1:N
d = lastd + d_step + d_noise.*randn(plant.nd)
y = evaloutput(plant, d) + y_step + y_noise.*randn(plant.ny)
ym = y[estim.i_ym]
Expand Down
6 changes: 6 additions & 0 deletions test/4_test_plot_sim.jl
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,8 @@
@test res.D_data[:, 1] ≈ model.dop
@test res.X_data[:, 1] ≈ zeros(model.nx)

@test_nowarn sim!(model, 15, progress=false)

res_man = SimResult(model, res.U_data, res.Y_data, res.D_data; X_data=res.X_data)
@test res_man.U_data ≈ res.U_data
@test res_man.Y_data ≈ res.Y_data
Expand Down Expand Up @@ -56,6 +58,8 @@ end
@test res.X_data[:, 1] ≈ zeros(estim.model.nx)
@test res.X̂_data[:, 1] ≈ zeros(estim.nx̂)

@test_nowarn sim!(estim, 15, progress=false)

res_man = SimResult(
estim, res.U_data, res.Y_data, res.D_data;
X_data=res.X_data, X̂_data=res.X̂_data
Expand Down Expand Up @@ -140,6 +144,8 @@ end
@test res.X_data[:, 1] ≈ zeros(mpc1.estim.model.nx)
@test res.X̂_data[:, 1] ≈ zeros(mpc1.estim.nx̂)

@test_nowarn sim!(mpc1, 15, progress=false)

mpc2 = ExplicitMPC(LinModel(sys, Ts, i_d=[3]))
res = sim!(mpc2, 15)
@test isa(res.obj, ExplicitMPC)
Expand Down