An improved approach of robust constrained model predictive tracking control for polytopic description systems

Abstract

A novel design of robust constrained model predictive tracking control is proposed for systems with polytopic description. Unlike the conventional robust model predictive tracking control, the proposed method adopts an improved state space model in which the process state variables and tracking error are combined such that they can be tuned in the cost function optimization separately. Based on the proposed new model, more degrees of freedom are provided for the subsequent controller design, which leads to improved control performance. The relevant feasibility and robust stability issues are further discussed, and the effectiveness of the proposed approach is tested on the control of a system which is open-loop unstable with dead time and reverse responses.

Introduction

Due to the advantages of convenience in dealing with various constraints, low requirements in model accuracy, and so on, model predictive control (MPC) has acquired a lot of progress in both theory and applications.^1,2 However, with the increasing control requirements, there is also strict requirement for the control performance improvement of MPC strategies.

Disturbances and uncertainties exist in industrial systems inevitably, which cause the control performance of relevant controller to deteriorate to be unacceptable. In order to cope with these situations, many researchers contributed to the development of MPC schemes. A combined strategy of MPC and least squares online parameter estimation for the control of the hypnotic depth, measured by the Bispectral Index, under uncertainty was presented by Krieger and Pistikopoulos. ³ Ding and Pan proposed the robust MPC for a linear polytopic uncertain system with bounded disturbance and immeasurable states. ⁴ Robust MPC for uncertain discrete-time Takagi-Sugeno fuzzy systems with input constraints and persistent disturbances was considered in Yang et al. ⁵ In Ojaghi et al., ⁶ a new robust MPC control strategy in which the linear matrix inequalities are employed to acquire the corresponding feedback control law was developed for uncertain nonlinear systems. A predictive control formulation for uncertain discrete-time non-linear uniformly continuous plant models, where the controller output data are transmitted over an unreliable communication channel, was studied by Quevedo and Nešić. ⁷ In order to cope with the chemical processes with parameter uncertainty, the real time optimization and economic MPC objectives were combined to form a novel control algorithm in Santander et al. ⁸

For the nominal MPC scheme, its closed-loop stability has been solved by imposing artificial constraints in the optimization problem. ⁹ On the basis of the dual problem of an MPC optimization problem, a distributed MPC scheme was proposed for linear systems with local and global constraints by Wang and Ong. ¹⁰ In order to improve the feasibility and optimality of the off-line MPC strategy, Bloemen et al., proposed an on-line linear model-based MPC algorithm. ¹¹

Among the existing MPC strategies, robust MPC is a significant branch. During the development history of robust MPC scheme, there are many representative results. As to the tracking problem in robust MPC subject to state and input constraints, some computationally efficient approaches were presented in Razi and Haeri. ¹² By applying parameter-dependent Lyapunov function (PDLF), Cuzzola, Geromel, and Morari improved the technology presented in Cuzzola et al. ¹³ Inspired by the approach in Ding et al. ⁹ and Bloemen et al. ¹¹ further improved the technique in Cuzzola et al. ¹³ Besides, the application of the robust MPC in the train regulation in underground railway transportation was investigated by Li et al. ¹⁴ In Rhouma and Bouani, ¹⁵ the robust MPC strategy was developed for the fractional-order systems with uncertain parameters.

With regard to the polytopic systems, the applications of robust MPC are also numerous. In Li and Xi, ¹⁶ the feedback MPC and an augmented polytopic uncertainty description were combined to design the robust MPC for the constrained polytopic systems. Adopting the polyhedral invariant sets, an off-line robust MPC was developed for the uncertain polytopic systems by Bumroongsri and Kheawhom. ¹⁷ For the polytopic systems under exogenous disturbance and constraints, the application of robust MPC strategy was discussed in Li et al. ¹⁸ In Muñoz-Carpintero et al., ¹⁹ the robust MPC scheme was presented for the polytopic systems with additive and multiplicative uncertainty. A fast ellipsoidal MPC approach was put forward for the polytopic systems under bounded disturbances by Casavola et al. ²⁰

In Zhang et al., ²¹ a novel state space model structure in which process state variables and tracking errors are combined and regulated separately was introduced to improve the control performance of MPC strategy, and there are some subsequent results.^22–24 Based on the improved state space model, there are more degrees of freedom for the relevant controller design and improved ensemble control performance can be obtained. In this paper, an improved robust constrained model predictive tracking control is proposed. The main contribution of this work include the followings: (1) A new state space model in which state variables and tracking error are combined and regulated separately is proposed. (2) The technique in Bloemen et al. ¹¹ is also employed to formulate an N free control moves plus a linear feedback law based robust constrained MPC strategy, which leads to improved control performance. (3) Based on the proposed model structure and the subsequent control law derivation, the feasibility and robust stability analysis of the proposed scheme are further discussed. The effectiveness of the proposed MPC approach is tested on the control of a system which is open-loop unstable with dead time and reverse responses under uncertainty.

The chapters of this paper are organized as follows. The problem needed to be resolved is introduced in section 2. In section 3, the conventional robust MPC strategy is briefly discussed. The improved state space model based robust MPC scheme is presented in section 4. The numerical examples are discussed in section 5, and the conclusion is given in section 6.

Problem statement

