# Data-Driven Optimization Framework for Nonlinear Model Predictive Control

Authors: Shiliang Zhang; Hui Cao; Yanbin Zhang; Lixin Jia; Zonglin Ye; Xiali Hei
Journal: Mathematical Problems in Engineering
 (2017)

---

## Abstract

The structure of the optimization procedure may affect the control quality of nonlinear model predictive control (MPC). In this paper, a data-driven optimization framework for nonlinear MPC is proposed, where the linguistic model is employed as the prediction model. The linguistic model consists of a series of fuzzy rules, whose antecedents are the membership functions of the input variables and the consequents are the predicted output represented by linear combinations of the input variables. The linear properties of the consequents lead to a quadratic optimization framework without online linearisation, which has analytical solution in the calculation of control sequence. Both the parameters in the antecedents and the consequents are calculated by a hybrid-learning algorithm based on plant data, and the data-driven determination of the parameters leads to an optimization framework with optimized controller parameters, which could provide higher control accuracy. Experiments are conducted in the process control of biochemical continuous sterilization, and the performance of the proposed method is compared with those of the methods of MPC based on linear model, the nonlinear MPC with neural network approximator, and MPC nonlinear with successive linearisations. The experimental results verify that the proposed framework could achieve higher control accuracy.


---

## Body

## 1. Introduction

Model predictive control (MPC) has been recognised as an efficient means and applied in many industrial applications successfully [1–5]. MPC refers to a class of computer algorithms that use a dynamic model to predict future behavior of the process [6–8]. In the family of MPC, the MPC based on linear model is more applicable as it brings about a quadratic programming optimization problem in online control [9], where analytical solutions exist. Since properties of many technological processes are nonlinear, the linear MPC algorithms are not likely to work properly all the time [10]. In the case of severe nonlinear process, linear MPC algorithms may result in poor control performance and even instability. Thus, nonlinear MPC is considered necessary.Nonlinear MPC schemes based on the direct use of nonlinear models have been proposed, like the nonlinear quadratic dynamic matrix control [11] and the modified nonlinear internal model control [12] as well as their extensions. Both the two schemes require a fundamental model. The derivation of fundamental models can be very time consuming and even elusive if the process is not well understood, which limits its application in practice. With the plant data, a nonlinear empirical model can be identified with less computational intensity. Some MPC algorithms employ the data-driven method for more effective implementation of control schemes. A nonlinear MPC scheme based on the second-order Volterra series model was proposed, where a successive substitution algorithm is employed [13]. Neural networks provide an alternative nonlinear modeling approach for MPC [14, 15]. Nonlinear MPC algorithms based on neural network [16, 17] imply the minimization of a cost function using computational methods to obtain the optimal control actions [18]. In the control procedure of MPC based on neural network model, gradients of the cost function are approximated numerically and the nonlinear optimization problem is solved online [19, 20]. The gradient-based optimization techniques may terminate in local minima, while global minima substantially increase the computational burden, yet they still give no guarantee that the global optimal solution is found [21].An efficient optimization framework is the key to successful implementation of nonlinear MPC. Some MPC algorithms with improved optimization structures have been proposed. Model predictive control with neural network approximator [22–24] (MPC-NNA) and its extension [25] are proposed to approximate the solution of the optimization procedure in nonlinear MPC. With an approximator trained by offline calculations [26, 27], the whole MPC algorithm in online control is replaced by the neural network approximator [28]. Though this method provides computational convenience in online control, the training process of the neural network approximator still cannot get rid of the nonconvex nonlinear optimization problem. Linearisation technologies could provide quadratic optimization framework in nonlinear MPC, and analytical solution exists in such case. In most practical applications, the accuracy of the linearisation method is sufficient compared with full nonlinear optimization in model predictive control [29, 30]. Nonlinear MPC based on multiple piecewise linear models [31] is proposed to calculate the manipulated variable based on a series of local linear MPC controllers, where quadratic optimization problem exists and analytical solution could be obtained. In this approach, the recurrent neurofuzzy [32] model is employed to represent the process, which is partitioned into several fuzzy operating regions. Within each region, a local linear model is used to represent the process. Nevertheless, the coefficients of the membership function in the local linear models are determined subjectively, which may result in undesirable control accuracy. MPC nonlinear with successive linearisation (MPC-NSL) has been described [33, 34], where online linearisation is performed at each sampling instant. With the linearisation procedure, a linear approximation model is obtained and a quadratic programming optimization can be derived. Some extensions have been reported like the MPC with nonlinear prediction and linearisation [35–37], which considers the free response of the nonlinear model. Though simple in implementation, the online linearisation approach may be susceptible to system disturbances and lower the control accuracy. A linguistic modeling method [38, 39] is reported where a series of linear models are used to represent the nonlinear model. A quadratic optimization problem without online linearisation could be obtained in the control procedure of MPC for the sake of the linear property of the consequents of the rules in the linguistic model [40]. The coefficients of the linguistic model are identified based on the plant data without subjective factors, resulting in more effective optimization procedure in MPC.In this paper, a data-driven optimization framework for nonlinear MPC is proposed, which employs the linguistic model as the prediction model. The linguistic model consists of a series of fuzzy rules, whose antecedents and consequents are the membership functions of the input variables and the predicted output calculated based on linear combinations of the input variables, respectively. Because of the linear properties of the consequent, a quadratic optimization framework can be derived where analytical solution exists. With the resulting quadratic optimization, the nonconvex nonlinear optimization or the susceptible-to-disturbance online linearisation procedure is avoided in nonlinear MPC. Based on the plant data, the antecedent parameters are initialized by the Gaussian kernel fuzzy clustering algorithm, and a hybrid-learning algorithm is employed to tune the coefficients of both the antecedents and the consequents of the model. The data-driven determination of the antecedent and consequent in the linguistic model provides optimized parameters free of subjective factors, which could bring better control accuracy. The process control of continuous sterilization is employed in the experiments to evaluate the effectiveness of the proposed framework, and the performance of the proposed framework is compared with those of the MPC based on linear model, MPC-NNA, and the MPC-NSL.The article is organized as follows: the MPC problem formulation is briefly shown in Section2; the proposed optimization framework is explained in detail in Section 3; Section 4 shows the experiments results of the process control of biochemical continuous sterilization; conclusions are given in Section 5.


## 2. MPC Problem Formulation

MPC algorithms aim to achieve accurate and fast control. The most important task is to minimize the difference between the controlled output variabley and its desired reference trajectory yref. Also, the changes of the value of the manipulated variable are expected not to be very big [24]. A set of control increments is calculated at each sampling instant t by the MPC algorithms as(1)Δut=Δut+1∣tΔut+2∣t⋯Δut+Hc-1∣tT,where Hc is the control horizon, and the increment in the equation is defined as(2)Δut+k∣t=ut+k∣t-ut+k-1∣t,k∈1,2,…,Hc-1.To minimize the difference between the predicted output y^t+k∣t and the reference trajectory yreft+k∣t over the prediction horizon N and to penalize excessive control increments, the following quadratic cost function is usually used in MPC algorithms:(3)Ct=∑k=1Ny^t+k∣t-yreft+k∣t2+∑k=0HcλΔut+k∣t2,where λ is the weighting coefficient. The control increments are calculated by minimizing the cost function online. Only the first value of the obtained control increments is applied to the control process, and the whole procedure is repeated at the next sampling instant.An explicit model is used in MPC algorithms to predict the value of future outputs of the process, based on which the MPC algorithms are derived. The selection of a suitable model is a critical task in constructing an effective MPC algorithm, as different models lead to different MPC algorithms and different optimization framework in minimizing the cost function illustrated in the Introduction. In this paper, the linguistic model is employed to improve control performance in MPC.


## 3. The Proposed Framework

