10 minute read

1. Predictive equations

This section derives an expression for the predicted state trajectory $x(k)$ in terms of the predicted input trajectory $u(k)$.

This expression allows the general quadratic cost:

$$ J(k) = \sum^{N-1}_{i=0} \left[ x^T(k + i \vert k) Q x(k + i \vert k) + u^T(k + i \vert k) R u(k + i \vert k) \right] + x^T(k + N \vert k) \bar{Q} x(k + N \vert k) $$

to be written as a quadratic function of $u(k)$ in following form.

$$ J(k) = \mathbf{u}^T(k)H\mathbf{u}(k) + 2f^T \mathbf{u}(k) + g $$

The predicted state sequence generated by the linear state-space model with input sequence $\mathbf{u}(k)$ can be written.

$$ x(k + i \vert k) = A x(k + i \vert k) + B u(k + i \vert k), \;\;\; i = 0, 1, \dots $$
$$ \begin{align*} x(k \vert k) &= x(k) \\ x(k + 1 \vert k) &= Ax(k) + Bu(k \vert k) \\ x(k + 2 \vert k) &= A^2x(k) + ABu(k \vert k) + Bu(k + 1 \vert k) \end{align*} $$

In conpact notation:

$$ x(k + i \vert k) = A^i x(k) + C_i \mathbf{u}(k), \;\;\; i = 0, \dots, N $$

or

$$ \mathbf{x}(k) = \mathcal{M} x(k) + \mathcal{C} \mathbf{u}(k), \;\;\; \text{where} \;\;\mathcal{M} = \begin{bmatrix} A \\ A^2 \\ \vdots \\ A^N \end{bmatrix}, \;\; \mathcal{C} = \begin{bmatrix} B & 0 & \dots & 0 \\ AB & B & \dots & 0 \\ \vdots & \vdots & \ddots \\ A^{N-1} & A^{N-2} & \dots & B \\ \end{bmatrix} $$

Substituting for $x(k + i \vert k)$ and collecting terms gives:

$$ J(k) = \mathbf{u}^T H \mathbf{u} + 2x^T(k) F^T \mathbf{u}(k) + x^T(k) G x(k) $$

where

$$ \begin{align*} H &= \mathcal{C}^T \tilde{Q} \mathcal{C} + \tilde{R} \\ F &= \mathcal{C}^T \tilde{Q} \mathcal{M} \\ G &= \mathcal{M}^T \tilde{Q} \mathcal{M} + Q \\ \end{align*} $$

with

$$ \tilde{Q} = \begin{bmatrix} Q & 0 & \dots & 0 \\ 0 & \ddots & & \vdots \\ \vdots & & Q & 0 \\ 0 & \dots & 0 & \bar{Q} \\ \end{bmatrix}, \;\;\; \tilde{R} = \begin{bmatrix} R & 0 & \dots & 0 \\ 0 & \ddots & & \vdots \\ \vdots & & R & 0 \\ 0 & \dots & 0 & R \\ \end{bmatrix} $$

Note that matrices $H, F$ and $G$ can be computed offline.

1.1 LTV prediction models

The above formulation applies to systems with linear time-varying models

\[x(k+1) = A(k) x(k) + B(k) u(k).\]

In this case the state predictions can be written

\[x(k + i \vert k) = \prod_{j=i-1}^0 A(k + j) x(k) + C_i (k) \mathbb{u}(k) \;\;\; i = 0, \dots, N\]

where

\[C_0(k) = 0 \\ C_i(k) = \begin{bmatrix} \prod_{j=i-1}^1 A(k + j)B(k) & \prod_{j=i-1}^2 A(k + j)B(k) & 0 & \dots & 0 \end{bmatrix}\]

2. Unconstrainted optimization

In the absense of constraints, the optimization $\mathbb{u}^*(k) = \underset{u}{\mathrm{argmin}} J(k)$ has a closed-form soluation which can be derived by considering the graident of $J$ with respect ot $\mathbb{u}$:

\[\nabla_{\mathbb{u}} J = 2 H \mathbb{u} + 2 F x.\]

If

  • $H$ is positive definite (or positive semidefinite)
  • $H$ is nonsingular

