 Open access peer-reviewed chapter

# Fuzzy–neural Model Predictive Control of Multivariable Processes

Written By

Michail Petrov, Sevil Ahmed, Alexander Ichtev and Albena Taneva

Submitted: October 15th, 2010 Reviewed: February 3rd, 2011 Published: July 5th, 2011

DOI: 10.5772/16828

From the Edited Volume

Edited by Tao Zheng

Chapter metrics overview

View Full Metrics

## 1. Introduction

Predictive control is a model-based strategy used to calculate the optimal control action, by solving an optimization problem at each sampling interval, in order to maintain the output of the controlled plant close to the desired reference. Model predictive control (MPC) based on linear models is an advanced control technique with many applications in the process industry (Rossiter, 2003). The next natural step is to extend the MPC concept to work with nonlinear models. The use of controllers that take into account the nonlinearities of the plant implies an improvement in the performance of the plant by reducing the impact of the disturbances and improving the tracking capabilities of the control system.

In this chapter, Nonlinear Model Predictive Control (NMPC) is studied as a more applicable approach for optimal control of multivariable processes. In general, a wide range of industrial processes are inherently nonlinear. For such nonlinear systems it is necessary to apply NMPC. Recently, several researchers have developed NMPC algorithms (Martinsen et al., 2004) that work with different types of nonlinear models. Some of these models use empirical data, such as artificial neural networks and fuzzy logic models. The model accuracy is very important in order to provide an efficient and adequate control action. Accurate nonlinear models based on soft computing (fuzzy and neural) techniques, are increasingly being used in model-based control (Mollov et al., 2004).

On the other hand, the mathematical model type, which the modelling algorithm relies on, should be selected. State-space models are usually preferred to transfer functions, because the number of coefficients is substantially reduced, which simplifies the computation; systems instability can be handled; there is no truncation error. Multi-input multi-output (MIMO) systems are modelled easily (Camacho et al., 2004) and numerical conditioning is less important.

A state-space representation of a Takagi-Sugeno type fuzzy-neural model (Ahmed et al., 2010; Petrov et al., 2008) is proposed in the Section 2. This type of models ensures easier description and direct computation of the gradient control vector during the optimization procedure. Identification procedure of the proposed model relies on a training algorithm, which is well-known in the field of artificial neural networks.

Obtaining an accurate model is the first stage of the of the NMPC predictive control strategy. The second stage involves the computation of a future control actions sequence. In order to obtain the control actions, a previously defined optimization problem has to be solved. Different types of objective and optimization algorithms (Fletcher, 2000) can be used in the optimization procedure. Two different approaches for NMPC are proposed in Section 3. They consider the unconstrained and constrained model predictive control problem. Both of the approaches use the proposed Takagi-Sugeno fuzzy-neural predictive model.

The proposed techniques of fuzzy-neural MPC are studied in Section 4, by experimental simulations in Matlab® environment in order to control the levels in a multi tank system (Inteco, 2009). The case study is capable to show how the proposed NMPC algorithms handle multivariable processes control problem.

## 2. Multivariable fuzzy-neural predictive model

The Takagi-Sugeno fuzzy-neural models are powerful modelling tools for a wide class of nonlinear systems. Fuzzy reasoning is capable of handling uncertain and imprecise information while neural networks can learn from samples. Fuzzy-neural networks combine the advantages of both artificial intelligent techniques and incorporate them in adaptive features. Those futures, based on a real time learning algorithm are the main advantage of the fuzzy-neural models.

The importance of the used in MPC strategy models and their adaptive characteristics is obvious. The accuracy of the model determines the accuracy of the control action. The proposed fuzzy-neural model is implemented in a classical NMPC scheme (Fig. 1) as a predictor (Camacho et al., 2004).

In this chapter a nonlinear discrete time state-space implementation is considered to represent the system dynamic:

x(k+1)=fx(x(k),u(k))y(k)=fy(x(k),u(k))E1

where x(k)n, u(k)mand y(k)qare state, control and output variables of the system, respectively. The unknown nonlinear functions fxand fycan be approximated by Takagi-Sugeno type fuzzy rules in the next form:

Rl:ifz1(k)isMl1andandzi(k)isMliandzp(k)isMlpthen{xl(k+1)=Axl(k)+Blu(k)yl(k)=Clx(k)+Dlu(k)E2

where Rlis the l-thrule of the rule base. Each rule is represented by an if-thenconception. The antecedent part of the rules has the following form “zi(k) is Mli” where zi(k)is an i-thlinguistic variable (i-thmodel input) and Mliis a membership function defined by a fuzzy set of the universe of discourse of the input zi.Note that the input regression vector z(k)pin this chapter contains the system states and inputs z(k)=[x(k) u(k)]T. The consequent part of the rules is a mathematical function of the model inputs and states. A state-space implementation is used in the consequent part of Rl, where Aln×n, Bln×m, Clq×nandDlq×mare the state-space matrices of the model (Ahmed et al., 2009).

The states in the next sampling time x^(k+1)and the system output y^(k)can be obtained by taking the weighted sum of the activated fuzzy rules, using

x^(k+1)=l=1Lμ¯yl(k)(Alx(k)+Blu(k))y^(k)=l=1Lμ¯yl(k)(Clx^(k)+Dlu(k))E3

On the other hand the state-space matrices A, B, C, and Dfor the global state-space plant model could be calculated as a weighted sum of the local matrices Al, Bl, Cl, and Dlfrom the activated fuzzy rules (2):

A(k)=l=1LAlμ¯yl(k)B(k)=l=1LBlμ¯yl(k)C(k)=l=1LClμ¯yl(k)D(k)=l=1LDlμ¯yl(k)E4

where μ¯yl=μyl/l=1Lμylis the normalized value of the membership function degree μylupon the l-thactivated fuzzy rule and Lis the number of the activated rules at the moment k.

Fuzzy implication in the l-thrule (2) can be realized by means of a product composition

μyl=i=1pμijE5

where μijspecifies the membership degree (Fig. 2) upon the activated j-thfuzzy set of the corresponded i-thinput signal and it is calculated according to the chosen here Gaussian membership function (6) for the l-thactivated rule:

μij(zi)=exp(zicGij)2σij22E6

where ziis the current input value of the i-thmodel input, cGijis the centre (position) and σijis the standard deviation (wide) of the j-thmembership function (j=1, 2,.., s) (Fig.2).

### 2.1. Identification procedure for the fuzzy–neural model

The proposed identification procedure determines the unknown parameters in the Takagi-Sugeno fuzzy model, i.e. the parameters of membership functions, according to their shape and the parameters of the functions fxand fyin the consequent part of the rules (2). It is realised by a five-layer fuzzy-neural network (Fig. 3). Each of the layers performs typical fuzzy logic strategy operations:

Layer 1. The first layer represents the model inputs through its own input nodes Z1, Z2, …, Zp. The network synaptic weights are set to one, so the model inputs are directly passed through the nodes to the next layer. Neurons here are represented by the elements of the regression vector z(k).

Layer 2. The fuzzification procedure of the input variables is performed in the second layer. The weights in this layer are the parameters of the chosen membership functions. Their number depends on the type and the number of the applied functions. All these parameters aijare adjustable and take part in the premise term of the Takagi-Sugeno type fuzzy rule base (2). In that section the membership functions for each model input variable are represented by Gaussian functions (Fig. 2). Hence, the adjustable parameters aijare the centres cGijand standard deviations σijof the Gaussian functions (6). The nodes in the second layer of the fuzzy-neural architecture represent the membership degrees μij(zi)of the activated membership functions for each model input zi(k)according to (6).The number of the neurons depends on the number of the model inputs pand the number of the membership functions sin corresponding fuzzy sets. It is calculated as.

Layer 3. The third layer of the network interprets the fuzzy rule base (2). Each neuron in the third layer has as many inputs as the input regression vector size p. They are the corresponding membership degrees for the activated membership functions calculated in the previous layer. Therefore, each node in the third layer represents a fuzzy rule Rl, defined by Takagi-Sugeno fuzzy model. The outputs of the neurons are the results of the applied fuzzy rule base.

Layer 4. The fourth layer implements the fuzzy implication (5). Weights in this layer are set to one, in case the rule Rlfrom the third layer is activated, otherwise weights are zeros.

Layer 5. The last layer (one node layer) represents the defuzzyfication procedure and forms the output of the fuzzy-neural network (3). This layer also contains a set of adjustable parameters – βl.These are the parameters in the consequent part of Takagi-Sugeno fuzzy model (2). The single node in this layer computes the overall model output signal as the summation of all signals coming from the previous layer.

I5=l=1LfylμylorI5=l=1LfxlμylorO5=l=1Lfylμyll=1LμylorO5=l=1Lfxlμyll=1LμylE7

where fxl=Alx(k)+Blu(k)andfyl=Clx(k)+Dlu(k).

### 2.2. Learning algorithm of the fuzzy–neural model

Two-step gradient learning procedure is used as a learning algorithm of the internal fuzzy-neural model. It is based on minimization of an instant error function EFNN. At time kthe function is obtained from the following equation

EFNN(k)=ε2(k)/2E8

where the error ε(k) is calculated as a difference between the controlled process output y(k)and the fuzzy-neural model output ŷ(k):

ε(k)=y(k)y^(k)E9

During step one of the procedure, the consequent parameters of Takagi-Sugeno fuzzy rules are calculated according to summary expression (10) (Petrov et al., 2002).

βl(k+1)=βl(k)+η(EFNNβl)E10

where ηis a learning rate and βlrepresents an adjustable coefficient aij, bij, cij, dij(11) for the activated fuzzy rule Rl(2). The coefficients take part in the state matrix Al, control matrix Bland output matrices Cland Dlof the l-thactivated rule (Ahmed et al., 2009). The matrices approximate the unknown nonlinear functions fxand fyaccording to defined fuzzy rule model (2). The matrix dimensions are specified by the system parameters – numbers of inputs m, outputs qand states nof the system.

Al=[a11a1nan1ann]Bl=[b11b1mbn1bnm]Cl=[c11c1ncq1cqn]Dl=[d11d1mdq1dqm]E11

In order to find a weight correction for the parameters in the last layer of the proposed fuzzy-neural network the derivative EFNNβlof the instant error should be determined. Following the chain rule, the derivative is calculated considering the expressions (7) and (8)

EFNNβl=EFNNy^y^I5I5βlE12

After the calculation of the partial derivatives, the matrix elements for each matrix of the state-space equations corresponding to the l-thactivated rule (2) are obtained according to the summary expression (12) (Petrov et al., 2002; Ahmed et al., 2010):

aij(k+1)=aij(k)+ηε(k)μ¯yl(k)xi(k)i=j=1÷nbij(k+1)=bij(k)+ηε(k)μ¯yl(k)uj(k)i=1÷n,j=1÷mcij(k+1)=cij(k)+ηε(k)μ¯yl(k)xj(k)i=1÷q,j=1÷ndij(k+1)=dij(k)+ηε(k)μ¯yl(k)uj(k)i=1÷q,j=1÷mE13

The proposed fuzzy-neural architecture allows the use of the previously calculated output error (8) in the next step of the parameters update procedure. The output error EFNNis propagated back directly to the second layer, where the second group of adjustable parameters are situated (Fig. 3). Depending on network architecture, the membership degrees calculated in the fourth and the second network layer are related as μylμij. Therefore, the learning rule for the second group adjustable parameters can be done in similar expression as (10):

αij(k+1)=αij(k)+η(EFNNαij)E14

where the derivative of the output error EFNNis calculated by the separate partial derivatives:

EFNNαij=EFNNy^y^μijμijαijE15

The adjustable premise parameters of the fuzzy-neural model are the centre cGijand the deviation σijof the Gaussian membership function (6). They are combined in the representative parameter aij, which corresponds to the i-thmodel input and its j-thactivated fuzzy set. Following the expressions (14) and (15) the parameters are calculated as follows (Petrov et al., 2002; Ahmed et al., 2010):

cGij(k+1)=cGij(k)+ηε(k)μ¯yl(k)[fyly(k)][zi(k)cGij(k)]σij2(k)E16
σij(k+1)=σij(k)+ηε(k)μ¯yl(k)[fyly(k)][zi(k)cGij(k)]2σij3(k)E17

The proposed identification procedure for the fuzzy-neural model could be summarized in the following steps (Table 1).

 Step 1. Initialize the membership functions – number, shape, parameters; Step 2. Assign initial values for the network inputs; Step 3. Start the algorithm at the current moment k; Step 4. Fuzzify the network inputs and calculate the membership degrees upon the activated fuzzy set of the membership functions according to (6); Step 5. Perform fuzzy implication according to (5); Step 6. Calculate the fuzzy-neural network output, which is represented by state-space description of the modelled system – (3) and (4); Step 7. Calculate the instant error according to (8) and (9); Step 8. Start training procedure for fuzzy-neural network; Step 9. Adjust the consequent parameters according to (13); Step 10. Adjust the premise parameters according to (16) and (17). Repeat the algorithm from Step 3 for each sampling time.

### Table 1.

Fuzzy-neural model identification procedure

## 3. Optimization algorithm of multivariable model predictive control strategy

The model provided by the Takagi-Sugeno type fuzzy-neural network is used to formulate the objective function for the optimization algorithm and to calculate the future control actions. The second stage of the predictive control strategy includes an optimization procedure. It utilizes the obtained results during the first (modelling) stage predictive model of the system. Using the Takagi-Sugeno fuzzy-neural model (3), the optimization algorithm computes the future control actions at each sampling period, by minimizing the typical for MPC strategy (Generalized Predictive Control – GPC) cost function (Akesson, 2006):

J(k)=i=HwHp+Hw1y^(k+i)r(k+i)Q˜2+i=0Hu1Δu(k+i)R˜2E18

where ŷ(k), r(k)and ∆u(k)are the predicted outputs, the reference trajectories, and the predicted control increments at time k, respectively. The length of the prediction horizon is Hp, and the first sample to be included in the horizon is Hw. The control horizon is given by Hu.Q˜≥0and R˜>0are weighting matrices representing the relative importance of each controlled and manipulated variable and they are assumed to be constant over the Hp.

The cost function (18) may be rewritten in a matrix form as follows

J(k)=Y(k)T(k)Q2+ΔU(k)R2E19

where Y(k), T(k), ∆U(k), Qand Rare predicted output, system reference, control variable increment and weighting matrices, respectively,

Y(k)=[y^(k|k)y^(k+Hp-1|k)],T(k)=[r(k|k)r(k+Hp-1|k)],ΔU(k)=[Δu(k|k)Δu(k+Hu-1|k)]E20
Q=[Q˜(1)00Q˜(Hp)]E21
R=[R˜(1)00R˜(Hu)]E22

The linear state-space model used for Takagi-Sugeno fuzzy rules (2) could be represented in the following form:

x^(k+1)=Ax(k)+Bu(k1)+BΔu(k)y^(k)=Cx^(k)+Du(k1)+DΔu(k)E23

Based on the state-space matrices A, B, Cand D(4), the future state variables are calculated sequentially using the set of future control parameters:

x^(k+1)=Ax(k)+Bu(k1)+BΔu(k)x^(k+2)=A2x(k)+(AB+B)u(k1)+(AB+B)Δu(k)+BΔu(k+1)x^(k+3)=A3x(k)+(A2B+AB+B)u(k1)+(A2B+AB+B)Δu(k)+(AB+B)Δu(k+1)+BΔu(k+2)x^(k+j)=Ajx(k)+i=0j1AiBu(k1)+i=0j1AiBm=0j1iΔu(k+m)x^(k+Hp)=AHpx(k)+i=0Hp1AiBu(k1)+i=0Hp1AiBΔu(k)+i=0Hp2AiBΔu(k+1)++AHpHuBΔu(k+Hu1)E24

The predictions of the output y^for jsteps ahead could be calculated as follows

y^(k+1)=Cx(k+1)+Du(k+1)=CAx(k)+(CB+D)u(k1)+(CB+D)Δu(k)+DΔu(k+1)y^(k+2)=CA2x(k)+(CAB+CB+D)u(k1)+(CAB+CB+D)Δu(k)+(CB+D)Δu(k+1)+DΔu(k+2)E25
y^(k+j)=CAjx(k)+(Ci=0j1AiB+D)u(k1)+(Ci=0j1AiB+D)Δu(k)+(CB+D)Δu(k+j1)+DΔu(k+jy^(k+Hp1)=CAHp1x(k)+(Ci=0Hp2AiB+D)u(k1)+(Ci=0Hp2AiB+D)Δu(k)++(Ci=0Hp3AiB+D)Δu(k+1)++(Ci=0HpHu1AiB+D)Δu(k+Hu1)E26

The recurrent equation for the output predictionsy^(k+jp), where jp= 1, 2,..., Hp–1, is in the next form:

.
y^(k+jp)=CAjpx(k)+(Ci=0jp1AiB+D)u(k1)+{i=0jp(Cj=0jp1AjB+D)Δu(k+i),jp<Hui=0Hu1(Cj=0jpi1AjB+D)Δu(k+i),jp>HuE27

The prediction model defined in (27) can be generalized by the following matrix equality

Y(k)=Ψx(k)+Γu(k-1)+ΘΔU(k)E28

where

Ψ=[CCACA2CAHp1]E29
Γ=[DCB+DCAB+CB+DCi=0Hp2AiB+D]E30
Θ=[D00CB+DDCAB+CB+DCB+D0Ci=0Hu2AiB+DDCi=0Hp2AiB+DCi=0HpHu1AiB+D]E31

All matrices, which take part in the equations above, are derived by the Takagi-Sugeno fuzzy-neural predictive model (4).

It is also possible to define the vector

E(k)=T(k)-Γu(k-1)-ΘΔU(k)E32

This vector can be thought as a tracking error, in the sense that it is the difference between the future target trajectory and the free responseof the system, namely the response that would occur over the prediction horizon if no input changes were made, i.e. ∆U(k)=0. Hence, the quantity of the so called freeresponse F(k)is defined as follows

F(k)=Ψx(k)+Γu(k-1)E33

### 3.1. Unconstrained model predictive control

In this section, the study is focused on the optimization problem of the unconstrained nonlinear predictive control with the quadratic cost function (18). The section presents an approximate solution of the problem where the information given by the obtained fuzzy-neural model is used to solve the problem.

The unconstrained optimization problem can be formulated in a matrix form. First, the predictor can be constructed as follows

Y(k)=ΘΔU(k)+F(k)E34

Second, the cost function (19) can be rewritten as

J(ΔU)=(TY)TQ(TY)+ΔUTRΔUE35

Hence, substituting the predictive model (25) into expression (26), the cost function of the model predictive optimization problem can be specified as follows:

J(ΔU)=ΔUT(ΘTQΘ+R)ΔU+2(FT)TQΘΔU+(TF)TQ(TF)E36

The minimum of the function J(∆U)can be obtained by calculating the input sequence ΔUso that ∂J/∂ΔU= 0:

ΔUJ(ΔU)=ΔU[ΔUT(ΘTQΘ+R)ΔU+2(FT)TQΘΔU+(FT)TQ(FT)]=0E37

Then the optimal sequence ΔU*is

ΔU*=(ΘTQΘ+R)1ΘTQ(TF)E38

The input applied to the controlled plant at time kis computed according to the receding horizon principle, i.e. the first element from the control sequence Δu*(k)of the vector ΔU*is taken. Then, control signal is calculated from:

u(k)=u(k-1)+Δu*(k)E39

It is evident that the expression given by the matrix equation (29) is the same as expression obtained for the generalized predictive control. However, in the GPC formulation the components involved in the calculation of the formula (29) are obtained from a linear model. In the present case the components introduced in this expression are generated by the designed nonlinear fuzzy-neural model. A more rigorous formulation of (29) will be representation of the components as time-variant matrices, as they are shown in the expression (22). In this case the matrix Θ(k)and the vectors Ψ(k), T(k)are being reconstructed at each sampling time. The vector Ψ(k)is obtained by simulating the fuzzy model with the current input u(k);the matrix Θ(k)is also rebuilt using a method described below.

ΔU*(k)=[ΘT(k)QΘ(k)+R]1ΘT(k)Q[T(k)F(k)]E40

The proposed method solves the problem of unconstrained MPC. A system of equations is solved at each sampling time k.The proposed approach decreases computational burden avoiding the necessity to inverse the gain matrix in (31) at each sampling time k.

Applying this method, minimization of the GPC criterion (18) is based on a calculation of the gradient vector of the criterion cost function Jat the moment ksubject to the predicted control actions:

J(k)=[J(k)Δu(k),J(k)Δu(k+1),,J(k)Δu(k+Hu1)]TE41

Each element of this gradient vector (32) can be calculated using the following derivative matrix equation:

J(k)ΔU(k)=[2[T(k)Y(k)]TQY(k)ΔU(k)+2ΔU(k)TRΔU(k)ΔU(k)]E42

From the above expression (33) it can be seen that it is necessary to obtain two groups of partial derivatives. The first one is[Y(k)ΔU(k)], and the second one is[ΔU(k)ΔU(k)]. The first partial derivatives in (33) have the following matrix form:

Y(k)ΔU(k)=[y^(k+Hw)Δu(k)y^(k+Hw)Δu(k+Hu1)y^(k+Hp+Hw1)Δu(k)y^(k+Hp+Hw1)Δu(k+Hu1)]E43

For computational simplicity assume that Hw=0 (18). Then each element of the matrix (34) is calculated by the expressed equations according to the Takagi-Sugeno rules consequents (2). For example the derivatives from first column of the matrix (34) have the following form:

y^(k)Δu(k)=l=1LDlμ¯yl(k)E44
y^(k+1)Δu(k)=l=1L(ClBl+Dl)μ¯yl(k+1)E45
y^(k+2)Δu(k)=l=1L(ClAlBl+ClBl+Dl)μ¯yl(k+2)E46
y^(k+Hp1)Δu(k)=l=1L(Clj=0Hp2AljBl+Dl)μ¯yl(k+Hp1)E47

The second group partial derivatives in (33) has the following matrix form:

ΔU(k)ΔU(k)=[Δu(k)Δu(k)....Δu(k)Δu(k+Hu1)Δu(k+Hu1)Δu(k)....Δu(k+Hu1)Δu(k+Hu1)]E48

SinceΔu(k)=u(k)u(k1), the matrix (39) has the following form:

ΔU(k)ΔU(k)=E49

Following this procedure it is possible to calculate the rest column elements of the matrix (34) which belongs to the next gradient vector elements (32). Finally, each element of the gradient-vector (32) could be obtained by the following system of equations:

J(k)Δu(k+1)=2e^(k+1)Q˜(1)y^(k+1)Δu(k)+...+2e^(k+Hp)Q˜(Hp)y^(k+Hp)Δu(k)++2R˜(1)Δu(k)2R˜(2)Δu(k+1)+...2R˜(Hu)Δu(k+Hu1)=0E50
J(k)Δu(k+2)=2e^(k+1)Q˜(1)y^(k+1)Δu(k+1)+...+2e^(k+Hp)Q˜(Hp)y^(k+Hp)Δu(k+1)++2R˜(2)Δu(k+1)2R˜(3)Δu(k+2)+...2R˜(Hu)Δu(k+Hu1)=0E51
J(k)Δu(k+Hu2)=2e^(k+1)Q˜(1)y^(k+1)Δu(k+Hu2)+...+2e^(k+Hp)Q˜(Hp)y^(k+Hp)Δu(k+Hu2)++2R˜(Hu1)Δu(k+Hu2)2R˜(Hu)Δu(k+Hu1)=0E52
J(k)Δu(k+Hu1)=2e^(k+1)Q˜(1)y^(k+1)Δu(k+Hu1)+...+2e^(k+Hp)Q˜(Hp)y^(k+Hp)Δu(k+Hu1)++2R˜(Hu)Δu(k+Hu1)=0E53

where e^(k+j)=r(k+j)y^(k+j),j=1,2,,Hpis the predicted system error.

The obtained system of equations (41)-(44) can be solved very easily, starting from the last equation (44) and calculating the last control action u(k+Hu-1)first. Then, the procedure can continue with finding the previous control action u(k+Hu-2)from (43). The calculations continue until the whole number of the control actions over the horizon Huis obtained. The calculation order of the control actions is very important, since the calculations should contain only known quantities. After that, only the first control action u(k)(30) will be used at the moment kas an input to the controlled process. The software implementation of the proposed algorithm is realized easily according to the following equations:

Δu(k+Hu1)=R˜(Hu)1[e^(k+1)Q˜(1)y^(k+1)Δu(k+Hu1)++e^(k+Hp)Q˜(Hp)y^(k+Hp)Δu(k+Hu1)]E54
Δu(k+Hu2)=Δu(k+Hu1)++R˜(Hu1)1[e^(k+1)Q˜(1)y^(k+1)Δu(k+Hu2)++e^(k+Hp)Q˜(Hp)y^(k+Hp)Δu(k+Hu2)]E55
Δu(k+1)=Δu(k+2)+R˜(2)1[e^(k+1)Q˜(1)y^(k+1)Δu(k+1)++e^(k+Hp)Q˜(Hp)y^(k+Hp)Δu(k+1)]E56
Δu(k)=Δu(k+1)+R˜(1)1[e^(k+1)Q˜(1)y^(k+1)Δu(k)++e^(k+Hp)Q˜(Hp)y^(k+Hp)Δu(k)]E57

The proposed unconstrained predictive control algorithm could be summarized in the following steps (Table 2).

 Step 1. Initial identification of the Takagi-Sugeno fuzzy-neural predictive model; Step 2. Start the algorithm at the sample k with the initial parameters; Step 3. Calculate the predicted model output ŷ(k+j) using the tuned fuzzy-neural model (2); Step 4. Calculate the derivatives for the matrix (34) according to the equations (35)-(38); Step 5. Calculate predicted control actions according to (45)-(48) and update the sequence; Step 6. Apply the first optimal control action u(k); Step 7. Modify the model parameters into the rule (3) and update them for the next step 3

### Table 2.

Basic fuzzy-neural model unconstrained predictive control algorithm

### 3.2. Constrained model predictive control

The constrained nonlinear predictive control problem can be described as a problem of finding the “optimal” input sequence to move a dynamic system to a desired state, taking into account the constraints on the inputs and the outputs of the control systems. This section reveals the formulation of the constrained control problem for MPC uses. Essentially, the problem becomes a quadratic programming problem with linear inequality constraints (LICQP). It follows by the nature of the operational constraints, which are usually described by linear inequalities of the control and plant variables.

The problem of nonlinear constrained predictive control is formulated as a nonlinear quadratic optimization problem. By means of local linearization (20) the problem can be solved using QP. That way the solution to the linear constrained predictive control problem is obtained. At each sampling time the LICQP is solved with new parameters, which are obtained by the Takagi-Sugeno fuzzy-neural model. An active set method is used for solving the constructed quadratic programming problem.

#### 3.2.1. Constraint types in model predictive control

The operational constraints may be classified in three major types according to the type of the system variables, which they are imposed on. The first two types of constraints deal with the control variable incremental variation ∆u(k)and control variable u(k). The third type is concerned with output y(k)or state variable x(k)constraints.

Related to the origin model predictive control problem, the constraints are expressed in a set of linear equations. All types of constraints are taken into consideration for each moving horizon window.

Umin(k)U(k)Umax(k)   ΔUmin(k)ΔU(k)ΔUmax(k)    Ymin(k)Y(k)Ymax(k)   E58

Where

Umax(k)=[umax(k)umax(k+1)umax(k+Nu1)]Umin(k)=[umin(k)umin(k+1)umin(k+Nu1)]E59
ΔUmax(k)=[Δumax(k)Δumax(k+1)Δumax(k+Nu1)]ΔUmin(k)=[Δumin(k)Δumin(k+1)Δumin(k+Nu1)]E60
Ymax(k)=[ymax(k)ymax(k+1)ymax(k+Np1)]Ymin(k)=[ymin(k)ymin(k+1)ymin(k+Np1)]E61

Therefore, for multi-input case the number of the constraints for the change of the control variable ∆u(k)is m×Nu.Similarly, the number of the constraints for the control variable amplitude is also m×Nuand for the output constraints it is q×Np.

#### 3.2.2. Quadratic programming in use of constrained MPC

Since the cost function J(k)(19) is quadratic and the constraints are linear inequalities, the problem of finding an optimal predictive control becomes one of finding an optimal solution to a standard quadratic programming problem with linear inequality constraints

minJ(x)=12xTHx+fTxsubjecttoAxbE62

where Hand fare the Hessian and the gradient of the Lagrange function, xis the decision variable. Constraints on the QP problem (50) are specified by Ax ≤ baccording to (49).

The Lagrange function is defined as follows

L(x,λ)=J(x)+i=1Nλiai,i=1,2,NE63
,

where λiare the Lagrange multipliers, aiare the constraints on the decision variable x, Nis the number of the constraints considered in the optimization problem.

Several algorithms for constrained optimization are described in (Fletcher, 2000). In this chapter a primal active set method is used. The idea of active set method is to define a set Sof constraints at each step of algorithm. The constraints in this active set are regarded as equalities whilst the rest are temporarily disregarded and the method adjusts the set in order to identify the correct active constraints on the solution to (52)

minJ(x)=12xTHx+fTxsubjecttoaix=biaixbiE64

At iteration ka feasible point x(k) is known which satisfies the active constraints as equalities. Each iteration attempts to locate a solution to an equality problem (EP) in which only the active constraints occur. This is most conveniently performed by shifting the origin to x(k) and looking for a correction δ(k) which solves

minδ{12δTHx+fTδ}subjecttoaiδ=0aiSE65

where f(k)is defined by f(k)=f + Hx(k) and is J(x(k))for the function defined by (52). If δ(k) is feasible with regard to the constraints not included in S, then the feasible point in next iteration is taken as x(k+ 1) = x(k) + δ(k). If not, a line search is made in the direction of δ(k) to find the best feasible point. A constraint is active if the Lagrange multipliers λi≥ 0,i.e. it is at the boundary of the feasible region defined by the constraints. On the other hand, if there exist λi< 0,the constraint is not active. In this case the constraint is relaxed from the active constraints set Sand the algorithm continues as before by solving the resulting equality constraint problem (53). If there is more than one constraint with correspondingλi< 0,then the miniSλi(k)is selected (Fletcher, 2000).

The QP, described in that way, is used to provide numerical solutions in constrained MPC problem.

#### 3.2.3. Design the constrained model predictive problem

The fuzzy-neural identification procedure from the Section 2 provides the state-space matrices, which are needed to construct the constrained model predictive control optimization problem.

Similarly to the unconstrained model predictive control approach, the cost function (18) can be specified by the prediction expressions (22) and (23).

J(k)=[Ψx(k)+Γu(k-1)+ΘΔU(k)-T(k)]TQ[Ψx(k)+Γu(k-1)+ΘΔU(k)-T(k)]+ΔUT(k)RΔU(k)==[ΘΔU(k)-E(k)]TQ[ΘΔU(k)-E(k)]+ΔUT(k)RΔU(k)==ΔUT(k)[ΘTQΘ+R]ΔU(k)+ET(k)QE(k)-2ΔUT(k)ΘTQE(k)E66

Assuming that

H=ΘTQΘ+RandΦ=2ΘTQE(k),E67

the cost function for the model predictive optimization problem can be specified as follow

J(k)=ΔUT(k)ΗΔU(k)-ΔUT(k)Φ+ET(k)QE(k)E68

The problem of minimizing the cost function (55) is a quadratic programming problem. If the Hessian matrix His positive definite, the problem is convex (Fletcher, 2000). Then the solution is given by the closed form

ΔU=12H1ΦE69

The constraints (49) on the cost function may be rewritten in terms of ∆U(k).

Umin(k)Iuu(k-1)+IΔuΔU(k)Umax(k)   ΔUmin(k)ΔU(k)ΔUmax(k)    Ymin(k)Ψx(k)+Γu(k-1)+ΘΔU(k)Ymax(k)   E70

where Imm×mis an identity matrix

Iu=[ImImIm]mNu×m,IΔu=[Im00ImIm0ImImIm]mNu×mNu.E71

All types of constraints are combined in one expression as follows

[IΔuIΔuIIΘΘ]ΔU[Umin+Iuu(k1)UmaxIuu(k1)ΔUminΔUmaxYmin+(Ψx(k)+Γu(k1))Ymax(Ψx(k)+Γu(k1))]E72

where ImNu×mNuis an identity matrix.

Finally, following the definition of the LIQP (50), the model predictive control in presence of constraints is proposed as finding the parameter vector ∆Uthat minimizes (55) subject to the inequality constraints (58).

minJ(k)=ΔUTHΔU-ΔUTΦ+ETQEsubjecttoΩΔUωE73

In (59) the constraints expression (58) has been denoted by Ω∆U ≤ω, where Ω is a matrix with number of rows equal to the dimension of ωand number of columns equal to the dimension of ∆U.In case that the constraints are fully imposed, the dimension of ωis equal to 4×m×Nu+ 2×q×Np,where mis the number of system inputs and qis the number of outputs. In general, the total number of constraints is greater than the dimension of the ∆U.The dimension of ωrepresents the number of constraints.

The proposed model predictive control algorithm can be summarized in the following steps (Table 3).

 Step 1. Read the current states, inputs and outputs of the system; Step 2. Start identification of the fuzzy-neural predictive model following Algorithm 1; Step 3. With A(k), B(k), C(k), D(k) from Step 2 calculate the predicted output Y(k) according to (17); Step 4. Obtain the prediction error E(k) according to (23); Step 5. Construct the cost function (55) and the constraints (58) of the QP problem; Step 6. Solve the QP problem according to (59); Step 7. Apply only the first control action u(k).

### Table 3.

State-space implementation of fuzzy-neural model predictive control strategy

At each sampling time, LIQP (59) is solved with new parameters. The Hessian and the Lagrangian are constructed by the state-space matrices A(k),B(k),C(k)and D(k)(4) obtained during the identification procedure (Table 1). The problem of nonlinear constrained predictive control is formulated as a nonlinear quadratic optimization problem. By means of local linearization a relaxation can be obtained and the problem can be solved using quadratic programming. This is the solution of the linear constrained predictive control problem (Espinosa et al., 2005).

## 4. Fuzzy-neural model predictive control of a multi tank system. Case study

The case study is implemented in MATLAB/Simulink® environment with Inteco® Multi tank system. The Inteco® Multi tank System (Fig. 4) comprises from three separate tanks fitted with drain valves (Inteco, 2009). The additional tank mounted in the base of the set-up acts as a water reservoir for the system. The top (first) tank has a constant cross section, while others are conical or spherical, so they are with variable cross sections. This causes the main nonlinearities in the system. A variable speed pump is used to fill the upper tank. The liquid outflows the tanks by the gravity. The tank valves act as flow resistors C1, C2, C3. The area ratio of the valves is controlled and can be used to vary the outflow characteristic. Each tank is equipped with a level sensor PS1, PS2, PS3based on hydraulic pressure measurement.

The linearized dynamical model of the triple tank system could be described by the linear state-space equations (2) where the matrices A, B, Cand Dare as follow (Petrov et al., 2008).

A=[α1awH11w100α1w(c+bH1/H1max)H11w1α2w(c+bH2/H2max)H21w200α2wR2(H3maxH3)2H21w2α2wR2(H3maxH3)2H31w3]E74
B=[1aw1awH11α100001w(c=bH2/H2max)H21α200001wR2(H3maxH3)2H31α3]E75
C=D=E76

The parameters a1, a2and a3 are flow coefficients for each tank of the model. The described linearized state-space model is used as an initial model for the training process of the fuzzy-neural model during the experiments.

### 4.1. Description of the multi tank system as a multivariable controlled process

Liquid levels Н1, Н2, Н3in the tanks are the state variables of the system (Fig. 4). The Inteco Multi Tank system has four controlled inputs: liquid inflow qand valves settings C1, C2, C3.Therefore, several models of the tanks system can be analyzed (Fig. 5), classified as pump-controlled system, valve-controlled system and pump/valve controlled system (Inteco, 2009). Figure 5.Model of the Multi Tank system as a pump and valve-controlled system

In this case study a multi-input multi-output (MIMO) configuration of the Inteco Multi Tank system is used (Fig. 5). This corresponds to the linearized state-space model (60). Several issues have been recognized as causes of additional nonlinearities in plant dynamics:

• nonlinearities (smooth and nonsmooth) caused by shapes of tanks;

• saturation-type nonlinearities, introduced by maximum or minimum level allowed in tanks;

• nonlinearities introduced by valve geometry and flow dynamics;

• nonlinearities introduced by pump and valves input/output characteristic curve.

The simulation results have been obtained with random generated set points and following initial conditions (Table 4):

 Model predictive controller parameters Prediction horizon Hp=10First included sample of the prediction horizon Hw=1Control horizon Hu=3 Inteco Multi tank system parameters Flow coefficients for each tanka1=0.29; a2=0.2256; a3=0.2487 Operational constraints on the system Constraints on valve cross section ratio 0 ≤ Ci ≤ 2e-04, i=1,2,3Constraint on liquid inflow 0 ≤ q ≤ 1e-04 m3/sConstraints on liquid level in each tank 0 ≤ Hi ≤ 0.35 m, i=1,2,3 Simulation parameters Time of simulation 600 sSample time Ts=1 s

### Table 4.

Simulation parameters for unconstrained and constrained fuzzy-neural MPC

Figures below show typical results for level control problem. The reference value for each tank is changed consequently in different time. The proposed fuzzy-neural identification procedure ensures the matrices for the optimization problem of model predictive control at each sampling time Ts. The plant modelling process during the unconstrained and constrained MPC experiments are shown in Fig. 6 and Fig. 9, respectively.

### 4.2. Experimental results with unconstrained model predictive control

The proposed unconstrained model predictive control algorithm (Table 2) with the Takagi-Sugeno fuzzy-neural model as a predictor has been applied to the level control problem. The experiments have been implemented with the parameters in Table 4. The weighting matrices are specified as follow: Q˜=0.01*diag(1,1,1)andR˜=10e4*diag(1,1,1,1). Note that the weighting matrix R˜is constant over all prediction horizon, which allows to avoid matrix inversion at each sampling time with one calculation of R˜1at time k=0. Figure 6.Fuzzy-neural model identification procedure of the multi tank system – unconstrained NMPC

The next two figures - Fig. 7 and Fig. 8, show typical results regarding level control, where the references for H1, H2and H3are changed consequently in different time. The change of every level reference behaves as a system disturbance for the other system outputs (levels). It is evident that the applied model predictive controller is capable to compensate these disturbances. Figure 7.Transient responses of multi tank system outputs – unconstrained NMPC Figure 8.Transient responses of multi tank system inputs – unconstrained NMPC

### 4.3. Experimental results with fuzzy-neural constrained predictive control

The experiments with the proposed constrained model predictive control algorithm (Table 3) have been made with level references close to the system outputs constraints. The weighting matrices in GPC cost function (19) are specified as Q˜=diag(1,1,1)andR˜=15e4*diag(1,1,1,1).System identification during the experiment is shown on Fig. 9. The proposed identification procedure uses the linearized model (60) of the Multi tank system as an initial condition. Figure 9.Fuzzy-neural model identification procedure of the multi tank system – constrained NMPC

The proposed constrained fuzzy-neural model predictive control algorithm provides an adequate system response as it can be seen on Fig. 10 and Fig. 11. The references are achieved Figure 10.Transient responses of the multi tank system outputs – constrained NMPC Figure 11.Transient responses of the multi tank system inputs – constrained NMPC

without violating the operational constraints specified in Table 4. Similarly to the unconstrained case, the Takagi-Sugeno type fuzzy-neural model provides the state-space matrices A, Band C(the system is strictly proper, i.e. D=0) for the optimization procedure of the model predictive control approach. Therefore, the LIQP problem is constructed with “fresh” parameters at each sampling time and improves the adaptive features of the applied model predictive controller. It can be seen on the next figures that the disturbances, which are consequences of a sudden change of the level references, are compensated in short time without violating the proper system work.

## 5. Conclusions

This chapter has presented an effective approach to fuzzy model-based control. The effective modelling and identification techniques, based on fuzzy structures, combined with model predictive control strategy result in effective control for nonlinear MIMO plants. The goal was to design a new control strategy - simple in realization for designer and simple in implementation for the end user of the control systems.

The idea of using fuzzy-neural models for nonlinear system identification is not new, although more applications are necessary to demonstrate its capabilities in nonlinear identification and prediction. By implementing this idea to state-space representation of control systems, it is possible to achieve a powerful model of nonlinear plants or processes. Such models can be embedded into a predictive control scheme. State-space model of the system allows constructing the optimization problem, as a quadratic programming problem. It is important to note that the model predictive control approach has one major advantage – the ability to solve the control problem taking into consideration the operational constraints on the system.

This chapter includes two simple control algorithms with their respective derivations. They represent control strategies, based on the estimated fuzzy-neural predictive model. The two- stage learning gradient procedure is the main advantage of the proposed identification procedure. It is capable to model nonlinearities in real-time and provides an accurate model for MPC optimization procedure at each sampling time.

The proposed consequent solution for unconstrained MPC problem is the main contribution for the predictive optimization task. On the other hand, extraction of a “local” linear model, obtained from the inference process of a Takagi–Sugeno fuzzy model allows treating the nonlinear optimization problem in presence of constraints as an LIQP.

The model predictive control scheme is employed to reduce structural response of the laboratory system - multi tank system. The inherent instability of the system makes it difficult for modelling and control. Model predictive control is successfully applied to the studied multi tank system, which represents a multivariable controlled process. Adaptation of the applied fuzzy-neural internal model is the most common way of dealing with plant’s nonlinearities. The results show that the controlled levels have a good performance, following closely the references and compensating the disturbances.

The contribution of the proposed approach using Takagi–Sugeno fuzzy model is the capacity to exploit the information given directly by the Takagi–Sugeno fuzzy model. This approach is very attractive for systems from high order, as no simulation is needed to obtain the parameters for solving the optimization task. The model’s state-space matrices can be generated directly from the inference of the fuzzy system. The use of this approach is very attractive to the industry for practical reasons related with the capacity of this model structure to combine local models identified in experiments around the different operating points.