SD-ARX modeling and robust MPC with variable feedback gain for nonlinear systems

As a generalized input-output model, the state-dependent exogenous variable autoregressive (SD-ARX) model has been intensively utilized to model complex nonlinear systems. Considering that more freedom can be provided by the state feedback control with variable feedback gain for constructing robust controllers, we propose a robust model predictive control (RMPC) algorithm with variable feedback gain on the basis of the SD-ARX model. First, the polytopic state space models (SSMs) of the system are constructed and the prediction accuracy of the SSMs is further improved by using the parameter variation rate information of the SD-ARX model. Then, an RMPC algorithm with variable feedback gain is synthesized for increasing the design freedom and enhancing the control performance. Two simulation examples, that is, the modeling and control of a continuous stirred tank reactor (CSTR) and a water tank system, are provided to demonstrate the feasibility and effectiveness of the proposed RMPC algorithm.


Introduction
Nonlinear systems with physical constraints can be controlled by a model-based control technique called model predictive control (MPC), which has been attracting widespread attention from both the academic sphere and the industrial domain, respectively. 1For the MPC, accurate system model is very important, because usually this kind of algorithm needs to solve an optimal control problem based on model prediction at each sampling interval.However, in actual industrial systems, it is difficult to obtain accurate explicit models due to various uncertainties. 2 So far, dealing with systemic uncertainties such as modeling errors and external disturbances is one of the main topics in MPC research. 3ince the dynamics of nonlinear and/or uncertain systems can be conveniently wrapped by the models in the polytopic form, the polytopic model has attracted much attention. 4Moreover, based on this type of model, robust model predictive control (RMPC) algorithms for the system with disturbances and/or uncertainties can also be conveniently constructed to achieve a reasonable performance.For a system with uncertain parameters, Kothare et al. 5 proposed a linear matrix inequality (LMI)-based RMPC algorithm for state feedback control optimization of the controller.For input-constrained linear parameter varying (LPV) systems, Lee et al. 6 proposed an RMPC algorithm based on novel sufficient conditions for the cost monotonicity.Lu and Arkun 7 proposed two quasi-min-max MPCs for polytopic linear parameter varying (PLPV) systems.For the input-constrained LPV system with unknown but bounded disturbances, He et al. 8 proposed an input-to-state stable RMPC algorithm.For the constrained PLPV systems with bounded disturbances, Casavola et al. 9 used a novel RMPC for addressing feedback regulation problems.Li et al. 10 proposed a new RMPC algorithm with constrained feedback for uncertain systems.Polytopic form model-based RMPCs have produced effective results, [5][6][7][8][9] but the majority of the research has aimed to deal with the problems of system's state or output tracking by assuming the given or accurately measurable steady states of a system. 11evertheless, the steady states of most systems remain unknown or beyond quantification in practicality, owing to extraneous disturbances or inescapable approximative disparities. 12o solve these problems, attention has been focused on system state estimator-based RMPC synthesis schemes. 13For polytopic systems, Ding et al. 13 proposed a novel feedback RMPC with a state estimatorbased synthesis scheme.For LPV systems, Park et al. 12 proposed a quasi-min-max RMPC algorithm.Song et al. 14 proposed a distributed RMPC with output feedback.For uncertain systems, Hu and Ding 15 proposed a dynamic output feedback RMPC.However, employing observers for the purpose of estimating the states may introduce additional complexity, and when the constraints of the variables such as the state, input, and output are activated, the closed-loop performance may be greatly degraded in terms of input disturbance rejection. 11,16t is important to directly construct an output feedback RMPC with an input-output type model, since modeling an actual system with an input-output equation is generally easier. 16As a generalized input-output model, the state-dependent exogenous variable autoregressive (SD-ARX) model is intensively applied in complex nonlinear system modeling, 17,18 and it belongs to a quasi-LPV model and its state-dependent coefficients and flexible model structure can be convenient for the design of subsequent robust controller. 19,202][23] The SD-ARX model-based RMPC for systems without disturbance was first studied in the work of Peng et al. 24 Based on this type of model, an RMPC with unknown steady-state information was proposed in Zhou et al. 25 To decrease the conservativeness of the LPV models in studies such as Peng et al, 24 Zhou et al. 25 Zhou et al. 26 proposed an RMPC synthesis method using the model's parameter boundary information, and based on this model, they 27 further proposed a one-stage scheduling RMPC.So far, many meaningful researches of SD-ARX model-based RMPC have been reported, [24][25][26][27] whereas the research has mainly focused on studying single-constant feedback RMPC algorithms.According to the research results in Lu et al., 2 Li et al., 10 the improvement of the performance and the enlargement of the feasible region with more freedom of the RMPC can be achieved if variable feedback gains are used, compared with the case of using constant single feedback gains.
Considering that more freedom can be provided by the state feedback control with variable feedback gain for constructing robust controllers, in this paper, we propose a SD-ARX model-based RMPC with variable feedback gain for nonlinear systems.The information of parameter variation rate in the process of constructing the system's state space models is further explored to both decrease the conservatism and increase the predication accuracy.3][14] The rest of the work is constructed as follows: the construction process of an identified SD-ARX model-based polytopic state space model (SSM) is presented in Section II, the design of the RMPC with variable feedback gain is introduced in Section III, the numerical example is illustrated in Section IV.Finally, the conclusions are described in Section V.