from $\nabla_{\mathbb{u}} \mathbf{J} = 0$

\[\mathbb{u}^*(k) = - H^{-1} F x(k)\]

If $H$ is singular, then optimal $\mathbf{u}^*$ is non-unique.

…

2.1 Horizon length and performance

The linear feedback form $\mathbb{u}^(k) = - H^{-1} F x(k)$ is to be expected since $\mathbf{u}^$ is the solution of a LQ-optimal control problem.

  • (Ideal) The Infinite horizon LQ-problem control problem : there is no difference between the optimal predicted input sequence and its receding horizon implementation in the absence disturbances and model error.
  • (Real) The horizon is reduced : There can be significant fidvtrpsnvy between the prediction s of $\mathbb{u}^*(k) = - H^{-1} F x(k)$ and closed-loop responses with the receding horizon controller

3. Infinite horizon cost

Because of the difference between predicted and closed-loop responsed, there is no guarantee that a receding horizon controller based on a finite-horizon cost will achieve the optimal predicted preformace in closed-loop operation.

In fact the closed-loop system may even be unstable.

However, this problem is avioded entirely if performance is evaluated over an infinite prediction horizon, i.e. if the cost is defined as

\[J(k) = \sum_{i=0}^{\infty} \left[ x^T (k + i \vert k) Q x(k + i \vert k) + u^T(u + i \vert k R) u(k + i \vert k) \right]\]

To make sure that the problem of minimizing ihis cost is tractable, it is then necessary to define the predicted input sequence over the infinite prediction horizon in such a way that the number of free vaiables in the MPC optimization remains finite.

Dual model predictions

\[u(k + i \vert k) = \begin{cases} \text{optimization vaiable} \;\;\; i = 0, 1, \dots, N-1 \;\;\; (\text{mode 1}) \\ K x(k + i \vert k) \;\;\; i = N, N + 1, \dots \;\;\; (\text{mode 2}) \\ \end{cases}\]
  • mode 1 : initial horizon of $N$ samples over which the predicted inputs are variables in the MPC optimization
  • mode 2 : inputs are defined by a stabilizing feedback law over the remaining infinite horizon

For the dual mode input predictions, an infinite horizon cost nedd only be evaluated explicitly over mode 1 since $J(k)$ can be re-written in the form of

\[J(k) = \sum^{N-1}_{i=0} \left[ x^T(k + i \vert k) Q x(k + i \vert k) + u^T(k + i \vert k) R u(k + i \vert k) \right] + x^T(k + N \vert k) \bar{Q} x(k + N \vert k)\]

This is done by choosing the terminal wighting matrix $\bar{Q}$ so that the term $x^T(k + N \vert k) \bar{Q} x(k + N \vert k)$ is equal to the cost over the mode 2 prediction horizon, which is achieved by specifying $\bar{Q}$ as the solution of the Lyapunov equation:

\[\bar{Q} - (A + BK)^T \bar{Q}(A + BK) = Q + K^T R K\]

Theorem 2.1 (Terminal weighting matrix).

Along closed-loop trajectories of the model under the feedback law $u(k) = Kx(k)$, the infinite horizon quadratic cost is given by

\[\sum_{i=0}^{\infty} \left[ x^T(i)Qx(i) + u^T R u(i) \right] = x^T \bar{Q} x(0)\]

where $\bar{Q}$ satisfies $\bar{Q} - (A + BK)^T \bar{Q}(A + BK) = Q + K^T R K$

증명 생략

Note

  1. The Lyapunov equation $\bar{Q} - (A + BK)^T \bar{Q}(A + BK) = Q + K^T R K$ has a unique solution for $\bar{Q}$ if and only if the eigenvalues of $A + BK$ lie within the unit circle.
  2. It is easy to show that $\bar{Q}$ is positive definite if either $Q + K^T R K > 0$ or $Q = C^TC$ where (A + BK, C) is observable.
  3. From theorem 2.1, it is clear that the finite horizon cost fomulation is equal to the infinite horizon formulation.