Figure1 shows the block diagram of the proposed framework. The linguistic model is employed as the prediction model to make predictions on system output. The antecedent of the linear fuzzy rules in the model is the membership functions of the input variables, and the consequent is the predicted output represented by system input variables. Ki and fi mean the forced response parameters and the free response parameters corresponding to the ith rule, respectively. Both of Ki and fi are calculated by the consequent. The final value of the forced response K and the final value of the free response f are calculated by the antecedent considering all the rules. Given the reference trajectory of the controlled variable, a quadratic optimization problem could be derived where the manipulated values within the control horizon can be calculated. Only the first value of the manipulated values is employed by the plant, and the whole optimization procedure is repeated in the following control process.Figure 1
Data-driven quadratic optimization framework for nonlinear MPC.


### 3.1. Linguistic Modeling Based on Plant Data

In the proposed framework, predictions on the system dynamic are made based on the linguistic model. The parameters of the linguistic model are identified based on the plant data. The linguistic model consists of a series of Tagaki-Sugeno- (TS-) type fuzzy rules, where both the antecedent parameters and the consequent parameters of the fuzzy rules are calculated in an optimized way.


#### 3.1.1. Initialization of the Fuzzy Rules

Given that there existx1,x2,…,xn as the input variables and y as the output variable, the centers of the clusters of the process input-output dataset can be obtained by the Gaussian kernel fuzzy clustering algorithm. The centers of the clusters are used to generate the initial membership functions and a set of rules, which are used to calculate the weighted parameters of the linear models. Gaussian kernel fuzzy clustering algorithm is based on FCM clustering algorithm, and it replaces the Euclidean distance with a Gaussian kernel-induced distance [41, 42]. Assuming that there are Cn clusters obtained, then the number of linguistic terms and the number of the total TS-type rules are both Cn. The Gaussian function is employed as the membership function and the ith rule’s firing strength μi is calculated by(4)μi=∏j=1nexj-cij2/2σij2,j∈1,2,…,n,i∈1,2,…,Cn,where cij is the jth dimension value of the ith cluster’s center and σij is the width of the membership function, which is usually set as [max(xj)-min(xj)]/8. With the initial parameters of the antecedent in the fuzzy obtained rules, the initial parameters of the consequent in linear models are estimated as(5)lij=yxj·∑i=1Cnμiμi,where lij is the jth input variable’s coefficient corresponding to the ith rule.With both the initial antecedent and the consequent parameters of the fuzzy rules obtained, the predicted system outputy^ could be calculated by a set of n-dimensional inputs:(6)y^=∑i=1Cnμiyi∗μi,where yi∗ is the predicted output variable under the ith fuzzy rule and yi∗=li1x1+li2x2+⋯+linxn.


#### 3.1.2. Optimization of the Parameters in the Fuzzy Rules

For better prediction of the system output, this paper employs an iterative hybrid-learning algorithm to adjust the initial parameters of the fuzzy rules.Define the matricesP and S as(7)P=l11,…,l1n,l21,…,l2n,…,lCn1,…,lCnT,S=m00⋯00m0⋯000m⋯0⋮⋮⋮⋱⋮000⋯mN×N,where N=Cn·n and m represents a positive infinite value, which is often assigned as 1e6. For the kth entry of the training data, {x1(k),x2(k),…,xn(k),y(k)}, where k∈1,2,…,Nt and Nt is the total number of the training data. Assuming that (8)αk=μ1kx1k,μ1kx2k,…,μ1kxnk,μ2kx1k,μ2kx2k,…,μ2kxnk,…,μCnkx1k,μCnkx2k,…,μCnkxnkT,where μi(k),  i∈{1,2,…,Cn} is the degree of membership of the ith rule corresponding to the kth entry of the training data, P and S could be approximated as(9)Pk+1=Pk+Sk·αk·yk-αkTPk,Sk+1=Sk-Sk·αkαkT·Sk1+αkT·Sk·αk.The initial value of S is defined in (7). Let the initial P be a zero matrix, and the new values of lij can be obtained in the iterative calculation. The iteration will not stop until the RMSE is no larger than TH, which is the user-defined threshold. The update of the antecedent parameters of the fuzzy rules, cij and σij, in the Gaussian membership functions is shown as follows. For the training data {x1(k),x2(k),…,xn(k),y(k)},  k∈{1,2,…,Nt}, dek can be calculated as(10)dek=-n·Cn∑i=1Cn-2μikyi∗k·yk-y^k∑i=1Cnμik,where y^(k) is the predicted value of the linguistic model corresponding to the input {x1(k),x2(k),…,xn(k),y(k)} and yi∗(k) means the predicted output under the ith fuzzy rule. Assuming that Tc and Tσ are both Cn×n matrices, the jth column elements of the ith row of Tc and Tσ in the (k+1)th iterative calculation are(11)Ti,jck+1=Ti,jck+∑l=1Ntdel·xjl-cijkσijk·exjl-cijk2/2σij2k,Ti,jσk+1=Ti,jσk+∑l=1Ntdel·xjl-cijkσijk·exjl-cijk2/2σij2k,where the initial Tc and Tσ are zero matrices. The values of cij and σij in the (k+1)th iterative calculation can be obtained by (12)cijk+1=cijk-ts·Ti,jck∑i=1Cn∑j=1nTi,jck2+Ti,jσk2,σijk+1=σijk-ts·Ti,jck∑i=1Cn∑j=1nTi,jck2+Ti,jσk2,where ts is the training step.The iterative calculation goes on until the value of RMSE is less than the defined thresholdTH or the number of iterations reaches the defined value. When the iteration stops, all the optimized parameters are obtained and the linguistic model is established. As a result, the system dynamic is presented by the nonlinear linguistic model consisting of a series of weighted linear models. The linguistic model will be employed in nonlinear MPC for a quadratic optimization framework, which will be shown in detail in the next section.


### 3.2. The Quadratic Programming Optimization Framework Based on Linguistic Model for Nonlinear Model Predictive Control