Design of polytopic state space models
Considering that the nonlinearity of many practical problems depends on one or more system variables, 17,18 the SD-ARX model 19 can be used to present the dynamics of the system.By approximating the functional coefficients with RBF networks, the following RBF-ARX model can be obtained 28 : where u(t) denotes the input of the system; y(t) denotes the output; j(t + 1) is the white noise item; n u is the maximal input lag; n y is the maximal output lag; the process variable w(t) = w(t), w(t À 1), Á Á Á , w(t À d + 1) ½ T is a measurable state vector that causes the operating point of the system to vary with time, in many cases, it contains maybe the output series and/or the input series 17,28 , p 0 (w(t)), p j, k (w(t))jj = y, u È É are the RBF network-type functional coefficients composed of multiple exp-functions;  1) can be conducted by using the variable projection algorithm, 20 the normal 17 or the regularized structured nonlinear parameter optimization method (SNPOM). 28f the bounded external disturbance of the system is considered, then the following more general system model can be obtained: where d is a known upper bound and z(t + 1) is the bounded external unknown disturbance, containing modeling errors.
To obtain the polytopic model of the system, we first define the input and output increments as follows: where i and j are integers less than or equal to zero, and y r denotes the expected output reference of the system.Then, the one-step forward prediction y(t + 1jt) can be derived as follows: a k+1,t y r +b k+1,t u(tÀkÀ1) ð Þ À y r +p 0,t +z(t+1jt) where z(t + 1jt) denotes the bounded unknown disturbance with z(t + 1jt) j j 4d.The term u t j j in (5) can be used to indicate the controlled system's steady state information, since it equals zero under the steady state for a perfect input u(t) and a stabilized output of y r . 26herefore, the future controlled system dynamics can be specified by designing model (6), where the ''perfect'' input increment set u(t), u(t + 1jt), . . .f gwill be constructed to achieve future items where j51; if j4k, u(t + j À kjt) = u(t + j À k), and y(t + j À kjt) = y(t + j À k).To establish the state space model of the system, first, the state vector of the system is defined as follows: where j50 and i 2 2, 3, . . ., k n f g .Thus, the state space models corresponding to the polynomial models ( 4) and ( 6) can be obtained as follows: and : Based on model (1) and state vector w(t) at the current time, it is clearly that X(tjt), A t , and B t in the state space model ( 8) can be obtained.From ( 8) and ( 5), one can see that X(t) is an unknown item.But it varies in a determined range.In this paper, a convex polytopic set O X is designed for wrapping the range of X(t), in which two vertices of O X are defined as follows: Remark 1: , which depend on the future state vectors.However, more often than not, the future state vectors w(t + jjt) cannot be obtained.Therefore, using the special SD-ARX structure of model ( 1), we can achieve the ranges of all exp-functions in the functional coefficients, and then we can construct convex polytopic sets for wrapping the variation range of the set A t + jjt , B t + jjt jj51 È É in model (9).
First, based on the current state matrix set A t , B t f g in ( 8) and the information of the parameter variation rates in model ( 1) calculated with the global identification data, we can construct two convex polytopic sets to wrap the future state matrices A t + jjt and B t + jjt in (9) as follows: l a a, j A j a , l a a, j 50, where L h = 2 h , h is the model order; and their vertices can be defined as follows:   are defined as follows: where where n p = 4