\[\sum^{N-1}_{i=0} \left[ x^T(k + i \vert k) Q x(k + i \vert k) + u^T(k + i \vert k) R u(k + i \vert k) \right] + x^T(k + N \vert k) \bar{Q} x(k + N \vert k) \\ = \sum_{i=0}^{\infty} \left[ x^T (k + i \vert k) Q x(k + i \vert k) + u^T(u + i \vert k R) u(k + i \vert k) \right]\]

3.1 The relationship between unconstrainted MPC and LQ optimal control

The obvious choice for the mode 2 feedback gain $K$ is the LQ-optimal gain for the cost.

Due to the optimality of predictions over both modes 1 and 2, the optimal predicted trajectory $\mathbb{u}^*$ is then necessarily identical to the infinite horizon optimal input sequence

\[\mathbb{u}^*(k) = \begin{bmatrix} K \\ K ( A + BK) \\ \vdots \\ K ( A + B K)^{N-1} \\ \end{bmatrix} x(k)\]

implying the the receding horizon control law

\[u(k) = u^*(k \vert k) = K_N x(k), \;\;\;, K_N = - \left[ I_{n_u} 0 \dots o \right] H^{-1}F\]

is in fact equal to the LQ-optimal feedback law $u = Kx$.

This result should bot be surprising since the model and cost are the same for MPC and LQ-optimal control, here MPC simply provides an alternative method determining the optimal control law

4. Incorporating constraints

The real advantage of MPC lies in its ability ot determine nonlinear feedback laws which are optimal for constrained systems through numerical calculations that are performed online (i.e. in between samples)

\[\underline{u} \leq u(k) \leq \overline{u} \\ \underline{x} \leq x(k) \leq \overline{x} \\\]

or

The input constraint are equivalent

\[u(k) \leq \overline{u}, \;\;\; -u(k) \leq - \underline{u}\]

As applying mode 1 prediction $u(k + i \vert k), i = 0, \dots N-1$ can therefore be expressed in terms of $u(k)$ as

\[\begin{bmatrix} I \\ -I \\ \end{bmatrix} u(k) \leq \begin{bmatrix} \mathbf{1} \overline{u} \\ -\mathbf{1} \underline{u} \\ \end{bmatrix}\]

where $\mathbf{1}$ is a vector ones for the single input case.

Similary

\[\begin{bmatrix} C_i \\ -C_i \end{bmatrix} u(k) \leq \begin{bmatrix} \overline{x} \\ \underline{x} \\ \end{bmatrix} + \begin{bmatrix} -A^i \\ A^i \end{bmatrix} x(k) \;\;\; i = 1, \dots, N\]

Therefore

\[A_c \mathbf{u}(k) \leq b_0 + B_x x(k)\]

where $A_c, b_0, B_x$ are constant matrices that can be determined offline.

5. Quadratic Programming

Combining the objective function and constraints dervied above, the optimization of the infinite horizon cost subject to constraints requites the solution of a QP problem:

\[\underset{u}{\mathrm{minimize}} \;\;\; \mathbf{u}^T H \mathbf{u} + 2x^T(k) F^T \mathbf{u} \\ \mathrm{subject \; to} \;\;\; A_c \mathbf{u} \leq b_0 + B_x x(k)\]

Since $H$ is positive semi-definite and constraint are linear, this is a convex optimization problem which has a unique solution.

Purpose of this section

  • outlines the 2 types of algorithm commonly used to solve QP problem : active set and interior point algorithms
  • general result of constrainted optimization theory

Theorem 2.2 (Optimization with equality constraints)

If $\mathbf{u}^*$ satisfies

\[\mathbf{u}^* = \underset{u}{\mathrm{argmin}} f(\mathbf{u}) \\ \mathrm{subject \; to} \;\;\; c_i(\mathbf{u}) = 0, \;\;\; i = 1, \dots, m \\\]

where $f$ and $c_i$ are smooth functions, then scalar $\lambda_i^*, i = 1, \dots, m$ exist satisfying

