# Model Predictive Control Based on Kalman Filter for Constrained Hammerstein-Wiener Systems

Authors: Man Hong; Shao Cheng
Journal: Mathematical Problems in Engineering
 (2013)

---

## Abstract

To precisely track the reactor temperature in the entire working condition,  the constrained Hammerstein-Wiener model describing nonlinear chemical processes such as in the continuous stirred tank reactor (CSTR) is proposed.  A predictive control algorithm based on the Kalman filter for constrained Hammerstein-Wiener systems is designed. An output feedback control  law regarding the linear subsystem is derived by state observation. The size of reaction heat produced and its influence on the output are evaluated by the Kalman filter. The observation and evaluation results are calculated by the multistep predictive approach. Actual control variables are computed while considering the constraints of the optimal control problem in a finite horizon through the receding horizon. The simulation example of the CSTR tester shows the effectiveness and feasibility of the proposed algorithm.


---

## Body

## 1. Introduction

Temperature control of the batch reactor during chemical processes [1] is the key to producing qualified chemical products. The control requirement is to change the temperature in the reactor with the given temperature process curve set up in advance for the entire batch reaction process. The accuracy of following the temperature change is the scale used to judge the performance of the designed controller. Two factors affect temperature change. First is the change of the interlayer cooling capacity, where the size reflects the input change of the control system. Second is the heat production issue of the chemical reaction in the reactor with a variation in size that cannot be measured generally because of uncertain, strong, and nonlinear characteristics, raising higher requirements for the design of the controller.The Hammerstein-Wiener model can describe the processes in the continuous stirred tank reactor (CSTR), pH neutralization reaction process, and other nonlinear chemical processes. Many research results were achieved in the nonlinear model predictive control (NMPC) [2–4]. Considering that the state of the actual system is usually unpredictable, the literature [5] presents a new two-step NMPC algorithm that constrains the output feedback of the Hammerstein system. The literature [6] designs the predictive controller of the batch reaction based on the Kalman filter; the DMC is used as the framework and is combined with the size of the interference evaluated by the Kalman filter. Its impact on the system output for the problems regarding uncertainty and strong nonlinear characteristics is reflected by the system during the heat production of the batch reaction itself. However, no reasonable method was found to resolve the complex relationship among the states of the practical nonlinear system, the constraints of the intermediate mass, and the impact of the interference signal on the output of the system.Considering the Hammerstein-Wiener system of state, output, and intermediate mass constraints, this paper proposes the predictive control algorithm of the batch reaction based on the Kalman filter for the constrained Hammerstein-Wiener system by combining the literature [5, 6]. The paper designs the output feedback optimal control law of the linear subsystem with the state observer, and it evaluates the size of interference and its impact on the system output with the linear Kalman filter. By predicting the solving error of the nonlinear algebraic equation and by roll optimizing the error performance indicator of the finite horizon, the relationship between the outputs calculates the actual control variables. The simulation study of the CSTR proves that the control method in this paper overcomes the impact of the heat production in the reaction in terms of the accuracy of temperature tracking and improves the control effect to a certain extent.


## 2. Description of Batch Reaction Process and the Hammerstein-Wiener System

### 2.1. Description of Batch Reaction Process [1]

The batch reactor process is generally operated under a closed reactor condition so that the reactants are added to the reactor once to react in the given condition. The intensifying reaction will cause the heat to accumulate and the temperature to rise and further intensify. Obviously, this process, as well as the interlayer temperature, is unsteady. As the control process is complex, the batch reactor control system is generally required to design the cascade structure, in which the main circuit of the outer ring takes the interlayer temperature as the control input and the reaction temperature as the control output, usingTj and Tr, respectively. The secondary circuit takes the cooling water flow rate as the control input and the interlayer temperature as the control output, using Qc and Tj, respectively. The value of the secondary circuit is set up by the main circuit and the secondary circuit, is generally controlled with a simple PID. The external structure diagram of the reactor is shown in Figure 1.Figure 1
CSTR.


### 2.2. Description of the Hammerstein-Wiener System