Assuming thaty(t) is the system output and u(t) is the manipulated input variable at the sampling instant t, the system dynamic can be represented by the linguistic model as follows:(13)yt=∑i=1Cnμi·yi∗∑i=1Cn,where(14)μi=∏p=1Peyt-p-cyip2/2σyip2·∏q=0Qeut-q-1-cuiq2/2σuiq2and the integrated Controller Autoregressive Moving-Average (CARIMA) model [43] is used to describe the consequent of the linguistic model as follows:(15)yi∗=∑p=1Paip·yt-p+∑q=0Qbiqut-q-1+eitΔ,where cyip means the center of the cluster of y(t-p) under the ith fuzzy rule, cuiq means the center of the cluster of u(t-q) under the ith rule, σyip means the membership function width of y(t-p) under the ith rule, σuiq means the membership function width of u(t-q-1) under the ith rule, and Δ=1-q-1 means the differencing operator, while q-1 means the backward shift operator. All the centers of the cluster, the membership functions width, and the linear model parameters aip,  p∈{1,2,…,P} and biq,  q∈{1,2,…,Q} are determined in the linguistic modeling process. Equation (15) can be rewritten as(16)Aq-1yi∗t=Biq-1ut-1+eitΔwith(17)Aiq-1=1-ai1q-1-ai2q-2-⋯-aiPq-P,Biq-1=bi0+bi1q-1+bi2q-2+⋯+biQq-Q.Consider the Diophantine equation: (18)1=Eijq-1A~iq-1+q-jFijq-1withA~iq-1=ΔAiq-1.The degree of the polynomials Eij and Fij is defined as j-1 and P+Q, respectively, and both the two polynomials could be recursively calculated [43]. Multiply (16) by ΔEijqj; it can be obtained that (19)A~iEijq-1yi∗t+j=Eijq-1Biq-1Δut+j-1+Eijq-1eit+j.Given (18) and (19), it can be obtained that (20)1-q-jFijq-1yi∗t+j=Eijq-1Biq-1Δut+j-1+Eijq-1eit+jwhich could be rewritten as (21)yi∗t+j=Fijq-1yt+EijBiq-1Δut+j-1+Eijq-1eit+j.Because the degree of Eij(q-1) is j-1, the noise terms are all in the future in (20). Thus, the best prediction of yi∗(t+j) is calculated by(22)yi∗t+j∣t=KijΔut+j-1+Fijq-1yt,where Kij(q-1)=Eij(q-1)Bi(q-1). For the degree of Bi(q-1) is Q, it can be gained that the degree of Kij(q-1) is j+Q-1. So it could be obtained that(23)Kijq-1=ki,0+ki,1q-1+ki,2q-2+⋯+ki,j+Q-1qj+Q-1.For the ith rule of the linguistic model, the N ahead optimal predictions of the system output can be calculated as(24)y^i∗t+1∣t=Ki1q-1Δut+Fi1q-1yt,y^i∗t+2∣t=Ki2q-1Δut+1+Fi2q-1yt,⋮y^i∗t+N∣t=KiNq-1Δut+N-1+FiNq-1ytwhich can be rewritten as(25)y^i∗=Kiu+Fiq-1yt+Ki′Δut-1,where(26)y^i∗=y^i∗t+1∣ty^i∗t+2∣t⋮y^i∗t+N∣t,u=ΔutΔut+1⋮Δut+N-1,Ki=ki,000⋯0ki,1ki,00⋯0ki,2ki,1ki,0⋯0⋮⋮⋮⋱⋮ki,N-1ki,N-2ki,N-3⋯ki,0,Ki′q-1=Ki1q-1-ki,0qKi2q-1-ki,0-ki,1q-1q2⋮KiNq-1-ki,0-ki,1q-1-⋯-ki,N-1qN-1qN,Fiq-1=Fi1q-1Fi2q-1⋮FiNq-1.The last two terms in (26) depend only on the variables of the past instants, assuming that(27)fiq-1=fi,1fi,2⋮fi,N=Fiq-1yt+Ki′Δut-1and (25) could be rewritten as(28)y^i∗=Kiu+fi.Considering all Cn rules of the linguistic model, the final optimal predictions N ahead could be calculated as follows:(29)y^=∑i=1Cnμi·y^i∗∑i=1Cnμiwhich could be represented by(30)y^=Ku+f,where(31)K=k000⋯0k1k00⋯0k2k1k0⋯0⋮⋮⋮⋱⋮kN-1kN-2kN-3⋯k0,kj=∑i=1Cnμiki,j∑i=1Cnμi,fq-1=f1f2⋮fN,fj=∑i=1Cnμifi,j∑i=1Cnμi.Based on the predictions, the optimal control sequence could be calculated by minimizing a cost function. Considering the constraint on the control effort, the cost function could be presented as(32)C=∑l=1Ny^t+l∣t-wt+l2+∑l=1NλΔut+l-12which could be rewritten as(33)C=Ku-f-wTKu-f-w+λuTu,where(34)w=wt+1,wt+2,…,wt+NTand w(t+i),  i∈{1,2,…,N} means the reference trajectory of the controlled output variable at the instant t+i, and  λ means the penalty coefficient on the control effort. Assuming that no more constraint is considered, the minimum of the cost function is calculated by making the gradient of C equal 0, resulting in(35)u=KTK+λI-1KTw-f.Only the first value in u is passed to the control process. In the following sampling instants, the procedure will be repeated.Now the data-driven optimization framework for nonlinear MPC has been described. Experiments are conducted in the next section to verify the effectiveness of the proposed framework.


## 3.1. Linguistic Modeling Based on Plant Data

In the proposed framework, predictions on the system dynamic are made based on the linguistic model. The parameters of the linguistic model are identified based on the plant data. The linguistic model consists of a series of Tagaki-Sugeno- (TS-) type fuzzy rules, where both the antecedent parameters and the consequent parameters of the fuzzy rules are calculated in an optimized way.


### 3.1.1. Initialization of the Fuzzy Rules

Given that there existx1,x2,…,xn as the input variables and y as the output variable, the centers of the clusters of the process input-output dataset can be obtained by the Gaussian kernel fuzzy clustering algorithm. The centers of the clusters are used to generate the initial membership functions and a set of rules, which are used to calculate the weighted parameters of the linear models. Gaussian kernel fuzzy clustering algorithm is based on FCM clustering algorithm, and it replaces the Euclidean distance with a Gaussian kernel-induced distance [41, 42]. Assuming that there are Cn clusters obtained, then the number of linguistic terms and the number of the total TS-type rules are both Cn. The Gaussian function is employed as the membership function and the ith rule’s firing strength μi is calculated by(4)μi=∏j=1nexj-cij2/2σij2,j∈1,2,…,n,i∈1,2,…,Cn,where cij is the jth dimension value of the ith cluster’s center and σij is the width of the membership function, which is usually set as [max(xj)-min(xj)]/8. With the initial parameters of the antecedent in the fuzzy obtained rules, the initial parameters of the consequent in linear models are estimated as(5)lij=yxj·∑i=1Cnμiμi,where lij is the jth input variable’s coefficient corresponding to the ith rule.With both the initial antecedent and the consequent parameters of the fuzzy rules obtained, the predicted system outputy^ could be calculated by a set of n-dimensional inputs:(6)y^=∑i=1Cnμiyi∗μi,where yi∗ is the predicted output variable under the ith fuzzy rule and yi∗=li1x1+li2x2+⋯+linxn.


### 3.1.2. Optimization of the Parameters in the Fuzzy Rules

For better prediction of the system output, this paper employs an iterative hybrid-learning algorithm to adjust the initial parameters of the fuzzy rules.Define the matricesP and S as(7)P=l11,…,l1n,l21,…,l2n,…,lCn1,…,lCnT,S=m00⋯00m0⋯000m⋯0⋮⋮⋮⋱⋮000⋯mN×N,where N=Cn·n and m represents a positive infinite value, which is often assigned as 1e6. For the kth entry of the training data, {x1(k),x2(k),…,xn(k),y(k)}, where k∈1,2,…,Nt and Nt is the total number of the training data. Assuming that (8)αk=μ1kx1k,μ1kx2k,…,μ1kxnk,μ2kx1k,μ2kx2k,…,μ2kxnk,…,μCnkx1k,μCnkx2k,…,μCnkxnkT,where μi(k),  i∈{1,2,…,Cn} is the degree of membership of the ith rule corresponding to the kth entry of the training data, P and S could be approximated as(9)Pk+1=Pk+Sk·αk·yk-αkTPk,Sk+1=Sk-Sk·αkαkT·Sk1+αkT·Sk·αk.The initial value of S is defined in (7). Let the initial P be a zero matrix, and the new values of lij can be obtained in the iterative calculation. The iteration will not stop until the RMSE is no larger than TH, which is the user-defined threshold. The update of the antecedent parameters of the fuzzy rules, cij and σij, in the Gaussian membership functions is shown as follows. For the training data {x1(k),x2(k),…,xn(k),y(k)},  k∈{1,2,…,Nt}, dek can be calculated as(10)dek=-n·Cn∑i=1Cn-2μikyi∗k·yk-y^k∑i=1Cnμik,where y^(k) is the predicted value of the linguistic model corresponding to the input {x1(k),x2(k),…,xn(k),y(k)} and yi∗(k) means the predicted output under the ith fuzzy rule. Assuming that Tc and Tσ are both Cn×n matrices, the jth column elements of the ith row of Tc and Tσ in the (k+1)th iterative calculation are(11)Ti,jck+1=Ti,jck+∑l=1Ntdel·xjl-cijkσijk·exjl-cijk2/2σij2k,Ti,jσk+1=Ti,jσk+∑l=1Ntdel·xjl-cijkσijk·exjl-cijk2/2σij2k,where the initial Tc and Tσ are zero matrices. The values of cij and σij in the (k+1)th iterative calculation can be obtained by (12)cijk+1=cijk-ts·Ti,jck∑i=1Cn∑j=1nTi,jck2+Ti,jσk2,σijk+1=σijk-ts·Ti,jck∑i=1Cn∑j=1nTi,jck2+Ti,jσk2,where ts is the training step.The iterative calculation goes on until the value of RMSE is less than the defined thresholdTH or the number of iterations reaches the defined value. When the iteration stops, all the optimized parameters are obtained and the linguistic model is established. As a result, the system dynamic is presented by the nonlinear linguistic model consisting of a series of weighted linear models. The linguistic model will be employed in nonlinear MPC for a quadratic optimization framework, which will be shown in detail in the next section.


