INTRODUCTION
Distillation accounts for approximately 95% of the separation systems used for refining and chemical industries (Humphrey et al., 1991). It needs to be controlled close to optimum operating conditions because of economic incentives. Close control of distillation column improves the product quality, minimizes energy usage and maximizes the plant throughput and its economy. It is a complex multivariable system and exhibit highly nonlinear dynamic behavior due to their nonlinear vaporliquid equilibrium relationships. A clear and thorough understanding of the nonlinear dynamic characteristics of the system is necessary for selecting an appropriate structure for the nonlinear model.
Development of nonlinear model is the critical step in the application of nonlinear model based control strategies. Nonlinear behavior is the rule, rather than the exception, in the dynamic behavior of physical systems. Most physical devices have nonlinear characteristics outside a limited linear range. In most chemical processes, understanding the nonlinear characteristics is important for designing controllers that regulate the process (Eskinat et al., 1991). Many authors (Leontaritis and Billings, 1985; Juditsky et al., 1995; Sjoberg et al., 1995; Qin and Badgewell, 1998; Pearson, 2003) have noted the difficulty of developing the models required for nonlinear modelbased control strategies. With carefully designed data collection experiments, the dominant behavior of plant can be fitted into one of the several possible structures. The main challenge in this task is to select a reasonable structure for the nonlinear model to capture the process nonlinearities. The nonlinear model used in control purposes should be as simple as possible, warranting less computational load and at the same time retain most of the nonlinear dynamic characteristics of the system (Henson, 1998).
The practical difficulty of nonlinear dynamic model development arises from several sources, of which the following two are fundamental. First fact is that model utility can be measured in several conflicting ways and second, the fact that the class of nonlinear models does not exhibit the unity that the class of linear models does. The four extremely important measures of model utility are approximation accuracy, physical interpretation, suitability for control and ease of development (Pearson, 2003). Fundamental models are generally far superior to empirical and semiempirical models with respect to the first two of these criteria, but they also suffer badly with respect to the last two. At the other extreme, blackbox models are purely empirical, based entirely on input/output data. Consequently, these models generally lack of physical interpretation possible for whitebox models, on the other hand, the option of choosing convenient model structures that facilitate formulating and solving the associated control problem.
In engineering dynamics, control engineering and many other areas, autoregressive with exogenous inputs (ARX) models are widely utilized for describing dynamic data regimes for linear and nonlinear systems (Erik, 2000; Mahfouf et al., 2002; Karny and Pavelkova, 2007). However, the performance of these linear models for prediction and control has been limited. In particular, the nonlinear nature and strong directionality of the process present problems during identification. NARX model has been considered as alternatives to linear models in a number of chemical process applications such as distillation column (Sriniwas et al., 1995; Verhaegen, 1998), CSTR (Lee and Lee, 2005; Srinivasarao et al., 2006), pH process (Proll and Karim, 1994; Bomberger and Seborg, 1998) etc. Many algorithms have been employed to identify the parameters of the NARX models. Among these algorithms, Prediction Error Minimization (PEM) method (Spinelli et al., 2006), Least Squares (LS) algorithm (Chen et al., 1989; Previdi and Lovera, 2004) and genetic algorithm (Li and Jeon, 1993; Chen et al., 2007) are the most often used.
Sriniwas et al. (1995) have used polynomial NARX model for identification and control of a highpurity distillation column and they didn’t use any systematic analysis to select the past inputs and past output in the model. Verhaegen (1998) made the comparative analysis of Wiener identification approach and NARX neural network blackbox identification method and concluded that former one is better in the identification of the temperatureproduct quality relationship in a multicomponent distillation column. None of the papers have explained about the regressor analysis which is important to capture the data relation in NARX modeling. Hence in this work, new wavenet based NARX model is used to capture the nonlinear dynamics of the distillation column and regressor analysis is carried out to ensure the proper selection of regressors.
DISTILLATION COLUMN CASE STUDY
The schematic of pilot plant distillation column utilized in this study is
shown in Fig. 1. The methanolwater binary mixture was used
as feed. The top and bottom product compositions are the controlled variables
in distillation column. The reflux flow rate and reboiler vapor boil up rate
are used as manipulated variables, whereas feed flow rate and feed composition
are considered as disturbances. The nominal operating conditions and the column
parameters are listed in Table 1. The reboiler and feed temperatures
were controlled using separate PID controllers. In all the experiments, tray
2, tray 6, tray 10, tray 14, distillate and bottom product temperatures are
measured. The top and bottom product compositions are determined using refractive
index analysis.
The experimentally validated first principle model is used as a process model
in this work. In the process model, calculation of activity and fugacity coefficients
are included in order to account for the nonideality.