The Hammerstein-Wiener model consists of two static, nonlinear modules and one dynamic linear module. Its model structure is shown in Figure2.Figure 2
Hammerstein-Wiener model structure.Consider the following discrete-time multivariable Hammerstein-Wiener absolute system:(1)x(t+1)=Ax(t)+Bv(t),v(t)=g(u(t),t),w(t)=Cx(t),y(t)=h(w(t),t),
where x∈Rn,v∈Rr, and u∈Rm are the state variable, intermediate variable, and the input variable of the system, respectively; w is the output variable in the linear part; (A,B,C,D) is the constant coefficient matrix with the corresponding dimension; g,h are the continuously differentiable nonlinear function vectors, indicating the static nonlinear link between the system input and the intermediate variable and between the linear part and the intermediate variable and meeting g(0,t)=0,h(0,t)=0.It is assumed that the system controllability, state observability, and the origin as the equilibrium point of the system, the parameters of the system states, intermediate variables, and the input constraints are defined as follows:(2)x(t)∈Mx⊆Rn,v(t)∈Mv⊆Rr,u(t)∈Mu⊆Rm,w(t)∈Mw⊆Rr,
where Mx,Mv,Mu, and Mw represent the compact convex set, including the origin as the interior point.


## 2.1. Description of Batch Reaction Process [1]

The batch reactor process is generally operated under a closed reactor condition so that the reactants are added to the reactor once to react in the given condition. The intensifying reaction will cause the heat to accumulate and the temperature to rise and further intensify. Obviously, this process, as well as the interlayer temperature, is unsteady. As the control process is complex, the batch reactor control system is generally required to design the cascade structure, in which the main circuit of the outer ring takes the interlayer temperature as the control input and the reaction temperature as the control output, usingTj and Tr, respectively. The secondary circuit takes the cooling water flow rate as the control input and the interlayer temperature as the control output, using Qc and Tj, respectively. The value of the secondary circuit is set up by the main circuit and the secondary circuit, is generally controlled with a simple PID. The external structure diagram of the reactor is shown in Figure 1.Figure 1
CSTR.


## 2.2. Description of the Hammerstein-Wiener System

The Hammerstein-Wiener model consists of two static, nonlinear modules and one dynamic linear module. Its model structure is shown in Figure2.Figure 2
Hammerstein-Wiener model structure.Consider the following discrete-time multivariable Hammerstein-Wiener absolute system:(1)x(t+1)=Ax(t)+Bv(t),v(t)=g(u(t),t),w(t)=Cx(t),y(t)=h(w(t),t),
where x∈Rn,v∈Rr, and u∈Rm are the state variable, intermediate variable, and the input variable of the system, respectively; w is the output variable in the linear part; (A,B,C,D) is the constant coefficient matrix with the corresponding dimension; g,h are the continuously differentiable nonlinear function vectors, indicating the static nonlinear link between the system input and the intermediate variable and between the linear part and the intermediate variable and meeting g(0,t)=0,h(0,t)=0.It is assumed that the system controllability, state observability, and the origin as the equilibrium point of the system, the parameters of the system states, intermediate variables, and the input constraints are defined as follows:(2)x(t)∈Mx⊆Rn,v(t)∈Mv⊆Rr,u(t)∈Mu⊆Rm,w(t)∈Mw⊆Rr,
where Mx,Mv,Mu, and Mw represent the compact convex set, including the origin as the interior point.


## 3. The Constrained Predictive Controller Design Based on the Kalman Filter

### 3.1. The Constrained NMPC Algorithm [5]