Synthesis method of the robust MPC with variable feedback gain
In this section, an output tracking RMPC with variable feedback gain are designed based on the SSMs of the system which is established in the previous section.First, the following cost function is utilized in the RMPC synthesis: where J ' 0 (t) = P ' denotes the positive constant; W denotes the positive definite diagonal matrix; x k k 2 R = x T Rx; X(t + jjt) denotes the state predication at time t; u(t + jjt) denotes the input increments constrained by u(t + jjt) j j 4 u max ; and U(t) = u(tjt), u(t + 1jt), ::: f g .Then, by dividing function (18) into two parts, one can obtain: where W and u(tjt) is the free decision variable.
In this paper, the control strategy for the RMPC algorithm with variable feedback gain is designed as follows: where F j is the variable feedback gain with the infinite control horizon of j, and N53.For j5N À 1, F j in the control strategy L is always set to be F NÀ1 .
Next, the Lyapunov function of the robust controller is designed as follows: where j51, and if j .N À 1, P j is always set to be P NÀ1 .
From the selected Lyapunov function V(t + j), the robust stability conditions of the RMPC can be designed as follows: Add inequality (22) from j = 1 to ', one can obtain the upper limit of J ' 1 (t): Therefore, the optimization problem ( 18) can be rewritten as follows: Thus, as shown in Theorem 1, the optimization problem ( 18) can be converted into a linear programing problem (LPP) for solution.
Theorem 1: The optimizing problem ( 18) with the control strategy as in ( 20) is identical with the following LMI constrained convex optimization problem: subject to where R and W denote weights;A t , B t , X(tjt), X t, 1 , X t, 2 , and A 1 m , B 1 m , A m , B m jm = 1, :::, n p n o can be computed based on the obtained historical system data and model (1).Therefore, the control input subjected to the control system at time t is u(t) = u(tjt) + u(t À 1), where u(tjt) is the solution of the LPP in Theorem 1.
Thus, we can convert the form of the minimization problem (24) into the following form: Considering that the vector X(t) in ( 8) belongs to the convex polytopic set O X and introducing ( 8) into (36), the inequality J 1 0 (t) + X(t + 1jt) k k 2 P 1 4g can be finally expressed as LMIs (28).
Moreover, referring to Zhou et al., 25 one can convert the input and input increment constraints imposed in (36) into the LMIs (29-30).
In the end, the min-max optimization problem shown in ( 18) is transformed into the minimization problem min g, u(tjt), Y j , G j , Q j , Z j jj = 1, ..., NÀ1 f g g that subjects to the LMIs (26-30).
The computation process of the proposed RMPC algorithm is summarized as follows.
Step 1: SD-ARX model ( 1) is used to model the nonlinear system, and the regularized-SNPOM is used to optimize the parameters of the model (1).
Step 2: Initialize controller parameters W and R, predictive feedback horizon N, upper limit of disturbance d and constraint value u max in (29-30).Step 3: At sampling time t, based on the identified model ( 1) and the current system state vector w(t), calculate A t , B t , X(tjt) according to (7-8) and calculate X t, 1 , X t, 2 using (10).
Step 4: Based on the identified model ( 1) and the historical state vectors obtained from the global identification data, calculate vertices (17) according to (11-16).
Step 5: Solve the optimization problem (25) to get the input increment u(tjt).
Step 6: Calculate control law u(t) = u(tjt) + u(t À 1) and implement it on the controlled system at the current time; Increase t by 1 and go to Step 2.
Theorem 2: When the robust MPC synthesis method is implemented with a rolling-horizon strategy, the stability can be guaranteed by the feasibility of the LPP (25).
Remark 4: Based on the uncertain SSMs (8-9) and Theorems 1 and 2, when there is a feasible solution for the LPP (25) at time t, one can achieve the outputtracking RMPC without steady-state offset to nonlinear system (2) because there is a u s satisfying u(t) !u s as t !' since y(t) !y r when X(t) !0 and u(t) !0 as t !'.

