Functions: Plant Models
The SimModel
types represents discrete state-space models that can be used to construct StateEstimator
s and PredictiveController
s, or as plant simulators by calling evaloutput
and updatestate!
methods on SimModel
instances (to test estimator/controller designs). For time simulations, the states $\mathbf{x}$ are stored inside SimModel
instances. Use setstate!
method to manually modify them.
SimModel
ModelPredictiveControl.SimModel
— TypeAbstract supertype of LinModel
and NonLinModel
types.
(model::SimModel)(d=[]) -> y
Functor allowing callable SimModel
object as an alias for evaloutput
.
Examples
julia> model = NonLinModel((x,u,_,_)->-x + u, (x,_,_)->x .+ 20, 4, 1, 1, 1, solver=nothing);
julia> y = model()
1-element Vector{Float64}:
20.0
LinModel
ModelPredictiveControl.LinModel
— TypeLinModel(sys::StateSpace[, Ts]; i_u=1:size(sys,2), i_d=Int[])
Construct a linear model from state-space model sys
with sampling time Ts
in second.
The system sys
can be continuous or discrete-time (Ts
can be omitted for the latter). For continuous dynamics, its state-space equations are (discrete case in Extended Help):
\[\begin{aligned} \mathbf{ẋ}(t) &= \mathbf{A x}(t) + \mathbf{B z}(t) \\ \mathbf{y}(t) &= \mathbf{C x}(t) + \mathbf{D z}(t) \end{aligned}\]
with the state $\mathbf{x}$ and output $\mathbf{y}$ vectors. The $\mathbf{z}$ vector comprises the manipulated inputs $\mathbf{u}$ and measured disturbances $\mathbf{d}$, in any order. i_u
provides the indices of $\mathbf{z}$ that are manipulated, and i_d
, the measured disturbances. The constructor automatically discretizes continuous systems, resamples discrete ones if Ts ≠ sys.Ts
, computes a new realization if not minimal, and separates the $\mathbf{z}$ terms in two parts (details in Extended Help). The rest of the documentation assumes discrete dynamics since all systems end up in this form.
See also ss
Examples
julia> model1 = LinModel(ss(-0.1, 1.0, 1.0, 0), 2.0) # continuous-time StateSpace
LinModel with a sample time Ts = 2.0 s and:
1 manipulated inputs u
1 states x
1 outputs y
0 measured disturbances d
julia> model2 = LinModel(ss(0.4, 0.2, 0.3, 0, 0.1)) # discrete-time StateSpace
LinModel with a sample time Ts = 0.1 s and:
1 manipulated inputs u
1 states x
1 outputs y
0 measured disturbances d
Extended Help
Extended Help
The state-space equations are similar if sys
is discrete-time:
\[\begin{aligned} \mathbf{x}(k+1) &= \mathbf{A x}(k) + \mathbf{B z}(k) \\ \mathbf{y}(k) &= \mathbf{C x}(k) + \mathbf{D z}(k) \end{aligned}\]
Continuous dynamics are internally discretized using c2d
and :zoh
for manipulated inputs, and :tustin
, for measured disturbances. Lastly, if sys
is discrete and the provided argument Ts ≠ sys.Ts
, the system is resampled by using the aforementioned discretization methods.
Note that the constructor transforms the system to its minimal realization using minreal
for controllability/observability. As a consequence, the final state-space representation may be different from the one provided in sys
. It is also converted into a more practical form ($\mathbf{D_u=0}$ because of the zero-order hold):
\[\begin{aligned} \mathbf{x}(k+1) &= \mathbf{A x}(k) + \mathbf{B_u u}(k) + \mathbf{B_d d}(k) \\ \mathbf{y}(k) &= \mathbf{C x}(k) + \mathbf{D_d d}(k) \end{aligned}\]
Use the syntax LinModel{NT}(A, Bu, C, Bd, Dd, Ts)
to force a specific state-space representation.
LinModel(sys::TransferFunction[, Ts]; i_u=1:size(sys,2), i_d=Int[])
Convert to minimal realization state-space when sys
is a transfer function.
sys
is equal to $\frac{\mathbf{y}(s)}{\mathbf{z}(s)}$ for continuous-time, and $\frac{\mathbf{y}(z)}{\mathbf{z}(z)}$, for discrete-time.
See also tf
Examples
julia> model = LinModel([tf(3, [30, 1]) tf(-2, [5, 1])], 0.5, i_d=[2])
LinModel with a sample time Ts = 0.5 s and:
1 manipulated inputs u
2 states x
1 outputs y
1 measured disturbances d
LinModel(sys::DelayLtiSystem, Ts; i_u=1:size(sys,2), i_d=Int[])
Discretize with zero-order hold when sys
is a continuous system with delays.
The delays must be multiples of the sample time Ts
.
LinModel{NT}(A, Bu, C, Bd, Dd, Ts)
Construct the model from the discrete state-space matrices A, Bu, C, Bd, Dd
directly.
See LinModel(::StateSpace)
Extended Help for the meaning of the matrices. The arguments Bd
and Dd
can be the scalar 0
if there is no measured disturbance. This syntax do not modify the state-space representation provided in argument (minreal
is not called). Care must be taken to ensure that the model is controllable and observable. The optional parameter NT
explicitly set the number type of vectors (default to Float64
).
LinModel(model::NonLinModel; x=model.x0+model.xop, u=model.uop, d=model.dop)
Call linearize(model; x, u, d)
and return the resulting linear model.
NonLinModel
ModelPredictiveControl.NonLinModel
— TypeNonLinModel{NT}(f::Function, h::Function, Ts, nu, nx, ny, nd=0; p=[], solver=RungeKutta(4))
NonLinModel{NT}(f!::Function, h!::Function, Ts, nu, nx, ny, nd=0; p=[], solver=RungeKutta(4))
Construct a nonlinear model from state-space functions f
/f!
and h
/h!
.
Both continuous and discrete-time models are supported. The default arguments assume continuous dynamics. Use solver=nothing
for the discrete case (see Extended Help). The functions are defined as:
\[\begin{aligned} \mathbf{ẋ}(t) &= \mathbf{f}\Big( \mathbf{x}(t), \mathbf{u}(t), \mathbf{d}(t), \mathbf{p} \Big) \\ \mathbf{y}(t) &= \mathbf{h}\Big( \mathbf{x}(t), \mathbf{d}(t), \mathbf{p} \Big) \end{aligned}\]
where $\mathbf{x}$, $\mathbf{y}$, $\mathbf{u}$, $\mathbf{d}$ and $\mathbf{p}$ are respectively the state, output, manipulated input, measured disturbance and parameter vectors. If the dynamics is a function of the time, simply add a measured disturbance defined as $d(t) = t$. The functions can be implemented in two possible ways:
- Non-mutating functions (out-of-place): define them as
f(x, u, d, p) -> ẋ
andh(x, d, p) -> y
. This syntax is simple and intuitive but it allocates more memory. - Mutating functions (in-place): define them as
f!(ẋ, x, u, d, p) -> nothing
andh!(y, x, d, p) -> nothing
. This syntax reduces the allocations and potentially the computational burden as well.
Ts
is the sampling time in second. nu
, nx
, ny
and nd
are the respective number of manipulated inputs, states, outputs and measured disturbances. The keyword argument p
is the parameters of the model passed to the two functions. It can be of any Julia object but use a mutable type if you want to change them later e.g.: a vector.
Replace the d
or p
argument with _
in your functions if not needed (see Examples below).
A 4th order RungeKutta
solver discretizes the differential equations by default. The rest of the documentation assumes discrete dynamics since all models end up in this form. The optional parameter NT
explicitly set the number type of vectors (default to Float64
).
The two functions must be in pure Julia to use the model in NonLinMPC
, ExtendedKalmanFilter
, MovingHorizonEstimator
and linearize
.
See also LinModel
.
Examples
julia> f!(ẋ, x, u, _ , p) = (ẋ .= p*x .+ u; nothing);
julia> h!(y, x, _ , _ ) = (y .= 0.1x; nothing);
julia> model1 = NonLinModel(f!, h!, 5.0, 1, 1, 1, p=-0.2) # continuous dynamics
NonLinModel with a sample time Ts = 5.0 s, RungeKutta solver and:
1 manipulated inputs u
1 states x
1 outputs y
0 measured disturbances d
julia> f(x, u, _ , _ ) = 0.1x + u;
julia> h(x, _ , _ ) = 2x;
julia> model2 = NonLinModel(f, h, 2.0, 1, 1, 1, solver=nothing) # discrete dynamics
NonLinModel with a sample time Ts = 2.0 s, empty solver and:
1 manipulated inputs u
1 states x
1 outputs y
0 measured disturbances d
Extended Help
Extended Help
State-space functions are similar for discrete dynamics:
\[\begin{aligned} \mathbf{x}(k+1) &= \mathbf{f}\Big( \mathbf{x}(k), \mathbf{u}(k), \mathbf{d}(k), \mathbf{p} \Big) \\ \mathbf{y}(k) &= \mathbf{h}\Big( \mathbf{x}(k), \mathbf{d}(k), \mathbf{p} \Big) \end{aligned}\]
with two possible implementations as well:
- Non-mutating functions: define them as
f(x, u, d, p) -> xnext
andh(x, d, p) -> y
. - Mutating functions: define them as
f!(xnext, x, u, d, p) -> nothing
andh!(y, x, d, p) -> nothing
.
Set Variable Names
ModelPredictiveControl.setname!
— Functionsetname!(model::SimModel; u=nothing, y=nothing, d=nothing, x=nothing) -> model
Set the names of model
inputs u
, outputs y
, disturbances d
, and states x
.
The keyword arguments u
, y
, d
, and x
must be vectors of strings. The strings are used in the plotting functions.
Examples
julia> model = setname!(LinModel(tf(3, [10, 1]), 2.0), u=["\$A\$ (%)"], y=["\$T\$ (∘C)"])
LinModel with a sample time Ts = 2.0 s and:
1 manipulated inputs u
1 states x
1 outputs y
0 measured disturbances d
Set Operating Points
ModelPredictiveControl.setop!
— Functionsetop!(model; uop=nothing, yop=nothing, dop=nothing, xop=nothing, fop=nothing) -> model
Set the operating points of model
(both LinModel
and NonLinModel
).
Introducing deviations vectors around manipulated input uop
, model output yop
, measured disturbance dop
, and model state xop
operating points (a.k.a. nominal values):
\[\begin{align*} \mathbf{u_0}(k) &= \mathbf{u}(k) - \mathbf{u_{op}} \\ \mathbf{d_0}(k) &= \mathbf{d}(k) - \mathbf{d_{op}} \\ \mathbf{y_0}(k) &= \mathbf{y}(k) - \mathbf{y_{op}} \\ \mathbf{x_0}(k) &= \mathbf{x}(k) - \mathbf{x_{op}} \\ \end{align*}\]
The state-space description of LinModel
around the operating points is:
\[\begin{aligned} \mathbf{x_0}(k+1) &= \mathbf{A x_0}(k) + \mathbf{B_u u_0}(k) + \mathbf{B_d d_0}(k) + \mathbf{f_{op}} - \mathbf{x_{op}} \\ \mathbf{y_0}(k) &= \mathbf{C x_0}(k) + \mathbf{D_d d_0}(k) \end{aligned}\]
and, for NonLinModel
:
\[\begin{aligned} \mathbf{x_0}(k+1) &= \mathbf{f}\Big(\mathbf{x_0}(k), \mathbf{u_0}(k), \mathbf{d_0}(k), \mathbf{p}\Big) + \mathbf{f_{op}} - \mathbf{x_{op}} \\ \mathbf{y_0}(k) &= \mathbf{h}\Big(\mathbf{x_0}(k), \mathbf{d_0}(k), \mathbf{p}\Big) \end{aligned}\]
The state xop
and the additional fop
operating points are frequently zero e.g.: when model
comes from system identification. The two vectors are internally used by linearize
for non-equilibrium points.
Examples
julia> model = setop!(LinModel(tf(3, [10, 1]), 2.0), uop=[50], yop=[20])
LinModel with a sample time Ts = 2.0 s and:
1 manipulated inputs u
1 states x
1 outputs y
0 measured disturbances d
julia> y = model()
1-element Vector{Float64}:
20.0
Linearize
ModelPredictiveControl.linearize
— Functionlinearize(model::SimModel; x=model.x0+model.xop, u=model.uop, d=model.dop) -> linmodel
Linearize model
at the operating points x
, u
, d
and return the LinModel
.
The arguments x
, u
and d
are the linearization points for the state $\mathbf{x}$, manipulated input $\mathbf{u}$ and measured disturbance $\mathbf{d}$, respectively (not necessarily an equilibrium, details in Extended Help). The Jacobians of $\mathbf{f}$ and $\mathbf{h}$ functions are automatically computed with ForwardDiff.jl
.
See Extended Help if you get an error like: MethodError: no method matching (::var"##")(::Vector{ForwardDiff.Dual})
.
Examples
julia> model = NonLinModel((x,u,_,_)->x.^3 + u, (x,_,_)->x, 0.1, 1, 1, 1, solver=nothing);
julia> linmodel = linearize(model, x=[10.0], u=[0.0]);
julia> linmodel.A
1×1 Matrix{Float64}:
300.0
Extended Help
Extended Help
With the nonlinear state-space model:
\[\begin{aligned} \mathbf{x}(k+1) &= \mathbf{f}\Big(\mathbf{x}(k), \mathbf{u}(k), \mathbf{d}(k), \mathbf{p}\Big) \\ \mathbf{y}(k) &= \mathbf{h}\Big(\mathbf{x}(k), \mathbf{d}(k), \mathbf{p}\Big) \end{aligned}\]
its linearization at the operating point $\mathbf{x_{op}, u_{op}, d_{op}}$ is:
\[\begin{aligned} \mathbf{x_0}(k+1) &≈ \mathbf{A x_0}(k) + \mathbf{B_u u_0}(k) + \mathbf{B_d d_0}(k) + \mathbf{f(x_{op}, u_{op}, d_{op}, p)} - \mathbf{x_{op}} \\ \mathbf{y_0}(k) &≈ \mathbf{C x_0}(k) + \mathbf{D_d d_0}(k) \end{aligned}\]
based on the deviation vectors $\mathbf{x_0, u_0, d_0, y_0}$ introduced in setop!
documentation, and the Jacobian matrices:
\[\begin{aligned} \mathbf{A} &= \left. \frac{∂\mathbf{f(x, u, d, p)}}{∂\mathbf{x}} \right|_{\mathbf{x=x_{op},\, u=u_{op},\, d=d_{op}}} \\ \mathbf{B_u} &= \left. \frac{∂\mathbf{f(x, u, d, p)}}{∂\mathbf{u}} \right|_{\mathbf{x=x_{op},\, u=u_{op},\, d=d_{op}}} \\ \mathbf{B_d} &= \left. \frac{∂\mathbf{f(x, u, d, p)}}{∂\mathbf{d}} \right|_{\mathbf{x=x_{op},\, u=u_{op},\, d=d_{op}}} \\ \mathbf{C} &= \left. \frac{∂\mathbf{h(x, d, p)}}{∂\mathbf{x}} \right|_{\mathbf{x=x_{op},\, d=d_{op}}} \\ \mathbf{D_d} &= \left. \frac{∂\mathbf{h(x, d, p)}}{∂\mathbf{d}} \right|_{\mathbf{x=x_{op},\, d=d_{op}}} \end{aligned}\]
Following setop!
notation, we find:
\[\begin{aligned} \mathbf{f_{op}} &= \mathbf{f(x_{op}, u_{op}, d_{op}, p)} \\ \mathbf{y_{op}} &= \mathbf{h(x_{op}, d_{op}, p)} \end{aligned}\]
Notice that $\mathbf{f_{op} - x_{op} = 0}$ if the point is an equilibrium. The equations are similar if the nonlinear model has nonzero operating points.
Automatic differentiation (AD) allows exact Jacobians. The NonLinModel
f
and h
functions must be compatible with this feature though. See Automatic differentiation for common mistakes when writing these functions.
ModelPredictiveControl.linearize!
— Functionlinearize!(linmodel::LinModel, model::SimModel; <keyword arguments>) -> linmodel
Linearize model
and store the result in linmodel
(in-place).
The keyword arguments are identical to linearize
.
Examples
julia> model = NonLinModel((x,u,_,_)->x.^3 + u, (x,_,_)->x, 0.1, 1, 1, 1, solver=nothing);
julia> linmodel = linearize(model, x=[10.0], u=[0.0]); linmodel.A
1×1 Matrix{Float64}:
300.0
julia> linearize!(linmodel, model, x=[20.0], u=[0.0]); linmodel.A
1×1 Matrix{Float64}:
1200.0
Differential Equation Solvers
DiffSolver
ModelPredictiveControl.DiffSolver
— TypeAbstract supertype of all differential equation solvers.
RungeKutta
ModelPredictiveControl.RungeKutta
— TypeRungeKutta(order=4; supersample=1)
Create a Runge-Kutta solver with optional super-sampling.
Only the 4th order Runge-Kutta is supported for now. The keyword argument supersample
provides the number of internal steps (default to 1 step).