## 3.1.1. Initialization of the Fuzzy Rules

Given that there existx1,x2,…,xn as the input variables and y as the output variable, the centers of the clusters of the process input-output dataset can be obtained by the Gaussian kernel fuzzy clustering algorithm. The centers of the clusters are used to generate the initial membership functions and a set of rules, which are used to calculate the weighted parameters of the linear models. Gaussian kernel fuzzy clustering algorithm is based on FCM clustering algorithm, and it replaces the Euclidean distance with a Gaussian kernel-induced distance [41, 42]. Assuming that there are Cn clusters obtained, then the number of linguistic terms and the number of the total TS-type rules are both Cn. The Gaussian function is employed as the membership function and the ith rule’s firing strength μi is calculated by(4)μi=∏j=1nexj-cij2/2σij2,j∈1,2,…,n,i∈1,2,…,Cn,where cij is the jth dimension value of the ith cluster’s center and σij is the width of the membership function, which is usually set as [max(xj)-min(xj)]/8. With the initial parameters of the antecedent in the fuzzy obtained rules, the initial parameters of the consequent in linear models are estimated as(5)lij=yxj·∑i=1Cnμiμi,where lij is the jth input variable’s coefficient corresponding to the ith rule.With both the initial antecedent and the consequent parameters of the fuzzy rules obtained, the predicted system outputy^ could be calculated by a set of n-dimensional inputs:(6)y^=∑i=1Cnμiyi∗μi,where yi∗ is the predicted output variable under the ith fuzzy rule and yi∗=li1x1+li2x2+⋯+linxn.


## 3.1.2. Optimization of the Parameters in the Fuzzy Rules

For better prediction of the system output, this paper employs an iterative hybrid-learning algorithm to adjust the initial parameters of the fuzzy rules.Define the matricesP and S as(7)P=l11,…,l1n,l21,…,l2n,…,lCn1,…,lCnT,S=m00⋯00m0⋯000m⋯0⋮⋮⋮⋱⋮000⋯mN×N,where N=Cn·n and m represents a positive infinite value, which is often assigned as 1e6. For the kth entry of the training data, {x1(k),x2(k),…,xn(k),y(k)}, where k∈1,2,…,Nt and Nt is the total number of the training data. Assuming that (8)αk=μ1kx1k,μ1kx2k,…,μ1kxnk,μ2kx1k,μ2kx2k,…,μ2kxnk,…,μCnkx1k,μCnkx2k,…,μCnkxnkT,where μi(k),  i∈{1,2,…,Cn} is the degree of membership of the ith rule corresponding to the kth entry of the training data, P and S could be approximated as(9)Pk+1=Pk+Sk·αk·yk-αkTPk,Sk+1=Sk-Sk·αkαkT·Sk1+αkT·Sk·αk.The initial value of S is defined in (7). Let the initial P be a zero matrix, and the new values of lij can be obtained in the iterative calculation. The iteration will not stop until the RMSE is no larger than TH, which is the user-defined threshold. The update of the antecedent parameters of the fuzzy rules, cij and σij, in the Gaussian membership functions is shown as follows. For the training data {x1(k),x2(k),…,xn(k),y(k)},  k∈{1,2,…,Nt}, dek can be calculated as(10)dek=-n·Cn∑i=1Cn-2μikyi∗k·yk-y^k∑i=1Cnμik,where y^(k) is the predicted value of the linguistic model corresponding to the input {x1(k),x2(k),…,xn(k),y(k)} and yi∗(k) means the predicted output under the ith fuzzy rule. Assuming that Tc and Tσ are both Cn×n matrices, the jth column elements of the ith row of Tc and Tσ in the (k+1)th iterative calculation are(11)Ti,jck+1=Ti,jck+∑l=1Ntdel·xjl-cijkσijk·exjl-cijk2/2σij2k,Ti,jσk+1=Ti,jσk+∑l=1Ntdel·xjl-cijkσijk·exjl-cijk2/2σij2k,where the initial Tc and Tσ are zero matrices. The values of cij and σij in the (k+1)th iterative calculation can be obtained by (12)cijk+1=cijk-ts·Ti,jck∑i=1Cn∑j=1nTi,jck2+Ti,jσk2,σijk+1=σijk-ts·Ti,jck∑i=1Cn∑j=1nTi,jck2+Ti,jσk2,where ts is the training step.The iterative calculation goes on until the value of RMSE is less than the defined thresholdTH or the number of iterations reaches the defined value. When the iteration stops, all the optimized parameters are obtained and the linguistic model is established. As a result, the system dynamic is presented by the nonlinear linguistic model consisting of a series of weighted linear models. The linguistic model will be employed in nonlinear MPC for a quadratic optimization framework, which will be shown in detail in the next section.


## 3.2. The Quadratic Programming Optimization Framework Based on Linguistic Model for Nonlinear Model Predictive Control

