Synthesis of model predictive control based on data-driven learning

TL;DR

This paper introduces a data-driven model predictive control method that iteratively learns system dynamics without prior models, enabling adaptive control of complex, time-varying, or nonlinear systems.

Contribution

It presents a novel MPC approach that uses data-driven learning to approximate system parameters and adapt to changing dynamics without requiring system matrices beforehand.

Findings

Successfully predicts and optimizes future system behaviors.

Adapts to time-varying and nonlinear plant dynamics.

Provides a data-driven identification tool for complex systems.

Abstract

For the application of MPC design in on-line regulation or tracking control problems, several studies have attempted to develop an accurate model, and realize adequate uncertainty description of linear or non-linear plants of the processes. In this study, we employ the data-driven learning technique to iteratively approximate the dynamical parameters, without requiring a priori knowledge of system matrices. The proposed MPC approach can predict and optimize the future behaviors using multiorder derivatives of control input as decision variables. Because the proposed algorithm can obtain a linear system model at each sampling, it can adapt to the actual dynamics of time-varying or nonlinear plants. This methodology can serve as a data-driven identification tool to study adaptive optimal control problems for unknown complex systems.

Equations

$$\dot{x}(t)=$$

$$\frac{\text{buildrel}△}{=}$$

$$J(x(t_k),{\hat{u}}_k(\cdot )))=$$

$$+\Phi (e(t_k+T)),$$

$$\mathcal{F}(t)=\mathcal{G}(t)\hat{\Theta },\forall t\in {\mathbb{R}}_{\ge 0},$$

$$\mathcal{F}(t)=$$

$$\mathcal{G}(t)=$$

$${\Xi }_x(t)=$$

$${\Xi }_u(t)=$$

$$A=\left[\begin{array}{cc}{0}_{q\times q}& I_q\\ A_1& A_2\end{array}\right],B=\left[\begin{array}{c}{0}_{q\times m}\\ B_1\end{array}\right].$$

$${\mathcal{F}}_1(t)={\mathcal{G}}_1(t){\hat{\Theta }}_1,\forall t\in {\mathbb{R}}_{\ge 0},$$

$${\mathcal{F}}_1(t)=$$

$${\mathcal{G}}_1(t)=$$

$$({\Xi }_u^1(t)\otimes I_n{)}^{\text{T}}],$$

$${\Xi }_p(t)=$$

$${\Xi }_v(t)=$$

$${\Xi }_u^1(t)=$$

$$\frac{{\Gamma }_k\text{buildrel}△}{=}$$

$$\frac{{\Psi }_k\text{buildrel}△}{=}$$

$${\Gamma }_k={\Psi }_k\hat{\Theta },$$

$$\hat{\Theta }={\left({\Psi }_k^{\mathrm{T}}{\Psi }_k\right)}^{-1}{\Psi }_k^{\mathrm{T}}{\Gamma }_k.$$

$$\displaystyle{\mathbf{u}}_{k}(t)\buildrel\vartriangle\over{=}\Big{[}\hat{u}_{k}^{\textrm{T}}(t),(\hat{u}_{k}^{[1]})^{\textrm{T}}(t),\ldots,(\hat{u}_{k}^{[r]})^{\textrm{T}}(t)\Big{]}^{\textrm{T}},$$

$$\frac{{\mathcal{A}}_1\text{buildrel}△}{=}$$

$$\frac{{\mathcal{A}}_2\text{buildrel}△}{=}$$

$$\frac{\mathcal{B}\text{buildrel}△}{=}$$

$${\mathbf{u}}_k^{\star }(\cdot )=\underset{{\mathbf{u}}_k(\cdot )}{\mathrm{argmin}}J(x(t_k),{\mathbf{u}}_k(\cdot ))$$

$$\text{s.t.}x(t_k+\tau )=\left[\begin{array}{cc}{T}_1(\tau )& T_2(\tau )\end{array}\right]\left[\begin{array}{c}{X}_1\\ X_2\end{array}\right],\tau \in [0,T]$$

$$X_1={\mathcal{A}}_{1}x,X_2={\mathcal{A}}_{2}x+\mathcal{B}{\mathbf{u}}_k,$$

$$\hat{x}(t_k)=x(t_k),{\hat{u}}_k(t)\in \mathcal{U},$$