The Hammerstein-Wiener model can be divided into the Hammerstein model [7] and the Wiener model [8]. Different approaches are adopted for these two models to design the controller of the system and to deal with the impact of system output and interference signal on the system. This section mainly describes the designs of the output feedback optimal control law of the Hammerstein linear subsystem based on the state observer.The considered output signaly(t) is assumed to have no direct control relationship with the input signal u(t) in formula (1). Therefore, in its corresponding linear subsystem (A,B,C), x(t+1)=Ax(t)+Bv(t), w(t)=Cx(t); if the observer status is x^(t), the optimal control law of the linear subsystem can be described as v(t)≈Kx^(t), where K=-(R+BTPB)-1BTPA, R=rI is the positive definite weighting matrix of the intermediate mass, P is a symmetric positive definite solution to the Riccati equation ATPA-P+ATPB(R+BTPB)-1BTPA+Q=0, and Q=qI is the semipositive definite state weighting matrix.The actual control variableu(t) can be calculated by solving the nonlinear equations g(u,t)-v(t)=0, inevitably introducing the calculation error, and thereby affecting the closed-loop stability and optimality. Next, the actual control is determined through the prediction and the scroll optimization of the calculation error.Considering the system in formula (1), the calculation error objective function is defined as follows:
(3)J1(t)=∑i=tt+Tp-1{[g(i)-v(i)]T[g(i)-v(i)]+u(i)TWu(i)},
where the integer Tp>0 predicts the step length and the positive definite matrix W is the input weighting matrix.Thus, the NMPC issue of the constrained Hammerstein system is defined as follows:(4)u(t,x^(t),Tp)*=minuJ1,s.t.x(i+1)=Ax(i)+Bg(u(i),i),v(i)=Kx^(t),x(t)∈Mx,v(t)∈Mv,u(t)∈Mu,w(t)∈Mw,x(t)=x^(t),(i=t,t+1,…,t+Tp-1).The sequential quadratic programming is applied to calculate online the optimization problem in formula (4) and obtain the optimal control solution:
(5)u(t,x^(t),Tp)*={u(t,x^(t))*,u(t+1,x^(t))*,…,u(t+Tp-1,x^(t))*}.According to the principle of rolling optimization, the NMPC law based on the observation state is defined as follows:(6)u(t)NMPC=u(t,x^(t))*,(t=0,1,…).In formula (6), the observation state x^(t) can be obtained through the following all-dimensional observer:
(7)x^(t+1)=(A-LC)x^(t)+Bg(u(t)NMPC,t)+Lw(t),(t=0,1,…),
where L is the gain matrix of the observer.Control law (6) observes the system status with the state observer (7) at each sampling time. It then updates the initial conditions of the optimal control problem (4) on this basis and scrolls the optimization to calculate the predictive control variables of the current.The observation error of the defined state ise=x-x^; thus, the NMPC closed-loop Hammerstein system is
(8)x(t+1)=Ax(t)+Bg(u(t)NMPC,t)x(t+1)=(A+BK)x(t)-BKe(t)+Bve(t),e(t+1)=(A-LC)e(t),(t=0,1,2,…),
where ve(t)=g(uNMPC,t)-v(t) is the calculation error of the control variables; the controller gain K is given by v(t)≈Kx^(t).The NMPC closed-loop structure of the whole system is shown in Figure3.Figure 3
NMPC closed-loop system structure of the constrained hammerstein system.


### 3.2. The Predictive Controller Design Based on the Kalman Filter