Assume that an industrial system is time varying and a set of models are the candidates for its dynamics. However, it is possible to adopt one of the set of models to design the controller. Without loss of generality, the following nominal model (the kth model in the model set) is adopted to design the controller for the time varying system

$$\begin{array}{c}x(k+1)=A_{k}x(k)+B_{k}u(k)\hfill \\ y(k+1)=\mathrm{Cx}(k+1)\hfill \end{array}$$

(1)

with incremental input and state constraints

$$\begin{array}{c}-y_{min}\le y(k+i+1)\le y_{max}\forall i\ge 0\hfill \\ -\mathrm{\Delta }{u}_{min}\le \mathrm{\Delta }u(k+i)\le \mathrm{\Delta }{u}_{max}\forall i\ge 0\hfill \end{array}$$

(2)

where $y(k),u(k),x(k),\mathrm{\Delta }u(k)$ are the output, control input, system state and incremental input at time instant $k$ , respectively. $A_k,B_k,C$ are the kth relevant system matrix, input matrix and output matrix, respectively. $y_{min},y_{max},\mathrm{\Delta }{u}_{min},\mathrm{\Delta }{u}_{max}$ are the lower bound, the upper bound for the output and the lower bound, the upper bound for incremental input, respectively.

Suppose that $[A_k|B_k]\in \mathrm{\Omega }$ , $\forall k\ge 0$ , where $\mathrm{\Omega }=\mathrm{Co}\{A_1|B_1,A_2|B_2,\cdots ,A_L|B_L\}$ ( $\mathrm{Co}$ denotes set and | represents the delimiter), that is, there exist $L$ nonnegative coefficients ${\alpha }_l(k)$ (with $l=1,2,\cdots ,L$ ) such that

$$\sum _{l=1}^L{\alpha }_l(k)=1,[A_k|B_k]=\sum _{l=1}^L{\alpha }_l(k)[A_l|B_l]$$

(3)

Remark 1. Practical production processes generally have the characteristics of large inertia, long time-delay, and nonlinearity, and it is difficult to get accurate mathematical models. However, for systems with input/output data sets at different operating points, the polytopic system with vertices given by the linear models will apply to the real system. Alternatively, the nonlinear system can be approximated by a polytopic description system. Thus, this polytopic description system is suitable for problems of production process and utilizing polytopic description can incorporate model-plant mismatch. For using polytopic description system, both practical and theoretical have been achieved (Kothare et al. ²⁵ ).

The subsequent robust MPC design is aimed at letting the process output track the reference value as closely as possible, and meanwhile maintaining the desired control performance under various uncertainties. The corresponding cost function is selected as

$$\begin{array}{c}\multicolumn{1}{c}{\underset{\mathrm{\Delta }U(k)}{min}\underset{[A_k|B_k]\in \mathrm{\Omega }}{max}{J}_{\infty }(k)=}\\ \multicolumn{1}{c}{\sum _{i=0}^{\infty }[\Vert y(k+i)-y_r(k+i){\Vert }_Q^2+\Vert \mathrm{\Delta }u(k+i){\Vert }_R^2]}\end{array}$$

(4)

subject to

$$\begin{array}{c}-y_{min}\le y(k+i+1)\le y_{max}\hfill \\ -\mathrm{\Delta }{u}_{min}\le \mathrm{\Delta }u(k+i)\le \mathrm{\Delta }{u}_{max}\hfill \end{array}$$

(5)

where $\mathrm{\Delta }U(k)$ is the set of the future control input increments, that is, $\mathrm{\Delta }u(k),\mathrm{\Delta }u(k+1),\cdots ,\mathrm{\Delta }u(k+\infty )$ . $Q,R$ are the weighting matrices for tracking error and control input increment respectively. $y_r(k)$ is the set point at time instant $k$ .

Here the receding horizon strategy is adopted, that is, $\mathrm{\Delta }u(k)$ is implemented, and the optimization problem in equation (4) will be repeated at time instant $k+1$ .

Traditional robust MPC strategy

Consider the state space model in equation (1), it can be rewritten as follows by adding the difference operator $\mathrm{\Delta }$ .

$$\begin{array}{c}\mathrm{\Delta }x(k+1)=A_k\mathrm{\Delta }x(k)+B_k\mathrm{\Delta }u(k)\hfill \\ \mathrm{\Delta }y(k+1)=C\mathrm{\Delta }x(k+1)\hfill \end{array}$$

(6)

Define the reference value as $y_r(k)$ , then the tracking error can be calculated as

$$e(k)=y(k)-y_r(k)$$

(7)

For the traditional robust model predictive tracking control, we can obtain the formula of $e(k+1)$ further

$$\begin{array}{c}\multicolumn{1}{c}{e(k+1)=e(k)+\mathrm{\Delta }y(k+1)=e(k)+C\mathrm{\Delta }x(k+1)}\\ \multicolumn{1}{c}{=e(k)+C{A}_k\mathrm{\Delta }x(k)+C{B}_k\mathrm{\Delta }u(k)}\end{array}$$

(8)

then the general formula of $e(k+i)$ is also acquired.