Case I -CSTR process
Considering that continuous stirred tank reactor (CSTR) 29 is very important in the core equipment for chemical production, a CSTR process is selected and utilized to verify the performance of the RMPC strategy.The selected CSTR process is an irreversible firstorder exothermic reaction (A !B), and its characteristics described by continuous-time nonlinear differential equations are as follows 30 : where T(t) is the reaction temperature; q c (t) is the coolant flow rate; Ca(t) is the product concentration of A; T c0 is the inlet coolant temperature; T 0 is the feed temperature; q is the process flow rate; Ca 0 is the feed concentration; k 0 is the reaction rate constant; v is the reactor volume; E=R is the activation energy term; k 1 = (DHk 0 ) (rC p ), k 2 = 1=v and k 3 = hA (rC p ) are the constant coefficients of the chemical reaction; DH is the heat of reaction; r is the liquid density; C p is the specific heats; hA is the heat transfer term.The CSTR parameters in (37) are shown in Table 1.In this process, the reaction temperature T(t) is measured in real time and controlled by the manipulated variable q c (t) optimized by the proposed robust MPC algorithm.
Polytopic SSM construction process.First, the basic model structure of the process can be selected as: ,j=y,u: where the system input and the corresponding output are the coolant flow rate q c (t) and the reaction temperature T(t), respectively; the measurable state vector is chosen as w(t) = T(t), T(t À 1), Á Á Á , T(t À d + 1) ½ T ; and p 0 (w(t)), f p j, k (w(t))jj = y, u É are the functional coefficients composed of multiple exp-functions.The 38) are optimized by the regularized-SNPOM. 28n this paper, to get a set of globally varying input/ output data, a set of sine signals and step signals are set to be the required outputs under the PID controller for the CSTR, which make the input and output of the system adequately vary within a large range, and a set of Gaussian white noise signal with small power is also added to the inputs in order to excite the systems' various dynamic modes.The global identification data of the CSTR process with a sampling period of 0.05 min are shown in Figure 1.We use the first and the last 2500 data points for training and testing of the model, respectively.In the modeling process, the step response performance of the identified model is evaluated as the criterion for determining the model orders in (38) and the final selected model orders in (38) are n y = 2, n u = 2, h = 1 and d = 1.
Based on the state vector defined in (7) and the global identification data of the process, the SSMs of the system can be achieved at each sampling instant.Their poles, varying with the variable w(t), are indicated by different colors, as shown in Figure 2. It is obvious that the characteristics of the process vary with variable w(t).Therefore, the parameter changing rate and the parameter boundary of the system can be obtained through the identification process.Then, we can obtain the previously designed polytopic SSMs (8-9) in Section II to wrap the dynamics of the system at each sampling time.Finally, based on these models, the proposed RMPC algorithm can be constructed and used for the process control.
Control results and analysis of the RMPC algorithms.The numerical experiments are conducted under the MATLAB R2017b running on a computer equipped with an Intel Core i7-3770 CPU and 16 GB RAM.In this CSTR process, the constraints of the control variable q c (t) and the input increment are selected as 104q c (t)4150 and u max = 60L=min, respectively.Just as in (18), the constraints of the control variable u(t) is not considered in the cost function of the proposed RMPC algorithm.So, the constraint on the input variable in this paper is the physical constraint of the system.After a series of trials and errors, the determined controller weights of our RMPC are W = diag(1, 1) and R = 0:005.
Figure 3 shows the control results of our RMPC for several predictive feedback horizons N (as in (20)).The pink dotted line in Figure 3 represents the expected output trajectory of the system.A square wave disturbance with a 2K amplitude and a 0.5 min duration is added into the system output at the fourth min and ninth min for testing the anti-jamming performance of the proposed RMPC.Therefore, the upper limit d, as in (2), is set to 2 for the uncertain system disturbance z(t + 1). Figure 4 shows the real-time computational burden corresponding to the control results in Figure 3.
From Figure 3, it can be seen that the control results of our RMPC are improved with the increase of the predictive feedback horizon N, especially the antiinterference performance of the system when an external interference is introduced.Moreover, when the feedback horizon N is .5, the improvement of the control effect becomes weak.From Figure 4, one can see that the online computational cost corresponding to the results in Figure 3 also increases obviously with increasing the predictive feedback horizon N.This paper comprehensively considers the control performance of the RMPC algorithm and the online computing burden of the algorithm, the predictive feedback horizon N is finally selected as 7 for the CSTR.
In this subsection, two RBF-ARX model-based RMPC algorithms suggested in Peng et al., 24 Zhou et al. 26 (named RMPC1 and RMPC2, respectively) are also performed for comparative research.To be fair, the same RBF-ARX model obtained previously is used in both RMPC algorithms.After a series of trials and errors, the well-adjusted parameters of the three RMPCs are presented in Table 2.
The control results of the 3 RMPCs are shown in Figure 5, and the corresponding online computational burden is shown in Figure 6.
From Figure 5, one can see that the control results of the RMPC1 are worse than those of the others.The control results of the RMPC2 and RMPC(N = 7) are not significantly different when the system has no external interference.However, when the system is interfered with by the square wave disturbance, the RMPC(N = 7) has a better anti-interference performance than the RMPC2.From Figure 6, it is also clearly that the RMPC(N = 7) shows a heavier online computing burden compared with the other two algorithms.However, its maximal value is still far smaller    than that of the control system with a 0.05 min sampling period.