Assuming thaty(t) is the system output and u(t) is the manipulated input variable at the sampling instant t, the system dynamic can be represented by the linguistic model as follows:(13)yt=∑i=1Cnμi·yi∗∑i=1Cn,where(14)μi=∏p=1Peyt-p-cyip2/2σyip2·∏q=0Qeut-q-1-cuiq2/2σuiq2and the integrated Controller Autoregressive Moving-Average (CARIMA) model [43] is used to describe the consequent of the linguistic model as follows:(15)yi∗=∑p=1Paip·yt-p+∑q=0Qbiqut-q-1+eitΔ,where cyip means the center of the cluster of y(t-p) under the ith fuzzy rule, cuiq means the center of the cluster of u(t-q) under the ith rule, σyip means the membership function width of y(t-p) under the ith rule, σuiq means the membership function width of u(t-q-1) under the ith rule, and Δ=1-q-1 means the differencing operator, while q-1 means the backward shift operator. All the centers of the cluster, the membership functions width, and the linear model parameters aip,  p∈{1,2,…,P} and biq,  q∈{1,2,…,Q} are determined in the linguistic modeling process. Equation (15) can be rewritten as(16)Aq-1yi∗t=Biq-1ut-1+eitΔwith(17)Aiq-1=1-ai1q-1-ai2q-2-⋯-aiPq-P,Biq-1=bi0+bi1q-1+bi2q-2+⋯+biQq-Q.Consider the Diophantine equation: (18)1=Eijq-1A~iq-1+q-jFijq-1withA~iq-1=ΔAiq-1.The degree of the polynomials Eij and Fij is defined as j-1 and P+Q, respectively, and both the two polynomials could be recursively calculated [43]. Multiply (16) by ΔEijqj; it can be obtained that (19)A~iEijq-1yi∗t+j=Eijq-1Biq-1Δut+j-1+Eijq-1eit+j.Given (18) and (19), it can be obtained that (20)1-q-jFijq-1yi∗t+j=Eijq-1Biq-1Δut+j-1+Eijq-1eit+jwhich could be rewritten as (21)yi∗t+j=Fijq-1yt+EijBiq-1Δut+j-1+Eijq-1eit+j.Because the degree of Eij(q-1) is j-1, the noise terms are all in the future in (20). Thus, the best prediction of yi∗(t+j) is calculated by(22)yi∗t+j∣t=KijΔut+j-1+Fijq-1yt,where Kij(q-1)=Eij(q-1)Bi(q-1). For the degree of Bi(q-1) is Q, it can be gained that the degree of Kij(q-1) is j+Q-1. So it could be obtained that(23)Kijq-1=ki,0+ki,1q-1+ki,2q-2+⋯+ki,j+Q-1qj+Q-1.For the ith rule of the linguistic model, the N ahead optimal predictions of the system output can be calculated as(24)y^i∗t+1∣t=Ki1q-1Δut+Fi1q-1yt,y^i∗t+2∣t=Ki2q-1Δut+1+Fi2q-1yt,⋮y^i∗t+N∣t=KiNq-1Δut+N-1+FiNq-1ytwhich can be rewritten as(25)y^i∗=Kiu+Fiq-1yt+Ki′Δut-1,where(26)y^i∗=y^i∗t+1∣ty^i∗t+2∣t⋮y^i∗t+N∣t,u=ΔutΔut+1⋮Δut+N-1,Ki=ki,000⋯0ki,1ki,00⋯0ki,2ki,1ki,0⋯0⋮⋮⋮⋱⋮ki,N-1ki,N-2ki,N-3⋯ki,0,Ki′q-1=Ki1q-1-ki,0qKi2q-1-ki,0-ki,1q-1q2⋮KiNq-1-ki,0-ki,1q-1-⋯-ki,N-1qN-1qN,Fiq-1=Fi1q-1Fi2q-1⋮FiNq-1.The last two terms in (26) depend only on the variables of the past instants, assuming that(27)fiq-1=fi,1fi,2⋮fi,N=Fiq-1yt+Ki′Δut-1and (25) could be rewritten as(28)y^i∗=Kiu+fi.Considering all Cn rules of the linguistic model, the final optimal predictions N ahead could be calculated as follows:(29)y^=∑i=1Cnμi·y^i∗∑i=1Cnμiwhich could be represented by(30)y^=Ku+f,where(31)K=k000⋯0k1k00⋯0k2k1k0⋯0⋮⋮⋮⋱⋮kN-1kN-2kN-3⋯k0,kj=∑i=1Cnμiki,j∑i=1Cnμi,fq-1=f1f2⋮fN,fj=∑i=1Cnμifi,j∑i=1Cnμi.Based on the predictions, the optimal control sequence could be calculated by minimizing a cost function. Considering the constraint on the control effort, the cost function could be presented as(32)C=∑l=1Ny^t+l∣t-wt+l2+∑l=1NλΔut+l-12which could be rewritten as(33)C=Ku-f-wTKu-f-w+λuTu,where(34)w=wt+1,wt+2,…,wt+NTand w(t+i),  i∈{1,2,…,N} means the reference trajectory of the controlled output variable at the instant t+i, and  λ means the penalty coefficient on the control effort. Assuming that no more constraint is considered, the minimum of the cost function is calculated by making the gradient of C equal 0, resulting in(35)u=KTK+λI-1KTw-f.Only the first value in u is passed to the control process. In the following sampling instants, the procedure will be repeated.Now the data-driven optimization framework for nonlinear MPC has been described. Experiments are conducted in the next section to verify the effectiveness of the proposed framework.


## 4. Experiment Implementation

The process control of continuous sterilization is employed in the experiments, and the performance of the proposed framework is compared with those of the MPC based on linear model (MPC-L), the MPC based on neural network approximator (MPC-NNA), and the MPC nonlinear with successive linearisation (MPC-NSL). The continuous sterilization system is described briefly, and the modeling and the controlling on the continuous sterilization system are conducted based on plant data.


### 4.1. The Continuous Sterilization Process

Continuous sterilization is conducted on the medium before fermentation in biochemical engineering, and the nonsterilized medium is heated to the desirable temperature by steam in the sterilization process, so the undesired microorganisms are removed in the medium. After the continuous sterilization, necessary pure culture could be provided for the objective microorganism. Figure2 illustrates the diagrammatic sketch of the continuous sterilization system.Figure 2
The continuous sterilization process.In the continuous sterilization system, the quantity of steam used to heat the unsterilized medium is changed by the opening degree of the valve on the steam tube (Ov), while ΔOv(i)=Ov(i)-Ov(i-1) is the manipulated variable. The controlled output variable is the temperature of the heated medium at the outlet of the steam injector (Th). There exist measurable disturbances as the flow rate of medium (Fm), the temperature of unsterilized medium (Tm), and the temperature of steam (Ts). The control procedure aims to provide suitable heating on the medium, so that the temperature of the heated medium can be kept around the set value.


### 4.2. Experiment Settings

To evaluate the effectiveness of the proposed framework in nonlinear MPC, the plant data in a real workshop of continuous sterilization are collected and experiments are conducted based on the obtained data. The plant data includes the input variableOv, which adjusts the quantity of steam for the sterilization, while Tm,  Fm, and Ts are the measurable system disturbances, and Th is the output variable. The value ranges of Ov,  Tm,  Fm,  Ts, and Th are 0–100%, 0–70°C, 0–25 ton/h, 150–210°C, and 90–140°C, respectively.Modeling and controlling are conducted based on the plant data, and MPC based on different models are employed, including the linear model, the neural network, the adaptive neurofuzzy inference system (ANFIS) [44–48] with online linearisation (ANFIS-OL), and the linguistic model. With those models, four MPC algorithms are derived: MPC-L, MPC-NNA, MPC-NSL, and the proposed optimization framework for nonlinear MPC based on linguistic model, respectively. The proposed method is expressed as MPC-LOF for short.


### 4.3. Experiments on Modeling

In the modeling process, 788 samples are treated as training data, and 394 samples are treated as test data. The performances of the linear model, the neural network, the ANFIS with online linearisation, and the linguistic model are compared. For there are different parameters contained in different models, the parameters are selected by scanning for the best result based on the training data. All the involved parameters in the modeling are listed in Table1.Table 1
Parameters involved in modeling corresponding to different control methods.
Modeling method
H
r
r
N
n
Linear model
√
√
×
Linguistic model
√
√
×
Neural network
√
×
√
ANFIS-OL
√
√
×In Table1, the receding horizon Hr is the sampling length of the variables employed in modeling of the prediction model. For the linguistic model in the proposed framework, the radius value r affects the result of the Gaussian kernel clustering algorithm and the number of rules obtained. For the MPC-NNA, the three-layer back-propagation neural network is employed, and Nn means the number of nodes in the hidden layers. For the MPC-NSL, the ANFIS model is adopted. For the ANFIS, the subtractive clustering algorithm constructs the initial fuzzy model, and the radius r would affect the number of clusters obtained. In the following experiments, the parameters involved in the linguistic model are discussed, and the performances of different modeling methods are compared.The first experiment shows the effect ofHr and r on the modeling of the linguistic model. With the training data, experiment is conducted and the modeling results under different couples of Hr and r are shown in Figure 3(a). From Figure 3(a), it can be seen that Hr=5 tends to bring about lower mean square error (MSE) compared with the other values, and the best result is obtained when Hr=5 and r=0.5.Figure 3
Results on the modeling.
(a)
The effect ofHr and r on the results of linguistic modeling
(b)
Learning curves of the linguistic model withHr=5
(c)
Modeling result on training data based on different modeling methods
(d)
Modeling result on testing data based on different modeling methodsThe second experiment is conducted with fixedHr=5 to observe in detail the effect of r on the modeling result of the linguistic model. The result is shown is in Figure 3(b), where r in the middle of its range tends to bring about higher accuracy.To make comparison between different model methods, the parameters of the linear model, the ANFIS with online linearisation, and the neural network are selected by scanning for the best results based on the training data. The selected parameters are shown in Table2.Table 2
Parameters involved in modeling corresponding to different control methods.
Modeling method
H
r
r
N
n
Linear model
3
×
×
Linguistic model
5
0.5
×
Neural network
4
×
6
ANFIS-OL
4
0.6
×The performances of different methods with the selected parameters are shown in Figures3(c) and 3(d). From Figure 3(c), it can be seen that the neural network achieves the best performance than the linear model, the ANFIS with online linearisation (ANFIS-OL), and the linguistic model in the training process. The ANFIS-OL shows the worst accuracy, possibly because of low noise rejection. In Figure 3(d), the linguistic model provides the best modeling accuracy in the test data, while the neural network shows slightly lower accuracy than the linguistic model and the linear model, which shows lower generalization than the latter two methods. Yet, the ANFIS-OL method still shows the worst performance. Table 3 shows the detailed comparison between different modeling methods.Table 3
Performances of different modeling method on both training and test data.
Modeling method
Final MSE
Training data
Test data
Linguistic model
0.1524
0.4430
Linear model
0.2266
0.5026
Neural network
0.0823
0.5584
ANFIS-OL
11.5272
23.0631