This section mainly introduces the interference model [6] based on the Kalman filter and evaluates the size of the interference signal as well as its impact on the system output based on the Wiener model. Meanwhile, with the above constrained NMPC output, the relationship between their outputs solves the error through the predictive nonlinear algebraic equation, and roll optimizes the performance indicators in the finite horizon to calculate the actual control variables.The input interference model is the key to achieving unbiased control and to improving the tracking performance. To improve the tracking performance and to realize the unbiased control, the system is characterized by introducing the interference model. The model is built in the input interference along with the influence of heat produced by the materials in the chemical reactor on the output. Davison and Smith [8] regarded this method as the standard in eliminating the steady-state error in the design of the quadratic regulator:
(9)Tr(k+1)=αTr(k)+β[Tj(k)+d(k)].After the establishment of the input interference model, the input interferenced is evaluated. If the interference is the state vector, the following augmented-state space model is obtained:
(10)[Tr(k+1)d(k+1)]=[αβ01][Tr(k)d(k)]+[η(k)ξ(k)][β0]TTj(k),y(k)=[10][Tr(k)d(k)]T,
where η and ξ are the process noises when the mean value is zero.According to the measured value of the current output, the system stateT and the interference d of the current can be determined using the following formula:
(11)[Tk^∣kdk^∣k]=[Tk^∣k-1dk^∣k-1]+[L1L2](yk-Tk^∣k-1),
where L1 and L2 represent the corresponding Kalman gains of the state variable T and the interference d, respectively.Then, the predictive value of the augmented state at the next moment is(12)[Tk^+1∣kdk^+1∣k]=[αβ01][Tk^∣kdk^∣k]+[β0]Tj(k).For the sake of discussion, the following are set:A1=[αβ01], B2=[β0], H=[10], and L=[L1L2]. The Kalman gains of each cycle can then be calculated by Lk=PkHT(HPkHT+R)-1, where Pk=(I-Kk-1H)APk-1AT+Q is the covariance of the estimation error of the current; Q,R are the process incentive noise covariance matrix and the observation noise covariance matrix, respectively, and their values directly affect the performance of the Kalman filter. Greg Welch and Gary Bishop (2006) discussed it in detail in their writings. Interference d has no direct relationship with input Tj, but the observability of interference d ensures passing the Kalman filter and then evaluates d. When designing the controller, the evaluation value of d is used to remove its influence on the system output. After evaluating the size of d using the Kalman filter in each sampling period, feedforward compensation is made on this part, and its influence on the system output is included in the prediction model. The real-time control variable is
(13)ΔTj(k)=DT[Tr_sp(k)-bΔd(k)-T~op(k)]ΔTj(k)=DT[Trsp′(k)-T~op(k)],
where Δd(k)=d(k)-d(k-1); D=(ATQA+R)-1ATQ; and Trsp′(k) deducts the influence of interference d on the output from the expected value to form the new one.The predictive control algorithm flowchart of the constrained Hammerstein-Wiener model with the feedforward compensation is shown in Figure4.Figure 4
The predictive control algorithm of the constrained Hammerstein-Wiener model with the feedforward compensation.


## 3.1. The Constrained NMPC Algorithm [5]

The Hammerstein-Wiener model can be divided into the Hammerstein model [7] and the Wiener model [8]. Different approaches are adopted for these two models to design the controller of the system and to deal with the impact of system output and interference signal on the system. This section mainly describes the designs of the output feedback optimal control law of the Hammerstein linear subsystem based on the state observer.The considered output signaly(t) is assumed to have no direct control relationship with the input signal u(t) in formula (1). Therefore, in its corresponding linear subsystem (A,B,C), x(t+1)=Ax(t)+Bv(t), w(t)=Cx(t); if the observer status is x^(t), the optimal control law of the linear subsystem can be described as v(t)≈Kx^(t), where K=-(R+BTPB)-1BTPA, R=rI is the positive definite weighting matrix of the intermediate mass, P is a symmetric positive definite solution to the Riccati equation ATPA-P+ATPB(R+BTPB)-1BTPA+Q=0, and Q=qI is the semipositive definite state weighting matrix.The actual control variableu(t) can be calculated by solving the nonlinear equations g(u,t)-v(t)=0, inevitably introducing the calculation error, and thereby affecting the closed-loop stability and optimality. Next, the actual control is determined through the prediction and the scroll optimization of the calculation error.Considering the system in formula (1), the calculation error objective function is defined as follows:
(3)J1(t)=∑i=tt+Tp-1{[g(i)-v(i)]T[g(i)-v(i)]+u(i)TWu(i)},
where the integer Tp>0 predicts the step length and the positive definite matrix W is the input weighting matrix.Thus, the NMPC issue of the constrained Hammerstein system is defined as follows:(4)u(t,x^(t),Tp)*=minuJ1,s.t.x(i+1)=Ax(i)+Bg(u(i),i),v(i)=Kx^(t),x(t)∈Mx,v(t)∈Mv,u(t)∈Mu,w(t)∈Mw,x(t)=x^(t),(i=t,t+1,…,t+Tp-1).The sequential quadratic programming is applied to calculate online the optimization problem in formula (4) and obtain the optimal control solution:
(5)u(t,x^(t),Tp)*={u(t,x^(t))*,u(t+1,x^(t))*,…,u(t+Tp-1,x^(t))*}.According to the principle of rolling optimization, the NMPC law based on the observation state is defined as follows:(6)u(t)NMPC=u(t,x^(t))*,(t=0,1,…).In formula (6), the observation state x^(t) can be obtained through the following all-dimensional observer:
(7)x^(t+1)=(A-LC)x^(t)+Bg(u(t)NMPC,t)+Lw(t),(t=0,1,…),
where L is the gain matrix of the observer.Control law (6) observes the system status with the state observer (7) at each sampling time. It then updates the initial conditions of the optimal control problem (4) on this basis and scrolls the optimization to calculate the predictive control variables of the current.The observation error of the defined state ise=x-x^; thus, the NMPC closed-loop Hammerstein system is
(8)x(t+1)=Ax(t)+Bg(u(t)NMPC,t)x(t+1)=(A+BK)x(t)-BKe(t)+Bve(t),e(t+1)=(A-LC)e(t),(t=0,1,2,…),
where ve(t)=g(uNMPC,t)-v(t) is the calculation error of the control variables; the controller gain K is given by v(t)≈Kx^(t).The NMPC closed-loop structure of the whole system is shown in Figure3.Figure 3
NMPC closed-loop system structure of the constrained hammerstein system.


