Ind. Eng. Chem. Res. 2008, 47, 4929–4936
4929
Integral Identification of Continuous-Time Delay Systems in the Presence of Unknown Initial Conditions and Disturbances from Step Tests Qing-Guo Wang,* Min Liu, and Chang Chieh Hang Department of Electrical and Computer Engineering, National UniVersity of Singapore, Singapore 119260
Yong Zhang GE Global Research, General Electric Company, 3050 Red Hill AVenue, Costa Mesa, California 92626
Yu Zhang Global Research Shanghai, GE China Technology Center, 1800 Cai Lun Road, Zhang Jiang Hi-tech Park, Pudong, Shanghai 201203, China
Wei Xing Zheng School of Computing and Mathematics, UniVersity of Western Sydney, Locked Bag 1797, Penrith South DC NSW 1797, Australia
In this paper, an integral identification method is proposed for continuous-time delay systems in presence of both unknown initial conditions and static disturbances from a step test. The integration limits are specifically chosen to make the resulting integral equation independent of the unknown initial conditions. This enables identification of the process model from a step test by one-stage least-squares algorithm without any iteration. The proposed identification method is demonstrated through numerical simulation and real time testing. 1. Introduction System identification plays a important role in system analysis, control, and optimization.1 Time delay is present in most industrial processes and identification of continuous time processes with delay is of engineering importance. A considerable number of identification methods have been reported in the literatures.2–9 A differential equation with time delay can be transformed to an integral equation by means of multiple integration.10 The resulting integral identification has proven robust against noise in measurements (see the work of Golubev and Horwowitz11 and references therein for details). In general, the presence of time delay causes a nonlinear estimation problem for identification. Identification of continuous-time delay processes has been an active area in the last decades. Wang and Zhang2 address the problem using step test. They took the advantage of the step input and devised a linear identification algorithm for all the model parameters including delay. They, like the previous work on continuous system identification, assumed that the initial conditions are zero and there is no disturbance. It is possible that the underlying process is operated to the constant steady state and kept there so that the above assumption is met. On the other hand, these limitations are the major concerns from application perspectives, as also raised by the reviewers of Wang and Zhang.2 It is definitely desirable to remove the assumption for easy practical applications under the nonsteady state condition. Hwang and Lai5 and Wang, Liu, and Hang6 address the problem using pulse test. Hwang and Lai5 used a pulse signal as the input and presented a two-stage identification algorithm. Essentially, two changes of the pulse signal could be used to establish two independent integral equations so that estimation * To whom all correspondence should be addressed. E-mail:
[email protected]. Tel.: (+65) 6516 2282. Fax: (+65) 6779 1103.
or elimination of nonzero initial conditions becomes possible. In the work of Wang, Liu, and Hang,6 Hwang’s method5 is simplified with reduced two regression equations. Wang, et al.8 address the problem using relay test. Basically, a relay test is regarded as a sequence of step tests. In the work of Ahmed, Huang, and Shah,7 an alternative identification method for a process with time delay under nonzero initial conditions is presented through the linear filtering method. However, this method needs an iterative procedure for the time delay estimation. A improved identification method is developed by Liu, et al.9 Their method applies to many types of test signals including open-loop tests such as PRBS, rectangular pulses, and rectan-
Figure 1. Step input and response for example 1.
10.1021/ie071532s CCC: $40.75 2008 American Chemical Society Published on Web 06/11/2008
4930 Ind. Eng. Chem. Res., Vol. 47, No. 14, 2008 Table 1. Identification Results for Example 1 NSR 3% 5% 10% 15% 25%
estimated model [-0.4729s [-0.4527s [-0.4106s [-0.3534s [-0.2075s
+ + + + +
2
0.5062]/[s 0.5096]/[s2 0.5166]/[s2 0.5252]/[s2 0.5423]/[s2
+ + + + +
1.502s 1.503s 1.506s 1.508s 1.514s
+ + + + +
-1.02s
0.4971]e 0.4953]e-1.03s 0.4917]e-1.07s 0.4872]e-1.12s 0.4782]e-1.28s
gular doublet pulses, as well as close-loop tests such as relay tests. Note that in the works of Hwang and Lai,5 Wang, Liu, and Hang,6 Wang, et al.,8 Ahmed, Huang, and Shah,7 and Liu et al.,9 the test signal has at least two changes of its magnitude. These methods all fail if the step test, the most popular one in process control applications, is used, because a step test only has one change of its magnitude. In this paper, a new integral identification method is proposed for a continuous-time process with time delay. The test input under consideration is a step signal. The process can have nonzero initial conditions and static disturbance which are unknown. The key idea in our method is to make both upper
estimated disturbance
error ε
0.9864 0.9783 0.9622 0.9421 0.9016
8.32 × 10-4 2.1 × 10-3 6.3 × 10-3 1.48 × 10-2 4.3 × 10-2
and lower limits of the inner integral dependent on the dummy variable of the outer integral so that the initial conditions do not appear in the resulting integral equation. The used multiple integration is different from those in the works of Hwang and Lai,5 Wang, Liu, and Hang,6 Ahmed, Huang, and Shah,7 and Liu et al.9 and similar to that in Wang et al.8 A novel regression equation different is proposed, so that the presented method can identify the delay process in one-stage with no iteration from a step test. The effectiveness of the proposed method is demonstrated through examples.
Figure 4. Step response for example 3. Figure 2. Nyquist plot for example 1.
Figure 3. Step response for example 2.
Figure 5. Nyquist plot for example 4.
Ind. Eng. Chem. Res., Vol. 47, No. 14, 2008 4931
input test is applied to the process with the initial conditions of y(0) and y(1)(0). The task is to estimate the model parameters, a1, a0, b1, b0, and d, from the input u(t) and output measurement y(t) in presence of unknown c and y(1)(0) which could be nonzero. To avoid the use of various time derivatives, which are too sensitive to noise, (1) is transformed to an integral equation by multiple integration. Normally, the integral interval is chosen from 0 to t.10 Thus, integrating (1) from 0 to t twice gives
∫ y(δ ) dδ - y(0)t] + dδ ) b ∫ ∫ u (δ - d) dδ dδ + ct (2) b ∫ ∫ u(δ - d) dδ dδ + 2 t
[y(t) - y(0) - y(1)(0)t] + a1[ a0
∫∫
δ1
t
0
0
0
0
y(δ0) dδ0
1
0
δ1
t
1 0
(1)
0
0
0
δ1
t 0 0
1
2
0
0
0
1
(1)
where y (0) is present but unknown. This is the first obstacle which makes the existing integral identification methods from step tests impossible to work in presence of unknown initial conditions, while Hwang and Lai5 use a pulse test whose two signal levels (like two tests) give rise to two independent equations so that the unknown initial conditions can be obtained or eliminated. Under u(t) ) 1(t), the unit step function, (2) can be rewritten as
Figure 6. NI temperature chamber set.
y(t) ) [-
[ ]
∫ y(δ ) dδ t
0
0
-
0
∫∫
δ1
t
0
0
y(δ0) dδ0 dδ1 1 t t2 ] ×
a1 a0
b0d2 - b1d y(0) + 2 (1) y (0) + a1y(0) - b0d + b1 b0 + c 2
:) φT(t)θ
Figure 7. Step test of the temperature chamber system.
This paper is organized as follows. In section 2, a common problem of the existing integral identification methods is revealed. In section 3, the method is presented for second-order modeling from step and pulse tests, respectively. The methods are further extended to high-order modeling in section 4. The proposed method is applied for a real temperature control system in section 5. Conclusions are drawn in section 6. 2. Review of Integral Identification In this section, we will use a 2nd-order model to show why the existing integral methods are unable to identify such a model from a step test under unknown nonzero initial conditions and static disturbance. Assume that a stable process is represented by
where there are five linear independent functions in φ(t), which enables estimation of five parameters in θ. But there are seven unknowns, a1, a0, b1, b0, d, c, and y(1)(0). Not all of them can be found from θ. The presence of y(1)(0) in the regression equation also increases the number of unknowns. This forms the second obstacle for the current integral identification. The essential cause which leads to these two obstacles and failure of the existing methods is that when a differential equation is transformed to an integral equation by multiple integration, the output derivative will inevitably appear in the resultant integral equation as long as one of integration limits is fixed. It should be pointed out that all the existing methods including Hwang and Lai5 have one integration limit fixed, indeed. In view of the above observation, the key idea in this paper is to make both upper and lower limit of any inner integral dependent on the dummy variable of the immediate next outer integral so that all the terms in the outcome of inner integral are functions of the outer dummy variable, but not fixed. 3. Proposed Method To get rid of the problem in the existing methods, we employ the following double-integral operation on f(t):
∫ [∫
y(2)(t) + a1y(1)(t) + a0y(t) ) b1u(1)(t - d) + b0u(t - d) + c
τ
0
(1) where y(t) and u(t) are the output and input of the process, respectively, d is the time delay, and c is the static disturbance or a bias value of the process. Suppose that at t ) 0, a step
(3)
t+δ1
t-δ1
]
f(δ0) dδ0 dδ1
For y(2)(t), one sees
∫ [∫ τ
0
t+δ1 (2)
t-δ1
]
y (δ0) dδ0 dδ1 )
∫ [y τ
0
(t + δ1) - y(1)(t -
(1)
(4)
4932 Ind. Eng. Chem. Res., Vol. 47, No. 14, 2008
∫1
δ1)] dδ1 ) y(t + τ) - 2y(t) + y(t - τ)
t
(1)
which depends on y(t) only but not on y (t). If, on the other t+δ1 (2) hand, any term in the outcome of ∫t-δ (δ0) dδ0 was 1 y independent of δ1, then when integrated with respect to δ1, there would be y(1)(0) terms in (5), which are not available. The double-integral operation is applied to y(1)(t) and y(t), respectively,
∫ [∫ τ
0
t+δ1 (1)
t-δ1
y (δ0) dδ0] dδ1 )
∫ [y(t + δ ) - y(t - δ )] dδ τ
1
0
1
1
(6)
∫ [∫ τ
0
t+δ1
t-δ1
]
y(δ0) dδ0 dδ1
(7)
which can both be numerically evaluated with knowledge of y(t). For the right-hand side of (1), consider the step test first since the step testing is the simplest and dominant in process control. Let u(t) ) h1(t). Then u(t - d) ) h1(t - d), the unit step function delayed by time of d. It is straightforward to verify that
Figure 8. Flowchart of the mixing procedure.
Figure 9. Step test of the flow control system.
(δ0 - d) dδ0 ) 1(t - d)
(1)
0
∫ 1(δ t
0
0 - d)
dδ0 ) (t - d)1(t - d)
∫ (δ t
0 - d)1(δ0 - d)
0
dδ0 )
(t - d)2 1(t - d) 2
It then follows that
∫ [∫ 1(δ - d) dδ ] dδ ) ∫ [(t + δ - d)1(t + δ τ
0
t+δ1
t-δ1
0
0
τ
0
1
1
1 - d) - (t - δ1 - d)1(t - δ1 -
d)] dδ1 1 ) [(t + τ - d)21(t + τ - d) - 2(t - d)21(t - d) + (t - τ 2 d)21(t - τ - d)] (8)
Ind. Eng. Chem. Res., Vol. 47, No. 14, 2008 4933
∫ [∫ τ
0
t+δ1
t-δ1
1(1)(δ0 - d) dδ0] dδ1 )
∫ [1(t + δ τ
1 - d) - 1(t - δ1 - d)]
0
dδ1
) (t + τ - d)1(t + τ - d) - 2(t - d)1(t - d) + (t - τ - d)1(t τ - d) (9) Let τ be fixed and t satisfy t - τ < d e t, which causes 1(t + τ - d) ) 1(t - d) ) 1,
det
(10)
1(t - τ - d) ) 0, t - τ < d (11) Integrating (1) in form of (4) and making use of (5)-(7) and (8)-(11) yields φT(t)θ ) γ(t),
t-τ Ts - d
(35)
Choose 2(n - 1)τ ≈ Ts to meet (35), which results in Ts τ≈ 2(n - 1) Equation (32) can be rearranged as
(36)
d + (n - 2)τ e t < d + (n - 1)τ Suppose d ∈ [dmin dmax]. Once τ is calculated from (36), ti are chosen as dmax + (n - 2)τ e ti < dmin + (n - 1)τ
The Nyquist plots of the process and model are given in Figure 5. 5. Real Time Testing
(30)
where Γ ) [γ(t1), γ(t2),..., γ(tN)]T and Φ ) [φ(t1), φ(t2),..., φ(tN)]T. The ordinary least-squares algorithm can be applied to (30) to find its solution θ ) (ΦTΦ)-1ΦTΓ
0.4207 e-2.56s s3 + 1.785s2 + 1.516s + 0.4243
(37)
N is such that t1, t2,..., and tN meet (37). Example 4. Consider a high-order process y(4)(t) + 4y(3)(t) + 6y(2)(t) + 4y(1)(t) + y(t) ) u(t - 2) + c
Laboratory Test. The proposed method was tested on a temperature chamber system. The experiment setup consists of two parts: a thermal chamber set (which is made by National Instruments Corp. (NI), as shown in Figure 6) and a personal computer with data acquisition cards and LabView software. The system input is a 12 V light with a 20 W halogen bulb, and the system output is the temperature of the chamber. An identification test was performed on the system and the recorded inputs and the output are given in Figure 7. At t ) 0, y(0) and y(1)(0) are nonzero. Applying the proposed identification method yields y(2)(t) + 18.55y(1)(t) + 55.34y(t) ) 632.6u(t - 0.106) The response for this model under the same input is also shown in Figure 7. The effectiveness of the proposed method is obvious. Field Test. Xi-Hua-Feng pulp and paper mills is located in Wuzhi, Henan Province, P. R. China. Three kinds of pulps are made by the mills: wood pulp, grass pulp, and recycled-paper pulp. These pulps are mixed together in the mixing tank. The flowchart of this process is given in Figure 8. It is required to stabilize the pulp concentration without large deviations from the given operation conditions. The flow rate of the pulp is often tuned to meet different manufacture requirements, and it is important to monitor and control the flow rate. One needs to identify a model for the pulp flow rate. The process considered consists of a valve and a pipe (DN100). The input is the position of the valve and the output is the flow rate (cubic meters per hour) in the pipe. A step test was applied by moving the valve from the fully close to a 1/6 open position. The resultant response of the flow rate is given in Figure 9. The proposed method was applied and one model obtained as y(2)(t) + 2.267y(1)(t) + 0.9351y(t) ) 214.3u(t - 3.2) The response for this model under the same input is also shown in Figure 9. 6. Conclusions In this paper, a new integral method has been proposed for identification of linear continuous-time delay system with unknown initial conditions and static disturbance from step tests. The integration limits are specially chosen to make the integral equation independent of the unknown initial conditions. The process model is obtained by one-stage least-squares algorithm with no iteration. The effectiveness of the proposed method is demonstrated by simulation, laboratory test, and field experiment. Literature Cited (1) Ljung, L. System Identification: Theory for The User, second ed.; Prentice Hall: Upper Saddle River, NJ, 1999. (2) Wang, Q. G.; Zhang, Y. Robust identification of continuous systems with dead-Time from step response. Automatica 2001, 37, 377–390.
4936 Ind. Eng. Chem. Res., Vol. 47, No. 14, 2008 (3) Wang, Q. G.; Guo, X.; Zhang, Y. Direct identification of continuous time delay systems from step responses. J. Process Control 2001, 11, 531– 542. (4) Hwang, S. H.; Wang, L. W. On the identification and model reduction of processes under unknown load disturbances. Ind. Eng. Chem. Res. 2003, 42, 1905–1913. (5) Hwang, S. H.; Lai, S. T. Use of two-stage least-squares algorithms for identification of continuous systems with time delay based on pulse responses. Automatica 2004, 40, 1561–1568. (6) Wang, Q. G.; Liu, M.; Hang, C. C. Simplified identification of time delay systems with non-zero initial conditions from pulse tests. Ind. Eng. Chem. Res. 2005, 44, 7591–7595. (7) Ahmed, S.; Huang, B.; Shah, S. L. Parameter and delay estimation of continuous-time models using a linear filter. J. Process Control 2006, 16, 323–331. (8) Wang, Q. G.; Liu, M.; Hang, C. C.; Tang, W. Robust Process Identification from Relay Tests in the Presence of Nonzero Initial Conditions and Disturbance. Ind. Eng. Chem. Res. 2006, 45, 4063–4070.
(9) Liu, M.; Wang, Q. G.; Huang, B.; Hang, C. C. Improved Identification of Continuous-time Delay Processes from Piecewise Step Tests. J. Process Control 2007, 17, 51–57. (10) Whitfield, A. H.; Messali, N. Integral-equation approach to system identification. Int. J. Control 1987, 45, 1431–1445. (11) Golubev, B.; Horowitz, I. Plant rational transfer function approximation from input output data. Int. J. Control 1982, 36, 711–723. (12) Astrom, K. J.; Hagglund, T. PID Controllers: Theory, Design, and Tuning, second ed.; Instrument Society of America: Research Triangle Park, NC, 1995. (13) So¨derstro¨m, T.; Sto¨ica, P. G. Instrumental Variable Methods for System Identification; Springer: Berlin, Germany, 1983.
ReceiVed for reView November 10, 2007 ReVised manuscript receiVed April 15, 2008 Accepted April 23, 2008 IE071532S