$$\begin{array}{c}\multicolumn{1}{c}{e(k+i)=e(k)+(C{A}_k+C{A_k}^2+\cdots +C{A_k}^i)\mathrm{\Delta }x(k)}\\ \multicolumn{1}{c}{+(C{B}_k+C{A}_k{B}_k+\cdots +C{A_k}^{i-1}{B}_k)\mathrm{\Delta }u(k)}\\ +(C{B}_k+C{A}_k{B}_k+\cdots +C{A_k}^{i-2}{B}_k)\mathrm{\Delta }u(k+1)\hfill \\ \multicolumn{1}{c}{+\cdots +C{B}_k\mathrm{\Delta }u(k+i-1)}\end{array}$$

(9)

At the same time, the cost function in equation (4) is equivalent to

$$\underset{\mathrm{\Delta }U(k)}{min}\underset{[A_k|B_k]\in \mathrm{\Omega }}{max}{J}_{\infty }(k)=\sum _{i=0}^{\infty }[\Vert e(k+i){\Vert }_Q^2+\Vert \mathrm{\Delta }u(k+i){\Vert }_R^2]$$

(10)

subject to

$$\begin{array}{c}-y_{min}\le y(k+i+1)\le y_{max}\hfill \\ -\mathrm{\Delta }{u}_{min}\le \mathrm{\Delta }u(k+i)\le \mathrm{\Delta }{u}_{max}\hfill \end{array}$$

(11)

then the control law of traditional robust model predictive tracking control can be obtained by solving the optimization problem in equation (10).

Robust MPC with N free control moves

Model derivation

For the proposed robust MPC strategy, a novel state space model is adopted. Firstly, the formula of the new state space vector that contains the state variables and tracking error is selected as follows

$$z(k)=\left[\begin{array}{c}\mathrm{\Delta }x(k)\hfill \\ \multicolumn{1}{c}{e(k)}\end{array}\right]$$

(12)

then the new formulated state space model is

$$\begin{array}{c}z(k+1)=A_{\mathrm{mk}}z(k)+B_{\mathrm{mk}}\mathrm{\Delta }u(k)\hfill \\ \mathrm{\Delta }y(k+1)=C_{m}z(k+1)\hfill \end{array}$$

(13)

where

$$A_{\mathrm{mk}}=\left[\begin{array}{cc}{A}_k\hfill & 0\hfill \\ \multicolumn{1}{c}{C{A}_k}& 1\hfill \end{array}\right];B_{\mathrm{mk}}=\left[\begin{array}{c}{B}_k\hfill \\ \multicolumn{1}{c}{C{B}_k}\end{array}\right];C_m=\left[\begin{array}{cc}C\hfill & 0\hfill \end{array}\right]$$

Here 0 in $A_{\mathrm{mk}}$ is a zero vector with appropriate dimension.

Note that $[A_{\mathrm{mk}}|B_{\mathrm{mk}}]$ can also be cast into the following polytopic description

$$[A_{\mathrm{mk}}|B_{\mathrm{mk}}]\in {\mathrm{\Omega }}_m=\mathrm{Co}\{[A_{m1}|B_{m1}],\cdots ,[A_{\mathrm{mL}}|B_{\mathrm{mL}}]\}$$

(14)

$$\begin{array}{c}[A_{\mathrm{mk}}|B_{\mathrm{mk}}]=\sum _{l=1}^L{\xi }_l(k)[A_{\mathrm{ml}}|B_{\mathrm{ml}}],\sum _{l=1}^L{\xi }_l(k)=1\hfill \\ \begin{array}{cc}{\xi }_l(k)\ge 0,\hfill & l=1,2,\cdots ,L\hfill \end{array}\hfill \end{array}$$

(15)

Decomposition of the cost function

Consider the following cost function

$$\underset{\mathrm{\Delta }U(k)}{min}\underset{[A_{\mathrm{mk}}|B_{\mathrm{mk}}]\in {\mathrm{\Omega }}_m}{max}{J}_{\infty }(k)=\sum _{i=0}^{\infty }[\Vert z(k+i){\Vert }_Q^2+\Vert \mathrm{\Delta }u(k+i){\Vert }_R^2]$$

(16)

subject to

$$\begin{array}{c}-y_{min}\le y(k+i+1)\le y_{max}\hfill \\ -\mathrm{\Delta }{u}_{min}\le \mathrm{\Delta }u(k+i)\le \mathrm{\Delta }{u}_{max}\hfill \end{array}$$

(17)

where $\mathrm{\Delta }U(k)$ is the set of the future control input increments. $y_{min},y_{max},\mathrm{\Delta }{u}_{min}$ , and $\mathrm{\Delta }{u}_{max}$ are the relevant limits. $Q,R$ are the weighting matrices for the new state space vector and control input increment, respectively.

Note that

$$\begin{array}{c}\multicolumn{1}{c}{y(k+i+1)=y_r(k+i+1)+e(k+i+1)}\\ \multicolumn{1}{c}{=y_r(k+i+1)+C_{e}z(k+i+1)}\end{array}$$

(18)

where $C_e=[\begin{array}{cc}0\hfill & 1\hfill \end{array}]$ , $0$ in $C_e$ is a zero vector with appropriate dimension.

thus the constraints in equation (17) can be transformed into the following formula further

$$\begin{array}{c}-y_{min}-y_r(k+i+1)\le C_{e}z(k+i+1)\hfill \\ \multicolumn{1}{c}{\le y_{max}-y_r(k+i+1)-\mathrm{\Delta }{u}_{min}\le \mathrm{\Delta }u(k+i)\le \mathrm{\Delta }{u}_{max}}\end{array}$$