## 3.2. The Predictive Controller Design Based on the Kalman Filter

This section mainly introduces the interference model [6] based on the Kalman filter and evaluates the size of the interference signal as well as its impact on the system output based on the Wiener model. Meanwhile, with the above constrained NMPC output, the relationship between their outputs solves the error through the predictive nonlinear algebraic equation, and roll optimizes the performance indicators in the finite horizon to calculate the actual control variables.The input interference model is the key to achieving unbiased control and to improving the tracking performance. To improve the tracking performance and to realize the unbiased control, the system is characterized by introducing the interference model. The model is built in the input interference along with the influence of heat produced by the materials in the chemical reactor on the output. Davison and Smith [8] regarded this method as the standard in eliminating the steady-state error in the design of the quadratic regulator:
(9)Tr(k+1)=αTr(k)+β[Tj(k)+d(k)].After the establishment of the input interference model, the input interferenced is evaluated. If the interference is the state vector, the following augmented-state space model is obtained:
(10)[Tr(k+1)d(k+1)]=[αβ01][Tr(k)d(k)]+[η(k)ξ(k)][β0]TTj(k),y(k)=[10][Tr(k)d(k)]T,
where η and ξ are the process noises when the mean value is zero.According to the measured value of the current output, the system stateT and the interference d of the current can be determined using the following formula:
(11)[Tk^∣kdk^∣k]=[Tk^∣k-1dk^∣k-1]+[L1L2](yk-Tk^∣k-1),
where L1 and L2 represent the corresponding Kalman gains of the state variable T and the interference d, respectively.Then, the predictive value of the augmented state at the next moment is(12)[Tk^+1∣kdk^+1∣k]=[αβ01][Tk^∣kdk^∣k]+[β0]Tj(k).For the sake of discussion, the following are set:A1=[αβ01], B2=[β0], H=[10], and L=[L1L2]. The Kalman gains of each cycle can then be calculated by Lk=PkHT(HPkHT+R)-1, where Pk=(I-Kk-1H)APk-1AT+Q is the covariance of the estimation error of the current; Q,R are the process incentive noise covariance matrix and the observation noise covariance matrix, respectively, and their values directly affect the performance of the Kalman filter. Greg Welch and Gary Bishop (2006) discussed it in detail in their writings. Interference d has no direct relationship with input Tj, but the observability of interference d ensures passing the Kalman filter and then evaluates d. When designing the controller, the evaluation value of d is used to remove its influence on the system output. After evaluating the size of d using the Kalman filter in each sampling period, feedforward compensation is made on this part, and its influence on the system output is included in the prediction model. The real-time control variable is
(13)ΔTj(k)=DT[Tr_sp(k)-bΔd(k)-T~op(k)]ΔTj(k)=DT[Trsp′(k)-T~op(k)],
where Δd(k)=d(k)-d(k-1); D=(ATQA+R)-1ATQ; and Trsp′(k) deducts the influence of interference d on the output from the expected value to form the new one.The predictive control algorithm flowchart of the constrained Hammerstein-Wiener model with the feedforward compensation is shown in Figure4.Figure 4
The predictive control algorithm of the constrained Hammerstein-Wiener model with the feedforward compensation.


