- Email: [email protected]
Multi-linear model set design based on the nonlinearity measure and H-gap metric
ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
Contents lists available at ScienceDirect
ISA Transactions journal homepage: www.elsevier.com/locate/isatrans
Research article
Multi-linear model set design based on the nonlinearity measure and H-gap metric Davood Shaghaghi, Alireza Fatehi n, Ali Khaki-Sedigh APAC Research Group, Industrial Control Center of Excellence, Faculty of Electrical Engineering, K. N. Toosi University of Technology, 16317-14191 Tehran, Iran
art ic l e i nf o
a b s t r a c t
Article history: Received 18 November 2015 Received in revised form 8 November 2016 Accepted 16 January 2017
This paper proposes a model bank selection method for a large class of nonlinear systems with wide operating ranges. In particular, nonlinearity measure and H-gap metric are used to provide an effective algorithm to design a model bank for the system. Then, the proposed model bank is accompanied with model predictive controllers to design a high performance advanced process controller. The advantage of this method is the reduction of excessive switch between models and also decrement of the computational complexity in the controller bank that can lead to performance improvement of the control system. The effectiveness of the method is verified by simulations as well as experimental studies on a pH neutralization laboratory apparatus which confirms the efficiency of the proposed algorithm. & 2017 ISA. Published by Elsevier Ltd. All rights reserved.
Keywords: Nonlinearity measure H-gap metric Multiple-model controller Model bank selection Generalized predictive control
1. Introduction Classical control design techniques have been matured for linear systems to incorporate robustness and performance requirements within a narrow operating range. Nevertheless, for nonlinear processes with severe nonlinearities and large set-point changes, a single linear controller may not provide satisfactory performance. This has resulted in the development of nonlinear model predictive control (NMPC) that uses a comprehensive and precise model for process prediction and optimization [1]. The fundamental step in the development of NMPC is to elect an appropriate structure of a nonlinear model for identification of system dynamics. Unfortunately, the theory of nonlinear system identification has been not developed very well [2]. Furthermore, the resulting dynamic model may be too complex to be useful for NMPC design. The main disadvantage of the nonlinear modeling approach is that using a nonlinear model led to change the control problem solution from a convex quadratic programming to a nonconvex nonlinear case. The solution of this problem is much more difficult and there is no guarantee to find the global optimum [3]. An approach that tries to retain the features of a linear design and is applied to nonlinear systems is the multi-model approach for modeling and control [2]. The multiple model approach with switching and supervisory control is best suited for processes that n
Corresponding author. E-mail addresses: [email protected] (D. Shaghaghi), [email protected] (A. Fatehi), [email protected] (A. Khaki-Sedigh).
have fast changes in operating range [4] and [5]. The key concept is to represent the nonlinear system as a combination of linear subsystems in which classical controller design techniques can be applied. Control of nonlinear systems based on the multiple model strategy has been a popular method in the last decades. In [2] a review of the various methods of multiple models based on the classical control theory, statistical methods and fuzzy structures are presented. Various control theories based on multiple models have been evaluated on a wide range of industrial systems. Li and Jilkov [6] provide a comprehensive review of methods and applications of multiple model systems. Model bank selection addresses the number of models and their constructions. Nevertheless, this issue has not been generally solved. In recent years, many efforts have been made to overcome this problem. Accordingly, many various techniques are proposed. One of the simplest methods is decomposition of the operating region into smaller spaces. This decomposition may be based on system parameters [7] or operating space between models [8], but a priori knowledge of the model is required. In [9] the number of local models according to the accuracy of the models are determined automatically. In [10] statistical methods are used to generate the model. In [11] local models are constructed using data clustering and least squares. In [12] the self-organizing map (SOM) is used to clustering of linear model parameter space. The gap metric concept was introduced by Zames and El-Sakkary and is applied to set up a method for determining the uncertainty range to ensure closed-loop stability [13,14].
http://dx.doi.org/10.1016/j.isatra.2017.01.021 0019-0578/& 2017 ISA. Published by Elsevier Ltd. All rights reserved.
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
2
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
El-Sakkary showed that the gap metric was much more suitable to measure the distance between two linear systems than other metrics [13]. Later, Vinnicombe introduces the v-gap metric that has the same topology as the gap metric, but is less conservative than the gap metric [15]. These two metrics uses H∞ norm that may cause inaccuracies in measuring the distance between the two systems in all range of frequencies [16]. To overcome this deficiency, Hosseini et al. [17] introduce the H-gap metric that uses the H2 norm to compare two models in some desired frequency range. Galán et al. [18] used gap-metric for the first time to obtain a reduced model bank. Du et al. [19] proposed a gap metric based method to decompose the operating range and obtained a reduced linear model bank., [20] and [16] used the gap metric and the H-gap metric respectively to compute the weights of local controllers and design a multi-model scheduling algorithm. Also, Du and Johansen (2014) employed the gap metric to compute weighting functions for the local controller combination. Nonlinearity measure is introduced by [21]. It measures nonlinearity based on comparison of the output of nonlinear plant and output of its approximated linearized model. Then, various indices are presented to measure the nonlinearity [22,23]. Nonlinearity measure is categorized in two groups: model-based and signalbased [24]. Model-based approach is based on the best linear approximation of a nonlinear system. In this method, the nonlinearity of the dynamic system is defined as the normalized largest difference between the nonlinear process and the best linear approximation of the process [25]. The advantage of the modelbased approach is that it requires fewer observations and the disadvantage is that it relies on the accuracy of the linear and nonlinear models [22]. Signal-based method, also named nonparametric nonlinearity measure, is based on direct analysis of the input/output or only output measurement of nonlinear plant. This method is studied both in frequency and time domain [22]. The proposed model set design is applied to the multiple model generalized predictive controller (MMGPC). Generalized predictive controller (GPC) is a type of model predictive control (MPC) that has been widely implemented as an advanced controller in the process industries. GPC gives the possibility of tuning of the control parameters and prediction horizons and adding of the control signal rate to cost function as well as explicitly considering constraints [26]. GPC was first proposed by [27]. The multiple model controller is widely used in MPC. For example, in [28] multiple model-based predictive control (MMPC) is used for pH process control. MMPC design is proposed in [29] for linear MIMO systems. The pH control process is used as the experimental case study. The pH process was modeled by Seborg for a batch-type process based on the Wiener model [30]. Ylén [31] presents a review on many different ways of pH process control. In several cases, multiple model-based controller is used for pH process control. In [7] a multiple adaptive PID controller is designed for this purpose. In [18] the Gap metric is used to multi-modeling of the pH process and then a cascade PI controller is designed. In this paper, a novel method is proposed to design the model set using the nonlinearity measure and gap metric. First, we develop a weighted partitioning method that selects a finite collection of operating point using the nonlinearity measure. Then, for each selected operating point a first order plus dead time (FOPDT) model is constructed. After that we introduce a new method to eliminate the similar models in model bank using clustering algorithm based on the nearest neighborhood and H-gap metric. As a result, a model bank with the smaller model number is provided. Then, we use this model bank to design a controller bank using GPC. Finally, the effectiveness of the proposed algorithm is tested on the pH neutralization process as a difficult nonlinear process in the industry.
The paper is organized as follows. In Section 2, a summary of the mathematical relations of the nonlinearity measure, the H-gap metric and a brief introductory of multiple model control are expressed. In Section 3, the algorithm of model bank selection is presented. Details of each step of the algorithm for model bank selection, controller design and switching supervisor design are studied through a simulation example in Section 4 and then the experimental results of MMGPC on a pH neutralization process are presented in Section 5. Section 6 contains the conclusions.
2. Preliminary 2.1. Nonlinearity measure Nonlinearity measure index (NLI) is used to quantify the nonlinearity degree of a system. One of the NLIs is Bi-spectrum analysis [22]. Based on the normalized higher order cross correlation function, the following indices can be defined [24]
(
)
2 NLI1 = max ϕNyuu (ω1, ω2) ω1, ω 2
2 2 − E⎡⎣ ϕNyuu (ω1, ω2)⎤⎦ + 2σ ⎡⎣ ϕNyuu (ω1, ω2)⎤⎦
(
)
(1)
where E[ . ] and σ[ . ] are expected value and variance, respectively, and ϕNyuu(ω1, ω2) is normalized higher order cross correlation between input u( t ) and output y( t ), which is defined as 2 (ω1, ω2) = ϕNyuu
E⎡⎣ U ( ω1)U ( ω2)Y *( ω1, ω2)⎤⎦ 2 2 2 E⎡⎣ U ( ω1) U ( ω2) ⎤⎦. E⎡⎣ Y ( ω1, ω2) ⎤⎦
(2)
where X ( ωi ) and X ( ωi ) is Fourier transform and power spectrum of x(t) at frequency ωi , respectively, and * stands for conjugate complex. This index has the following properties 1. The value of the index lies in the interval [0,1] where “0” implies linear and “1” refers to sever nonlinear behavior. 2. Is frequency independent in the case of linear time invariant (LTI) systems. 3. Shows the second order high frequency nonlinear behavior. 4. Shows robustness to the level of measurement noise. 5. Requires a high volume of data. 6. It is only applicable for Gaussian inputs. 2.2. H-gap metric There are several ways for measuring the distance between two linear models. In [17] a new metric is presented to measure the similarity of models. By definition, the H-gap metric for two single input single output (SISO) dynamic models P1 and P2 is calculated as follows 1
δhv(P1, P2) =
⎛ 1 ωH ⎞2 ⎜ ∫ Dis(ω)dω⎟ ⎝ 2π ωL ⎠ ⎛ 1 ωH ( 1 + ⎜ ∫ ⎝ 2π ωL
1
) d ⎞⎟ 2
max( Dis(ω)) 2
ω ⎠
(3)
where
⎛ P2 − P1 Dis(ω) = ⎜⎜ ⎝ 1 + P1 2 . 1 + P2 2
⎞2 ⎟⎟ ⎠
(4)
is the distance of models P1 and P2 at every frequency ω .
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
3
After decomposition of the nonlinear region into local linear regions, a FOPDT model is constructed for each region. The FOPDT model is described as
G(z −1) =
Fig. 1. The structure of the multiple model.
2.3. Multiple model method The basic idea of control based on the multiple model is constructing a model bank that its models are determined based on different operating conditions of the system. A general structure of the multiple model controller is presented in Fig. 1. This structure consists of two main loops: outer loop that contains the supervisor and model bank and the inner loop that contains the controller and process. The controller could be any classical types. For each model, an appropriate controller can be designed off-line. Supervisor observes the process behavior continuously and compares it with the behavior of pre-determined related models based on an appropriate cost function.
3. Model bank selection This section introduces the proposed algorithm for model bank selection. The structure of models is selected in the form of FOPDT that is a simple yet popular model in process control. In this algorithm, first a nonlinearity index is employed to decompose the nonlinear operating region of the process to linear local operating points. Then for each of these operating points a FOPDT model is constructed. A process could have similar behavior in different and separated operating points. So it is most probably that some of these linear models are similar. The H-Gap metric is used to eliminate similar models using a clustering algorithm. The remaining linear models constitute the minimum number of models in the model bank. More details are proposed in the following subsections.
b1z −1 z −d 1 + a1z −1
(5)
where b1 and a1 are process model transfer function parameters and z is discrete transform operator. d is delay step of the plant and assumed to be integer constant. The model parameters are estimated using input-output data of the process. The last step is to select the proper minimum number of models. This is achieved by means of the clustering and H-gap metric. Clustering means partitioning of a collection of data into disjoint subsets or clusters so that maximum resemblance of members within each cluster is achieved. In this step, using clustering algorithm based on the nearest neighborhood clustering algorithm [32] and the H-gap metric, local models are divided into different clusters. Then, between these clusters, centers are selected for the model bank and other members are eliminated. 3.2. Model bank selection algorithm In this section, we propose the details of the algorithm for model bank selection. As stated in the previous section, this algorithm measures input-output data and uses the nonlinearity measure and H-gap metric tools to select a model bank for a nonlinear system. The algorithm comprises of 5 main steps as follows: Step 1. Enter the start and end of the operating region as Region=[Rstart , R end], compute the first operating point and it's NLI as (6) and (7) using the Eq. (1).
⎧ (R + R end ) ⎫ ⎬ Op = ⎨ start ⎩ ⎭ 2
(6)
{ NLI( Op1)}
(7)
NLI =
Step 2. Decompose the region Op to multi-linear regions through the sub-algorithm1. Step 3. Construct FOPDT models using least square (LS) for operating points of vector Op and put them into vector G, so named the primary model bank.
G = {g1, g2 , ... , gn}
(8)
3.1. The general structure of the presented algorithm The whole nonlinear operating space of the process is divided into some linear regions by means of the nonlinearity index. The start and end of operating point space are required as the prior knowledge. Then a middle point is selected as the operating point. The process is excited by an amplitude modulated random binary signal (APRBS) input with the amplitude equal to the half of the operating range. Afterward the nonlinearity measure is calculated for this operating point and is compared with a pre-set threshold. If it is smaller than the threshold, algorithm is stopped; otherwise the operating region is broken into smaller areas with respect to the NLI of every operating point. Details of the procedure to select the broken points are stated in Step 2 of section 3.2. Again, the NLI is calculated for these regions. If the calculated value for the index of either is greater than the given threshold, then the region is decomposed into smaller regions. This procedure is continued until all of the NLI values of operating points are less than the threshold value.
Step 4. Using Eq. (3) and measure the distance between each pair of linear models, construct the symmetric pairwise distance matrix Gap as (9):
Gap = Hgapmetric (gi , gj );
i, j = 1, ... , n
(9)
Step 5. Cluster the primary Model bank (8) using sub-algorithm2 and find the reduced final model bank. Details of the sub-algorithm1 of step 2 are as follows: Sub-algorithm1: Decomposition of nonlinear region to multilinear region
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
4
Table 1 Model parameters. Parameter
Notation
Value
Acid flow rate Stream flow rate Acid concentration Base concentration Acid equilibrium constant Base equilibrium constant Height of Solution in stirred tank Cross section of stirred tank
Fa Fo Ca Cb Ka Kb H A
5 ml/s 30 ml/s 0.3 mol 0.3 mol 4.47 * 10–14 5.62 * 10–14 10 cm 90 cm2
Table 2 Applied thresholds for the model set selection (simulation).
Fig. 2. pH process with continuous stirred tank reactor.
Step 2-1. Set i = 0 and j = 0. Step 2-2. Set n as the number of current operating points: n = length(Op). Step 2-3. Set i = i + 1. If i > n go to Step2-4, else go to Step2-5. Step 2-4. If j = 0 go to Step2-11 else put i = 1 , j = 0. Step 2-5. If NLIi > ThNLI go to Step2-6 else go back to Step2-3. Step 2-6. Add two new operating points before and after the current operating point as (10) and (11). Also, compute SNR as (12).
Parameter (Threshold)
Notation
Value
Nonlinearity measure SNR Distance of two adjacent operating points
ThNLI ThSNR ThOp
0.15 0.5 0.1
H-gap metric
ThGap
0.25
NLI = {NLI1, ... , NLIi − 1 ,
before NLInew,
NLIi , NLIi + 1 , ... ,
NLIn}
Opi − 1 . NLIi − 1 − Opi . NLIi = NLIi − 1 + NLIi
before Opnew
after Opnew
SNRnew =
(10)
where ThOp is a threshold for the minimum distance between two consecutive operating points and ThSNR is a threshold for the bound of measurement noise as defined in the Table 2. Step 2-8. If j = j + 1 and
Opi + 1 . NLIi + 1 − Opi . NLIi = NLIi + 1 + NLIi
(11)
2 σNoise
− Opi > ThOp and
Op = {Op1, ... , Opi − 1 , Opi ,
after Opnew ,
SNRnew > ThSNR put
Opi + 1 , ... , (16)
(12) after NLInew
before Opnew
Op = {Op1, ... , Opi − 1 ,
before NLInew
after Opnew
Opn }
2 σSignal
Step 2-7. If j = j + 1 and
(15)
− Opi > ThOp and
before Opnew ,
= NLI( before Opnew )
SNRnew > ThSNR put
Opi , Opi + 1 , ... , Opn }
(13)
(14)
= NLI( after Opnew )
NLI = {NLI1, ... , NLIi − 1 , NLIi ,
(17) after NLInew,
NLIi + 1 , ... ,
NLIn}
(18)
Step 2-9. Due to adding a new operating point(s) and change in the span of other operating points, compute the NLI again for Opm where m = {i − 2 , i − 1 , i , i + 1 , i + 2}. Step 2-10. Go back to Step2-3. Step 2-11. Stop, read Op and return to Step 3. Details of the sub-algorithm 2 of step 5 are as follows:
Fig. 3. The control structure of the pH process.
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
5
Titration curve of pH model
pH
10
8
6
0
0.5
1
1.5 2 2.5 Base feed rate/Acid feed rate Derivative of Titration curve
3
3.5
4
10
11
12
150 X: 7.499 Y: 142.1
Gain
100
50 X: 4.349 Y: 1.196
0
4
5
6
7
8 pH
9
Fig. 4. a (upper): The titration curve of pH process model, b (lower): derivative of the titration curve.
Fig. 5. Some steps of the region decomposition procedure.
Sub-algorithm2:Model clustering and selection Step 5-1. Set k = 1, v = 1, z = 1. Set C as the index of models that are chosen as the center of clusters. Take XC as a matrix that it's Ck th row is the index of models in kth cluster in C. Create a cluster with center G1 ( G1 = g1) and set C1 = { 1 }. Step 5-2. Set k = k + 1, z = z + 1. If z > n go to Step 5-8 else
go to Step 5-3. Step 5-3. Measure the distance between Gz and each cluster center as follows:
d(w ) = Gap( Gz , GCw ); w = 1, 2, ... , v
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
(19)
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
6
Define ind as the index of the cluster that has minimum distance from Gz and min _d as the value of minimum distance. Step 5-5. If min _d > Thgap , take Gz as the center of the new cluster: Ck = { z}. Set v = v + 1. Step 5-6. If min _d < Thgap , put gz in indth cluster: Xind = z . Step 5-7. Go back to Step 5-2. Step 5-8. Stop and read C as the index of the reduced model bank.
4. Simulation results In this section, the strategy described in the previous section is applied to the model of a pH neutralization process. Nonlinear model of this process is presented in [33]. Results are presented in terms of MPC comparison designed by the proposed algorithm and two other controller design algorithms: 1-Multiple PID controller and 2-Multiple model controller presented by [19]. 4.1. pH neutralization process This process consists of a continuous stirred tank reactor (CSTR) with a constant volume as shown in Fig. 2. Acid and base enter the tank with the flow rate Fa and Fb with concentration of Ca , Cb , respectively, and water has the flow of Fw . Mixed product discharges from the tank with the flow rate of Fo. The control objective is to tune the pH value of the outlet stream by manipulating the base flow rate. The value of the acid flow rate is considered as the disturbance. Also, tank level is controlled by manipulating the water flow in another decoupled loop control using a PI controller as illustrated in Fig. 3. The acid flow rate is constant and has been considered to be 5 ml/s . Other parameters of the system are presented in Table 1. One of the important characteristic of the pH process is a static curve, also called the titration curve. The titration curve is obtained from the steady state variation of the output (pH) to the steady state variation of the input (base flow). The titration curve of this process is shown in Fig. 4a. Fig. 4b is the derivative of the titration curve that represents the gain versus pH value. As can be seen, range of gain variation is from around 2 to around 142, which cause difficulty in control of this process.
Table 3 The primary model bank obtained from the proposed algorithm (simulation). #
Op. (pH)
NLI
Gain
1 2 3 4 5 6 7 8 9 10 11 12
5.63 6.25 6.65 6.91 7.37 7.50 7.57 7.86 8.06 8.22 8.36 8.75
0.04 0.04 0.10 0.12 0.09 0.04 0.09 0.05 0.09 0.04 0.08 0.03
4.87 17.4 44.1 75.48 142.25 142.15 135.94 97.38 69.4 50 36.56 15.7
Table 4 Identified models and corresponding steady state gain (simulation). #
Op. (pH)
1
5.63
0.1z−1
6.25
1 − 0.977z−1 0.27z−1
6.65
1 − 0.977z−1 0.806z−1
7.37
1 − 0.975z−1 2.88z−1
2 3 4
1 − 0.974z−1
Gain 4.87
z −5
17.4
z −5
44.1
z −5
142.25
z −5
Table 5 Control and prediction horizon (simulation). #
Op. (pH)
N2
Nu
1 2 3 4
5.63 6.25 6.65 7.37
101 100 92 91
86 85 77 76
Table 6 The optimal PID controller coefficient value. #
Op. (pH)
4.2. Model bank selection Sampling time is selected to be 3 s. In the simulations a noise signal with variance of 0.001 is added at the pH output. The upper and lower limits of the operating region are equal to pH ¼5 and pH ¼10, respectively. Other parameters and thresholds are given in the Table 2. In the first stage of the algorithm, we obtain local linear regions of the overall operating space through the implementation of sub algorithm 1 of Section 3.2. The results outcome of some steps are shown in Fig. 5. Finally, local operating points are obtained as shown in the Table 3. In the next step, FOPDT model is constructed for the operating points presented in the Table 3. Then, models are classified into the clusters with the similar distances by using clustering algorithm and the H-gap metric, as stated in sub algorithm 2 of Section 3.2. Table 4 presents the transfer function of the final model bank. In Fig. 4, the final model bank is marked by red points on the titration curve. As can be seen, models are selected on one side of the curve. It means that the algorithm detected the symmetry on the titration curve.
FOPDT Model
1 2 3 4
5.63 6.25 6.65 7.37
Cof. kp
ki
kd
0.38 0.11 0.03 0.01
4*10-3 1*10-3 5*10-4 1*10-1
3.08 0.86 0.34 0.11
4.3. Controller bank design GPC is the selected controller [34]. In the GPC, the model of the system is assumed to be as follows
A(z −1)y(t ) = B(z −1)z −du(t − 1) + C (z −1)
e(t ) Δ
(20)
where y(t ), u(t ) and e(t ) are the model output, model input and zero-mean white noise, respectively. Also A, B and C are the output dynamic, input dynamic and noise dynamic polynomial, respectively. Δ is the difference operator. Each GPC is designed to optimize the following cost function N2
J (t ) =
∑ j = N1
⎡ y^ (t + d + j ) − w (t + j )⎤2 + λ ⎣ ⎦
Nu− 1
∑ j=1
⎡⎣ Δu(t + j − 1)⎤⎦2
(21)
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
7
Multiple Model GPC in compare with Multiple Model PID, with Small reference changes 10 y
MMGPC
pH
8
MMPID
6 2000
4000
6000
8000
10000
12000
14000
16000
18000
12000
14000
16000
18000
12000
14000
16000
18000
Control Signal
Base feed
3 2 1 0
2000
4000
6000
8000
10000
Switching signal Model Index(pH=)
7.5 7 6.5 6 5.5 2000
4000
6000
8000
10000
time (second)
Fig. 6. a) Output, b) control signal & c) switching signal for MMGPC and MMPID for small change of set-point. Multiple Model GPC in compare with Multiple Model PID, with Medium reference changes 10
pH
8 y MMGPC
6
MMPID
2000
4000
6000
8000
10000
12000
14000
16000
18000
12000
14000
16000
18000
12000
14000
16000
18000
Base feed
Control Signal
2 1 0
2000
4000
6000
8000
10000
Model Index(pH=)
Switching signal 7.5 7 6.5 6 5.5
2000
4000
6000
8000
10000
time (second)
Fig. 7. a) Output, b) control signal & c) switching signal for MMGPC and MMPID for medium change of set-point.
where w (t ) is the desired reference trajectory. λ is the suppressing coefficient of the control signal. N1 is the lower bound of the prediction horizon and is usually set to d + 1. The upper bounds of the prediction and control horizons are denoted by Nu and N2, respectively. Eq. (21) can be written as the vector form as follows
J = (Gu + f−w)T (Gu + f−w) + λ uT u
(22)
J=
1 T u Hu + bT u + f0 2
(24)
where H = 2(GT G + λI ), bT = 2(f−w)T G and f0 = (f−w)T (f−w). With optimization of cost function J with respect to input signal u , control signal is reached
where w is the vector of the feature reference trajectory given as follows
u = − H−1b = (GT G + λI)−1GT (w−f)
T w = ⎡⎣ w (t + d + 1) w (t + d + 2) ⋯ w (t + d + N )⎤⎦
Eq. (25) is the optimized vector of the control signal that is computed at each sample but only the first array of this vector is applied to the system. The GPC bank is designed based on the model bank in the Table 4. The prediction and control horizon parameters of the controllers are given in the Table 5.
(23)
f is the free response vector and G is a matrix extracted from the step response of the plant. Details derivation of these matrices can be found in [34]. Sorting of Eq. (22) in terms of u , we have a quadratic cost function as the form of Eq. (24).
(25)
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
8
Multiple Model GPC in compare with Multiple Model PID, with High reference changes 10
y MMGPC
pH
8
MMPID
6 2000
4000
6000
8000
10000
12000
14000
16000
18000
12000
14000
16000
18000
12000
14000
16000
18000
Control Signal
Base feed
4 3 2 1 0
2000
4000
6000
8000
10000
Model Index(pH=)
Switching signal 7.5 7 6.5 6 5.5 2000
4000
6000
8000
10000
time (second)
Fig. 8. a) Output, b) control signal & c) switching signal for MMGPC and MMPID for high change of set-point.
Titration curve and Model bank based on Du algorithm
pH
10
8
6
0
0.5
1
1.5 2 2.5 Base feed rate/Acid feed rate
3
3.5
4
10
11
12
Derivative of Titration curve 150
Gain
100
50
0
4
5
6
7
8 pH
9
Fig. 9. Final model bank (marked by red points) obtained by Du's algorithm.
4.4. Switching supervisor As it was mentioned in Section 2.3, the multiple model controller at each sampling time requires a switching supervisor to detect the operating region of the process and then select the appropriate controller in the controller loop. Switching supervisor uses a cost function as follows in order to select the suitable model [35]
cost function is selected as the best model. In order to reduce the chattering in the switching process, we use switching with hysteresis logic [36]. In this method, switching is not performed unless the cost-function of the ith model ( Ji (t )) is some percent less than the current selected model. 4.5. Comparison study
(26)
Simulation results are presented in two parts as follows in order to compare the proposed controller and some other multiple model controllers.
where a and ρ are weighting parameters and ei( . ) is the difference between the process output and the ith linear model. The first term is instance squared error and the second is the sum of weighted squared errors. The cost function (26) is calculated for each model at every sampling time. The model that has the lowest
4.5.1. Comparison with the multiple model PID controller(MMPID) In this section, the proposed MMGPC is compared with the MMPID controller presented in [37]. In this design, three constraints are considered; phase margin on cutoff frequency, restriction on gain crossover frequency and derivative of phase
t
Ji (t ) = a. ei2(t ) +
∑ ρt − k ei2(k) k=0
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
diagram at the gain crossover frequency. Parameters of the PID controllers are shown in Table 6. Simulation results are given in Figs. 6–8 for different step change. As can be seen, in all of three cases, there are some operating points that MMPID controller is unstable while MMGPC has more appropriate response. 4.5.2. Comparison with the model bank selection of Du [19] In this comparison, the proposed MMGPC is compared with another MMMPC that its model bank is designed according to the algorithm presented by [19]. According to [19], the operating region is uniformly divided into a large number (50 or more) of operating points. Then similar models are removed based on the gap metric criterion. In Fig. 9, the final operating points obtained by this algorithm are shown, when the whole operating space was initially modeled by 50 models. Note that the threshold of gap
9
metric is selected so that the number of models obtained from Du's algorithm is the same as the proposed algorithm. Simulation results are given in Figs. 10–12 for different step change. As can be seen, due to improper selection of initial operating points in the model bank, Du's algorithm has not desirable tracking performances. This illustrates the effectiveness of the nonlinearity index to appropriate selection of operating points. In fact, NLI introduces the appropriate set of linear behavior subspaces, to the clustering algorithm to select the model bank.
5. Implementation results The experimental apparatus, shown in Fig. 13, is a continuous pH neutralization pilot plant in the Process Control laboratory at K. N. Toosi University of Technology. A schematic diagram of the
Multiple Model GPC with two different algorithm, with Small reference changes
pH
10 Prop. alg.
8
Du alg.
6 2000
4000
6000
8000
10000
12000
14000
16000
18000
12000
14000
16000
18000
12000
14000
16000
18000
Control Signal
Base feed
3 2 1 0
2000
4000
6000
8000
10000
Model Index(pH=)
Switching signal 9 8 7 6 2000
4000
6000
8000 10000 time (second)
Fig. 10. a) Output, b) control signal & c) switching signal of MMGPC obtained by proposed algorithm and Du's algorithm for small change of set-point.
Fig. 11. a) Output, b) control signal & c) switching signal of MMGPC obtained by proposed algorithm and Du's algorithm for medium change of set-point.
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
10
Multiple Model GPC with two different algorithm, with High reference changes 10
pH
8 Prop. alg.
6
Du alg.
2000
4000
6000
8000
10000
12000
14000
16000
18000
12000
14000
16000
18000
12000
14000
16000
18000
Control Signal
Base feed
2 1.5 1 0.5 0
2000
4000
6000
8000
10000
Model Index(pH=)
Switching signal 9 8 7 6 2000
4000
6000
8000 10000 time (second)
Fig. 12. a) Output, b) control signal & c) switching signal of MMGPC obtained by proposed algorithm and Du's algorithm for high change of set-point.
Fig. 14. The schematic diagram of the pH pilot plant [38].
Table 7 Tunable parameters of pH process pilot plant. Fig. 13. pH process pilot plant.
system is shown in Fig. 14. It contains a continuous stirred tank reactor (CSTR) with a motorized mixer in order to have a wellmixed tank. The effluent stream flow rate is regulated by a manual valve which affects transportation delay and time constant of the process. A pH sensor is located after the outlet valve to measure the pH of the effluent stream. Acid, base, and water are pumped into the CSTR by the corresponding dosing pumps. The control objective is to achieve the desirable pH values by manipulating the base flow rate in the presence of changes in acid stream flow rate as a disturbance to evaluate ability of different control strategies. The advanced pH controller and level controller are performed by
Parameter
Notation
Value
Acid flow rate
fa
0.45 ml/s
Acid concentration pH of acid solution Base concentration pH of base solution Stream flow rate
Ca – Cb – fo
1.5 mlit./lit. 3 1 gr./lit. 12 0.45 mlit/s
Height of Solution in stirred tank Sample time
H Ts
15 cm 3s
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎ Table 8 Applied thresholds for the model set selection (implementation).
11
Table 10 The final model bank in experimental study.
Parameter (Threshold)
Notation
Value
#
Op. (pH)
Time Constant (s)
Gain
Delay (s)
Nonlinearity measure SNR Distance of two adjacent operating points
ThNLI ThSNR ThOp
0.35 0.5 0.5
H-gap metric
ThGap
0.25
1 2 3 4 5
4.32 4.94 5.45 6.41 7.55
462 453 384 462 390
0.16 0.29 0.42 0.73 2.52
66 87 48 42 36
the MATLAB™ software and SIMULINK™; in addition, a PLC is installed for data acquisition which communicates to the computer through a PROFIBUS cable. Some of the tunable parameters of the pH pilot plant are shown in Table 7. This apparatus is severely nonlinear with ample source of disturbances and uncertainty. The static gain changes by a ratio of more than 10 even for a very simplified experiment. One of the main sources of disturbances is the effect of the minerals in the water tab, which is used instead of buffer in our experience. In addition, acid and base tanks are filled using some materials that prepared in the lab manually. So, the concentration of materials changes every time that they are prepared. Adding to that neutralization of acid in the source tank in connection with air, change in the temperature of the room and change in noise variance makes this apparatus one of the most challenging problem in process control study. 5.1. Model bank selection In this section, an appropriate model bank is designed using the proposed algorithm. Threshold parameters are assumed as presented in Table 8. NLI1 needs a considerable amount of data. Due to the inherent time invariance of the process, for instance, when the base and acid tanks are filled by the new source of materials, it is not possible to acquire enough experimental data for the nonlinearity index from this process. So, in the experimental study, we use the cross correlation index to compute the nonlinearity measure. This
Table 11. Control and prediction horizon. #
Op. (pH)
N2
Nu
1 2 3 4 5
4.32 4.94 5.45 6.41 7.55
44 34 27 25 22
15 12 11 11 10
Table 12 Weighted coefficient of control signal ( λ ). # Op. (pH)
1 2 3 4 5
4.32 4.94 5.45 6.41 7.55
λ Small change set-point
Medium change set-point
High change set-point
0.05 0.17 .35 2.66 6.35
0.38 1.26 3.35 5.33 76.2
0.38 1.26 3.35 5.33 50.8
index is defined as follows [22] 2 ⎛ ⎞ ϕyu(jω) ⎜ ⎟ NLI2 = 1 − mean⎜ ω ϕuu(jω)ϕyy(jω) ⎟ ⎝ ⎠
(27)
Table 9 The initial linear local regions. pH NLI2
4.11 0.16
4.32 0.15
4.94 0.11
5.14 0.15
5.45 0.18
5.8 0.14
6.33 0.11
6.41 0.16
6.58 0.16
7.1 0.23
7.44 0.25
7.55 0.33
7.9 0.12
Fig. 15. a) Output, b) control signal & c) switching signal of MMGPC for the small change of set-point.
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
12
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
Fig. 16. a) Output, b) control signal & c) switching signal of MMGPC for the medium change of set-point.
Fig. 17. a) Output, b) control signal & c) switching signal of MMGPC for the high change of set-point.
where ϕyu( jω) = Y ( jω)U *( jω) is the cross-correlation function. Using the proposed algorithm in Section 3 and consideration of threshold parameters as presented in Table 8, we obtained 13 initial linear regions. The operating points of these regions and their NLI indices are shown in Table 9. Finally, after applying the clustering algorithm, reduced model bank is obtained as presented in Table 10.
Table 13 Time domain properties of response of the implemented MMGPC on pH pilot plant for high change of set-point presented in Figs. 15-17. Setpoint Situation
MSE
Av. Settling time (s)
Av. Overshoot (%)
small change medium change high change
0.014 0.064 0.364
1007.1 3691 3582
26.1 123.9 78.38
5.2. Controller bank design and switching supervisor For each model of the model bank, a GPC is designed. Parameters of these controllers are shown in Table 11. Weighted coefficients of control signal λ with respect to the change in setpoint value (low, medium and high change) have been achieved experimentally as Table 12. Parameters of cost function (26) are chosen to be α = 0.1 and η = 0.9. Three modes of low, medium and high change in the reference signal are studied experimentally. Hysteresis constant (h) for these three modes are 0.7, 0.05 and 0.05, respectively. Results are shown in Figs. 15–17. Time domain properties of set points, are shown in Table 13. As can be seen, the performance of
the control system is satisfactory for small and medium step change except when there is a step from a medium steady state gain operating point to a high gain operating point. This is predictable, since even a step change in pH value from around 6.5 to 7.5 means 3 times changes in the steady state gain. The high value of over-shoot is avoidable by reducing the speed of the desired closed loop system. Fig. 17 shows that the control system provides a stable closed loop system even in large step test. However, a step change from pH of around 5 to 7.5 means a steady state gain change by a factor of around 10, which is a very large value of uncertainty for every control system. Again this large overshoot is
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
avoidable sacrificing the speed of the closed loop system.
6. Conclusion In the context of the controller design based on the multiple model technique, the model bank selection problem is one of the main issues that is still under study for a general yet reliable solution. This problem is even more challenging in the case that there is not any first principle model of the processes. In this paper, an algorithm is proposed that can provide an appropriate model bank using the nonlinearity measure and the H-gap metric indices, using only the input/output data of the process. The operating space of the process is divided into some small linear behavior sub-spaces using the proposed nonlinearity measure. This may leads to some excessive local models. So, the next step a clustering algorithm based on H-gap metric is employed to prune the model bank. The simulation and implementation results show the functionality of the algorithm. The proposed method gives a systematic solution for the modeling and control of a large class of nonlinear stable systems with wide operating ranges. For the highly nonlinear pH neutralization process with plenty of source of disturbances and uncertainty, this method provides a stable experimental control system in all operating points at step changes with proper performance.
References [1] Camacho EF, Bordons C. Model predictive control. 2nd ed. London: SpringerVerlag; 2004. [2] Murray-Smith R, Johansen TA. Multiple model approaches to modeling and control. London: Taylor & Francis; 1997. [3] Scokaert PO, Mayne DQ, Rawlings JB. Suboptimal model predictive control (feasibility implies stability). IEEE Trans Autom Control 1999;44:648–54. [4] Johansen TA, Foss BA. Operating regime based process modeling and identification. Comput Chem Eng 1997;21:159–76. [5] Rodriguez J, Romagnoli J, Goodwin G. Supervisory multiple regime control. J Process Control 2003;13:177–91. [6] Li XR, Jilkov VP. Survey of maneuvering target tracking. Part V. Mult-model methods Aerosp Electron Syst, IEEE Trans 2005;41:1255–321. [7] Böling JM, Seborg DE, Hespanha JP. Multi-model adaptive control of a simulated pH neutralization process. Control Eng Pract 2007;15:663–72. [8] Aufderheide B, Bequette BW. A variably tuned multiple model predictive controller based on minimal process knowledge. In: Proceedings of the 2001 American control conference. IEEE 2001. p. 3490–3495; 2001. [9] Peymani E, Fatehi A, Sedigh AK. Automatic learning in multiple model adaptive control. International Conference on Control, UKACC, Manchester, UK; 2008. [10] Liang Y, Li X, Han C, Duan Z. A general systematic method for model-set design. Aerosp Electron Syst, IEEE Trans 2012;48:1505–20. [11] Novák J, Chalupa P. , Bobál V. Modelling and predictive control of a nonlinear system using local model network. In: Proceedings of the World Congress. p.
13
1622–7; 2011. [12] Fatehi A, Abe K. Flexible structure multiple modeling using irregular self-organizing maps neural network. Int J Neural Syst 2008;18:233–56. [13] El-Sakkary A. The gap metric: robustness of stabilization of feedback systems. Autom Control, IEEE Trans 1985;30:240–7. [14] Zames G, El-Sakkary A. Unstable systems and feedback: the gap metric. Proc Allerton Conf 1980:380–5. [15] Vinnicombe G. Frequency domain uncertainty and the graph topology. Autom Control, IEEE Trans 1993;38:1371–83. [16] Hosseini S, Fatehi A, Johansen TA, Sedigh AK. Multiple model bank selection based on nonlinearity measure and H-gap metric. Journal of Process Control; 2012. [17] Hosseini S, Fatehi A, Khakisedigh A A new metric for measuring the distance between two linear systems. Control, Automation and Systems (ICCAS). In: Proceedings of 11th international conference on. IEEE; 2011. p. 1552–1557. [18] Galán O, Romagnoli JA, Palazoglu A, Arkun Y. Gap metric concept and implications for multilinear model-based controller design. Ind Eng Chem Res 2003;42:2189–97. [19] Du J, Song C, Li P. Application of gap metric to model bank determination in multilinear model approach. J Process Control 2009;19:231–40. [20] Arslan E, Çamurdan MC, Palazoglu A, Arkun Y. Multimodel scheduling control of nonlinear systems using gap metric. Ind Eng Chem Res 2004;43:8275–83. [21] Desoer C, Wang Y-T. Foundations of feedback theory for nonlinear dynamical systems. Circuits Syst, IEEE Trans 1980;27:104–23. [22] Hosseini M, Johansen T, Fatehi A. Comparison of nonlinearity measures based on time series analysis for nonlinearity detection. Model Identif Control 2011;32:123–40. [23] Le Caillec J-M, Garello R. Comparison of statistical indices using third order statistics for nonlinearity detection. Signal Process 2004;84:499–525. [24] MAAS Choudhury, Shah SL, Thornhill NF. Diagnosis of Process Nonlinearities and Valve Stiction: data Driven Approaches. Springer; 2008. [25] Helbig A, Marquardt W, Allgöwer F. Nonlinearity measures: definition, computation and applications. J Process Control 2000;10:113–23. [26] Qin SJ, Badgwell TA. A survey of industrial model predictive control technology. Control Eng Pract 2003;11:733–64. [27] Clarke DW, Mohtadi C, Tuffs P. Generalized predictive control – Part I. The basic algorithm. Automatica. 1987;23:137–48. [28] Novák J, Chalupa P, Bobál V. Multiple model modeling and predictive control of the pH neutralization process. Int J Math Model Methods Appl Sci 2011. [29] Bagheri P, Mardanlou V, Fatehi A. Multiple model predictive control of multivariable pH process using adaptive weighting matrices.In: Proceedings of The 18th IFAC World congress2011; 2011. p. 12366–71. [30] Hall RC, Seborg DE. Modelling and self-tuning control of a multivariable pH neutralization process part I: modelling and multiloop control. In: Proceedings of the American control conference. IEEE1989; 1989. p. 1822–1827. [31] Ylén J-P. Measuring, modelling and controlling the pH value and the dynamic chemical state. Helsinki University of Technology; 2001. [32] Theodoridis S, koutroumbas L. Pattern Recognit 1999. [33] Sadeghassadi M, Fatehi A, Khaki Sedigh A, Hosseini S. On the structural optimization of a neural network model predictive controller. Ind Eng Chem Res 2012;52:394–407. [34] Camacho EF, Bordons C. Model Predictive Control.London: Springer; 2004. [35] Narendra KS, Balakrishnan J. Adaptive control using multiple models. Autom Control, IEEE Trans 1997;42:171–87. [36] Liberzon D. Switching In Systems and Control. Springer; 2003. [37] Li H, Luo Y, Chen Y. A fractional order proportional and derivative (FOPD) motion controller: tuning rule and experiments. IEEE Trans Control Syst Technol 2010;18:516–20. [38] Peymani E, Fatehi A, Bashivan P, Sedigh AK. An experimental comparison of adaptive controllers on a pH neutralization pilot plant. In Proceedings of the India conference, 2008 INDICON Annual IEEE. IEEE2008; 2008. p. 377-382.
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
Contents lists available at ScienceDirect
ISA Transactions journal homepage: www.elsevier.com/locate/isatrans
Research article
Multi-linear model set design based on the nonlinearity measure and H-gap metric Davood Shaghaghi, Alireza Fatehi n, Ali Khaki-Sedigh APAC Research Group, Industrial Control Center of Excellence, Faculty of Electrical Engineering, K. N. Toosi University of Technology, 16317-14191 Tehran, Iran
art ic l e i nf o
a b s t r a c t
Article history: Received 18 November 2015 Received in revised form 8 November 2016 Accepted 16 January 2017
This paper proposes a model bank selection method for a large class of nonlinear systems with wide operating ranges. In particular, nonlinearity measure and H-gap metric are used to provide an effective algorithm to design a model bank for the system. Then, the proposed model bank is accompanied with model predictive controllers to design a high performance advanced process controller. The advantage of this method is the reduction of excessive switch between models and also decrement of the computational complexity in the controller bank that can lead to performance improvement of the control system. The effectiveness of the method is verified by simulations as well as experimental studies on a pH neutralization laboratory apparatus which confirms the efficiency of the proposed algorithm. & 2017 ISA. Published by Elsevier Ltd. All rights reserved.
Keywords: Nonlinearity measure H-gap metric Multiple-model controller Model bank selection Generalized predictive control
1. Introduction Classical control design techniques have been matured for linear systems to incorporate robustness and performance requirements within a narrow operating range. Nevertheless, for nonlinear processes with severe nonlinearities and large set-point changes, a single linear controller may not provide satisfactory performance. This has resulted in the development of nonlinear model predictive control (NMPC) that uses a comprehensive and precise model for process prediction and optimization [1]. The fundamental step in the development of NMPC is to elect an appropriate structure of a nonlinear model for identification of system dynamics. Unfortunately, the theory of nonlinear system identification has been not developed very well [2]. Furthermore, the resulting dynamic model may be too complex to be useful for NMPC design. The main disadvantage of the nonlinear modeling approach is that using a nonlinear model led to change the control problem solution from a convex quadratic programming to a nonconvex nonlinear case. The solution of this problem is much more difficult and there is no guarantee to find the global optimum [3]. An approach that tries to retain the features of a linear design and is applied to nonlinear systems is the multi-model approach for modeling and control [2]. The multiple model approach with switching and supervisory control is best suited for processes that n
Corresponding author. E-mail addresses: [email protected] (D. Shaghaghi), [email protected] (A. Fatehi), [email protected] (A. Khaki-Sedigh).
have fast changes in operating range [4] and [5]. The key concept is to represent the nonlinear system as a combination of linear subsystems in which classical controller design techniques can be applied. Control of nonlinear systems based on the multiple model strategy has been a popular method in the last decades. In [2] a review of the various methods of multiple models based on the classical control theory, statistical methods and fuzzy structures are presented. Various control theories based on multiple models have been evaluated on a wide range of industrial systems. Li and Jilkov [6] provide a comprehensive review of methods and applications of multiple model systems. Model bank selection addresses the number of models and their constructions. Nevertheless, this issue has not been generally solved. In recent years, many efforts have been made to overcome this problem. Accordingly, many various techniques are proposed. One of the simplest methods is decomposition of the operating region into smaller spaces. This decomposition may be based on system parameters [7] or operating space between models [8], but a priori knowledge of the model is required. In [9] the number of local models according to the accuracy of the models are determined automatically. In [10] statistical methods are used to generate the model. In [11] local models are constructed using data clustering and least squares. In [12] the self-organizing map (SOM) is used to clustering of linear model parameter space. The gap metric concept was introduced by Zames and El-Sakkary and is applied to set up a method for determining the uncertainty range to ensure closed-loop stability [13,14].
http://dx.doi.org/10.1016/j.isatra.2017.01.021 0019-0578/& 2017 ISA. Published by Elsevier Ltd. All rights reserved.
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
2
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
El-Sakkary showed that the gap metric was much more suitable to measure the distance between two linear systems than other metrics [13]. Later, Vinnicombe introduces the v-gap metric that has the same topology as the gap metric, but is less conservative than the gap metric [15]. These two metrics uses H∞ norm that may cause inaccuracies in measuring the distance between the two systems in all range of frequencies [16]. To overcome this deficiency, Hosseini et al. [17] introduce the H-gap metric that uses the H2 norm to compare two models in some desired frequency range. Galán et al. [18] used gap-metric for the first time to obtain a reduced model bank. Du et al. [19] proposed a gap metric based method to decompose the operating range and obtained a reduced linear model bank., [20] and [16] used the gap metric and the H-gap metric respectively to compute the weights of local controllers and design a multi-model scheduling algorithm. Also, Du and Johansen (2014) employed the gap metric to compute weighting functions for the local controller combination. Nonlinearity measure is introduced by [21]. It measures nonlinearity based on comparison of the output of nonlinear plant and output of its approximated linearized model. Then, various indices are presented to measure the nonlinearity [22,23]. Nonlinearity measure is categorized in two groups: model-based and signalbased [24]. Model-based approach is based on the best linear approximation of a nonlinear system. In this method, the nonlinearity of the dynamic system is defined as the normalized largest difference between the nonlinear process and the best linear approximation of the process [25]. The advantage of the modelbased approach is that it requires fewer observations and the disadvantage is that it relies on the accuracy of the linear and nonlinear models [22]. Signal-based method, also named nonparametric nonlinearity measure, is based on direct analysis of the input/output or only output measurement of nonlinear plant. This method is studied both in frequency and time domain [22]. The proposed model set design is applied to the multiple model generalized predictive controller (MMGPC). Generalized predictive controller (GPC) is a type of model predictive control (MPC) that has been widely implemented as an advanced controller in the process industries. GPC gives the possibility of tuning of the control parameters and prediction horizons and adding of the control signal rate to cost function as well as explicitly considering constraints [26]. GPC was first proposed by [27]. The multiple model controller is widely used in MPC. For example, in [28] multiple model-based predictive control (MMPC) is used for pH process control. MMPC design is proposed in [29] for linear MIMO systems. The pH control process is used as the experimental case study. The pH process was modeled by Seborg for a batch-type process based on the Wiener model [30]. Ylén [31] presents a review on many different ways of pH process control. In several cases, multiple model-based controller is used for pH process control. In [7] a multiple adaptive PID controller is designed for this purpose. In [18] the Gap metric is used to multi-modeling of the pH process and then a cascade PI controller is designed. In this paper, a novel method is proposed to design the model set using the nonlinearity measure and gap metric. First, we develop a weighted partitioning method that selects a finite collection of operating point using the nonlinearity measure. Then, for each selected operating point a first order plus dead time (FOPDT) model is constructed. After that we introduce a new method to eliminate the similar models in model bank using clustering algorithm based on the nearest neighborhood and H-gap metric. As a result, a model bank with the smaller model number is provided. Then, we use this model bank to design a controller bank using GPC. Finally, the effectiveness of the proposed algorithm is tested on the pH neutralization process as a difficult nonlinear process in the industry.
The paper is organized as follows. In Section 2, a summary of the mathematical relations of the nonlinearity measure, the H-gap metric and a brief introductory of multiple model control are expressed. In Section 3, the algorithm of model bank selection is presented. Details of each step of the algorithm for model bank selection, controller design and switching supervisor design are studied through a simulation example in Section 4 and then the experimental results of MMGPC on a pH neutralization process are presented in Section 5. Section 6 contains the conclusions.
2. Preliminary 2.1. Nonlinearity measure Nonlinearity measure index (NLI) is used to quantify the nonlinearity degree of a system. One of the NLIs is Bi-spectrum analysis [22]. Based on the normalized higher order cross correlation function, the following indices can be defined [24]
(
)
2 NLI1 = max ϕNyuu (ω1, ω2) ω1, ω 2
2 2 − E⎡⎣ ϕNyuu (ω1, ω2)⎤⎦ + 2σ ⎡⎣ ϕNyuu (ω1, ω2)⎤⎦
(
)
(1)
where E[ . ] and σ[ . ] are expected value and variance, respectively, and ϕNyuu(ω1, ω2) is normalized higher order cross correlation between input u( t ) and output y( t ), which is defined as 2 (ω1, ω2) = ϕNyuu
E⎡⎣ U ( ω1)U ( ω2)Y *( ω1, ω2)⎤⎦ 2 2 2 E⎡⎣ U ( ω1) U ( ω2) ⎤⎦. E⎡⎣ Y ( ω1, ω2) ⎤⎦
(2)
where X ( ωi ) and X ( ωi ) is Fourier transform and power spectrum of x(t) at frequency ωi , respectively, and * stands for conjugate complex. This index has the following properties 1. The value of the index lies in the interval [0,1] where “0” implies linear and “1” refers to sever nonlinear behavior. 2. Is frequency independent in the case of linear time invariant (LTI) systems. 3. Shows the second order high frequency nonlinear behavior. 4. Shows robustness to the level of measurement noise. 5. Requires a high volume of data. 6. It is only applicable for Gaussian inputs. 2.2. H-gap metric There are several ways for measuring the distance between two linear models. In [17] a new metric is presented to measure the similarity of models. By definition, the H-gap metric for two single input single output (SISO) dynamic models P1 and P2 is calculated as follows 1
δhv(P1, P2) =
⎛ 1 ωH ⎞2 ⎜ ∫ Dis(ω)dω⎟ ⎝ 2π ωL ⎠ ⎛ 1 ωH ( 1 + ⎜ ∫ ⎝ 2π ωL
1
) d ⎞⎟ 2
max( Dis(ω)) 2
ω ⎠
(3)
where
⎛ P2 − P1 Dis(ω) = ⎜⎜ ⎝ 1 + P1 2 . 1 + P2 2
⎞2 ⎟⎟ ⎠
(4)
is the distance of models P1 and P2 at every frequency ω .
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
3
After decomposition of the nonlinear region into local linear regions, a FOPDT model is constructed for each region. The FOPDT model is described as
G(z −1) =
Fig. 1. The structure of the multiple model.
2.3. Multiple model method The basic idea of control based on the multiple model is constructing a model bank that its models are determined based on different operating conditions of the system. A general structure of the multiple model controller is presented in Fig. 1. This structure consists of two main loops: outer loop that contains the supervisor and model bank and the inner loop that contains the controller and process. The controller could be any classical types. For each model, an appropriate controller can be designed off-line. Supervisor observes the process behavior continuously and compares it with the behavior of pre-determined related models based on an appropriate cost function.
3. Model bank selection This section introduces the proposed algorithm for model bank selection. The structure of models is selected in the form of FOPDT that is a simple yet popular model in process control. In this algorithm, first a nonlinearity index is employed to decompose the nonlinear operating region of the process to linear local operating points. Then for each of these operating points a FOPDT model is constructed. A process could have similar behavior in different and separated operating points. So it is most probably that some of these linear models are similar. The H-Gap metric is used to eliminate similar models using a clustering algorithm. The remaining linear models constitute the minimum number of models in the model bank. More details are proposed in the following subsections.
b1z −1 z −d 1 + a1z −1
(5)
where b1 and a1 are process model transfer function parameters and z is discrete transform operator. d is delay step of the plant and assumed to be integer constant. The model parameters are estimated using input-output data of the process. The last step is to select the proper minimum number of models. This is achieved by means of the clustering and H-gap metric. Clustering means partitioning of a collection of data into disjoint subsets or clusters so that maximum resemblance of members within each cluster is achieved. In this step, using clustering algorithm based on the nearest neighborhood clustering algorithm [32] and the H-gap metric, local models are divided into different clusters. Then, between these clusters, centers are selected for the model bank and other members are eliminated. 3.2. Model bank selection algorithm In this section, we propose the details of the algorithm for model bank selection. As stated in the previous section, this algorithm measures input-output data and uses the nonlinearity measure and H-gap metric tools to select a model bank for a nonlinear system. The algorithm comprises of 5 main steps as follows: Step 1. Enter the start and end of the operating region as Region=[Rstart , R end], compute the first operating point and it's NLI as (6) and (7) using the Eq. (1).
⎧ (R + R end ) ⎫ ⎬ Op = ⎨ start ⎩ ⎭ 2
(6)
{ NLI( Op1)}
(7)
NLI =
Step 2. Decompose the region Op to multi-linear regions through the sub-algorithm1. Step 3. Construct FOPDT models using least square (LS) for operating points of vector Op and put them into vector G, so named the primary model bank.
G = {g1, g2 , ... , gn}
(8)
3.1. The general structure of the presented algorithm The whole nonlinear operating space of the process is divided into some linear regions by means of the nonlinearity index. The start and end of operating point space are required as the prior knowledge. Then a middle point is selected as the operating point. The process is excited by an amplitude modulated random binary signal (APRBS) input with the amplitude equal to the half of the operating range. Afterward the nonlinearity measure is calculated for this operating point and is compared with a pre-set threshold. If it is smaller than the threshold, algorithm is stopped; otherwise the operating region is broken into smaller areas with respect to the NLI of every operating point. Details of the procedure to select the broken points are stated in Step 2 of section 3.2. Again, the NLI is calculated for these regions. If the calculated value for the index of either is greater than the given threshold, then the region is decomposed into smaller regions. This procedure is continued until all of the NLI values of operating points are less than the threshold value.
Step 4. Using Eq. (3) and measure the distance between each pair of linear models, construct the symmetric pairwise distance matrix Gap as (9):
Gap = Hgapmetric (gi , gj );
i, j = 1, ... , n
(9)
Step 5. Cluster the primary Model bank (8) using sub-algorithm2 and find the reduced final model bank. Details of the sub-algorithm1 of step 2 are as follows: Sub-algorithm1: Decomposition of nonlinear region to multilinear region
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
4
Table 1 Model parameters. Parameter
Notation
Value
Acid flow rate Stream flow rate Acid concentration Base concentration Acid equilibrium constant Base equilibrium constant Height of Solution in stirred tank Cross section of stirred tank
Fa Fo Ca Cb Ka Kb H A
5 ml/s 30 ml/s 0.3 mol 0.3 mol 4.47 * 10–14 5.62 * 10–14 10 cm 90 cm2
Table 2 Applied thresholds for the model set selection (simulation).
Fig. 2. pH process with continuous stirred tank reactor.
Step 2-1. Set i = 0 and j = 0. Step 2-2. Set n as the number of current operating points: n = length(Op). Step 2-3. Set i = i + 1. If i > n go to Step2-4, else go to Step2-5. Step 2-4. If j = 0 go to Step2-11 else put i = 1 , j = 0. Step 2-5. If NLIi > ThNLI go to Step2-6 else go back to Step2-3. Step 2-6. Add two new operating points before and after the current operating point as (10) and (11). Also, compute SNR as (12).
Parameter (Threshold)
Notation
Value
Nonlinearity measure SNR Distance of two adjacent operating points
ThNLI ThSNR ThOp
0.15 0.5 0.1
H-gap metric
ThGap
0.25
NLI = {NLI1, ... , NLIi − 1 ,
before NLInew,
NLIi , NLIi + 1 , ... ,
NLIn}
Opi − 1 . NLIi − 1 − Opi . NLIi = NLIi − 1 + NLIi
before Opnew
after Opnew
SNRnew =
(10)
where ThOp is a threshold for the minimum distance between two consecutive operating points and ThSNR is a threshold for the bound of measurement noise as defined in the Table 2. Step 2-8. If j = j + 1 and
Opi + 1 . NLIi + 1 − Opi . NLIi = NLIi + 1 + NLIi
(11)
2 σNoise
− Opi > ThOp and
Op = {Op1, ... , Opi − 1 , Opi ,
after Opnew ,
SNRnew > ThSNR put
Opi + 1 , ... , (16)
(12) after NLInew
before Opnew
Op = {Op1, ... , Opi − 1 ,
before NLInew
after Opnew
Opn }
2 σSignal
Step 2-7. If j = j + 1 and
(15)
− Opi > ThOp and
before Opnew ,
= NLI( before Opnew )
SNRnew > ThSNR put
Opi , Opi + 1 , ... , Opn }
(13)
(14)
= NLI( after Opnew )
NLI = {NLI1, ... , NLIi − 1 , NLIi ,
(17) after NLInew,
NLIi + 1 , ... ,
NLIn}
(18)
Step 2-9. Due to adding a new operating point(s) and change in the span of other operating points, compute the NLI again for Opm where m = {i − 2 , i − 1 , i , i + 1 , i + 2}. Step 2-10. Go back to Step2-3. Step 2-11. Stop, read Op and return to Step 3. Details of the sub-algorithm 2 of step 5 are as follows:
Fig. 3. The control structure of the pH process.
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
5
Titration curve of pH model
pH
10
8
6
0
0.5
1
1.5 2 2.5 Base feed rate/Acid feed rate Derivative of Titration curve
3
3.5
4
10
11
12
150 X: 7.499 Y: 142.1
Gain
100
50 X: 4.349 Y: 1.196
0
4
5
6
7
8 pH
9
Fig. 4. a (upper): The titration curve of pH process model, b (lower): derivative of the titration curve.
Fig. 5. Some steps of the region decomposition procedure.
Sub-algorithm2:Model clustering and selection Step 5-1. Set k = 1, v = 1, z = 1. Set C as the index of models that are chosen as the center of clusters. Take XC as a matrix that it's Ck th row is the index of models in kth cluster in C. Create a cluster with center G1 ( G1 = g1) and set C1 = { 1 }. Step 5-2. Set k = k + 1, z = z + 1. If z > n go to Step 5-8 else
go to Step 5-3. Step 5-3. Measure the distance between Gz and each cluster center as follows:
d(w ) = Gap( Gz , GCw ); w = 1, 2, ... , v
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
(19)
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
6
Define ind as the index of the cluster that has minimum distance from Gz and min _d as the value of minimum distance. Step 5-5. If min _d > Thgap , take Gz as the center of the new cluster: Ck = { z}. Set v = v + 1. Step 5-6. If min _d < Thgap , put gz in indth cluster: Xind = z . Step 5-7. Go back to Step 5-2. Step 5-8. Stop and read C as the index of the reduced model bank.
4. Simulation results In this section, the strategy described in the previous section is applied to the model of a pH neutralization process. Nonlinear model of this process is presented in [33]. Results are presented in terms of MPC comparison designed by the proposed algorithm and two other controller design algorithms: 1-Multiple PID controller and 2-Multiple model controller presented by [19]. 4.1. pH neutralization process This process consists of a continuous stirred tank reactor (CSTR) with a constant volume as shown in Fig. 2. Acid and base enter the tank with the flow rate Fa and Fb with concentration of Ca , Cb , respectively, and water has the flow of Fw . Mixed product discharges from the tank with the flow rate of Fo. The control objective is to tune the pH value of the outlet stream by manipulating the base flow rate. The value of the acid flow rate is considered as the disturbance. Also, tank level is controlled by manipulating the water flow in another decoupled loop control using a PI controller as illustrated in Fig. 3. The acid flow rate is constant and has been considered to be 5 ml/s . Other parameters of the system are presented in Table 1. One of the important characteristic of the pH process is a static curve, also called the titration curve. The titration curve is obtained from the steady state variation of the output (pH) to the steady state variation of the input (base flow). The titration curve of this process is shown in Fig. 4a. Fig. 4b is the derivative of the titration curve that represents the gain versus pH value. As can be seen, range of gain variation is from around 2 to around 142, which cause difficulty in control of this process.
Table 3 The primary model bank obtained from the proposed algorithm (simulation). #
Op. (pH)
NLI
Gain
1 2 3 4 5 6 7 8 9 10 11 12
5.63 6.25 6.65 6.91 7.37 7.50 7.57 7.86 8.06 8.22 8.36 8.75
0.04 0.04 0.10 0.12 0.09 0.04 0.09 0.05 0.09 0.04 0.08 0.03
4.87 17.4 44.1 75.48 142.25 142.15 135.94 97.38 69.4 50 36.56 15.7
Table 4 Identified models and corresponding steady state gain (simulation). #
Op. (pH)
1
5.63
0.1z−1
6.25
1 − 0.977z−1 0.27z−1
6.65
1 − 0.977z−1 0.806z−1
7.37
1 − 0.975z−1 2.88z−1
2 3 4
1 − 0.974z−1
Gain 4.87
z −5
17.4
z −5
44.1
z −5
142.25
z −5
Table 5 Control and prediction horizon (simulation). #
Op. (pH)
N2
Nu
1 2 3 4
5.63 6.25 6.65 7.37
101 100 92 91
86 85 77 76
Table 6 The optimal PID controller coefficient value. #
Op. (pH)
4.2. Model bank selection Sampling time is selected to be 3 s. In the simulations a noise signal with variance of 0.001 is added at the pH output. The upper and lower limits of the operating region are equal to pH ¼5 and pH ¼10, respectively. Other parameters and thresholds are given in the Table 2. In the first stage of the algorithm, we obtain local linear regions of the overall operating space through the implementation of sub algorithm 1 of Section 3.2. The results outcome of some steps are shown in Fig. 5. Finally, local operating points are obtained as shown in the Table 3. In the next step, FOPDT model is constructed for the operating points presented in the Table 3. Then, models are classified into the clusters with the similar distances by using clustering algorithm and the H-gap metric, as stated in sub algorithm 2 of Section 3.2. Table 4 presents the transfer function of the final model bank. In Fig. 4, the final model bank is marked by red points on the titration curve. As can be seen, models are selected on one side of the curve. It means that the algorithm detected the symmetry on the titration curve.
FOPDT Model
1 2 3 4
5.63 6.25 6.65 7.37
Cof. kp
ki
kd
0.38 0.11 0.03 0.01
4*10-3 1*10-3 5*10-4 1*10-1
3.08 0.86 0.34 0.11
4.3. Controller bank design GPC is the selected controller [34]. In the GPC, the model of the system is assumed to be as follows
A(z −1)y(t ) = B(z −1)z −du(t − 1) + C (z −1)
e(t ) Δ
(20)
where y(t ), u(t ) and e(t ) are the model output, model input and zero-mean white noise, respectively. Also A, B and C are the output dynamic, input dynamic and noise dynamic polynomial, respectively. Δ is the difference operator. Each GPC is designed to optimize the following cost function N2
J (t ) =
∑ j = N1
⎡ y^ (t + d + j ) − w (t + j )⎤2 + λ ⎣ ⎦
Nu− 1
∑ j=1
⎡⎣ Δu(t + j − 1)⎤⎦2
(21)
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
7
Multiple Model GPC in compare with Multiple Model PID, with Small reference changes 10 y
MMGPC
pH
8
MMPID
6 2000
4000
6000
8000
10000
12000
14000
16000
18000
12000
14000
16000
18000
12000
14000
16000
18000
Control Signal
Base feed
3 2 1 0
2000
4000
6000
8000
10000
Switching signal Model Index(pH=)
7.5 7 6.5 6 5.5 2000
4000
6000
8000
10000
time (second)
Fig. 6. a) Output, b) control signal & c) switching signal for MMGPC and MMPID for small change of set-point. Multiple Model GPC in compare with Multiple Model PID, with Medium reference changes 10
pH
8 y MMGPC
6
MMPID
2000
4000
6000
8000
10000
12000
14000
16000
18000
12000
14000
16000
18000
12000
14000
16000
18000
Base feed
Control Signal
2 1 0
2000
4000
6000
8000
10000
Model Index(pH=)
Switching signal 7.5 7 6.5 6 5.5
2000
4000
6000
8000
10000
time (second)
Fig. 7. a) Output, b) control signal & c) switching signal for MMGPC and MMPID for medium change of set-point.
where w (t ) is the desired reference trajectory. λ is the suppressing coefficient of the control signal. N1 is the lower bound of the prediction horizon and is usually set to d + 1. The upper bounds of the prediction and control horizons are denoted by Nu and N2, respectively. Eq. (21) can be written as the vector form as follows
J = (Gu + f−w)T (Gu + f−w) + λ uT u
(22)
J=
1 T u Hu + bT u + f0 2
(24)
where H = 2(GT G + λI ), bT = 2(f−w)T G and f0 = (f−w)T (f−w). With optimization of cost function J with respect to input signal u , control signal is reached
where w is the vector of the feature reference trajectory given as follows
u = − H−1b = (GT G + λI)−1GT (w−f)
T w = ⎡⎣ w (t + d + 1) w (t + d + 2) ⋯ w (t + d + N )⎤⎦
Eq. (25) is the optimized vector of the control signal that is computed at each sample but only the first array of this vector is applied to the system. The GPC bank is designed based on the model bank in the Table 4. The prediction and control horizon parameters of the controllers are given in the Table 5.
(23)
f is the free response vector and G is a matrix extracted from the step response of the plant. Details derivation of these matrices can be found in [34]. Sorting of Eq. (22) in terms of u , we have a quadratic cost function as the form of Eq. (24).
(25)
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
8
Multiple Model GPC in compare with Multiple Model PID, with High reference changes 10
y MMGPC
pH
8
MMPID
6 2000
4000
6000
8000
10000
12000
14000
16000
18000
12000
14000
16000
18000
12000
14000
16000
18000
Control Signal
Base feed
4 3 2 1 0
2000
4000
6000
8000
10000
Model Index(pH=)
Switching signal 7.5 7 6.5 6 5.5 2000
4000
6000
8000
10000
time (second)
Fig. 8. a) Output, b) control signal & c) switching signal for MMGPC and MMPID for high change of set-point.
Titration curve and Model bank based on Du algorithm
pH
10
8
6
0
0.5
1
1.5 2 2.5 Base feed rate/Acid feed rate
3
3.5
4
10
11
12
Derivative of Titration curve 150
Gain
100
50
0
4
5
6
7
8 pH
9
Fig. 9. Final model bank (marked by red points) obtained by Du's algorithm.
4.4. Switching supervisor As it was mentioned in Section 2.3, the multiple model controller at each sampling time requires a switching supervisor to detect the operating region of the process and then select the appropriate controller in the controller loop. Switching supervisor uses a cost function as follows in order to select the suitable model [35]
cost function is selected as the best model. In order to reduce the chattering in the switching process, we use switching with hysteresis logic [36]. In this method, switching is not performed unless the cost-function of the ith model ( Ji (t )) is some percent less than the current selected model. 4.5. Comparison study
(26)
Simulation results are presented in two parts as follows in order to compare the proposed controller and some other multiple model controllers.
where a and ρ are weighting parameters and ei( . ) is the difference between the process output and the ith linear model. The first term is instance squared error and the second is the sum of weighted squared errors. The cost function (26) is calculated for each model at every sampling time. The model that has the lowest
4.5.1. Comparison with the multiple model PID controller(MMPID) In this section, the proposed MMGPC is compared with the MMPID controller presented in [37]. In this design, three constraints are considered; phase margin on cutoff frequency, restriction on gain crossover frequency and derivative of phase
t
Ji (t ) = a. ei2(t ) +
∑ ρt − k ei2(k) k=0
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
diagram at the gain crossover frequency. Parameters of the PID controllers are shown in Table 6. Simulation results are given in Figs. 6–8 for different step change. As can be seen, in all of three cases, there are some operating points that MMPID controller is unstable while MMGPC has more appropriate response. 4.5.2. Comparison with the model bank selection of Du [19] In this comparison, the proposed MMGPC is compared with another MMMPC that its model bank is designed according to the algorithm presented by [19]. According to [19], the operating region is uniformly divided into a large number (50 or more) of operating points. Then similar models are removed based on the gap metric criterion. In Fig. 9, the final operating points obtained by this algorithm are shown, when the whole operating space was initially modeled by 50 models. Note that the threshold of gap
9
metric is selected so that the number of models obtained from Du's algorithm is the same as the proposed algorithm. Simulation results are given in Figs. 10–12 for different step change. As can be seen, due to improper selection of initial operating points in the model bank, Du's algorithm has not desirable tracking performances. This illustrates the effectiveness of the nonlinearity index to appropriate selection of operating points. In fact, NLI introduces the appropriate set of linear behavior subspaces, to the clustering algorithm to select the model bank.
5. Implementation results The experimental apparatus, shown in Fig. 13, is a continuous pH neutralization pilot plant in the Process Control laboratory at K. N. Toosi University of Technology. A schematic diagram of the
Multiple Model GPC with two different algorithm, with Small reference changes
pH
10 Prop. alg.
8
Du alg.
6 2000
4000
6000
8000
10000
12000
14000
16000
18000
12000
14000
16000
18000
12000
14000
16000
18000
Control Signal
Base feed
3 2 1 0
2000
4000
6000
8000
10000
Model Index(pH=)
Switching signal 9 8 7 6 2000
4000
6000
8000 10000 time (second)
Fig. 10. a) Output, b) control signal & c) switching signal of MMGPC obtained by proposed algorithm and Du's algorithm for small change of set-point.
Fig. 11. a) Output, b) control signal & c) switching signal of MMGPC obtained by proposed algorithm and Du's algorithm for medium change of set-point.
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
10
Multiple Model GPC with two different algorithm, with High reference changes 10
pH
8 Prop. alg.
6
Du alg.
2000
4000
6000
8000
10000
12000
14000
16000
18000
12000
14000
16000
18000
12000
14000
16000
18000
Control Signal
Base feed
2 1.5 1 0.5 0
2000
4000
6000
8000
10000
Model Index(pH=)
Switching signal 9 8 7 6 2000
4000
6000
8000 10000 time (second)
Fig. 12. a) Output, b) control signal & c) switching signal of MMGPC obtained by proposed algorithm and Du's algorithm for high change of set-point.
Fig. 14. The schematic diagram of the pH pilot plant [38].
Table 7 Tunable parameters of pH process pilot plant. Fig. 13. pH process pilot plant.
system is shown in Fig. 14. It contains a continuous stirred tank reactor (CSTR) with a motorized mixer in order to have a wellmixed tank. The effluent stream flow rate is regulated by a manual valve which affects transportation delay and time constant of the process. A pH sensor is located after the outlet valve to measure the pH of the effluent stream. Acid, base, and water are pumped into the CSTR by the corresponding dosing pumps. The control objective is to achieve the desirable pH values by manipulating the base flow rate in the presence of changes in acid stream flow rate as a disturbance to evaluate ability of different control strategies. The advanced pH controller and level controller are performed by
Parameter
Notation
Value
Acid flow rate
fa
0.45 ml/s
Acid concentration pH of acid solution Base concentration pH of base solution Stream flow rate
Ca – Cb – fo
1.5 mlit./lit. 3 1 gr./lit. 12 0.45 mlit/s
Height of Solution in stirred tank Sample time
H Ts
15 cm 3s
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎ Table 8 Applied thresholds for the model set selection (implementation).
11
Table 10 The final model bank in experimental study.
Parameter (Threshold)
Notation
Value
#
Op. (pH)
Time Constant (s)
Gain
Delay (s)
Nonlinearity measure SNR Distance of two adjacent operating points
ThNLI ThSNR ThOp
0.35 0.5 0.5
H-gap metric
ThGap
0.25
1 2 3 4 5
4.32 4.94 5.45 6.41 7.55
462 453 384 462 390
0.16 0.29 0.42 0.73 2.52
66 87 48 42 36
the MATLAB™ software and SIMULINK™; in addition, a PLC is installed for data acquisition which communicates to the computer through a PROFIBUS cable. Some of the tunable parameters of the pH pilot plant are shown in Table 7. This apparatus is severely nonlinear with ample source of disturbances and uncertainty. The static gain changes by a ratio of more than 10 even for a very simplified experiment. One of the main sources of disturbances is the effect of the minerals in the water tab, which is used instead of buffer in our experience. In addition, acid and base tanks are filled using some materials that prepared in the lab manually. So, the concentration of materials changes every time that they are prepared. Adding to that neutralization of acid in the source tank in connection with air, change in the temperature of the room and change in noise variance makes this apparatus one of the most challenging problem in process control study. 5.1. Model bank selection In this section, an appropriate model bank is designed using the proposed algorithm. Threshold parameters are assumed as presented in Table 8. NLI1 needs a considerable amount of data. Due to the inherent time invariance of the process, for instance, when the base and acid tanks are filled by the new source of materials, it is not possible to acquire enough experimental data for the nonlinearity index from this process. So, in the experimental study, we use the cross correlation index to compute the nonlinearity measure. This
Table 11. Control and prediction horizon. #
Op. (pH)
N2
Nu
1 2 3 4 5
4.32 4.94 5.45 6.41 7.55
44 34 27 25 22
15 12 11 11 10
Table 12 Weighted coefficient of control signal ( λ ). # Op. (pH)
1 2 3 4 5
4.32 4.94 5.45 6.41 7.55
λ Small change set-point
Medium change set-point
High change set-point
0.05 0.17 .35 2.66 6.35
0.38 1.26 3.35 5.33 76.2
0.38 1.26 3.35 5.33 50.8
index is defined as follows [22] 2 ⎛ ⎞ ϕyu(jω) ⎜ ⎟ NLI2 = 1 − mean⎜ ω ϕuu(jω)ϕyy(jω) ⎟ ⎝ ⎠
(27)
Table 9 The initial linear local regions. pH NLI2
4.11 0.16
4.32 0.15
4.94 0.11
5.14 0.15
5.45 0.18
5.8 0.14
6.33 0.11
6.41 0.16
6.58 0.16
7.1 0.23
7.44 0.25
7.55 0.33
7.9 0.12
Fig. 15. a) Output, b) control signal & c) switching signal of MMGPC for the small change of set-point.
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
12
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
Fig. 16. a) Output, b) control signal & c) switching signal of MMGPC for the medium change of set-point.
Fig. 17. a) Output, b) control signal & c) switching signal of MMGPC for the high change of set-point.
where ϕyu( jω) = Y ( jω)U *( jω) is the cross-correlation function. Using the proposed algorithm in Section 3 and consideration of threshold parameters as presented in Table 8, we obtained 13 initial linear regions. The operating points of these regions and their NLI indices are shown in Table 9. Finally, after applying the clustering algorithm, reduced model bank is obtained as presented in Table 10.
Table 13 Time domain properties of response of the implemented MMGPC on pH pilot plant for high change of set-point presented in Figs. 15-17. Setpoint Situation
MSE
Av. Settling time (s)
Av. Overshoot (%)
small change medium change high change
0.014 0.064 0.364
1007.1 3691 3582
26.1 123.9 78.38
5.2. Controller bank design and switching supervisor For each model of the model bank, a GPC is designed. Parameters of these controllers are shown in Table 11. Weighted coefficients of control signal λ with respect to the change in setpoint value (low, medium and high change) have been achieved experimentally as Table 12. Parameters of cost function (26) are chosen to be α = 0.1 and η = 0.9. Three modes of low, medium and high change in the reference signal are studied experimentally. Hysteresis constant (h) for these three modes are 0.7, 0.05 and 0.05, respectively. Results are shown in Figs. 15–17. Time domain properties of set points, are shown in Table 13. As can be seen, the performance of
the control system is satisfactory for small and medium step change except when there is a step from a medium steady state gain operating point to a high gain operating point. This is predictable, since even a step change in pH value from around 6.5 to 7.5 means 3 times changes in the steady state gain. The high value of over-shoot is avoidable by reducing the speed of the desired closed loop system. Fig. 17 shows that the control system provides a stable closed loop system even in large step test. However, a step change from pH of around 5 to 7.5 means a steady state gain change by a factor of around 10, which is a very large value of uncertainty for every control system. Again this large overshoot is
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i
D. Shaghaghi et al. / ISA Transactions ∎ (∎∎∎∎) ∎∎∎–∎∎∎
avoidable sacrificing the speed of the closed loop system.
6. Conclusion In the context of the controller design based on the multiple model technique, the model bank selection problem is one of the main issues that is still under study for a general yet reliable solution. This problem is even more challenging in the case that there is not any first principle model of the processes. In this paper, an algorithm is proposed that can provide an appropriate model bank using the nonlinearity measure and the H-gap metric indices, using only the input/output data of the process. The operating space of the process is divided into some small linear behavior sub-spaces using the proposed nonlinearity measure. This may leads to some excessive local models. So, the next step a clustering algorithm based on H-gap metric is employed to prune the model bank. The simulation and implementation results show the functionality of the algorithm. The proposed method gives a systematic solution for the modeling and control of a large class of nonlinear stable systems with wide operating ranges. For the highly nonlinear pH neutralization process with plenty of source of disturbances and uncertainty, this method provides a stable experimental control system in all operating points at step changes with proper performance.
References [1] Camacho EF, Bordons C. Model predictive control. 2nd ed. London: SpringerVerlag; 2004. [2] Murray-Smith R, Johansen TA. Multiple model approaches to modeling and control. London: Taylor & Francis; 1997. [3] Scokaert PO, Mayne DQ, Rawlings JB. Suboptimal model predictive control (feasibility implies stability). IEEE Trans Autom Control 1999;44:648–54. [4] Johansen TA, Foss BA. Operating regime based process modeling and identification. Comput Chem Eng 1997;21:159–76. [5] Rodriguez J, Romagnoli J, Goodwin G. Supervisory multiple regime control. J Process Control 2003;13:177–91. [6] Li XR, Jilkov VP. Survey of maneuvering target tracking. Part V. Mult-model methods Aerosp Electron Syst, IEEE Trans 2005;41:1255–321. [7] Böling JM, Seborg DE, Hespanha JP. Multi-model adaptive control of a simulated pH neutralization process. Control Eng Pract 2007;15:663–72. [8] Aufderheide B, Bequette BW. A variably tuned multiple model predictive controller based on minimal process knowledge. In: Proceedings of the 2001 American control conference. IEEE 2001. p. 3490–3495; 2001. [9] Peymani E, Fatehi A, Sedigh AK. Automatic learning in multiple model adaptive control. International Conference on Control, UKACC, Manchester, UK; 2008. [10] Liang Y, Li X, Han C, Duan Z. A general systematic method for model-set design. Aerosp Electron Syst, IEEE Trans 2012;48:1505–20. [11] Novák J, Chalupa P. , Bobál V. Modelling and predictive control of a nonlinear system using local model network. In: Proceedings of the World Congress. p.
13
1622–7; 2011. [12] Fatehi A, Abe K. Flexible structure multiple modeling using irregular self-organizing maps neural network. Int J Neural Syst 2008;18:233–56. [13] El-Sakkary A. The gap metric: robustness of stabilization of feedback systems. Autom Control, IEEE Trans 1985;30:240–7. [14] Zames G, El-Sakkary A. Unstable systems and feedback: the gap metric. Proc Allerton Conf 1980:380–5. [15] Vinnicombe G. Frequency domain uncertainty and the graph topology. Autom Control, IEEE Trans 1993;38:1371–83. [16] Hosseini S, Fatehi A, Johansen TA, Sedigh AK. Multiple model bank selection based on nonlinearity measure and H-gap metric. Journal of Process Control; 2012. [17] Hosseini S, Fatehi A, Khakisedigh A A new metric for measuring the distance between two linear systems. Control, Automation and Systems (ICCAS). In: Proceedings of 11th international conference on. IEEE; 2011. p. 1552–1557. [18] Galán O, Romagnoli JA, Palazoglu A, Arkun Y. Gap metric concept and implications for multilinear model-based controller design. Ind Eng Chem Res 2003;42:2189–97. [19] Du J, Song C, Li P. Application of gap metric to model bank determination in multilinear model approach. J Process Control 2009;19:231–40. [20] Arslan E, Çamurdan MC, Palazoglu A, Arkun Y. Multimodel scheduling control of nonlinear systems using gap metric. Ind Eng Chem Res 2004;43:8275–83. [21] Desoer C, Wang Y-T. Foundations of feedback theory for nonlinear dynamical systems. Circuits Syst, IEEE Trans 1980;27:104–23. [22] Hosseini M, Johansen T, Fatehi A. Comparison of nonlinearity measures based on time series analysis for nonlinearity detection. Model Identif Control 2011;32:123–40. [23] Le Caillec J-M, Garello R. Comparison of statistical indices using third order statistics for nonlinearity detection. Signal Process 2004;84:499–525. [24] MAAS Choudhury, Shah SL, Thornhill NF. Diagnosis of Process Nonlinearities and Valve Stiction: data Driven Approaches. Springer; 2008. [25] Helbig A, Marquardt W, Allgöwer F. Nonlinearity measures: definition, computation and applications. J Process Control 2000;10:113–23. [26] Qin SJ, Badgwell TA. A survey of industrial model predictive control technology. Control Eng Pract 2003;11:733–64. [27] Clarke DW, Mohtadi C, Tuffs P. Generalized predictive control – Part I. The basic algorithm. Automatica. 1987;23:137–48. [28] Novák J, Chalupa P, Bobál V. Multiple model modeling and predictive control of the pH neutralization process. Int J Math Model Methods Appl Sci 2011. [29] Bagheri P, Mardanlou V, Fatehi A. Multiple model predictive control of multivariable pH process using adaptive weighting matrices.In: Proceedings of The 18th IFAC World congress2011; 2011. p. 12366–71. [30] Hall RC, Seborg DE. Modelling and self-tuning control of a multivariable pH neutralization process part I: modelling and multiloop control. In: Proceedings of the American control conference. IEEE1989; 1989. p. 1822–1827. [31] Ylén J-P. Measuring, modelling and controlling the pH value and the dynamic chemical state. Helsinki University of Technology; 2001. [32] Theodoridis S, koutroumbas L. Pattern Recognit 1999. [33] Sadeghassadi M, Fatehi A, Khaki Sedigh A, Hosseini S. On the structural optimization of a neural network model predictive controller. Ind Eng Chem Res 2012;52:394–407. [34] Camacho EF, Bordons C. Model Predictive Control.London: Springer; 2004. [35] Narendra KS, Balakrishnan J. Adaptive control using multiple models. Autom Control, IEEE Trans 1997;42:171–87. [36] Liberzon D. Switching In Systems and Control. Springer; 2003. [37] Li H, Luo Y, Chen Y. A fractional order proportional and derivative (FOPD) motion controller: tuning rule and experiments. IEEE Trans Control Syst Technol 2010;18:516–20. [38] Peymani E, Fatehi A, Bashivan P, Sedigh AK. An experimental comparison of adaptive controllers on a pH neutralization pilot plant. In Proceedings of the India conference, 2008 INDICON Annual IEEE. IEEE2008; 2008. p. 377-382.
Please cite this article as: Shaghaghi D, et al. Multi-linear model set design based on the nonlinearity measure and H-gap metric. ISA Transactions (2017), http://dx.doi.org/10.1016/j.isatra.2017.01.021i