$$J(x(t_k),$$

$$+{\tilde{X}}_2^{\text{T}}{\mathcal{T}}_{2,2}{\tilde{X}}_2+{{\mathbf{u}}_k}^{\text{T}}\mathcal{T}{\mathbf{u}}_k+\Phi ({\tilde{X}}_i(t_k+T)),$$

$$\displaystyle w(t)=\mathcal{H}\left(x(t),\hat{u}_{k}^{\star}(t)\right)\Big{(}\hat{\Theta}-\Theta\Big{)},$$

Full text

\ArticleType

LETTER \Year2018 \Month \Vol61 \No \DOI \ArtNo \ReceiveDate \ReviseDate \AcceptDate \OnlineDate

Synthesis of model predictive control based on data-driven learning

dwli@sjtu.edu.cn

\AuthorMark

Zhou Y

\AuthorCitation

Zhou Y, Li D, Xi Y, Gan Z

Synthesis of model predictive control based on data-driven learning

Yuanqiang ZHOU

  

Dewei LI

  

Yugeng XI

  

Zhongxue GAN

ENN Science & Technology Development Co.,Ltd, Langfang, Hebei, 065001, China

Department of Automation, Shanghai Jiao Tong University, Shanghai 200240, China

Key Laboratory of System Control and Information Processing, Ministry of Education of China, Shanghai 200240, China

\deareditor

Model predictive control (MPC) is a practically effective and attractive approach in the field of industrial processes [1] owing to its excellent ability to handle constraints, nonlinearity, and performance/cost trade-offs. The core of all model-based predictive algorithms is to use “open-loop optimal control” instead of “closed-loop optimal control” within a moving horizon [2]. It is assumed in this letter that the reader is familiar with MPC as a control design methodology.

Because the dynamic model of a system predicting its evolution is usually inaccurate, the actual behaviors may deviate significantly from the predicted ones. Thus, acquiring accurate knowledge of the physical model is essential to ensure satisfactory performance of MPC controllers. Owing to the well-developed information technology, copious amounts of measurable process data can be easily collected, and such data can then be employed to predict and assess system behaviors and make control decisions, especially for the establishment and development of learning MPC.

For the application of MPC design in on-line regulation or tracking control problems, several studies have attempted to develop an accurate model, and realize adequate uncertainty description of linear or non-linear plants of the processes [3, 4, 5]. In this work, we employ the data-driven learning technique specified in [6] to iteratively approximate the dynamical parameters, without requiring a priori knowledge of system matrices. The proposed MPC approach can predict and optimize the future behaviors using multi-order derivatives of control input as decision variables. Because the proposed algorithm can obtain a linear system model at each sampling, it can adapt to the actual dynamics of time-varying or nonlinear plants. This methodology can serve as a data-driven identification tool to study adaptive optimal control problems for unknown complex systems.

\lettersection

Problem Formulation In this work, we consider a continuous-time industrial process given by

[TABLE]

where $t\geq{t_{0}}$, ${x}\in\mathbb{R}^{n}$, and $u\in\mathbb{R}^{m}$ are the system states and input, respectively. $\mathcal{H}(\cdot,\cdot):\mathbb{R}^{n}\times\mathbb{R}^{m}\rightarrow\mathbb{R}^{n\times({n^{2}+mn})}$ is defined as $\mathcal{H}(x,u)\buildrel\vartriangle\over{=}\left[(x\otimes I_{n})^{\textrm{T}}~{}(u\otimes I_{n})^{\textrm{T}}\right]$, where $\otimes$ denotes the Kronecker product. $\Theta$ denotes the vector of the system parameters given by $\Theta\buildrel\vartriangle\over{=}\left[\textrm{vec}(A)^{\textrm{T}}~{}\textrm{vec}(B)^{\textrm{T}}\right]^{\textrm{T}}\in\mathbb{R}^{n^{2}+nm}$, where $A\in\mathbb{R}^{n\times n}$ is the system matrix, $B\in\mathbb{R}^{n\times m}$ is the input matrix, and $\mathrm{vec}(\cdot)$ denotes the vectorization operator, that is, $\textrm{vec}(P)=\left[p_{1}^{\textrm{T}},\ldots,p_{m}^{\textrm{T}}\right]^{\textrm{T}}$, where $p_{i}\in\mathbb{R}^{n}$ is the $i$th column of a matrix $P\in\mathbb{R}^{n\times m}$. We assume that $(A,B)$ is controllable and $(A,C)$ is observable.

In this work, we consider the following input constraint: $u\in\mathcal{U}\subset\mathbb{R}^{m}$, where $\mathcal{U}$ denotes a nonempty compact convex set and contains the origin as its interior point. In this case, since $\Theta$ is unknown, the primary objective of this work is to design a data-driven MPC formulation to obtain an open-loop optimal control policy that tracks a given reference $x_{\textrm{d}}$ and, at each sampling time $t_{k},k=1,2,\ldots$, minimizes the following cost function

[TABLE]

where $e(\cdot)=x(\cdot)-x_{\textrm{d}}(\cdot)$ denotes the error, $\varPhi(\cdot)$ denotes the terminal cost, and $Q=Q^{T}\!\succ\!0$ and $R=R^{T}\!\succeq\!0$ are the symmetric weighting matrices.

\lettersection

Methodology To facilitate MPC design, at time $t=t_{k}$, the states $x$ and the parameter $\Theta$ of the predicted model over the moving horizon $[t,t+T]$ are both learned from the input–output measurements, using a data-driven learning technique. In this work, we consider two situations.

• •

All the states and input information are available to us. Then, by rearranging (1), we have the linear error system in the form

[TABLE]

where $\hat{\Theta}$ is an estimate of the unknown parameter $\Theta$; the matrices $\mathcal{F}(\cdot):\mathbb{R}_{\geq 0}\rightarrow\mathbb{R}^{n}$ and $\mathcal{G}(\cdot):\mathbb{R}_{\geq 0}\rightarrow\mathbb{R}^{n\times(n^{2}+mn)}$ are defined as

[TABLE]

where $\delta$ denotes the sampling period; the vectors $\Xi_{x}(\cdot):\mathbb{R}_{\geq 0}\rightarrow\mathbb{R}^{n}$ and $\Xi_{u}(\cdot):\mathbb{R}_{\geq 0}\rightarrow\mathbb{R}^{m}$ are defined as

[TABLE]

• •

Only partial states and input information are available; we assume the available states as the first $q=n/2<n$ components of the states and denote them as $\xi\in\mathbb{R}^{q}$. We assume that the pair $(A,B)$ has the form

[TABLE]

Then, the linear error system is given by

[TABLE]

where $\hat{\Theta}_{1}$ is an estimate of the unknown parameter $\Theta_{1}=\left[\textrm{vec}(A_{1})^{\textrm{T}}~{}\textrm{vec}(A_{2})^{\textrm{T}}~{}\textrm{vec}(B_{1})^{\textrm{T}}\right]^{\textrm{T}}\in\mathbb{R}^{2q^{2}+mq}$; $\mathcal{F}_{1}(\cdot):\mathbb{R}_{\geq 0}\rightarrow\mathbb{R}^{q}$ and $\mathcal{G}_{1}(\cdot):\mathbb{R}_{\geq 0}\rightarrow\mathbb{R}^{q\times(2q^{2}+mq)}$ are defined as

[TABLE]

where $\delta_{1}$ and $\delta_{2}$ ($\delta_{1}\neq\delta_{2}$) denote the different periods; $\Xi_{p}(\cdot):\mathbb{R}_{\geq 0}\rightarrow\mathbb{R}^{q}$, $\Xi_{v}(\cdot):\mathbb{R}_{\geq 0}\rightarrow\mathbb{R}^{q}$, and $\Xi_{u}^{1}(\cdot):\mathbb{R}_{\geq 0}\rightarrow\mathbb{R}^{m}$ are defined as

[TABLE]

Further, from (3) and using the measurements, for a positive integer $l\leq k$, we define the vector $\varGamma_{k}\in\mathbb{R}^{ln}$ and matrix $\Psi_{k}\in\mathbb{R}^{ln\times(n^{2}+mn)}$ such that

[TABLE]

where $0\leq t_{0}<t_{1}<\cdots<t_{l}$ and $t_{i}=i\delta,i=0,1,\ldots,l$. Then, (3) implies the linear equation

[TABLE]

Notice that if $\Psi_{k}$ has full column rank, (6) can be directly solved as

[TABLE]

Similarly, for (5), we let $\delta_{1}=\delta$ and $\delta_{2}=2\delta$; thus, we obtain the same results for $\hat{\Theta}_{1}$. To guarantee $\text{rank}(\Psi_{k}^{\mathrm{T}}\Psi_{k})=n^{2}+nm$, we let the states and inputs collected over a sufficiently large number of data samples be $l\gg n^{2}+nm$. In practice, we assume that there exists a nominal control input $u=-K_{0}x$, where $K_{0}$ denotes a stabilizing feedback gain matrix, such that $\varGamma_{k}$ and $\Psi_{k}$ in (6) can be implemented using $2l$ integrators to collect information about the states and inputs. Using $\mathfrak{D}_{k}\buildrel\vartriangle\over{=}\bigcup\limits_{i=0}^{l}{\{\mathcal{F}(t_{i}),\mathcal{G}(t_{i})\}}$ by (7), we have $\hat{\Theta}=\left[\mathrm{vec}(\hat{A})^{\mathrm{T}}~{}\mathrm{vec}(\hat{B})^{\mathrm{T}}\right]^{\mathrm{T}}$, which can be used to predict and optimize the future behaviors over a finite horizon $[0,T]$.

First, let us consider the MPC formulation. As mentioned in [4], the majority (if not all) of existing formulations consider only $u(\cdot)$ as the decision variable. For (1), we can extend it to a higher order derivative of $u(\cdot)$, that is,

[TABLE]

with some control order $r\in\mathbb{N}_{+}$ larger than $\rho\geq 1$, where $\rho$ denotes the input relative degree of (1). This will improve the efficacy of our learning MPC, and the first term of $\mathbf{u}_{k}(\cdot)$ in (8) is the control input $\hat{u}_{k}(\cdot)$ that is to be optimized in (2). Then, we let $\tilde{B}=\textrm{vec}^{-1}(\hat{B})$, $\tilde{A}=\textrm{vec}^{-1}(\hat{A})$, where $\textrm{vec}^{-1}(\cdot)$ denotes the inverse operation of $\textrm{vec}(\cdot)$, and define the following matrices:

[TABLE]

At time instant $t=t_{k}$, the MPC formulation can be given by Eq. (9) (see Appendix A for details).

[TABLE]

where ${X}_{1}=\left[x^{\textrm{T}},(x^{[1]})^{\textrm{T}},\ldots,(x^{[\rho-1]})^{\textrm{T}}\right]^{\textrm{T}}$, ${X}_{2}=\left[(x^{[\rho]})^{\textrm{T}},(x^{[\rho+1]})^{\textrm{T}},\ldots,(x^{[r]})^{\textrm{T}}\right]^{\textrm{T}}$, $T_{1}(\tau)=\left[1,\tau,\ldots,\frac{\tau^{\rho-1}}{(\rho-1)!}\right]$, $T_{2}(\tau)=\left[\frac{\tau^{\rho}}{\rho!},\ldots,\frac{\tau^{r}}{r!}\right]$, $T_{3}(\tau)=\left[T_{1}(\tau),T_{2}(\tau)\right]$, $X_{1,\textrm{d}}=\left[x_{\textrm{d}}^{\textrm{T}},(x_{\textrm{d}}^{[1]})^{\textrm{T}},\ldots,(x_{\textrm{d}}^{[\rho-1]})^{\textrm{T}}\right]$, $X_{2,\textrm{d}}=\left[(x_{\textrm{d}}^{[\rho]})^{\textrm{T}},(x_{\textrm{d}}^{[\rho+1]})^{\textrm{T}},\ldots,(x_{\textrm{d}}^{[r]})^{\textrm{T}}\right]$, and

[TABLE]

with $\tilde{X}_{i}=X_{i}-X_{i,\textrm{d}}$, $\Xi_{i}(\tau)=\sqrt{Q}T_{i}(\tau)$, ${\mathcal{T}}_{i,j}=\int_{0}^{T}\Xi_{i}^{\textrm{T}}\Xi_{j}\textrm{d}\tau$, $i,j\in\{1,2\}$, and ${\mathcal{T}}=\int_{0}^{T}T_{3}^{\textrm{T}}RT_{3}\textrm{d}\tau$. This is a standard quadratic programming (QP) problem that can be solved by many available tools. In particular, for the case with box constraints, with consideration of the possible model error with the data-driven method, we can handle the input constraints by the sub-optimal method in Appendix B. We thus have the optimal control policy $\hat{u}_{k}^{\star}(t)=I_{u}{\mathbf{u}}_{k}^{\star}$, where $I_{u}=[1,0,\ldots,0]_{1\times(r+1)}$. We summarize this proposed approach as Algorithm 1 in Appendix C.

Furthermore, we consider the linear error system with the control policy $\hat{u}_{k}^{\star}(t)$ applied to (1). For the actual state trajectory, we have the continuous error as

[TABLE]

where $t\in[t_{k},t_{k}+T]$. In Appendix D, we show that $w(t)$ is bounded and has the upper bounded rate of change with time $t$; if $\tilde{A}=A$ and $\tilde{B}=B$, then $w(t)=0,t\geq 0$; with the updated control policy $\hat{u}_{k}^{\star}(t)$ at each time $t=t_{k}$, $\lim_{t\rightarrow\infty}w(t)=0$, which implies the asymptotic stability of the closed-loop system.

\lettersection

Simulation results An illustrative numerical example for two continuous stirred tank reactor (CSTR) systems is provided to validate the performance of the proposed approach. More details and discussions are presented in Appendix E.

\Acknowledgements

This work was supported by the National Key Basic Research Special Foundation of China (2014CB249200) and the authors would like to thank Prof. Zhong-ping Jiang and his CAN Lab at Tandon School of Engineering, New York University, Brooklyn, NY, USA, for many inspirations and the help of this work.

\Supplements

Appendix A, B, C, D and E.