### 4.4. Experiments on Controlling

To evaluate the effectiveness of the proposed control method, the control methods of MPC-L, MPC-NNA, MPC-NSL, and MPC-LOF are conducted to make comparison based on the plant data, which contain 1187 samples as training data and anther 1097 samples as test data. The process model used in those MPC methods is employed as the same neural network model generated by the plant data withHr=7. In the experiments, the performance indicator is the control accuracy, which represents the approximation of the controlled variable to the expected reference trajectory. The control accuracy is also the most important indicator in a practical continuous sterilization. The indicator of the control accuracy Ic is calculated as(36)Ic=1Nv∑t=1NvTht-Ret2,where Nv means the number of data samples, while Re(t) means the reference trajectory at the time instant t. Re(t) is set to be 123°C in the control process. Different MPC methods are conducted, with the involved parameters discussed. All the involved parameters in the controlling experiments are listed in Table 4, where λ is the penalty coefficients on the control effort and Hc is the control horizon.Table 4
Parameters involved in MPC methods.
Control method
Controller parameters
λ
H
c
H
r
r
N
n
MPC-L
√
√
√
×
×
MPC-NNA
√
×
√
×
√
MPC-NSL
√
√
√
√
×
MPC-LOF
√
√
√
√
×In the following experiments, the parameters involved in the MPC-LOF are discussed, and the performances of MPC-L, MPC-NNA, MPC-NSL, and the MPC-LOF are compared.The first experiment investigates the effects of parameters in MPC-LOF on the controlling performance, and the best combination of parameters in MPC-LOF is selected by scanning for the best result. With the training data, the optimal parameters combination is selected asHr=5,  Hc=5,  r=0.5, and  λ=3.78. The parameter selection in MPC-LOF is shown in Figures 4–7 for illustration. The effect of Hr and r on the performance of MPC-LOF is shown for illustration in Figure 4, where the indicator of TestingIc is employed to show the control performance in detail. The ith TestingIc is calculated as ∑t=1iTht-Ret2/i. It can be seen from Figure 4 that, when Hr=5 and λ=3.78, desirable result could be found when Hr=4or5, and the best result is obtained when Hr=5 and r=0.5. The final Ic corresponding to different combinations of Hr and r with fixed Hc=5 and λ=3.78 is detailed in Table 5.Table 5
The finalIc in MPC-LOF corresponding to different couples of Hr and r when Hc=5 and λ=3.78.
I
c
r
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
H
r
2
107.9560
265.6859
257.2303
252.4744
246.9592
245.1016
274.8260
254.3659
164.3792
3
224.1424
251.8802
282.2562
263.9272
6.3428
262.2148
265.0530
244.6227
2.0565
4
282.3905
280.4251
274.4064
258.4806
3.6191
212.4163
262.9407
662.4718
706.5165
5
9.3565
88.5872
29.2493
30.5592
1.6593
161.6644
291.7744
2.0702
87.0250
6
197.6085
187.1392
262.6352
259.0817
258.1882
258.5736
7.0774
39.0084
17.5822Figure 4
The effect ofHr and r on performance of MPC-LOF when Hc=5 and λ=3.78.Figure 5
The effect ofHc and r on performance of MPC-LOF when Hr=5 and λ=3.78, part 1.Figure 6
The effect ofHc and r on performance of MPC-LOF when Hr=5 and λ=3.78, part 2.Figure 7
The effect ofλ and r on performance of MPC-LOF when Hc=5 and Hr=5.When there isHr and λ=3.78, the best performance is found when Hc=5 and r=0.5, which performs better than the other couples of parameters as shown in Figures 5-6 for illustration. The details of the final Ic with different couples of Hc and r with fixed Hr=5 and λ=3.78 are shown in Table 6.Table 6
The finalIc in MPC-LOF corresponding to different couples of Hc and r when Hr=5 and λ=3.78.
I
c
r
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
H
c
2
593.6593
140.8252
283.2953
261.1969
35.1593
303.4396
267.9231
396.3824
439.1347
3
283.0946
137.2669
274.4064
258.3611
3.6191
212.4163
262.9395
662.4733
706.5170
4
297.0199
133.8319
263.1655
255.9483
1.7700
38.8709
259.8053
692.2658
738.5095
5
9.3565
88.5872
29.2493
30.5592
1.6593
161.6644
291.7744
2.0702
87.0250
6
338.4374
129.4590
262.1812
250.5345
31.3147
31.4507
250.3436
737.1146
766.6540
7
184.7563
125.0401
261.9619
248.3564
32.4505
33.1497
245.2064
753.7656
771.5174
8
187.6874
55.7913
261.8151
246.7770
32.7710
35.2300
241.9463
378.1504
733.8148
9
203.9504
7.2457
261.5909
245.5387
32.6305
35.8943
240.0754
373.3230
561.3975
10
186.2029
4.7204
261.3996
244.4238
32.9118
33.3927
238.0571
371.1604
570.3603
11
162.1932
3.9902
261.0869
243.2266
33.9692
34.4692
235.5092
334.7633
587.1572
12
142.0463
3.8393
260.5998
241.9458
34.0403
37.4656
231.6528
311.7689
548.6458
13
120.2090
3.6480
259.9733
240.5188
35.1123
37.9268
183.0111
174.5661
556.9916
14
80.8538
3.3041
258.9066
238.9687
35.6674
38.3702
169.9955
154.8572
541.8535
15
56.8366
2.8585
256.4588
237.1937
35.9254
38.7456
22.6742
150.4693
510.2067
16
61.1071
2.4136
254.7232
235.3396
36.1617
39.2147
9.7830
152.9248
354.7270Figure7 shows the effect of λ and r on the MPC-LOF, and the final Ic under different couples of λ and r is shown for comparison. When Hc=5 and Hr=5, the best result is found when λ=3.78 and r=0.5.The second experiment is conducted to show the performances comparison between the MPC-L, MPC-NNA, MPC-NSL, and the proposed MPC-LOF. For there are different parameters contained in different MPC methods, the parameters are selected by scanning for the best result based on the training data, whose selection results are shown in Table7.Table 7
Optimal parameters for four MPC methods.
Control method
Controller parameters
λ
H
c
H
r
r
N
n
MPC-L
3.16
5
3
×
×
MPC-NNA
4.26
×
3
×
6
MPC-NSL
0.50
5
4
0.2
×
MPC-LOF
3.78
5
5
0.5
×The plot of the manipulated variable (ΔOv) and controlled variable (Th) corresponding to the controller parameters in Table 7 are shown in Figure 8. Given the reference trajectory of 123°C, the controlled variable Th can be kept within the interval [122,124.5] under the MPC-NNA method, and there exists deviation from the reference trajectory in the controlled variable; for example, Hh is kept around 122°C from the 90th to the 140th iteration and rises to about 123.5°C in the following 200 iterations. There are no significant fluctuations in the trajectory of the controlled variable in MPC-NNA, but the controlled variable could not be held around the reference trajectory stably. As for the MPC-NSL method, the controlled variable shows better approximation to the reference compared with MPC-NNA, and Hh is kept within [121,125] generally. Yet there exists fluctuation during the control process; for example, the controlled variable fluctuation appears from the 700th to the 860th iteration, whose corresponding manipulated variable is also in severe fluctuation as shown in Figure 8. The fluctuations show low resistance of the online linearisation approach in MPC-NSL to the disturbances. Besides, there are deviations in the control process of MPC-NSL; for example, the controlled variable rises from 120.5°C to 127°C at the 20th iteration and rises from 123°C to 130°C at the 690th iteration. In MPC-L, the controlled output variable is generally kept within [120.5,125], and both severe deviation and fluctuation exist in the control process. For example, fluctuations appear from the 200th to the 320th iteration, and deviations exist at the 825th iteration with a decrease from 130°C to 117°C and the increase from 123°C to 128°C at the 690th iteration. Though the modeling accuracy of linear model is not the worst in the experiments on modeling, the bad control performance compared with the other methods shows the limit of the pure linear optimization framework in nonlinear system. As for the performance of the proposed MPC-LOF, Th is kept within [122.5,124]. There exist fluctuations in the controlled variable, yet they are much less than those in MPC-NSL and MPC-L, and the value of the controlled variable is kept closely around the reference compared with MPC-NNA. The value of indicator Ic of the MPC-LOF is 0.3829 and is less than that of the MPC-NNA (0.6452), MPC-NSL (1.1689), and MPC-L (2.2194), indicating that the proposed method could provide higher control accuracy than MPC-NNA, MPC-NSL, and MPC-L.Figure 8
Manipulated variable and controlled variable under four control methods.