Fig. 1: 
Pilot plant distillation column 
Table 1: 
Nominal operating conditions and the column parameters 

The activity coefficients are calculated using UNIFAC model and the fugacity
coefficients are calculated using virial equation of state. The detailed first
principle model equations and the experimental validation of model were discussed
in somewhere else (Ramesh et al., 2007).
NONLINEARITY STUDIES
The nonlinear dynamic characteristic of the system is necessary to select an
appropriate structure of the nonlinear model along with an input and output
variables for the nonlinear model. The dynamic behavior of the process is studied
by giving positive and negative step changes in reflux flow rate, vapor boil
up rate, feed flow rate and feed composition.

Fig. 2: 
Responses in top product composition for ±10% change
in R 
The response in top product composition change for ±10% change in reflux
flow rate (R) is shown in Fig. 2. It was observed that the
magnitude of % change in top product composition for +10% R is different compared
to 10% change in R, i.e., asymmetric responses to symmetric input changes and
it clearly indicated the violation of odd symmetry of the linear systems. The
similar types of responses were obtained for top product composition for ±10%
change in vapor boil up rate (VB), feed flow rate (F) and feed composition.
All the dynamic simulation studies have indicated that the violation of odd
symmetry of the linear systems. The asymmetric responses to symmetric input
changes in the dynamic behavior specified the necessity of nonlinear model (Pearson,
2003) for distillation column.
NARX MODEL
The nonlinear ARX structure models dynamic systems using a parallel combination
of nonlinear and linear blocks, as shown in the Fig. 3.
The nonlinear and linear functions are expressed in terms of variables called
regressors, which are functions of measured inputoutput data. The predicted
output Ŷ(t) of
a nonlinear model at time t is given by the following general equation.
where, u(t) represents the regressors. F is a nonlinear regression function,
which is approximated by the nonlinearity estimators. The function F can include
both linear and nonlinear functions of u(t), as shown in Fig.
3.

Fig. 3: 
NARX model structure 
The past output regressors are represented as:
The past input regressors are represented as:
Where:
na 
= 
Number of past output terms used to predict the current output. 
nb 
= 
Number of past input terms used to predict the current output. 
nk 
= 
Delay from input to the output in terms of the number of samples. This
value defines the least delayed input regressor. 
The meaning of na and nb is similar to the linearARX model parameters in the sense that na represents the number of output terms, nb represents the number of input terms and, nk represents the minimum input delay from an input to an output.
Wavenet structure based nonlinear function x = F(u) is used to represent the
static nonlinearity of the Hammerstein model.
where, f(u) is the scaling function given by
g(u) is the wavelet function given by
Parameters with the s are scaling parameters and with the w are wavelet parameters.
In this study, two past outputs (na = 2), current input and one past input (nb = 2) without delay (nk = 0) were used as regressors in the fourth order nonlinear function. The parameters of the NARX model were estimated using standard prediction error method similar to that one described by Eskinat et al. (1991). The model parameters are given as follows.
The values of nonlinear subspace P and linear subspace Q associated with wavenet function model are same for this case.
The values of other nonlinear parameters are given by
The output offset d 
= 
4.9427x10¯^{4} 
Regressor mean r 
= 
[5.0717x10¯^{4} 4.9979x10¯^{4}2.5621x10¯^{5}
3.119x10¯^{5} ] 
The different sets of data generated from validated first principle model were
used for parameter estimation and model validation. The data used in the NARX
model validation is the one which is not used in the parameter estimation. The
multilevel changes were made in reflux flow rate to generate u(t) which is the
input to the NARX model. u(t) is the difference between the steadystatevalue
and the value of the reflux flow rate in the particular time instant.

Fig. 4: 
Input u(t) for multilevel changes in reflux flow rate 

