Abstract
This paper is taking actual control demand of rotary kiln as background and builds a calcining belt state space model using PO-Moesp subspace method. A novel order-delay double parameters error criterion (ODC) is presented to reduce the modeling order. The proposed subspace order identification method takes into account the influence of order and delay on model error criterion simultaneously. For the introduction of the delay factors, the order is reduced dramatically in the system modeling. Also, in the data processing part sliding-window method is adopted for stripping delay factor from historical data. For this, the parameters can be changed flexibly. Some practical problems in industrial kiln process modeling are also solved. Finally, it is applied to an industrial kiln case.
1. Introduction
Rotary kiln is the key equipment in metallurgy, building materials, and many other industries. A lot of rotary kiln control theories and methodologies have been studied, for the reason that kiln belongs to the typical complex process, which has the characteristics of multivariables, time delay, serious coupling, and difficult modelling [1–5].
In the field of rotary kiln modeling, researchers started studying internal status of calcinations process as early as 1960s. The initial researchers established control model based on solid mechanics, heat transfer mechanism, and dry calcined kinetics. Imber and Paschkis calculated the optimum kiln length for heating process in 1961 [6, 7]. These early approaches are all based upon the rmodynamic equilibrium and analysis of the kiln structure, and the establishment of these equations harshly demands material flow balance. So, it is not a good way to put in use widely. In 1980s dynamic system identification technology was used in rotary kiln modeling, but the identification model can only characterize part states of kiln. So this technology was restricted in application and exploitation. In the 1990s, after the emergence of the new intelligent modeling and control methods, a large number of intelligent modeling ways started to use in building model of kiln successfully [8–10]. After a long-term development many mature experiences have already been obtained, but there are still many problems such as imprecise model, excess restrictions, and difficulty in promotion. Consequently, building model has become a bottleneck on the development road of kiln process control.
Subspace identification method as an effective identification modelling way for multi-input and multi-output system has drawn much attention recently [11–13]. In comparison with the conventional methods, it has obvious advantage. This is because the model can be directly got from the input/output data without nonlinear optimization and iteration. Furthermore, the computations are based on robust tools such as QR-factorization and singular value decomposition (SVD), for which numerically reliable algorithms are available [14, 15].
Order and delay are the most important construction parameters of industrial process model. There exists an extensive literature for order estimation algorithms of linear, dynamical, state space systems [16–20]. Among all of these methods, the most famous one is AIC (Akaike information criterion), which is introduced by Akaike [21]. And also, great deal of subsequent researches were done about the properties effects of the penalty term [22, 23]. And then Baure introduces another criterion for the Larimore type of procedures, which is similar to AIC [24, 25]. From above literature, an obvious problem can be found that the delay parameter is always neglected in systems order selection. The information about the delay of any process is valuable for industrial process model. However, there exist only few references dealing with the estimation of the order in delay system [26–28]. When the system has input delay, many methods attempt to increase the model order as the cost to improve model accuracy. This enhances the system complexity.
In this paper, the calcining belt state space model of rotary kiln is built by using PO-Moesp subspace method. After analyzing the main components technology and calcinations reaction mechanism in detail, the input and output variables of modeling are selected according to industrial field operational procedures. And more, a novel order-selecting method is put forward. The original error criterion is replaced by the new one, which includes the two parameters—the order and the delay. This improvement effectively solves the problem of model order too high when modelling object with time delay, and it also increases the precision of model.
The outline of this paper is as follows: in Section 2, the basic process of rotary kiln is described; in Section 3 the general subspace identification algorithm is introduced and the order selection problem is addressed. Section 4 presents the main results of this paper. The method to identify the order and the delay is deduced in detail. In Section 5, we discuss the performance of the proposed method by means of some simulation studies on an industrial rotary kiln illustration.
2. Lime Rotary Kiln Process Description
The rotary lime kiln is a long cylinder that has a 3–5°d slope, and it rotates about its axis at roughly 0.5–1.5 rpm. Lime mud is fed from the elevated cold end and moves down the kiln due to rotation and gravity. The hot end, in which the burner operates, is typically maintained at about 1200°C by burning fuel [1]. Active lime rotary kiln system has many technological processes, such as long kiln technology system, chain grate cooling technology system, and vertical preheater technology system. The process is mainly divided into material system and air flow system. Raw material becomes product after preheating, high temperature calcinations, and cooling. Primary air, secondary air, and coke oven gas constitute the air flow system. The whole technological process is shown in Figure 1. In this paper, the lime kiln under study is 600 t/d, 4.0 m × 60 m.
The main calcinations process parameters are listed in Table 1.
3. Problem Formulation and Assumptions
A discrete-time model in state space form is described by the equations: where input and output , and are input and output dimensions, respectively; is the system state at time , and is the system order to be identified; process noise vector and measurement noise vector are zero mean value white noise series.
We introduce the following hypotheses (van Overschee and de Moor, 1995) [29].(a)The system is asymptotically stable.(b)The system is observable and reachable.(c)All the process noise and measurement noise are statistically independent of the input , that is mean: where .
The goal of subspace methods consist in the online estimating order , delay , and system matrices , , , and .
3.1. Matrices Construct
Suppose at time the past input and output data samples are available; then construct past/future output vectors and the Hankel output matrices as follows:where , , represent the past information, and represent the future information, , .
The state matrix and is defined as follows: , .
The augmented observation matrix
Block-triangular matrix Toeplitz matrix and
The augmented observation inverse matrix of is . And the augmented observation inverse matrix of is .
Suppose the input and state matrices are all row full rank matrices, and they are row space orthogonal. Then the augmented input/output matrices can be written as follows:
The key step of Moesp method is to estimate the augmented observability matrix through the projection future data onto past data. The augmented observability matrix and state space vector estimation can be got through SVD decomposition.
3.2. General Algorithm
In this section, we introduce the general subspace method framework which is proffered by Favoreel et al. [30]. Subspace identification algorithm consists of two main steps. The first step always performs a weighted projection of the row space of the previously defined data Hankel matrices. From this projection, the augmented observability matrix and estimate of the state sequence can be retrieved. Then the order is got from SVD decomposition. In the second step, the system matrices , , , and , , are determined through least square method. The concrete Moesp algorithm is as follows.(1)Project the row space into the orthogonal complement of the row space: Since it is assumed that the noise is uncorrelated with the inputs, so and ; therefore .(2)Select the weighting matrices and : The weighting matrices can be chosen appropriately according to different subspace methods, including N4SID, MOESP, CVA, basic-4SID, and IV-4SID [31–34].
Then, we can get (3)Carry SVD decomposition: And then take the number of nonzero eigenvalue as the system order .(4)The augmented observability matrix or is derived from the third step.(5)Extract estimate , , , from or .
Remark. By reference to [30], the weighting matrices and should satisfy the following three conditions:(1);(2);(3).
The first two conditions guarantee that the rank- property of is preserved after projection onto and weighting by and . The third condition expresses that should be uncorrelated with the noise sequences and . By choosing the appropriate weighting matrices and , all subspace algorithms for LTI systems can be interpreted in the above framework, including N4SID, MOESP, CVA, Basic-4SID, and IV-4SID.
4. The Proposed Method
In the classical system identification theory, the actual model structure is usually assumed to be known. However, in practical it is always not clear. Subspace system identification method determines the order of the system by the nonzero eigenvalue of the augmented observability matrix. However, the system nonzero singular values may be very small. This may lead to the wrong system order and large identification error.
4.1. The Order-Delay Double Parameters Error Criterion
The most directly order-selection method is based on the error performance criterion. This idea is to choose the smallest possible order that keeps the error below a certain level. Then the MRSE (mean relative squared error) index is introduced by model error as follows: where is the model prediction error and is the sample number. In [35] use the AIC, which was originally developed by Akaike and then adapted by Larimore for SMI. Given a set of samples, for a sequence of system order , for example, , the order of the model will be the one which makes the following AIC index minimum: where
For calculating the criterion, we first suppose the upper bound of the system order and then calculate the sequence; the appropriate system order is the one which decrees the AIC index obviously, and the order should be as small as possible. The MRSE and AIC index can be analyzed in the same manner.
However, the performance index based on a single order parameter cannot provide an effective solution to the delay system, which is shown in (1). This led to the problem that the original identification method had to increase the order as cost to improve model accuracy. Here, we introduce an order-delay double parameters error criterion, which identifies the two key structural parameters at the same time. That means that the index is changed into the form.
For each given individual , a state space model can be identified using the Moesp algorithm described in Section 2. Then, its model error can be deserved as , then the with respect to individual as
Also, the MRSE criterion has the similar form:
Other performance index such as SVC IVC, NIC criteria mentioned in [25] can also be modified as this method.
The original performance index just identifies the order . Suppose ; then are calculated, respectively; then the inflection point is the best order and the corresponding system matrices , , , are the best model. After improvement, we add the delay as other optimal parameters. So the system order and input delay are all embedded in . Then calculate , , respectively. By searching the minimum point of surface , the best order, delay, and system matrices can all be got. Thus, taking the delay as another parameter in the modelling methods, it can effectively avoid high order results in the delay system.
4.2. The Delay Factor Stripping from Historical Data
The introduction of delay parameters in performance criterion has resulted to a notable problem. The modelling historical data matrices have already included delay information. It is difficult to change in the performance criterion artificially. To solve this problem, the sliding-window method is adopted here. Sliding-window principle is shown in Figure 2. When new samples are added to the window, the oldest data inside the window will be discarded.
We use sliding window to change delay, which is shown in Figure 3. Suppose the output data length is ; select the input data region . Then the input data can be moved according to different delay from to .
4.3. The Algorithm Description
The detailed procedure of subspace identification based on ODC algorithm can be expressed as Algorithm 1.
| ||||||||||||||||||||||||
Then the best combination , and the optimum matching model parameters are all obtained.
5. Simulation Results
To demonstrate the superiority of the proposed order selection method in this paper over the conventional method, their performance is evaluated through a numerical example and an industrial illustrate.
5.1. Example 1: MIMO Process with Input Delay
The first model (15) is a MIMO system:
The time delay is and the sample is . We use the zero mean white noise as the input and sequence to excite the system. The input and output curves are shown in Figure 4.
(a) Input signal
(b) Input signal
(c) Output signal
(d) Output signal
Firstly, we use conventional subspace model selection method, which is dependant on the number of eigenvalues. Figure 5(a) is the situation of delay = 0 and Figure 5(b) is the situation of delay = 5. According to this strategy, the order of model increases from 4 to 13 when the system has delay.
(a) situation
(b) situation
Also the MSE (mean squared error), MRSE (mean relative squared error), and AIC criteria are all tested as shown in Figures 6 and 7. These methods all have the distinct inflection point at 13, which is for the 4-order system with input delay 5. Obviously, for obtaining the ideal model of a delay system, these conventional methods have to increase the order.
Next, the performances of the proposed ODC method in this paper are presented. As the two search parameters are available, so the index performance is shown as a surface. The of output is shown in Figure 8 and the of output is shown in Figure 9. Three axes are the order , delay , and , respectively. The AIC index surface is shown in Figure 10. The minimum value of these indexes can be got easily; they are also the inflection points. The order and delay results are . This is the same with the actual system model.
The corresponding identified model matrices
For verifying the model, the model output error is drawn as and in Figure 11.
(a)
(b)
5.2. Example 2: The Kiln Industrial Illustration
To demonstrate the superiority of the proposed order selection method in this paper over the conventional method, their performance is evaluated through an industrial illustration. The data come from actual kiln production data of an enterprise. Here, the gas flow and the second air flow are selected as the control input, and the calcination temperature and kiln tail temperature are taken as output variables. The sampling time is ; then use Moesp method modelling the kiln based on the input/output data after preprocessing. Here, three main practical problems are solved.
5.2.1. The First Problem
The problem is that these two input variables have different delay. They need to be identified, respectively; that is to say, a triple loop about , , and should be carried for solving .
In order to obtain the inflection point information more directly, we identify the order firstly. According to the industrial field situation, choose the possibility maximum values which are , . At first, travel all the possible , and compute the smallest , , and corresponding to the different order . These curves are shown in Figure 12. As can be seen from the curves, all of these three error criteria achieve the inflection point at , so the order is got. Then specify ; the triple loop is reduced to double loop, which is to compute and corresponding to and between .
5.2.2. The Second Problem
From Figure 13 we notice another problem that the outputs and generate different inflection point.
In order to solve this problem, the criterion should choose MRSE and AIC, which take into account and both together. Then get Figures 14 and 15.
It can be conducted that .
5.2.3. The Third Problem
However, there is still a problem which needs to be resolved. In general, the model with large time delay will increase the difficulty of controller designing. We know that when the sampling frequency is higher than the actual needed frequency, there will be lots of redundant data. And this will raise the model order and the delay. Therefore, the delay can be reduced by properly decreasing sampling frequency. Considering the kiln is a slow time-varying process, changing the sampling time from 1 s to 10 s will not affect the model accuracy.
From Table 2 it can be seen that when set , the order and inflexion point is still . We can also change the delay to , . The surface can be got as in Figure 16; the results show that . This is in agreement with the analysis before.
The corresponding rotary kiln calcining zone temperature model is and the order ; delay ;
Compare the measured value and the predicted value generated by the model in Figures 17 and 18.
6. Conclusion
In this paper, the calcining belt state space model of rotary kiln is built using PO-Moesp subspace method. And a novel double parameters error performance criterion for the order choosing in subspace modelling is introduced. Since the presented method considering the order and delay simultaneously, it can reduce the model order of the delay system effectively. And also, a strategy for stripping the delay factors from the historical data is also proposed. The algorithm is verified in identifying an industrial lime kiln. In this example, we solve the practical problem of industrial process with multidelay and also reduce the order by adjusting sampling time. Further research could shed more light on the issue of applying the model online. The problem in industrial field is more complex than simulation environment. How to extract problems from industrial practice and guide the direction of modeling research has become the study focus.
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.