Case II -water tank system
Level control of the water tank system is often used in literatures 30,31 for comparison study.The mass and energy balance equations for the water tank system are given in Zhou et al. 30 In this case study, the control problem is to control the level H 2 (t) of the lower tank by adjusting the input voltage U(t) of the pump.The detailed RBF-ARX modeling process and results of the water tank system have been given in Zhou et al., 30 that is, the estimated RBF-ARX model orders are n y = 5, n u = 4, h = 1 and d = 2. Based on this identified model, the proposed RMPC algorithm can be obtained and used for the level control of the water tank system.
The numerical experiments are conducted under the same software and hardware environment as in the previous subsection.The sampling period of the control system is 1.5 s.In this water tank system, the constraints of the control variable H 2 (t) and the input increment U(t) are selected as 24H 2 (t)4100 and u max = 150, respectively.After a series of trials and errors, the determined controller weights of our RMPC are W = diag (1, 1, 1, 1, 1) and R = 0:001.Similar to the previous section, the RMPC1 and RMPC2 based on the same RBF-ARX model (n y = 5, n u = 4, h = 1 and d = 2) are also performed for comparative research.After a series of trials and errors, the well-adjusted parameters of the three RMPCs are presented in Table 3.
The control results of the three RMPCs for the water tank system are shown in Figure 7, and the corresponding online computational burden is shown in Figure 8.
From Figure 7, one can see that the control results of the RMPC(N = 3) are better than those of the   Figure 8. Online computational burden corresponding to the results in Figure 7.
others, especially in terms of the step response performance.From Figure 8, it is also clearly that the RMPC(N = 3) shows a heavier online computing burden compared with the others.However, its maximal value is still a little smaller than that of the system sampling period 1.5 s.It is worth noting that when N is .3, the calculation time of the RMPC has exceeded the sampling time of the water tank system, so the predictive feedback horizon N is finally selected as 3 for the water tank system.In general, for the actual system, it is recommended to choose the maximum value of N that can complete the RMPC algorithm execution within a sample time of the system as the infinite control horizon.

Conclusions
In this paper, a SD-ARX model-based RMPC algorithm with variable feedback gain for nonlinear systems were proposed.First, the construction process of the polytopic SSMs of the system was introduced.Then, the synthesis method of the proposed RMPC was designed.In general, the proposed approach has two advantages: by using the information of the variation rate of the SD-ARX model parameters, our scheme has the ability to shrink the variation regions of the state matrices in the system's polytopic SSMs, and by using the control policy with variable feedback gain, it also presents the ability to enlarge the robust controller design freedom.The feasibility and effectiveness of the algorithm were verified by two numerical examples.However, the limitation of the RMPC is its relatively large computation burden, especially as the predictive feedback horizon N increases.In order to reduce the online computation burden, we shall try to study a partially offline computing and partially online synthesizing RMPC in the future research.

Declaration of conflicting interests
The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.
of the RBF networks; and z j m and l j m are the RBF network centers and scaling parameters, respectively.The optimization of the linear parameters u L = c 0 where e u m (t), e u, m , e u, m , D e u, m , De u, m ju = y, u; m = 1, :::, h n o 1, :::, n p )

Figure 1 .
Figure 1.Global identification data of the process.

Figure 3 .
Figure 3.Control results of our RMPC for selected predictive horizons.

Figure 4 .
Figure 4. Online computational burden corresponding to the results in Figure3.

2 Figure 5 .
Figure 5.Control results of the three RMPCs for CSTR.
= e u, m for all j52.So, the corresponding vertices A m , B m of the convex polytopic set defined as O 1 D , to which the set A t + jjt , B t + jjt D that wraps the state matrix set A t + jjt , B t + jjt jj52 È É .Then their corresponding polytopic LPV models can be obtained to wrap the dynamics of the SSMs(8-9).
h and the state matrices A j m , B j