(19)

Remark 2. Note that the new state space vector contains state variables and tracking error simultaneously, and there are more degrees of freedom for the subsequent controller design because by adjusting the relevant weighting coefficients, the closed-loop system performance will be tuned.

The cost function in equation (16) can be divided into two parts

$$\underset{\mathrm{\Delta }{U}_1(k)}{min}\underset{[A_{\mathrm{mk}}|B_{\mathrm{mk}}]\in {\mathrm{\Omega }}_m}{max}{J}_1(k)=\sum _{i=0}^{N-1}[\Vert z(k+i){\Vert }_Q^2+\Vert \mathrm{\Delta }u(k+i){\Vert }_R^2]$$

(20)

subject to (19) for all $i=0,1,\cdots ,N-1$ .

and

$$\underset{\mathrm{\Delta }{U}_2(k)}{min}\underset{[A_{\mathrm{mk}}|B_{\mathrm{mk}}]\in {\mathrm{\Omega }}_m}{max}{J}_2(k)=\sum _{i=N}^{\infty }[\Vert z(k+i){\Vert }_Q^2+\Vert \mathrm{\Delta }u(k+i){\Vert }_R^2]$$

(21)

subject to (19) for all $i\ge N$ .

where $\mathrm{\Delta }{U}_1(k)$ is the set of $\mathrm{\Delta }u(k),\cdots ,\mathrm{\Delta }u(k+N-1)$ , and $\mathrm{\Delta }{U}_2(k)$ is the set of $\mathrm{\Delta }u(k+N),\cdots ,\mathrm{\Delta }u(k+\infty )$ . $N$ is the switching horizon.

Remark 3. Equation (20) in fact considers the tracking because the tracking error is included in it. Moreover, since the state changes are considered, improved control performance may be anticipated, and the control law is composed of the first N step freedom control moves and linear state feedback control law.

As for the infinite horizon optimization problem in equation (21), the linear state feedback control law is adopted.

$$\begin{array}{c}\mathrm{\Delta }u(k+i)=F(k)z(k+i),i\ge N\hfill \end{array}$$

(22)

Define the following quadratic function

$$\begin{array}{c}V(i,k)=z(k+i{)}^{\mathrm{T}}P(i,k)z(k+i),i\ge N\hfill \end{array}$$

(23)

where $P(i,k)>0$ for $\forall k$ and $i\ge N$ .

meanwhile, suppose that for $\forall [A_{\mathrm{mk}}|B_{\mathrm{mk}}]\in {\mathrm{\Omega }}_m$ and $i\ge N$ , $V(i,k)$ satisfies the following robust stability constraint

$$V(i+1,k)-V(i,k)\le -\Vert z(k+i){\Vert }_Q^2-\Vert \mathrm{\Delta }u(k+i){\Vert }_R^2$$

(24)

Summing equation (24) from $i=N$ to $\infty$ , we can obtain

$$\underset{[A_{\mathrm{mk}}|B_{\mathrm{mk}}]\in {\mathrm{\Omega }}_m,i\ge N}{max}{J}_2(k)\le V(N,k)=z(k+N{)}^{\mathrm{T}}P(N,k)z(k+N)$$

(25)

hence the optimization problem in equation (21) is turned into the minimization of $V(N,k)$ .

Finally, the optimization problem in equation (16) is equivalent to

$$\overline{J}(k)=\sum _{i=0}^{N-1}[\Vert z(k+i){\Vert }_Q^2+\Vert \mathrm{\Delta }u(k+i){\Vert }_R^2]+\Vert z(k+N){\Vert }_{P(N,k)}^2$$

(26)

with respect to $\mathrm{\Delta }{U}_1(k),F(k)$ , and $P(N,k)$ .

Controller design

At first, the prediction of the process state is obtained on the basis of equation (13)

$$\begin{array}{c}\multicolumn{1}{c}{\left[\begin{array}{c}z(k+1)\hfill \\ \multicolumn{1}{c}{z(k+2)}\\ \multicolumn{1}{c}{⋮}\\ \multicolumn{1}{c}{z(k+N)}\end{array}\right]=\left[\begin{array}{c}{A}_{\mathrm{mk}}\hfill \\ \multicolumn{1}{c}{{A_{\mathrm{mk}}}^2}\\ \multicolumn{1}{c}{⋮}\\ \multicolumn{1}{c}{{A_{\mathrm{mk}}}^N}\end{array}\right]z(k)}\\ \multicolumn{1}{c}{+\left[\begin{array}{cccc}{B}_{\mathrm{mk}}\hfill & 0\hfill & \cdots \hfill & 0\hfill \\ \multicolumn{1}{c}{A_{\mathrm{mk}}{B}_{\mathrm{mk}}}& B_{\mathrm{mk}}\hfill & \cdots \hfill & ⋮\hfill \\ \multicolumn{1}{c}{⋮}& ⋮\hfill & \ddots \hfill & 0\hfill \\ \multicolumn{1}{c}{{A_{\mathrm{mk}}}^{N-1}{B}_{\mathrm{mk}}}& \cdots \hfill & \cdots \hfill & B_{\mathrm{mk}}\hfill \end{array}\right]\times \left[\begin{array}{c}\mathrm{\Delta }u(k)\hfill \\ \multicolumn{1}{c}{\mathrm{\Delta }u(k+1)}\\ \multicolumn{1}{c}{⋮}\\ \multicolumn{1}{c}{\mathrm{\Delta }u(k+N-1)}\end{array}\right]}\end{array}$$