Fig. 5: 
Comparison of the output of the NARX model and process model 
The input signal u(t) used for model validation was shown in Fig.
4.
y(t) is the output of the NARX model which is the difference between
the steadystate value and the current value of the top product composition
in the particular time instant. The comparison of the NARX model output and
the experimentally validated first principle model (process model) output was
shown in Fig. 5. It has been found that the NARX model results
have shown 92.46% agreement with the process model results. The results proved
that the wavenet based NARX model was capable of capturing the nonlinear dynamics
of distillation column.
REGRESSOR ANALYSIS
The regressor analysis plot displays the characteristics of model nonlinearities
as a function of one or two regressors and this will helpful in identify which
regressors have the strongest effect on the model output.

Fig. 6: 
Effect of regressor y(t1) on nonlinearity 

Fig. 7: 
Effect of regressor y(t2) on nonlinearity 
Understanding the relative importance of the regressors on the output is needed to know which regressors should be included in the nonlinear function. The linear relationship between the regressor and nonlinearity reveal better performance of the regressor in capturing the nonlinearity.
In this research four regressors namely current input regressor u(t), past input regressor u(t1), past out put regressors y(t1) and y(t2) were used. The other regressors like u(t2) and y(t3) were studied and their inclusion doesn’t show any significant effect in the model results, but increased the complexity of the model. So only above mentioned four regressors were used in this NARX model and their effect on the nonlinearity of the process was discussed below.
The effect of regressor y(t1) on nonlinearity of the process was shown in
Fig. 6.

Fig. 8: 
Effect of regressor u(t) on nonlinearity 

Fig. 9: 
Effect of regressor u(t1) on nonlinearity 
It has been seen from the graph that the relationship between the nonlinearity
of the process and the value of regressor y(t1) is linear and proportional
to each other. So the nonlinearity of the process can be better captured by
the immediate past output of the process.
The effect of regressor y(t2) on nonlinearity of the distillation column was
shown in Fig. 7. It has been observed that the nonlinearity
decreases with increase in regressor y(t2) although the relationship is linear.
The effect of current input regressor u(t) on nonlinearity of the process was
shown in Fig. 8. It can be evident from the graph that the
relationship is linear, however the effect of regressor u(t) on nonlinearity
is not significant compare to regressor y(t1). The effect of input regressor
u(t1) on nonlinearity of the process was shown in Fig. 9.
It can be seen from the graph that the relation between regressor u(t1) is
slightly nonlinear compare to other three regressors. The study of effect of
each regressor in the NARX model on nonlinearity of the process has concluded
that the regressors of past outputs having more effect on nonlinearity of the
process compare to input regressors.
CONCLUSION
The studies on nonlinear dynamic characteristics of the distillation column using the experimentally validated first principle model have proven the necessity of nonlinear model for the distillation column. The wavenet based NARX model using input and output regressors was developed and applied to distillation column. The model validation results and the regressor analysis have concluded that the selected regressors have the potential of capturing the nonlinearity of the process.
ACKNOWLEDGMENT
The authors gratefully acknowledge the financial support by Ministry of Science, Technology and Innovation (MOSTI), Malaysia through the IRPA8 project No: 03024279EA019.
NOMENCLATURE
as_{k} 
: 
:Scaling coefficient 
AT 
: 
Analysis transmitter 
aw_{k} 
: 
Wavelet coefficient 
bs_{k} 
: 
Scaling dilation 
bw_{k} 
: 
Wavelet dilation 
cs_{k} 
: 
Scaling translation 
cw_{k} 
: 
Wavelet translation 
D 
: 
Distillate flow rate (kmol h¯^{1}) 
d 
: 
Output offset 
f(u) 
: 
Scaling function 
FL 
: 
Feed flow rate (kmol h¯^{1}) 
g(u) 
: 
Wavelet function 
LC 
: 
Level controller 
na 
: 
Number of past outputs 
nb 
: 
Number of past inputs 
nk 
: 
Delay from input to the output 
P 
: 
Linear subspace matrix 
PT 
: 
Pressure transmitter 
Q 
: 
Nonlinear subspace matrix 
R 
: 
Reflux flow rate (kmol h¯^{1}) 
R 
: 
Regressor mean vector 
Tc 
: 
Temperature controller 
TT 
: 
Temperature transmitter 
u(t) 
: 
Input to the NARX model 
VB 
: 
Vapor boilup rate (kmol h¯^{1}) 
XB 
: 
Bottom product composition 
XD 
: 
Top product composition 
XF 
: 
Feed composition 
y(t) 
: 
Output of the NARX model 