Copyright © IFAC Dynamics and Control of Process Systems. Corfu , Greece, 1998
CONSTRAINED PARAMETER ESTIMATION WITH APPLICATIONS TO BLENDING OPERATIONS Kenya Murakami t and Dale E. Seborg Dept. of Chemical Engineering, University afCalifornia, Santa Barbara. CA 93106 USA
[email protected] and
[email protected]
Abstract: The classical least squares approach ignores a priori information about the feasible values of the estimated parameters. But in many practical problems, such information is available in the form of upper and lower limits. In this paper, two alternative techniques are evaluated for this important class of constrained parameter estimation problems for linear. dynamic systems. Simulation results for two blending problems illustrate that more accurate parameter estimates and better predictions can be obtained using a quadratic programming approach. Copyright © 1998 IFAC Keywords: Blending Systems, Least Squares, Constrained Parameter Estimation, Inequality Constraints, Quadratic Programming
I. INTRODUCTION
constrained estimation has been derived (Gorman and Hero, 1990). Several authors have analyzed the effects of uncertainties in the assumed inequality constraints (O'Hagan and Leonard, 1976). A simple computational approach for recursive least squares is to orthogonally project the parameter estimates onto the constraint set after a violation occurs (Goodwin and Sin 1984).
Classical least squares estimation techniques have been widely used in engineering problems and provide optimal estimates under ideal conditions. However, the classical approach ignores a priori information about feasible values of the estimated parameters. Such information is widely available in practical problems, often in the form of inequality constraints. For example, if a mole fraction or mass fraction is to be estimated, it is bounded by zero and one. Physical considerations and process knowledge often allow upper and lower limits to be placed on parameters such as rate constants, heat transfer coefficients, steady-state gains and time constants.
In this paper we compare two alternative methods for solving constrained parameter estimation problems for dynamic systems, based on linear or non linear models. The first method is a simple variation of the standard least squares (LS) approach, namely, that an estimated parameter is "clipped" by setting it equal to an upper or lower limit whenever it exceeds its allowable range . This modification will be referred to as "clipped least squares (CLS)". The second method uses a quadratic programming CQP) formulation where inequality constraints for the estimated parameters are included explicitly in the problem statement. The performance of these two methods is compared to standard LS for a simulated blending process.
Despite their widespread occurrence, constrained estimation problems have received surprisingly little attention in the engineering literature. For example, system identification textbooks usually ignore this topic. Norton (1976, 1979) treated inequality violations as "psuedo-observations" in a recursive estimation problem. A Cramer-Rao type bound for tVisiting researcher from Fuji Electric Co., Ltd.
229
2.
BLENDING PROBLEM
~ ru)
a o =ex p( -
Blending operations are widely used in industrial plants including oil refineries. chemical plants. and cement plants. However. it can be difficult to achieve the desired blend characteristics if the compositions of the inlet streams vary significantly in an unknown manner. But if the unknown compositions could be estimated accurately on-line. the blending control problem would be much easier.
(5)
bj={I_exp(_~ru)}(C,;c)
(6)
~ {l-exp( - ~ ~t)}
(7)
eo =
c
~t
is the sampling period and value of c, etc.
Consider the stirred tank blending system with n inlet streams shown in Fig. I. A dynamic model is developed based on the following assumptions: perfect mixing. isothermal operation. each process stream has the same constant density p. and the liquid volume in the tank V. is kept constant with overflow line. The model based on mass balance concerning the liquid and component has the form.
denotes the nominal
Euler integration for (2) gives the estimation model:
nonlinear
c(k + 1) =
~ -hoq(k) }:(k) + ho! c (k)qj (k)+ hoS(k) j
(8)
;=1
ru/V.
where. ho =
2.2 Least Squares Estimation: Linear Model (LS-L)
(1)
It will be shown that the problem of estimating the inlet compositions can be put into a standard form for LS estimation based on a linear model. First. it should be noted that the estimation of Ci is equivalent to the estimation of in view of (6). The LS-L formulation for the blending problem can be written as:
;=1
dc n (2) V-= Lcjqj -cq dt j=1 where C is the concentration of exit flow. c j is the concentration of i-th inlet flow. q is the exit flow rate. qj is the i-th inlet flow rate.
bj,
(9) where,
1
c'(k)-aoc'(k-l)
y= [
$=[q;(~-l)
v <-
_J
(l0)
c'(k-w+l)~aoc'(k-w)
c.q
e =[b Figure I. Schematic diagram of the blending system.
(12)
is the residual, and w is the data window size for estimation. Denote the solution to (9) as,
e
L5 L -
Three models, a process model. a linear estimation model and a nonlinear estimation model are used . Process model is used for simulating actual process. Combining (1) and (2) and introducing a stochastic disturbance Syields the process model:
-V
bolT
•••
l
E
2.1 Process Models
dc(t) 1 ~{ } 1 -;;- = ~ Cl (t) - c(t) qj (t) +-VS(t)
(11 )
q;(k - w)
(k) =
[b lL5 - L (k) . .. bL5 - L (k)r
2.3 Clipped Least Squares Model (CLS-L)
Estimation:
Denote the estimated parameters as
e
CL5 L -
(3)
(k) =
(13)
n
[b lCL5 - L (k) . .. bnCL5 - L (k)
r
Linear
(14)
Then,
Linearizing (3) introducing deviation variables, and discretizing gives the linear estimation model:
b
CL5
j
n
c'(k + 1) = aoc'(k) + ~)jq:(k) + eos
where, (4)
-
L
(k) =
CIL
bounds for
;=1
min(max(b
j
L5 L -
(k),
b/), bu) j
(15)
U
and cj denote the lower and upper Cj •
and
where. (16)
230
I
Symmetric Lower limit Upper limit 0 .8 1.2 SI 0.5 1.5 S2 0.0 2.0 S3 Asymmetric 0.8 1.5 Al 0.5 A2 1.2 0.5 2.0 A3 0.0 1.2 A4 0.0 1.5 AS 0 .8 2.0 A6
2.4 Quadratic Programming Formulation Based on Linear Model (QP-L)
I
Minimization of the mean squared residual (or equation error) can be expressed as, min
J.. eTHe + le
bl
bn
(18)
2
. ..
subject to
b/ $b; $b;u (i = I, ·· ·, n) where
(19)
Table 2. Inequality Constraints
e = [bl ...b rand
In order to compare the effectiveness of the six estimation techniques (LS-L, CLS-L, QP-L, LS-NL, CLS-NL, QP-NL), different sets of constraints were considered and are shown in Table 2. The constraint values in Table 2 have been normalized by the true values of the parameters. For example, for the two stream problem, constraint SI is expressed as,
n
H =cpTcp
(20)
r =- cpTy
(21 )
y is the vector in (10), and cp is the matrix in (11) .
2.5 Estimation Based on Nonlinear Model LS-NL, CLS-NL and QP-NL are formulated based on the nonlinear estimation model in analogy to LSL, CLS-L and QP-L, respectively (the detailed derivation is omitted).
O.8c l ~ Cl ~ 1.2c I ' O.8c 1 ~
The simulation study was performed in the following manner. Gaussian excitation and a noise sequence (either PI, P2 or M) were used to generate data for 2000 sampling instants based on process model (3) with Runge-Kutta integration method assuming a 0.3 min sampling period. A moving window of the w most recent samples was used to estimate Cl and C2 at each sampling instant. The parameter estimates were then used to predict c(k + p) , the future exit
3. SIMULATION RESULTS Representative simulation results are included for two cases, n = 2 and n = 5, where n is the number of streams to be blended. For each case, all of the inlet stream compositions are estimated and all flow rates are assumed to be known.
composition, p sampling instants ahead . (In a subsequent study, the future predictions will be used in a predictive control strategy.) A preliminary investigation indicated that values of w 10 and p I provided satisfactory results.
3.1 Two stream blending (n = 2). The nominal steady-state conditions are shown Table I .
v
Cl
= 30 moUm 3
= 0.2 m /min
c2
= 10 moUm 3
=0. 1 m3/min
C
=23.3 moUm 3
= I m3 3
lJl
lJ2
In
=
=
in the estimated Representative variatIons parameters for a simulation with noise P I and constraint set S I are shown in Fig. 2 (LS-UCLSUQP-L, LS-NUCLS-NUQP-NL). It can be seen that the differences between CLS-L and QP-L parameter variations are quite small. Estimates by CLS-NL and QP-NL are different at all, and the variations of estimates by LS-NL and QP-NL are strongly correlated. In both approaches, estimated values Cl and 2 vary in a quite symmetric way,
Table I. Nominal Conditions for the two stream problem. In order to provide the necessary excitation for the parameter estimation, the inlet flows were varied in a Gaussian manner, N(O, 10-4). Two types of process disturbances were considered: a zero mean Gaussian disturbance, "P I ", and a non-stationary disturbance with periodic and drifting components, "P2" : ~(k)=O.OIO{ 1.5+(1/150)k+sin(kJ13.6)+11(k)}
c2 ~ 1.2c 2
c
though the amplitude of them are different.
(22)
Especially, in the case of QP-NL (the right bottom figure in Fig. 2) , the estimate Cl does not hit the
where 11 is a Gaussian noise sequence. PI and P2 give a signal-to-noise ratio of one. As a measurement noise, Gaussian noise denoted as "M", giving a signal-to-noise ratio of five is considered.
corresponding constraints.
This is because Cl is
restricted by the constraints for
c
2
due to the
correlation between them, before it reaches to constraints for Cl .
231
Linear
Two inlet stream, nonlinear model
Nonlinear 0 .14
0 .12
o.1
...,
-~~------------~
..:
~ o.oe
~r-------------,
c: o
....130
~~I
~ O.06
1ok_·. .... ~ ,--
~\
t ...
#
J....." _' .\
f -
~
-
CJ LS CJ CLS OP
I
~
- ,,"
r
0 .0<
~r-------------,
r
0 .02
o
I SI
S2
S3
102
Al
A3
A4
AS
AS
Fig. 3. Prediction error J for PI and nonlinear model.
It can be seen from Fig. 3, that QP-NL is slightly better than or almost same as LS-NL, but CLS-NL is much worse than QP-NL and even LS-NL.
Fig. 2. Estimated parameters, Pl/S I case. (Left : linear model, right: nonlinear model,
-- :c,' - - - :Cl ).
This result can be explained in the following way. The LS-NL or QP-NL methods yield strongly correlated estimates, which was demonstrated in Fig. 2. This correlation minimizes the (one step) prediction error based on nonlinear model. Estimates by CLS-NL, on the other hand, clipped from the LS-NL estimates, and it breaks this correlation.
Two performance measures, J and Jp, are used to evaluate the estimation techniques. The mean squared prediction error J defined by (23) is calculated for each case, N-p
L
{c(k+ p)-c(k+ p)Y
k=w+1
J=
N-w-p
(23)
Figure 3 illustrates that similar trends occur for all of the constraints in table 2, consequently, S2 and A2 are chosen as representative of nine pairs of constraints in Table 2. A preliminary simulation showed that these pairs of constraints are representative of other pairs. J and Jp values for LSUCLS-UQP-L, LS-NUCLS-NUQP-NL are plotted in Figs. 4-6.
where N is the number of the data points, w is the data window size for estimation, p is the prediction horizon and c(k + p) is the predicted value of c(k + p) based on the estimated parameters at time, k + p.
In the case for linear model, c'(k + p)
(prediction for the deviation variable c'(k + p) ) is calculated, then adding the nominal value c gives
Two inlet stream . linear model 0.04
c(k+ p) .
I~
0 .035
The mean squared parameter estimation error Jp defined by (24) is also calculated,
l-
0.0 3
LS CLS OP
..., ~ 0. 025
-:
W § ~ '0
0 .02
I
!:? 0 .0 15
i
where C; is the true value of i-th parameter, cj(k) is the estimated value at time, k. In our study, the parameter values are, w = 10, p = I, N = 2000, and !1t =0.3 min.
0. 0 1 0 .005
o
P lIS2
PlIlo2
1
P2IS2
P2I102
MlS2
MlIo2
Fig. 4. Prediction error J for linear model.
For each set of conditions, the simulations were repeated for 21 different sets of random number seeds for the Gaussian excitation and random disturbances. The mean values for these 21 runs were used for purposes of comparison. The result J for LS-NL, CLS-NL and QP-NL for all pairs of constraints and the PI case are shown in Fig. 3.
From Fig. 4. it can be seen that both the CLS-L and QP-L methods provide significant improvements over the LS-L method. However. the QP-L values are only slightly lower than the CLS-L values. For disturbances P I and M, similar trends can be seen. i.e., QP-NL gives slightly better results than those of LS-NL. but the CLS-NL results are much worse. On the other hand, for disturbance P2, the results of LS-NL, CLS-NL and QP-NL are almost
232
values Cl and cz , respectively. On the other hand, LS-NL gives biased estimates, but correctly detects the changes in the parameters.
same. When criterion Jp is used , CLS-NL and QPNL provide significant improvements over LS-NL, though QP-NL gives only slight improvements over CLS-NL.
3.3 Five stream blending (n
Two inlet stream . nonlinear model 0 .2
-
CJ CLS
0 .16
QP
.., 0 . 14
eO.12
w 5
Next, the problem of five inlet streams and five compositions to be estimated is considered. The nominal operating conditions are shown in Table 3.
CJ LS
0.16
i
= 5 moUm3 3 C = 21 moUm 3 if; =0.06 m /min (i = 1,00 .,5)
v = 1m3
Cs
=40moUm
Cl
o. 1
U
3
=30 moUm3 3 c3 = 20 moUm c = 10 moUm3
Cz
~0. 06
.t
= 5).
0 .06 o.~
4
0.0 2 0
Table 3. Nominal Conditions for the five stream problem.
i PlIS2
PlIA2
P2IS2
P2IA2
MlS2
MlA2
Nonlinear
Linear
Fig. 5. Prediction error J for nonlinear model. Two inlet stream . nonlinaar model
o•
_20 L-___
i-_~~~
en
6O~
U
0 - - - - _, - _ - , . - . - " , , '
-'
Z40 I
.
-'20 ~
60
~0 .4
E
~
-'40 I
'" 0 .2
Cl.
PlIA2
P2I52
P2IA2
M/52
MlA2
LS-L
LS-NL
.s 35
'"
.
v..",
1Sr----------.
,.
Q)
' of---__.....
Q)
0
- '0 -20
L --::-20-----:~----:!SO
0
bmt> ( InIn )
320
340
~- _ .-
360
lime (mln)
Once again, the CLS-L and QP-L methods provide significant improvements over LS-L, for either kind of disturbance. The QP-L is consistently better than CLS-L, but the amount of improvement is modest. It can be seen that, based on criterion J, for disturbance PI and M, QP-NL gives significantly better result than that of LS-NL, and CLS-NL gives much worse results even than LS-NL.
LS-NL
1020
v '"
0
300
The mean values of J and Jp for LS-UCLS-UQP-L and LS-NUCLS-NUQP-NL for 21 sets of random number seeds are shown in Figs. 9-12 .
30+----117«1 ~ 25
30
, - , - \ " , , _ - , - \ . ..
360
The estimated parameters for LS-NL and QP-NL in Fig. 8 are not highly correlated, in contrast to Fig. 2.
~~-------
.~
~r---------.
trT'lln )
The same numerical values were used for parameters w, p, N, and 111, and for disturbances PI, P2 and M, as for the two stream case. The estimated parameters for 200 samples, disturbance PI, and constraint set S I are shown in Fig. 8. The estimated parameters for CLS-L and QP-L vary quite differently, in contrast to Fig. 2.
40
LS-L
E '
Fig. 8. Estimated parameters for five stream, Pl/S I case. (Left: linear model, right : nonlinear model. ---- : c,' - - - :Cs ).
To investigate the influence of plant-model mismatch, a case is studied in which the parameters change from nominal values in Table 1. No noise or disturbance is assumed . The parameter changes occur at t = 20 and t = 40. Figure 7 shows that LS-L gives quite poor estimates after the parameters Cl and Cz change from nominal 4'~-------
30&0
320
om.
3.2 Influence of plant-model mismatch.
'"
I
Cl. 020
, , - - - , ... ... - _ - .. r \ , r - ' . -
300
P1I52
Fig. 6. Parameter estimation error Jp for nonlinear model.
§
60
~4O
. .
Cl. 020
o
o
.
90
20
40
60
I.".... (mlnl
Fig. 7. Estimates with parameters changed.
233
Five inlet stream. linear model 0 .14
-
0 .12
[=:J
LS
[=:J
CLS
OP
.., o.1
.n5
-
r-
, ,i
Five inlet stream. nonlinear model
,.... I
0. 35
i
,
-
0. 3
I i
0 .25
§ W
0 .08
J~
:
[=:J
LS
[=:J
CLS
OP
0 .2
.~
~O. 06
~O.' 5
e
Cl.
Cl.
0 .04
o.1
0 .0 2
0.05
0
I
PlIS2
PlIA2
P21S2
11
P2JA2
MlS2
0
MlA2
Fig. 9. Prediction error J for five stream case and linear model.
P1 /S2
Five inlet stream. linear model
B
P2JA2
MlS2
MlA2
i
Five inlet stream. non linear model
LS
[=:J
B
CJ CLS -
P21S2
Fig. 11 . Prediction error J for five stream case and nonlinear model.
9'~--------~====~------~--1 [=:J
PlIA2
LS
CJ CLS
OP
-
o
PlIS2
P lIA2
OP
P2JS2
P2/A2
Fig. 12. Parameter estimation error Jp for five stream case and nonlinear model. For disturbance P2, LS-NL, CLS-NL and QP-NL give very similar results. Based on Jp , CLS-NL gives significantly bener results than LS-NL and QP-NL gives moderately better results than CLS-NL. 4.
ACKNOWLEDGMENT Financial support for Kenya Murakami from his employer, Fuji Electric Co., Ltd., is gratefully acknowledged.
CONCLUDING REMARKS
REFERENCES
In order to incorporate a priori information into online parameter estimation, two alternatives to the standard least squares (LS) approach have been considered: clipped least squares (CLS) and a quadratic programming (QP) formulation . Simulation results show that more accurate parameter estimates and better predictions can be obtained when a QP approach is used to estimate inlet compositions for a simulated blending system. Results for the CLS approach depend whether a linear or nonlinear model used. The nonlinear process model provided better parameter estimates even in the presence of plant/model mismatch.
Goodwin, G. C. and K. S. Sin (1984), Adaptive Filtering Prediction and Control, Prentice-Hall, Inc., Englewood Cliffs, NJ, 07632 Gorman, J. D. and A. O. Hero (1990), Lower Bounds for Parametric Estimation with Constraints, IEEE Trans. Inform. Theory , 26, 1285- I 30 1. Norton, J. P. (1976), Error Reduction by Linear Constraints in Recursive Identification, Electronics Letters, 12, 524-525 . Norton, J. P. (1979), Recursive Linear Estimation of Linearly Constrained Parameters, Electronics Letters, 15,715-716.
Recursive versions for these methods also can be derived . This approach is equivalent to those proposed in Goodwin and Sin (1984), if the inequality constraints are given in linear form.
O'Hagan, A., and T. Leonard (1976), Bayes Estimation Subject to Uncertainty about Parameter Constraints, Biometrika , 63, 201203 .
234