(27)

Further, equation (27) can be simplified as

$$\left[\begin{array}{c}\tilde{z}(k+1)\hfill \\ \multicolumn{1}{c}{z(k+N)}\end{array}\right]=\left[\begin{array}{c}{\tilde{A}}_{\mathrm{mk}}\hfill \\ \multicolumn{1}{c}{{\tilde{A}}_{\mathrm{mNk}}}\end{array}\right]z(k)+\left[\begin{array}{c}{\tilde{B}}_{\mathrm{mk}}\hfill \\ \multicolumn{1}{c}{{\tilde{B}}_{\mathrm{mNk}}}\end{array}\right]\mathrm{\Delta }\tilde{u}(k)$$

(28)

where $\tilde{z}(k+1),\mathrm{\Delta }\tilde{u}(k),{\tilde{A}}_{\mathrm{mk}},{\tilde{A}}_{\mathrm{mNk}},{\tilde{B}}_{\mathrm{mk}}$ and ${\tilde{B}}_{\mathrm{mNk}}$ can be acquired easily from equation (27).

then equation (26) can be rewritten as

$$\begin{array}{c}\multicolumn{1}{c}{\overline{J}(k)=\Vert z(k){\Vert }_Q^2+\Vert {\tilde{A}}_{\mathrm{mk}}z(k)+{\tilde{B}}_{\mathrm{mk}}\mathrm{\Delta }\tilde{u}(k){\Vert }_{\tilde{Q}}^2}\\ \multicolumn{1}{c}{+\Vert \mathrm{\Delta }\tilde{u}(k){\Vert }_{\tilde{R}}^2+\Vert {\tilde{A}}_{\mathrm{mNk}}z(k)+{\tilde{B}}_{\mathrm{mNk}}\mathrm{\Delta }\tilde{u}(k){\Vert }_{P(N,k)}^2}\end{array}$$

(29)