## 4.1. The Continuous Sterilization Process

Continuous sterilization is conducted on the medium before fermentation in biochemical engineering, and the nonsterilized medium is heated to the desirable temperature by steam in the sterilization process, so the undesired microorganisms are removed in the medium. After the continuous sterilization, necessary pure culture could be provided for the objective microorganism. Figure2 illustrates the diagrammatic sketch of the continuous sterilization system.Figure 2
The continuous sterilization process.In the continuous sterilization system, the quantity of steam used to heat the unsterilized medium is changed by the opening degree of the valve on the steam tube (Ov), while ΔOv(i)=Ov(i)-Ov(i-1) is the manipulated variable. The controlled output variable is the temperature of the heated medium at the outlet of the steam injector (Th). There exist measurable disturbances as the flow rate of medium (Fm), the temperature of unsterilized medium (Tm), and the temperature of steam (Ts). The control procedure aims to provide suitable heating on the medium, so that the temperature of the heated medium can be kept around the set value.


## 4.2. Experiment Settings

To evaluate the effectiveness of the proposed framework in nonlinear MPC, the plant data in a real workshop of continuous sterilization are collected and experiments are conducted based on the obtained data. The plant data includes the input variableOv, which adjusts the quantity of steam for the sterilization, while Tm,  Fm, and Ts are the measurable system disturbances, and Th is the output variable. The value ranges of Ov,  Tm,  Fm,  Ts, and Th are 0–100%, 0–70°C, 0–25 ton/h, 150–210°C, and 90–140°C, respectively.Modeling and controlling are conducted based on the plant data, and MPC based on different models are employed, including the linear model, the neural network, the adaptive neurofuzzy inference system (ANFIS) [44–48] with online linearisation (ANFIS-OL), and the linguistic model. With those models, four MPC algorithms are derived: MPC-L, MPC-NNA, MPC-NSL, and the proposed optimization framework for nonlinear MPC based on linguistic model, respectively. The proposed method is expressed as MPC-LOF for short.


## 4.3. Experiments on Modeling

In the modeling process, 788 samples are treated as training data, and 394 samples are treated as test data. The performances of the linear model, the neural network, the ANFIS with online linearisation, and the linguistic model are compared. For there are different parameters contained in different models, the parameters are selected by scanning for the best result based on the training data. All the involved parameters in the modeling are listed in Table1.Table 1
Parameters involved in modeling corresponding to different control methods.
Modeling method
H
r
r
N
n
Linear model
√
√
×
Linguistic model
√
√
×
Neural network
√
×
√
ANFIS-OL
√
√
×In Table1, the receding horizon Hr is the sampling length of the variables employed in modeling of the prediction model. For the linguistic model in the proposed framework, the radius value r affects the result of the Gaussian kernel clustering algorithm and the number of rules obtained. For the MPC-NNA, the three-layer back-propagation neural network is employed, and Nn means the number of nodes in the hidden layers. For the MPC-NSL, the ANFIS model is adopted. For the ANFIS, the subtractive clustering algorithm constructs the initial fuzzy model, and the radius r would affect the number of clusters obtained. In the following experiments, the parameters involved in the linguistic model are discussed, and the performances of different modeling methods are compared.The first experiment shows the effect ofHr and r on the modeling of the linguistic model. With the training data, experiment is conducted and the modeling results under different couples of Hr and r are shown in Figure 3(a). From Figure 3(a), it can be seen that Hr=5 tends to bring about lower mean square error (MSE) compared with the other values, and the best result is obtained when Hr=5 and r=0.5.Figure 3
Results on the modeling.
(a)
The effect ofHr and r on the results of linguistic modeling
(b)
Learning curves of the linguistic model withHr=5
(c)
Modeling result on training data based on different modeling methods
(d)
Modeling result on testing data based on different modeling methodsThe second experiment is conducted with fixedHr=5 to observe in detail the effect of r on the modeling result of the linguistic model. The result is shown is in Figure 3(b), where r in the middle of its range tends to bring about higher accuracy.To make comparison between different model methods, the parameters of the linear model, the ANFIS with online linearisation, and the neural network are selected by scanning for the best results based on the training data. The selected parameters are shown in Table2.Table 2
Parameters involved in modeling corresponding to different control methods.
Modeling method
H
r
r
N
n
Linear model
3
×
×
Linguistic model
5
0.5
×
Neural network
4
×
6
ANFIS-OL
4
0.6
×The performances of different methods with the selected parameters are shown in Figures3(c) and 3(d). From Figure 3(c), it can be seen that the neural network achieves the best performance than the linear model, the ANFIS with online linearisation (ANFIS-OL), and the linguistic model in the training process. The ANFIS-OL shows the worst accuracy, possibly because of low noise rejection. In Figure 3(d), the linguistic model provides the best modeling accuracy in the test data, while the neural network shows slightly lower accuracy than the linguistic model and the linear model, which shows lower generalization than the latter two methods. Yet, the ANFIS-OL method still shows the worst performance. Table 3 shows the detailed comparison between different modeling methods.Table 3
Performances of different modeling method on both training and test data.
Modeling method
Final MSE
Training data
Test data
Linguistic model
0.1524
0.4430
Linear model
0.2266
0.5026
Neural network
0.0823
0.5584
ANFIS-OL
11.5272
23.0631


## 4.4. Experiments on Controlling