## 4. Simulation Analysis

The experimental platform is the CSTR tester that can simulate in real time the entire process of the batch reaction. The main parameters selected include the volume of the reactants: 200 L, the spin speed of the stirrer: 90r/min, the power of the stirrer: 1.5 Kw, the feeding of the reactant m1: 50 kg, and the specific heat: 1.2kJ/kgK. The valve operation of the heating steam comprises manual and automatic modes. The initial stage of the reaction mainly involves the evoked response within the reactor. In other words, the temperature of the reactor is improved to 60°C with steam. The steam valve is manually operated in this process. The temperature in the reactor continues to rise because of the heat production in the reaction. The steam valve is switched to automatic mode to implement the control strategies. The process curve of the reaction temperature is set as follows: rising at 0.2°C/s, remaining constant at 120°C, and the overshoot of the temperature being less than 1.5°C. The result is compared with that of the traditional PID control, as shown in Figure 5.The process temperature curve rising at 0.2°C/s and remaining constant at 120°C.
(a)
The predictive control algorithm control result diagram referred to in the paper
(b)
The traditional PID control algorithm control result diagramCompared with the PID control, the PID control method cannot make the reactor temperature follow the process curve well because of the inherent nonlinearity and time lag characteristics of the object. However, the method in this paper achieves a good control effect.In practical production, different products are produced through different processes, requiring corresponding control tasks to be changed, that is, selecting different process curves while keeping the algorithm unchanged. In the same simulation platform, the process curve is set to rise at 0.3°C/s and be constant at 120°C, with other parameters unchanged. The control result is shown in Figure6(a), whereas the traditional PID algorithm control result is shown in Figure 6(b).The process temperature curve rising at 0.3°C/s and remaining constant at 120°C.
(a)
The predictive control algorithm control result diagram referred to in the paper
(b)
The traditional PID control algorithm control result diagramFor batch production, the different natures of the reactants result in different exothermic characteristics in the entire reaction process, requiring the adopted control algorithm to have good robustness. In other words, the original control result remains even when the material property changes. Now, the reactantm2 is used with the specific heat of 1.5kJ/kgK. The reaction temperature process curve is set as follows: rising at 0.2°C/s and remaining constant at 120°C. The result is compared with that of the traditional PID control, as shown in Figure 7.Comparison of the result after changing the special heat of the material.
(a)
The predictive control algorithm control result diagram referred to in the paper after changing the special heat of the material
(b)
The traditional PID control algorithm control result diagram after changing the special heat of the materialThe traditional PID control algorithm is difficult to use in controlling the reaction temperature in a satisfactory range because of the specific heat change of the reactant. However, the method proposed in this paper shows good robustness, indicating practical significance to the complex and changeable strong, nonlinear, chemical process system.


## 5. Conclusions

With the temperature tracking issue in the batch reaction process as the study object, as well as with the condition that the value of the heat production cannot be measured and its influence law on the output is unknown in the whole process of the reactant reaction, this paper proposes the predictive control algorithm of the batch reaction based on the Kalman filter for the constrained Hammerstein-Wiener system. Combined with state observer, the paper designs the output feedback optimal control law of the linear subsystem and characterizes the unmodeled dynamic system with the output interference model, which compensates for the impact of the heat production from the reaction on the system output and improves the accuracy of temperature tracking. The simulation experiment is conducted to verify the control results of this algorithm on different process curves and determine its robustness. The proposed algorithm, when compared with the traditional PID algorithm, fully demonstrates effectiveness and feasibility.



---

*Source: 104702-2013-10-22.xml*