\[\nabla_{\mathbf{u}} f(\mathbf{u}^*) + \sum_{i=1}^m \lambda_i \nabla_{\mathbf{u}} c_i (\mathbf{u}^*) = 0 \\ c_i (\mathbf{u}^*) = 0, \;\;\; i = 1, \dots, m\]
  • Condition 1 : extension of the gradient condition, $\nabla_{\mathbf{u}} f(\mathbf{u}^)$, which must be satisfied at a minimum point of $f(u)$ if $\mathbf{u}$ is unconstrained. The addition of the second term simply ensures that $f(\mathbf{u})$ cannot be reduced by perturbing $\mathbf{u}^$ by incremental distance in any direction for which the constraint $c_i(\mathbf{u}) = 0$ remains satisfied.

Scalar $\lambda_i$ : Lagrange multipliers

Note that

  1. Condition 1 applies to the case of inequality constraints $c_i (\mathbf{u}) \leq 0$.
  2. Problems with inequality constraints are generally harder to solve than similar equality constrainted problems since only a subset of the constraints may be active at the solution.

5.1 Active set algorithms

$a_i^T$ : $i$th row of $A_c$

$b_i$ : $i$th element of $b_0 + B_x x(k), \;\;\; i = 1, \dots, ,$

individual constraint $a_i^T \mathbf{u} \leq b_i$ is active at the solution if $a_i^T \mathbf{u}^* = b_i$, and inactive if $a_i^T \mathbf{u}^* < b_i$.

Clearly the inactive constraints can be removed from the problem without affecting the solution, and it follow that $\mathbf{u}^*$ is also the solution of the equality constrained problem:

\[\underset{u}{\mathrm{minimize}} \;\;\; \mathbf{u}^T H \mathbf{u} + 2x^T(k) F^T \mathbf{u} \\ \mathrm{subject \; to} \;\;\; a_i^T \mathbf{u} = b_i, \;\;\; i \in \mathcal{A}^*\]

where

\[\mathcal{A}^* = \left[ i : a_i^T \mathbf{u}^* = b_i \right]\]

is the set of active constraints at the solution.

The solution of

\[\underset{u}{\mathrm{minimize}} \;\;\; \mathbf{u}^T H \mathbf{u} + 2x^T(k) F^T \mathbf{u} \\ \mathrm{subject \; to} \;\;\; A_c \mathbf{u} \leq b_0 + B_x x(k)\]

can be found by iteratively

  1. selecting a possible combination of active constraints
  2. solving the corresponding equality constrainted problem (thm 2.2)
  3. testing the optimality of the solution ot the closest point to the solution at which all of the inequality constraint satisfied
  • The optimality of a trial solution can be determined from the associated Lagrange multifpliers
  • The successive active sets are chosed so as to make the objective function decrease at each iteratino
  • Wll-design active set solvers manage to avoid testing large numers of possible combinations of active constraints

5.2 Interior point QP algorithms

This approach solves an unconstrainted problem:

\[\underset{u}{\mathrm{minimize}} \;\;\; \mu \left[ \mathbf{u}^T H \mathbf{u} + 2 F^T \mathbf{u} \right] + \phi (\mathbf{u})\]

at each iteration, where

  • $\mu$ is a scalar, and
  • $\phi (\mathbf{u})$ is a barrier function which is finite whenever $\mathbf{u}$ satisfies constraints but which tends to infinity as $\mathbf{u}$ approaches a constraint doundary

For constraints $a_i^T \mathbf{u} \leq b_i$, the barrier function is typically defined as

\[\phi (\mathbf{u}) = \sum_{i=1}^{m} - \log (b_i - a_i^T \mathbf{u})\]

The solution of unconstrainted problem for any given $\mu$ satisfies constraints.

It can be shown the solution of unconstrainted problem tends to $\mathbf{u}^*$ as $\mu \rightarrow \infty$.

Interior point methods increase $\mu$ over successive iterations until constraint are met to within a given tolerance.

The interior point approach has a lower computational load for large-scale problems involving hendreds of optimization vairblaes.

Unlike active set solvers, it is not possible to initialize the iteration at an estimate close to the actual solution because the corresponding value of $\mu$ is unknown.

This can be a big disadvantage in predictive control, where a good estimate of the current optimal control sequence can usually be determined from the solution of the opimization problem computed at the preivous sampling instant.

You may also enjoy