where $\tilde{Q}=\mathrm{diag}\{Q,Q,\cdots ,Q\}$ and $\tilde{R}=\mathrm{diag}\{R,R,\cdots ,R.$

If and only if there exist $L$ symmetric positive matrices $P_l$ ( $l=1,2,\cdots ,L$ ) satisfying

$$P(i,k)=\sum _{l=1}^L{\xi }_l(k+i)P_l,\begin{array}{c}i\ge N\hfill \end{array}$$

and

$$\begin{array}{c}[A_{\mathrm{ml}}+B_{\mathrm{ml}}F(k){]}^{\mathrm{T}}{P}_t[A_{\mathrm{ml}}+B_{\mathrm{ml}}F(k)]-P_l+Q\hfill \\ \multicolumn{1}{c}{+F(k{)}^T\mathrm{RF}(k)\le 0,}\\ t=1,2,\cdots ,L;\begin{array}{c}l=1,2,\cdots ,L\hfill \end{array}.\hfill \end{array}$$

(30)

then equations (22) and (24) will be tenable.

Define

$${\gamma }_1\ge \Vert {\tilde{A}}_{\mathrm{mk}}z(k)+{\tilde{B}}_{\mathrm{mk}}\mathrm{\Delta }\tilde{u}(k){\Vert }_{\tilde{Q}}^2+\Vert \mathrm{\Delta }\tilde{u}(k){\Vert }_{\tilde{R}}^2$$

(31)

and

$${\gamma }_2\ge \Vert {\tilde{A}}_{\mathrm{mNk}}z(k)+{\tilde{B}}_{\mathrm{mNk}}\mathrm{\Delta }\tilde{u}(k){\Vert }_{P(N,k)}^2$$

(32)

then the optimization problem in equation (29) is transformed as

$$\underset{{\gamma }_1,{\gamma }_2,\mathrm{\Delta }\tilde{u}(k),F(k),P_l}{min}\Vert z(k){\Vert }_Q^2+{\gamma }_1+{\gamma }_2$$

(33)

subject to equation (19), equations (30)–(32).

Equations (30)–(32) can be converted into linear matrix inequalities (LMIs) by applying the Schur complement.

Define $S_l={\gamma }_2{P_l}^{-1}$ , $l=1,2,\cdots ,L$ , $F(k)=Y{G}^{-1}$ , then equation (30) is equivalent to the following linear matrix inequality (LMI)

$$\begin{array}{c}\left[\begin{array}{cccc}G+G^{\mathrm{T}}-S_l\hfill & *\hfill & *\hfill & *\hfill \\ \multicolumn{1}{c}{A_{\mathrm{ml}}G+B_{\mathrm{ml}}Y}& S_t\hfill & *\hfill & *\hfill \\ \multicolumn{1}{c}{Q^{1/2}G}& 0\hfill & {\gamma }_{2}I\hfill & *\hfill \\ \multicolumn{1}{c}{R^{1/2}Y}& 0\hfill & 0\hfill & {\gamma }_{2}I\hfill \end{array}\right]\ge 0\hfill \\ t=1,2,\cdots ,L;l=1,2,\cdots ,L.\hfill \end{array}$$

(34)

Similar with equation (3), $[{\tilde{A}}_{\mathrm{mk}}|{\tilde{B}}_{\mathrm{mk}}]$ can be cast into the following polytopic description

$$[{\tilde{A}}_{\mathrm{mk}}|{\tilde{B}}_{\mathrm{mk}}]\in {\tilde{\mathrm{\Omega }}}_m=\mathrm{Co}\{[{\tilde{A}}_{m1}|{\tilde{B}}_{m1}],\cdots ,[{\tilde{A}}_{\mathrm{mL}}|{\tilde{B}}_{\mathrm{mL}}]\}$$

(35)

$$\begin{array}{c}[{\tilde{A}}_{\mathrm{mk}}|{\tilde{B}}_{\mathrm{mk}}]=\sum _{l=1}^L{\varsigma }_l(k)[{\tilde{A}}_{\mathrm{ml}}|{\tilde{B}}_{\mathrm{ml}}],\sum _{l=1}^L{\varsigma }_l(k)=1\hfill \\ \begin{array}{cc}{\varsigma }_l(k)\ge 0,\hfill & l=1,2,\cdots ,L\hfill \end{array}\hfill \end{array}$$

(36)

then LMI of equation (31) can be described as

$$\left[\begin{array}{ccc}{\tilde{Q}}^{-1}\hfill & *\hfill & *\hfill \\ \multicolumn{1}{c}{0}& {\tilde{R}}^{-1}\hfill & *\hfill \\ \multicolumn{1}{c}{\mathrm{\Delta }\tilde{u}{(k)}^{\mathrm{T}}{{\tilde{B}}^{\mathrm{T}}}_{\mathrm{ml}}+z{(k)}^{\mathrm{T}}{{\tilde{A}}^{\mathrm{T}}}_{\mathrm{ml}}}& \mathrm{\Delta }\tilde{u}{(k)}^{\mathrm{T}}\hfill & {\gamma }_1\hfill \end{array}\right]\ge 0l=1,2,\cdots ,L.$$

(37)

The same as equations (35)–(37), $[{\tilde{A}}_{\mathrm{mNk}}|{\tilde{B}}_{\mathrm{mNk}}]$ can be cast into the following polytopic description

$$[{\tilde{A}}_{\mathrm{mNk}}|{\tilde{B}}_{\mathrm{mNk}}]\in {\tilde{\mathrm{\Omega }}}_{\mathrm{mN}}=\mathrm{Co}\{[{\tilde{A}}_{\mathrm{mN}1}|{\tilde{B}}_{\mathrm{mN}1}],\cdots ,[{\tilde{A}}_{\mathrm{mNL}}|{\tilde{B}}_{\mathrm{mNL}}]\}$$

(38)

$$\begin{array}{c}[{\tilde{A}}_{\mathrm{mNk}}|{\tilde{B}}_{\mathrm{mNk}}]=\sum _{l=1}^L{\tau }_l(k)[{\tilde{A}}_{\mathrm{mNl}}|{\tilde{B}}_{\mathrm{mNl}}],\sum _{l=1}^L{\tau }_l(k)=1\hfill \\ \begin{array}{cc}{\tau }_l(k)\ge 0,\hfill & l=1,2,\cdots ,L\hfill \end{array}\hfill \end{array}$$

(39)

and equation (32) can be rewritten as

$$\left[\begin{array}{cc}1\hfill & *\hfill \\ \multicolumn{1}{c}{{\tilde{A}}_{\mathrm{mNl}}z(k)+{\tilde{B}}_{\mathrm{mNl}}\mathrm{\Delta }\tilde{u}(k)}& S_l\hfill \end{array}\right]\ge 0l=1,2,\cdots ,L.$$

(40)

Finally, the optimization problem in equation (33) can be transformed as

$$\underset{{\gamma }_1,{\gamma }_2,\mathrm{\Delta }\tilde{u}(k),Y,G,S_l}{min}\Vert z(k){\Vert }_Q^2+{\gamma }_1+{\gamma }_2$$

(41)

subject to equation (19), equation (34), equation (37) and equation (40).

It is clearly that the control inputs before horizon $N$ are parameterized by $\mathrm{\Delta }\tilde{u}(k)$ , then the following constraints are acquired.

$$\begin{array}{c}-y_{min}-y_r(k+i+1)\le C_{e}z(k+i+1)\le y_{max}\hfill \\ \multicolumn{1}{c}{-y_r(k+i+1)-\mathrm{\Delta }{\tilde{u}}_{min}\le \mathrm{\Delta }\tilde{u}(k+i)\le \mathrm{\Delta }{\tilde{u}}_{max}}\\ i=0,1,\cdots ,N-1.\hfill \end{array}$$

(42)

where $\mathrm{\Delta }{\tilde{u}}_{min}$ and $\mathrm{\Delta }{\tilde{u}}_{max}$ are vectors constructed by $\mathrm{\Delta }{u}_{min}$ and $\mathrm{\Delta }{u}_{max}$ properly.

Meanwhile, the control inputs beyond horizon $N$ are parameterized by the feedback control law in equation (22). The following result is obtained referring to Cuzzola et al. ¹³

Lemma 1. Suppose that there exist $L$ symmetric matrices $S_l={P_l}^{-1}>0$ , two values $\{\mathrm{\Phi },\mathrm{\Gamma }\}$ and a set of matrices $\{G,Y\}$ satisfying equation (34) and

$$\left[\begin{array}{cc}\mathrm{\Phi }\hfill & Y\hfill \\ \multicolumn{1}{c}{Y^{\mathrm{T}}}& G+G-S_l\hfill \end{array}\right]\ge 0,\begin{array}{c}\mathrm{\Phi }\le {{\phi }^2}_{inf}\hfill \end{array}$$

(43)

$$\left[\begin{array}{cc}G+G^{-1}-S_l\hfill & *\hfill \\ \multicolumn{1}{c}{C_e(A_{\mathrm{ml}}G+B_{\mathrm{ml}}Y)}& \mathrm{\Gamma }\hfill \end{array}\right]\begin{array}{cc}\ge 0,\hfill & \mathrm{\Gamma }\le {{\psi }^2}_{inf}\hfill \end{array},l=1,2,\cdots ,L.$$

(44)

where ${\phi }_{inf}=min\{\mathrm{\Delta }{u}_{min},\mathrm{\Delta }{u}_{max}\}$ and ${\psi }_{inf}=min\{(y_{min}+y_r(k+i+1)),(y_{max}-y_r(k+i+1))\}$ . Then the state feedback control law in equation (22) satisfies the constraint in equation (19) for all $i\ge N$ , and the stability of the system is guaranteed.

Hence, the whole optimization problem of the proposed robust constrained model predictive tracking control can be described as

$$\underset{{\gamma }_1,{\gamma }_2,\mathrm{\Delta }\tilde{u}(k),Y,G,S_l}{min}\Vert z(k){\Vert }_Q^2+{\gamma }_1+{\gamma }_2$$

(45)

subject to equation (19), equation (34), equation (37), equation (40) and equations (42)–(44).

At time instant $k$ , only the first item of $\mathrm{\Delta }\tilde{u}(k)$ , that is, $\mathrm{\Delta }u(k)$ is implemented, and the optimization problem in equation (45) will be repeated at next cycle.

Remark 4. Note that the infinite horizon MPC problem with input, and output constraints can be formulated as a convex optimization problem involving LMIs, and this on-line optimization can be solved in polynomial time, which means that it has low computational complexity.

Robust stability analysis

Theorem 1. Once a feasible solution of the proposed robust MPC strategy whose objective function can be expressed as equation (45) is found, then exponential closed-loop stability of the system in equation (13) can be guaranteed by the receding horizon implementation of the robust MPC scheme.

Proof. Denote the feasible solution at time instant $k$ as $\{\mathrm{\Delta }{\tilde{u}}^{*}(k),F^{*}(k),P^{*}(i,k)\}$ , where $P^{*}(i,k)=\sum _{l=1}^L{\xi }_l(k+i)P_l^{*}$ with ${\xi }_l(k+i)$ , $i\ge N$ usually unknown. According to Lemma 1, the following solution is feasible at time instant $k+1$ .

$$\begin{array}{cc}\mathrm{\Delta }u(k+i+1|k+1)=\mathrm{\Delta }{u}^{*}(k+i+1|k),\hfill & i=0,\cdots ,N-2\hfill \end{array}.$$

(46)

$$\begin{array}{cc}\mathrm{\Delta }u(k+i+N|k+1)=F^{*}(k)z^{*}(k+i+N|k),\hfill & i\ge 0\hfill \end{array}$$

(47)

where $\mathrm{\Delta }u(k+i+1|k+1)$ is the future control input increment from time instant $k+1$ . $z^{*}(k+i+N|k)$ , $i\ge 0$ are usually unknown.

Combing equation (26) and equations (46)–(47), the following formula is obtained.

$$\begin{array}{l}\overline{J}(k)=\sum _{i=0}^{N-1}[{\Vert z(k+i+1|k+1)\Vert }_Q^2+{\Vert \mathrm{\Delta }u(k+i+1|k+1)\Vert }_R^2]\\ +{\Vert z(k+N+1|k+1)\Vert }_{P(N,k+1)}^2\\ =\sum _{i=1}^{N-1}[{\Vert z^{\ast }(k+i|k)\Vert }_Q^2+{\Vert \mathrm{\Delta }{u}^{\ast }(k+i|k)\Vert }_R^2]\\ +{\Vert z^{\ast }(k+N|k)\Vert }_Q^2+{\Vert \mathrm{\Delta }{u}^{\ast }(k+N|k)\Vert }_R^2\\ +{\Vert z^{\ast }(k+N+1|k)\Vert }_{P^{\ast }(N+1,k)}^2\end{array}$$

(48)

where $\mathrm{\Delta }{u}^{*}(k+N|k)=F^{*}(k)z^{*}(k+N|k)$ .

Based on equations (23)–(24), the following holds.

$$\begin{array}{c}\multicolumn{1}{c}{\Vert z^{*}(k+N+1|k){\Vert }_{P^{*}(N+1,k)}^2\le \Vert z^{*}(k+N|k){\Vert }_{P^{*}(N,k)}^2}\\ \multicolumn{1}{c}{-\Vert z^{*}(k+N|k){\Vert }_Q^2-\Vert \mathrm{\Delta }{u}^{*}(k+N|k){\Vert }_R^2}\end{array}$$

(49)

Synthesizing equations (48)–(49), we can obtain

$$\begin{array}{c}\multicolumn{1}{c}{\overline{J}(k+1)\le \sum _{i=1}^{N-1}[\Vert z^{*}(k+i|k){\Vert }_Q^2+\Vert \mathrm{\Delta }{u}^{*}(k+i|k){\Vert }_R^2]}\\ \multicolumn{1}{c}{+\Vert z^{*}(k+N|k){\Vert }_{P^{*}(N,k)}^2\le {\overline{J}}^{*}(k)}\end{array}$$

(50)

At time instant $k+1$ , when $\overline{J}(k+1)$ is optimized, the optimum ${\overline{J}}^{*}(k+1)\le \overline{J}(k+1)$ . Thus, ${\overline{J}}^{*}(k)$ is the Lyapunov function.

This completes the proof.

Numerical example

Consider a nominal system which is open-loop unstable with dead time and inverse-response

$$\begin{array}{c}x(k+1)=\mathrm{Ax}(k)+\mathrm{Bu}(k-2)\hfill \\ y(k+1)=\mathrm{Cx}(k+1)\hfill \end{array}$$

(51)

where

$$A=\left[\begin{array}{cc}1.1053\hfill & 0\hfill \\ \multicolumn{1}{c}{-0.03}& 0.8186\hfill \end{array}\right];B=\left[\begin{array}{c}1\hfill \\ \multicolumn{1}{c}{0.0858}\end{array}\right];C=\left[\begin{array}{cc}0\hfill & 1\hfill \end{array}\right]$$

(52)

Since uncertainties exist inevitably in practice, and it is more meaningful to investigate the control performance of the proposed robust MPC strategy under uncertainties. In addition, equation (51) serves as the nominal model for the controller design of the industrial system, while in fact, the real industrial system dynamics will not be as perfect as the dynamics shown in equation (51). In view of this and to mimic this situation, two cases of the real industrial system, that is, model/plant mismatch cases, are generated through Monte Carlo algorithm as

• Case 1:

$$A=\left[\begin{array}{cc}1.1526\hfill & 0\hfill \\ \multicolumn{1}{c}{-0.0318}& 0.8471\hfill \end{array}\right];B=\left[\begin{array}{c}0.91\hfill \\ \multicolumn{1}{c}{0.0834}\end{array}\right];C=\left[\begin{array}{cc}0\hfill & 1\hfill \end{array}\right]$$

(53)

• Case 2:

$$A=\left[\begin{array}{cc}1.0041\hfill & 0\hfill \\ \multicolumn{1}{c}{-0.0285}& 0.7931\hfill \end{array}\right];B=\left[\begin{array}{c}1.06\hfill \\ \multicolumn{1}{c}{0.0879}\end{array}\right];C=\left[\begin{array}{cc}0\hfill & 1\hfill \end{array}\right]$$

(54)

Now, the controller will be designed using the model shown in equation (51) to control the real industrial systems shown in equations (53) and (54). Meanwhile, the traditional robust MPC strategy is introduced as the comparison. The set point of the controlled system is 1, and the control parameters of both MPC strategies are listed in Table 1.

Table 1.

Tuning parameters for both methods.

Parameters	Proposed	Conventional
$N$	20	20
$Q$	$\mathrm{diag}(0.5,0.1,0.01,0.01,1)$	1
$R$	0.01	0.01

In equation (16), Q, R are the weighting matrix for the state variables and control inputs respectively. The output error and the state variables can both be regulated by tuning Q, so more degrees of freedom are provided and improved control performance can be achieved. Thus, the control parameters for the two MPCs are the same with N = 20, R = 0.01, and the output error weighting factors are both 1. For the proposed MPC the weighting matrix is $Q=\mathrm{diag}(0.5,0.1,0.01,0.01,1)$ and for conventional MPC it can only adopt $Q=1$ for the output error weighting.

The relevant constraints are selected as

$$\{\begin{array}{c}-1\le \mathrm{\Delta }u(k)\le 1\hfill \\ \multicolumn{1}{c}{-1.5\le y(k+1)\le 1.5}\end{array}$$

(55)

Figures 1 and 2 show the closed-loop responses under the two cases. In Figure 1, both methods achieve the set-point tracking successfully. The traditional approach shows bigger overshoot and oscillation. From an overall perspective, the proposed method provides better ensemble control performance with better control efforts. In Figure 2, set-point tracking under constraints was still achieved for the proposed methods, and it is obvious that responses of the proposed strategy are smoother with smaller overshoot, but the oscillation of the traditional approach is drastic. It is easy to find that the incremental control input is kept in the range of the given constraints, which verified the superiority of the proposed approach further. In a word, the proposed method shows better control performance.

Figure 1.

Closed-loop responses under case 1.

Figure 2.

Closed-loop responses under case 2.

Conclusion

An improved robust constrained model predictive tracking control is proposed for the systems described by polytopic description in this paper. In the framework of the proposed method, a novel state space model in which state variables and tracking error are combined and regulated separately is adopted such that more degrees of freedom are provided for the subsequent controller design. The effectiveness of the proposed strategy is tested on a open-loop unstable system with dead time and reverse responses.