Copyright @ IF AC Advanced Control of Chemical Processes, Pisa, Italy, 2000
ISSUES IN DESIGN OF INPUT SIGNALS FOR THE IDENTIFICATION OF NONLINEAR MODELS OF PROCESS SYSTEMS Raghunathan Rengasamy·· Robert S. Parker· Francis J. Doyle III . ,1
• Department of Chemical Engineering, University of Delaware, U.S.A •• Department of Chemical Engineering, Indian Institute of Technology - Bombay, India
Abstract: This paper addresses the design of input signals for the identification of Volterra models for nonlinear process systems. In this work, we review and compare deterministic and stochastic signals for the identification of Volterra models. We also discuss the various issues involved in the design of input signals. Accurate identification of Volterra kernels through the suppression of higher order nonlinearities is also discussed . Copyright © 2000 IFAC Keywords: empirical model-based control, input sequence design, nonlinear identification, Volterra series
1. INTRODUCTION
In et al., 1999) , is the Volterra functional series model.
Traditionally, model-based control schemes have been implemented using discrete-time linear convolution models. One of the very popular control schemes which employs linear models is the Model Predictive Controller (MPC) . However, since most physical systems are nonlinear to some extent, these models often perform poorly in practice. Thus, to obtain improved control of a nonlinear process, an alternative model structure that incorporates nonlinearities is needed. There are a number of model structures that include nonlinearities, however, the artificial neural network appears to be gaining acceptance with the vendors (Qin and Badgwell, 1999). Unfortunately, there is little formal theory for input sequence design for such models. Another popular model class, which has received considerable attention in process identification as well as nonlinear control (Doyle et al., 1995; Zheng and Zafiriou, 1995; Genceli and Nikolaou, 1995; Ling and Rivera, 1998; Doyle 1
An Nth-order discrete-time Volterra model of the form :
is capable of capturing significant nonlinearity. Particularly important is the ability of these models to exhibit asymmetric responses to symmetric inputs, a phenomenon often observed in chemical processes. The parameters ha and hn(i}, ... , in) are referred to as the Volterra kernels. The identification problem is one of determining the values of these coefficients. There is a large body of literature that addresses the problem of identification of Volterra kernels . Provided the input sequence is white, a crosscorrelation technique may be used to determine the Volterra kernels. This method originates from the work of Wiener (1958). Using a Gram-Schmidt procedure, he converted the Volterra series into a sum of functionals that were mutually orthogonal
Corresponding author:
[email protected]
839
with respect to a Gaussian white noise input. The resulting series is called a Wiener series and its coefficients are referred to as the Wiener kernels. Lee and Schetzen (1965) developed a simple crosscorrelation formula for the computation of these terms. Once all relevant Wiener kernels have been measured, they can be rearranged into a Volterra model that gives the least-squares fit for the particular Gaussian input employed. From a chemical engineering perspective, there are two major problems with this approach. First, Gaussian white noise is continuously distributed. As such, at any given time instant k, the input is going to switch to a new value. This constant movement would result in significant actuator wear. Second, there is considerable estimation error in the kernels resulting from the use of finite length input/output data records. This is a consequence of the incomplete formation of white autocorrelation properties of the input stimulus and the fact that the cross-correlation estimates converge at a rate proportional to the square-root of the record length (Marmarelis, 1993) . Thus, long data records are required for reasonably accurate parameter estimates. Since the collection of this data is done off-line and may entail production of off-specification product, a lengthy experiment may have serious economic repercussions. Multi-level (m-level) sequences have also been used in the identification of Volterra kernels. A Pseudo Random Binary Sequence (PRBS) is one example of a maximum length binary sequence that has been used in the identification of the linear portion of the Volterra kernel. The advantages with these signals are: they are easily implementable; deterministic; and need much shorter experimentation time than the random signals. The downside though, is that there might still be considerable valve movement during experimentation. There has been a significant amount of research in the theory and applications of PRBS signals. PRBS signals were tried on diverse range of industrial applications. A detailed review of the theory and applications is provided in (Godfrey, 1969a; Godfrey, 1969b) . Time domain techniques for identifying the second-order kernel have been discussed in (Barker et al., 1972; Barker and Davy, 1978) . In (Barker et al., 1972) they consider three and five level sequences for identification of the Volterra kernels, whereas in (Barker and Davy, 1978) a detailed analysis of the three level sequence is provided. More work can be found in (Parker and Moore, 1980; Parker and Moore, 1982; Baheti et al., 1980) .
the m-level sequences. We also address two other important questions in the design of input signals. First, we address the multi-objective nature of the problem of identification of fixed structure Volterra models (for example, a second- or thirdorder model). Various design issues are discussed in detail. Second, the problem of identification of accurate Volterra kernels is discussed. We discuss the suppression of the effect of higher order kernels in the identification process. As a specific example, we show how even-order and odd-order kernels could be cancelled to suppress the effect of higher order nonlinearities in the identification process. 2. ISSUES IN INPUT SIGNAL DESIGN FOR IDENTIFICATION In the design of input signals for identifying Volterra kernels , the most important requirement is one of persistence of excitation. The Volterra kernel identification problem can be posed as least squares estimation problem for a suitable 'X' matrix in the following equation, where /3 represents the Volterra coefficients:
y=X/3+t As Volterra models are linear in the parameters, persistence of excitation would require that the XT X matrix be non-singular and positive definite. It can be shown that for estimating nth-order kernels a signal of at least n +1 levels is needed (Nowak and Veen, 1994). There are many kinds of signals that could be used for identification of Volterra kernels. Two of these, as mentioned earlier, are the CSRS and m-level sequences. It can be easily shown that the persistence of excitation criteria is satisfied by both these signals. Other than the rank criteria for persistence of excitation, other measures to evaluate the input signals have also been discussed in the literature. Plant friendliness, characterized by a friendliness index, has been discussed in (Doyle III et al., 1999). IT the total length of the signal is N and if there are nT transitions then the friendliness factor is defined as
f
= 100 x (1 -
Nn~ 1)
There are other factors that one has to consider in the design of input signals. The length of the input signal (duration of experimentation) is an important factor that should be minimized. Input signals could be evaluated on this basis. Other objectives are the amplitude of the input signal needed for sufficient excitation of the process and the variability in the process variable during identification. Clearly there is a trade-off between the input variability in terms of the friendliness index and output variability. Since most of the chemical processes act as high frequency filters,
In this work, we review and compare a broad class of white inputs called Constant-switchingpace Symmetric Random Signals (CSRS's) and
840
the process variability is minimized when inputs of higher frequencies are used. As long as the persistence of excitation condition is satisfied, these high frequency signals having low friendliness index in terms of input, and high friendliness index in terms of output, could be used for identifying Volterra models. Hence, the input signal design problem is clearly a multi-objective problem with competing objectives. Another important goal of the identification signal would be reduction of the effect of noise and outliers on the estimates. In the least squares estimation, for deterministic signals, it can be shown that var({3) = (XT X)-lvar(e) . One of the objectives of the identification process should be the minimization of the variance of {3 and it can be seen that it would be directly related to the signal designed. We are in the process of developing an optimal framework for handling these competing objectives in the design of input signal and the results will be communicated in future publications. In this work, though, we contrast two signals for identifying Volterra kernels in terms of some of the objectives mentioned above. 3. IDENTIFICATION USING M-LEVEL AND CSRS SIGNALS
= ho + L(k) + D(k) + O(k)
A 5-level sequence with a shift register length of 5 was sufficiently exciting to identify the 230 parameters. The length of the designed sequence was 624. An irreducible primitive polynomial of the form X = D3 X + D4 X (implies feedback from registers 3 and 4 with a multiplicative factor of one) was used in the design of the 5-level sequence. Following (Godfrey, 1993), the sequence value logic was converted to signal level by mapping [0,1,2,3,4] to [0,0,20,-20,-0], where 0 is the amplitude of the signal, and this generates a signal which has zero mean with the inverse repeat property. The levels 0, 20, -20, -0 occur 5n - 1 times and zero occurs 5n - 1 -1 times in the signal.
(1)
where:
= ~t!lhl(i)u(k - i) = ~t!l h2(i , i)u 2(k - i) O(k) = 2~t!I~~:~h2(i , j)u(k -
Many of the m-level signals used in literature are usually based on maximum length sequences of multiple levels. There are other kinds of binary signals also, examples being, Hall (Hall Jr. , 1956) and Baker (Donce, 1966) sequences. In this work, we show the utility of a 5-level maximum length sequence for the identification of Volterra models. m-level signals are developed using a combination of shift registers and modulo addition. The theory behind m-level sequence generation is well known and more information on this can be found in (Godfrey, 1993). The maximum length of the repeating signal that can be generated is qn_1 where q is the number of levels and n is the length of the shift register. The elements of the sequence can take values of the integers 0,1 , .. ., q - 1 and the elements can be generated by a shift register with the feedback to the first stage comprising of the modulo q sum of the outputs of the other stages multiplied by coefficients Cl , C2, . .. , en which are also integers 0, 1, ... , q - 1 (Godfrey, 1993). Primitive irreducible polynomial feedback generates sequences that repeat after a maximum length of elements. Such polynomials for generating different m-level signals can be found in (Godfrey, 1993) . There are a number of m-sequences that can be generated through the choice of different polynomial functions. For example, when one tries to identify a particular kind of nonlinearity, some signals might give the "best" results. An example is the preference for certain binary msequences for identification of linear systems with an additive quadratic drift term at the output (Godfrey, 1993). In this work, we generated a 5-level sequence for the identification of secondorder Volterra models. The number of parameters to be identified based on a memory of 20 is 230 (20 linear and 210 quadratic). A number of 3-level sequences with a shift register length of up to 6 (length of sequence = 728) were tried and it was not possible to find a sequence with persistence of excitation (rank of XT X matrix equals 230) for estimating all the required parameters.
In this section we show the design of two different kinds of signals for the identification of the Volterra kernels. The starting point for the second-order Volterra model that serves as a real plant is a first-principles description of the isothermal free-radical polymerization of methyl methacrylate using azo-bis-isobutryonitrile as initiator and toluene as solvent (Congalidis et al., 1989; Doyle et al., 1995) . In our SISO case study, the number average molecular weight is the output variable and the initiator flow rate is the input variable. For more details on the model, operating point and normalization, the interested reader is referred to (Doyle et al., 1995). Assuming a memory M = 20, a second-order Volterra model acting as the plant can be derived through Carlemann linearization. The following model decomposition will act as a reference for further discussion of the input signals (Doyle III et al., 1999) : y(k)
3.1 Identification Using m-level Signal
L(k)
D(k)
i)u(k - j)
represent the linear, diagonal, and off-diagonal model components, respectively. Here, we assume h2 (i, j) = h2 (j, i) without loss of generality. 841
This property creates a low input friendliness index which will be discussed later in this paper. This signal was generated and applied to the case study discussed above. All the parameters (linear and diagonal) were identified together using least squares identification. The identified kernels were then compared with the true kernels and the sum of squared errors is reported in Table 1. It is clearly seen that this persistently exciting signal was able to identify the second-order Volterra kernels to a high degree of accuracy. Table 1. Sum of squares error in identification using m-level sequences. Term Linear Diagonal Off-diagonal
Coefficient Error 2.113 x 10 ., 4.156 x 10- 7 9.715 x 10- 7
amplitude. The second sub-sequence consisted of a 789 point binary sequence (levels equal to ±56 in scaled variables) . Based on the prediction error variance expression, the theoretical kurtosis of -2 guaranteed the diagonal coefficients would not be excited by this signal (see (Doyle III et al., 1999)). Pre-whitening was used to generate a residual sequence from which the off-diagonal coefficients were calculated using cross-correlation. To guarantee that the results were not corrupted by a bias based on the specific characteristics of a particular random binary sequence, the identification procedure was performed 100 times, and the results were averaged . It should be noted that the deterministic nature of the first sub-sequence leads to identical linear and diagonal coefficients for noise-free simulations. As above, the identification error was examined term-by-term assuming a noise-free process output signal. The results are shown in Table 2.
3.2 Identification using CSRS signal For comparison purposes, Algorithm 2 from (Doyle III et al., 1999) was used to identify a second-order Volterra model from the Carlemann linearization of the polymerization reactor. Algorithm 2 was composed of 2 sub-sequences: the first was a 211 point sequence containing a deterministic binary pulse replicated 5 times (see Figure 1) . The bias, 60r----T----~----~--~----~--~~
Table 2. Sum of squares error in identification using CSRS sequences (Algorithm 2).
Term Linear Diagonal Off-Diagonal
Coefficient Error 4.11 x lO- b 1.47 x 10- 7 8.61 x 10- 7
4. DISCUSSION
40
It is clearly seen that both the signals were able to provide good identification of Volterra kernels. Deterministic signals, like the m-level signals, have the advantage that a block of experiments can be repeated and identification results verified, whereas stochastic CSRS signals cannot be repeated. Clearly, the friendliness factor of the mlevel sequence is not very high as the number of times each level is attained is equally distributed (property of a zero mean input signal) and the signal level switches quite often in the switching time period. The friendliness factor for the m-level signal that was used to identify the second-order Volterra kernel is 20%. The friendliness index of the CSRS like signal that is used to identify the linear and diagonal terms is 91%, while the theoretical friendliness index of the CSRS signal that is used to identify the non diagonal terms is 50%.
20
-20
-40 -60~--~--------~----~--~----~ o 2 4 6 8 10 12 Time (hours)
Fig. 1. Binary pulse CSRS used for linear and diagonal identification in Algorithm 2. linear, and diagonal second-order coefficients were calculated without pre-whitening using the following equations: ho = hdk)
The predicted variance due to measurement noise on the estimated parameters can be calculated explicitly. If the noise variance is var(e) , then the variance in the parameter estimation in the CSRS case is of the order of 1O-5var(e) for the linear terms and is approximately of the order of 1O-8var(e) for the quadratic terms . In the 5-level sequence case, for the linear terms, the variance is of the order of lO-Svar(e) and for the quadratic
---.L " Nr - I {y(2j[M+l])+y([2j+l][M+l])}, 'IN; ~J-O
= ut;:x L;.:; I{!,t(2j[M+l]+k) -y([2j+l][M+l]+k)},
h2(k,k)
= ~2N "N -l{y(2j[M+l]+k)-y(2j[M+l]) "x L.Jp_o r
+y([2j+l][M+l]+k)-y([2j+l][M+l])} .
Here N r is the number of repetitions of the binary pulse sequence and X is the scaled input sequence 842
terms it is of the order of lO-lOvar(e) . Hence, the 5-level sequence is a more robust signal for identification of Volterra kernels in the presence of random noise. Input and output variances are also important factors in the identification experiment. Because of the high friendliness factor, the input variance of the CSRS signal is quite low, whereas, the m-level sequence has a high input variance in comparison. In terms of the output variance, if the switching time for the m-sequence is chosen to be high, the output variance in the identification experiment can be decreased considerably. Another interesting aspect of the m-level signal designed is the suppression of the effect of the higher order nonlinearities in the identification experiment. These issues become important when one tries to identify reduced order Volterra kernels for higher order nonlinear plants. It can be shown that signals that are of the inverse repeat form can be used to selectively suppress the odd and even order nonlinearities. An example of a Binary Inverse Repeat Sequence (IRS) can be constructed as follows. The IRS signal is generated by doubling the PRBS sequence and toggling every other digit of the doubled sequence. This is illustrated in Table 3. Over a common period of 2N~t (where N is the length of the signal and ~t is the switching time), the auto-correlation function (ACF), 4>uu(7) is given by 4>uu(7) = 0 2 7=0 =_0 2 7=N 02 7 = 2,4,6 ... N 02 7=1,3,5 ... = N
(2)
as shown in Figure 2. As the IRS is a repeating signal with period T (T = 2N~t) and antisymmetric, that is (as can be seen from Table 3):
u(t)
= -u(t+Tj2)
(3)
Fig. 2. ACF of inverse repeat sequence
=~ =
1 '1'
+t
loT u(t)u(t +
J. JT
T 2 0 /
T/2
7 -
Ar}U(t + 7
-
A2)dt
u(t)u(t+r-A,)u(t+r-A2)dt u(t)u(Hr-A,)u(Hr-A2)dt
(5)
The two components on the right hand side of equation (5) are exactly equal and opposite. This could be easily seen from the antisymmetric property of the IRS - refer to equation (3). For the oddorder ACF of u(t), the two components become exactly equal, but positive. Because of this, the even-order kernels vanish and odd-order kernels remain. It can be shown that the 5-level sequence that was used in identification of Volterra kernels is also of the inverse repeat form by construction. Hence, in the estimation of the second-order kernel, the effects of the odd-order nonlinearities of order higher than 1 are minimized. Similarly, if one were to estimate the third-order Volterra kernel using a least squares algorithm, the effect of even-order nonlinearities of order higher than 2 will be minimized. Hence this would lead to better estimation of the Volterra kernels. Further resolution and suppression of nonlinearities other than even- and odd-order are extremely difficult and needs further investigation. These issues have to be investigated from the perspective of short record length.
In this paper, the identification of Volterra kernels using stochastic and deterministic input signals was presented. A number of issues in the design of input signal for identification were discussed and the design signals used in this work were on the basis of important identification issues. Further, better estimation of Volterra kernels through the suppression of higher order nonlinearities was also discussed.
= LimK .... oc-k foK u(t)u(t+r-A,)x u(t+r-A2) ... u(Hr-A,,)dt
·1
5. CONCLUSIONS
we have all the even order kernels (2,4,6, ... ) of the ACF ofu(t) going to zero. This could be explained by considering the 2nd order kernel. The pth order ACF of the input signal u is defined as given below:
1111...
.liN ·
(4)
6. ACKNOWLEDGMENTS
From equation (4) and the fact that IRS signal is repeating, the 2nd order ACF kernel could be written as:
This research was supported by the VD Process Control and Monitoring Consortium (VD PCMC) . 843
In: Proc. European Control Conf. . Karlsruhe, Germany. Genceli, H. and M. Nikolaou (1995) . Design ofrobust constrained model-predictive controllers with Volterra series. AIChE J. 41(9) , 20982107. Godfrey, K. (1993) . Introduction to perturbation signals for time-domain system identification. In: Perturbation signals for system identification (K. Godfrey, Ed.). pp. 1- 59. Prentice Hall International (UK) Limited. Godfrey, K. R (1969a). The theory of the correlation method of dynamic analysis and its application to industrial processes and nuclear power plant. Meas. and Control 2, T65-72. Godfrey, K. R (1969b). Dynamic analysis of an oil refinery under normal operating conditions. Proceedings lEE 116, 879--888. Hall Jr., M. (1956). A survey of difference sets. Proc. Am. Math. Society 1 , 975-. Lee, Y.W. and M. Schetzen (1965) . Measurement of the Wiener Kernels of a Non-linear System by Cross-correlation. Int. J. Control 2, 237254. Ling, W.-M. and D.E. llivera (1998). Control relevant model reduction of volterra series models. J. Proc. Cont. 8, 79--88. Marmarelis, V. Z. (1993). Identification of Nonlinear Biological Systems Using Laguerre Expansions of Kernels. Ann. Biomed. 21 , 573589. Nowak, R. D. and B. D. Van Veen (1994) . Random and pseudorandom inputs for Volterra filter identification. IEEE Trans . Signal Processing 42(8), 2124- 2135. Parker, G. A. and E . L. Moore (1980). A modified volterra series representation for a class of single-valued, continuous nonlinear systems. J. Dynam. Syst. Meas. Contr. 102, 163-167. Parker, G. A. and E. L. Moore (1982). Practical nonlinear system identification using a modified volterra series approach. Automatica 18, 85--91. Qin, S. J. and T . A. Badgwell (1999) . An overview of nonlinear MPC applications. In: Nonlinear Model Predictive Control: Assessment and Future Directions (F. Allgower and A. Zheng, Eds.). Birkhauser. in press. Wiener, N. (1958) . Nonlinear Problems in Random Theory. Wiley, New York. Zheng, Q. and E . Zafiriou (1995). Nonlinear system identification for control using VolterraLaguerre expansion. In: Proc. American Control Conf.. Seattle, WA. pp. 2195-2199.
Table 3. Generation of Inverse Repeat Sequence - Methodology Number of clock pulses
msequence
1
2 3 4 5 6 7 8 9 10 11
o o o
o 1
o
12 13 14 15 16 17 18 19
1
o o o o o
20 21 22 23 24
o
25 26
o
1
27
28
1
29
o o
30
PRBS ut{t)
+0 -0 -0 -0 +0 +0 +0 +0 -0 +0 -0 +0 +0 -0 -0 +0 -0 -0 -0 +0 +0 +0 +0 -0 +0 -0 +0 +0 -0 -0
Inverse repeat sequence
1 0 1 1 0 1 0 0 0 0 0
Inverse repeat signal U2(t)
+0 +0 -0
+0 +0 -0
+0 -0 -0
-0 -0 -0
+0 +0 0 0 0 1 0 0
-0
-0 -0
+0 -0 -0
+0 0 1 1
1 0 0
-0
+0 +0 +0 +0 +0 -0 -0
+0
7. REFERENCES Baheti, R S., R R Mohler and H. A. Spang (1980). Second-order correlation method for bilinear system identification. IEEE Trans. Autom. Contr. AC-25, 1141- 1146. Barker, H. A. and R W. Davy (1978) . Secondorder volterra kernel measurement using pseudorandom ternary signals and discrete fourier transforms. Int. J. of Control 21, 277291. Barker, H. A., S. N. Obidegwu and T . Pradisthayon (1972) . Performance of antisymmetric pseudorandom signals in the measurement of volterra kernels by cross-correlation. Proceedings lEE 119, 353-362. Congalidis, J. P., J. R llichards and W. H. Ray (1989). Feedforward and feedback control of a solution copolymerization reactor. AIChE J. 35(6), 891-907. Donce, J. 1. (1966) . Barker sequences. Electron. Letters 2(4) , 159--. Doyle, F .J ., B.A. Ogunnaike and RK. Pearson (1995). Nonlinear Model-Based Control Using Second Order Volterra Models. Automatica (31) , 697-714. Doyle Ill, F . J ., R S. Parker, R K. Pearson and B. A. Ogunnaike (1999). Plant-friendly identification of second-order Volterra models. 844