To evaluate the effectiveness of the proposed control method, the control methods of MPC-L, MPC-NNA, MPC-NSL, and MPC-LOF are conducted to make comparison based on the plant data, which contain 1187 samples as training data and anther 1097 samples as test data. The process model used in those MPC methods is employed as the same neural network model generated by the plant data withHr=7. In the experiments, the performance indicator is the control accuracy, which represents the approximation of the controlled variable to the expected reference trajectory. The control accuracy is also the most important indicator in a practical continuous sterilization. The indicator of the control accuracy Ic is calculated as(36)Ic=1Nv∑t=1NvTht-Ret2,where Nv means the number of data samples, while Re(t) means the reference trajectory at the time instant t. Re(t) is set to be 123°C in the control process. Different MPC methods are conducted, with the involved parameters discussed. All the involved parameters in the controlling experiments are listed in Table 4, where λ is the penalty coefficients on the control effort and Hc is the control horizon.Table 4
Parameters involved in MPC methods.
Control method
Controller parameters
λ
H
c
H
r
r
N
n
MPC-L
√
√
√
×
×
MPC-NNA
√
×
√
×
√
MPC-NSL
√
√
√
√
×
MPC-LOF
√
√
√
√
×In the following experiments, the parameters involved in the MPC-LOF are discussed, and the performances of MPC-L, MPC-NNA, MPC-NSL, and the MPC-LOF are compared.The first experiment investigates the effects of parameters in MPC-LOF on the controlling performance, and the best combination of parameters in MPC-LOF is selected by scanning for the best result. With the training data, the optimal parameters combination is selected asHr=5,  Hc=5,  r=0.5, and  λ=3.78. The parameter selection in MPC-LOF is shown in Figures 4–7 for illustration. The effect of Hr and r on the performance of MPC-LOF is shown for illustration in Figure 4, where the indicator of TestingIc is employed to show the control performance in detail. The ith TestingIc is calculated as ∑t=1iTht-Ret2/i. It can be seen from Figure 4 that, when Hr=5 and λ=3.78, desirable result could be found when Hr=4or5, and the best result is obtained when Hr=5 and r=0.5. The final Ic corresponding to different combinations of Hr and r with fixed Hc=5 and λ=3.78 is detailed in Table 5.Table 5
The finalIc in MPC-LOF corresponding to different couples of Hr and r when Hc=5 and λ=3.78.
I
c
r
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
H
r
2
107.9560
265.6859
257.2303
252.4744
246.9592
245.1016
274.8260
254.3659
164.3792
3
224.1424
251.8802
282.2562
263.9272
6.3428
262.2148
265.0530
244.6227
2.0565
4
282.3905
280.4251
274.4064
258.4806
3.6191
212.4163
262.9407
662.4718
706.5165
5
9.3565
88.5872
29.2493
30.5592
1.6593
161.6644
291.7744
2.0702
87.0250
6
197.6085
187.1392
262.6352
259.0817
258.1882
258.5736
7.0774
39.0084
17.5822Figure 4
The effect ofHr and r on performance of MPC-LOF when Hc=5 and λ=3.78.Figure 5
The effect ofHc and r on performance of MPC-LOF when Hr=5 and λ=3.78, part 1.Figure 6
The effect ofHc and r on performance of MPC-LOF when Hr=5 and λ=3.78, part 2.Figure 7
The effect ofλ and r on performance of MPC-LOF when Hc=5 and Hr=5.When there isHr and λ=3.78, the best performance is found when Hc=5 and r=0.5, which performs better than the other couples of parameters as shown in Figures 5-6 for illustration. The details of the final Ic with different couples of Hc and r with fixed Hr=5 and λ=3.78 are shown in Table 6.Table 6
The finalIc in MPC-LOF corresponding to different couples of Hc and r when Hr=5 and λ=3.78.
I
c
r
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
H
c
2
593.6593
140.8252
283.2953
261.1969
35.1593
303.4396
267.9231
396.3824
439.1347
3
283.0946
137.2669
274.4064
258.3611
3.6191
212.4163
262.9395
662.4733
706.5170
4
297.0199
133.8319
263.1655
255.9483
1.7700
38.8709
259.8053
692.2658
738.5095
5
9.3565
88.5872
29.2493
30.5592
1.6593
161.6644
291.7744
2.0702
87.0250
6
338.4374
129.4590
262.1812
250.5345
31.3147
31.4507
250.3436
737.1146
766.6540
7
184.7563
125.0401
261.9619
248.3564
32.4505
33.1497
245.2064
753.7656
771.5174
8
187.6874
55.7913
261.8151
246.7770
32.7710
35.2300
241.9463
378.1504
733.8148
9
203.9504
7.2457
261.5909
245.5387
32.6305
35.8943
240.0754
373.3230
561.3975
10
186.2029
4.7204
261.3996
244.4238
32.9118
33.3927
238.0571
371.1604
570.3603
11
162.1932
3.9902
261.0869
243.2266
33.9692
34.4692
235.5092
334.7633
587.1572
12
142.0463
3.8393
260.5998
241.9458
34.0403
37.4656
231.6528
311.7689
548.6458
13
120.2090
3.6480
259.9733
240.5188
35.1123
37.9268
183.0111
174.5661
556.9916
14
80.8538
3.3041
258.9066
238.9687
35.6674
38.3702
169.9955
154.8572
541.8535
15
56.8366
2.8585
256.4588
237.1937
35.9254
38.7456
22.6742
150.4693
510.2067
16
61.1071
2.4136
254.7232
235.3396
36.1617
39.2147
9.7830
152.9248
354.7270Figure7 shows the effect of λ and r on the MPC-LOF, and the final Ic under different couples of λ and r is shown for comparison. When Hc=5 and Hr=5, the best result is found when λ=3.78 and r=0.5.The second experiment is conducted to show the performances comparison between the MPC-L, MPC-NNA, MPC-NSL, and the proposed MPC-LOF. For there are different parameters contained in different MPC methods, the parameters are selected by scanning for the best result based on the training data, whose selection results are shown in Table7.Table 7
Optimal parameters for four MPC methods.
Control method
Controller parameters
λ
H
c
H
r
r
N
n
MPC-L
3.16
5
3
×
×
MPC-NNA
4.26
×
3
×
6
MPC-NSL
0.50
5
4
0.2
×
MPC-LOF
3.78
5
5
0.5
×The plot of the manipulated variable (ΔOv) and controlled variable (Th) corresponding to the controller parameters in Table 7 are shown in Figure 8. Given the reference trajectory of 123°C, the controlled variable Th can be kept within the interval [122,124.5] under the MPC-NNA method, and there exists deviation from the reference trajectory in the controlled variable; for example, Hh is kept around 122°C from the 90th to the 140th iteration and rises to about 123.5°C in the following 200 iterations. There are no significant fluctuations in the trajectory of the controlled variable in MPC-NNA, but the controlled variable could not be held around the reference trajectory stably. As for the MPC-NSL method, the controlled variable shows better approximation to the reference compared with MPC-NNA, and Hh is kept within [121,125] generally. Yet there exists fluctuation during the control process; for example, the controlled variable fluctuation appears from the 700th to the 860th iteration, whose corresponding manipulated variable is also in severe fluctuation as shown in Figure 8. The fluctuations show low resistance of the online linearisation approach in MPC-NSL to the disturbances. Besides, there are deviations in the control process of MPC-NSL; for example, the controlled variable rises from 120.5°C to 127°C at the 20th iteration and rises from 123°C to 130°C at the 690th iteration. In MPC-L, the controlled output variable is generally kept within [120.5,125], and both severe deviation and fluctuation exist in the control process. For example, fluctuations appear from the 200th to the 320th iteration, and deviations exist at the 825th iteration with a decrease from 130°C to 117°C and the increase from 123°C to 128°C at the 690th iteration. Though the modeling accuracy of linear model is not the worst in the experiments on modeling, the bad control performance compared with the other methods shows the limit of the pure linear optimization framework in nonlinear system. As for the performance of the proposed MPC-LOF, Th is kept within [122.5,124]. There exist fluctuations in the controlled variable, yet they are much less than those in MPC-NSL and MPC-L, and the value of the controlled variable is kept closely around the reference compared with MPC-NNA. The value of indicator Ic of the MPC-LOF is 0.3829 and is less than that of the MPC-NNA (0.6452), MPC-NSL (1.1689), and MPC-L (2.2194), indicating that the proposed method could provide higher control accuracy than MPC-NNA, MPC-NSL, and MPC-L.Figure 8
Manipulated variable and controlled variable under four control methods.


## 5. Conclusions

In this paper, the data-driven optimization framework for nonlinear MPC is proposed, which employs the linguistic model as the prediction model. The proposed framework has some advantages as follows. First, the proposed framework solves a quadratic optimization problem in online control without the calculation of piecewise linearisation or online linearisation procedure for local linear models. Second, the controller parameters are tuned based on the plant data, which is free of subjective factors, and the optimal parameters result in effective optimization framework that could bring about high control accuracy. Third, the data-driven determination of controller parameters makes the proposed method applicable in practice and easy to implement. The effectiveness is evaluated in comparison with MPC-L, MPC-NNA, and MPC-NSL in the control of continuous sterilization, and the experimental results show that the proposed framework has higher control accuracy. In the further research, we intend to apply the proposed method into other process control practices like the fermentation process in biochemical industry.



---

*Source: 9402684-2017-08-28.xml*