Functions: Plant Models

Functions: Plant Models

The SimModel types represents discrete state-space models that can be used to construct StateEstimators and PredictiveControllers, 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.

Info

The nomenclature in this page introduces the model manipulated input $\mathbf{u}$, measured disturbances $\mathbf{d}$, state $\mathbf{x}$ and output $\mathbf{y}$, four vectors of nu, nd, nx and ny elements, respectively. The $\mathbf{z}$ vector combines the elements of $\mathbf{u}$ and $\mathbf{d}$.

SimModel

ModelPredictiveControl.SimModel — Type

Abstract 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

source

LinModel

ModelPredictiveControl.LinModel — Type

LinModel(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{ẋ}(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.

NonLinModel

ModelPredictiveControl.NonLinModel — Type

NonLinModel{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{ẋ}(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) -> ẋ and h(x, d, p) -> y. This syntax is simple and intuitive but it allocates more memory.
Mutating functions (in-place): define them as f!(ẋ, x, u, d, p) -> nothing and h!(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 any Julia objects but use a mutable type if you want to change them later e.g.: a vector.

Tip

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).

Warning

The two functions must be in pure Julia to use the model in NonLinMPC, ExtendedKalmanFilter, MovingHorizonEstimator and linearize.

Set Variable Names

ModelPredictiveControl.setname! — Function

setname!(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

source

Set Operating Points

ModelPredictiveControl.setop! — Function

setop!(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

source

Linearize

ModelPredictiveControl.linearize — Function

linearize(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.

Warning

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.

source

ModelPredictiveControl.linearize! — Function

linearize!(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

source

Differential Equation Solvers

DiffSolver

ModelPredictiveControl.DiffSolver — Type

Abstract supertype of all differential equation solvers.

source

RungeKutta

ModelPredictiveControl.RungeKutta — Type

RungeKutta(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).

source