Feedforward and Feedback Controller Performance Assessment of

Dec 20, 2003 - ... and Materials Engineering, University of Alberta, Edmonton,. Alberta, Canada T6G 2G6, and Syncrude Canada Ltd., Fort McMurray, Albe...
0 downloads 0 Views 119KB Size
Ind. Eng. Chem. Res. 2004, 43, 589-596

589

Feedforward and Feedback Controller Performance Assessment of Linear Time-Variant Processes Folake B. Olaleye,† Biao Huang,*,† and Edgar Tamayo‡ Department of Chemical and Materials Engineering, University of Alberta, Edmonton, Alberta, Canada T6G 2G6, and Syncrude Canada Ltd., Fort McMurray, Alberta, Canada

This paper is concerned with the feedforward and feedback control performance assessment of linear time-variant processes. The developed algorithms provide a way to calculate the timevarying minimum-variance feedforward and feedback control benchmark from routine operating data. The proposed techniques are illustrated through a simulated stirred-tank heater example and an industrial application. 1. Introduction Minimum-variance control and performance assessment techniques for linear time-variant (LTV) feedback control loops have been discussed by Li and Evans1 and Huang,2,3 respectively. Minimum-variance control provides a measure of the lower bound on the process variance that serves as an appropriate benchmark against which the control-loop performance can be evaluated.4 The key to the control performance assessment is to find a benchmark, which is the minimum-variance feedforward and feedback control in this work. Studies and discussions on the design of minimum-variance feedback-only controllers can be found in work by Astrom,5 Box and Jenkins,6 and other references. Box and Jenkins6 and Sternad and Stoderstrom7 have also discussed the design of linear time-invariant (LTI) minimum-variance feedforward and feedback controllers. Stanfelj et al.8 presented a hierarchical method for monitoring and diagnosing the cause of poor performance of LTI feedforward/feedback control systems using autocorrelation and cross-correlation functions. Desborough and Harris9 developed a performance assessment algorithm for LTI feedback plus feedforward control loops. Huang and Shah10 and Huang et al.11 developed methods for the performance assessment of multivariate LTI feedback plus feedforward control systems using minimum-variance control as the benchmark. One of the main differences between the LTI and LTV systems is the noncommutativity of the LTV transfer function multiplications. This paper is an extension of the performance assessment of LTV feedback control loops of Huang3 to that of LTV feedforward plus feedback control systems. The remainder of this paper is organized as follows: As a prerequisite of the LTV feedforward and feedback performance assessment, the performance assessment of LTI systems is discussed in section 2. In section 3, we present several general properties of LTV feedforward and feedback control systems. In considering the performance assessment of feedforward and feedback * To whom correspondence should be addressed. Tel.: (780)492-2971 or 9016. Fax: (780) 492-2881. E-mail: Biao.Huang@ ualberta.ca. † University of Alberta. E-mail: [email protected]. ‡ Syncrude Canada Ltd. E-mail: [email protected].

control loops, two scenarios are possible: mutually correlated measured disturbances and mutually uncorrelated measured disturbances. Because the latter is the special case of the former, we will only present the more general one in section 4. A simulation example is presented in section 5, followed by concluding remarks in section 6. 2. LTI Feedforward and Feedback Control System Performance Assessment The performance assessment of a LTI feedforward and feedback control scheme is the foundation to performance analysis for LTV processes. Desborough and Harris,9 Huang and Shah,10 Kadali et al.,12 and Huang et al.11 have reported that the closed-loop response to both unmeasured and measured disturbances can be modeled as

yt ) Ga(q-1) at +

n

Gi,b(q-1) bi,t ∑ i)1

(1)

where n is the number of measured disturbance variables. yt is the process output variable, and Ga and Gi,b are proper and rational transfer functions. The white noise at is the driving force of the unmeasured disturbances, while bi,t represents the driving force of the ith measured disturbance. The driving force bi,t, in general, is not measurable but can be calculated by applying time-series analysis to the ith measured disturbance variable Di,t to yield an ARMA model:

Ai(q-1) Di,t ) Bi(q-1) bi,t

(2)

Now write eq 1 in an impulse response form (a) (a) yt ) (f(a) 0 at + ... + fd-1 at-(d-1)) + (fd at-d + ...) + n

) (b ) (b ) [(f(b ∑ 0 bi,t + ... + fd-1 bi,t-(d-1)) + (fd bi,t-d + ...)] i)1 i

i

i

(3)

where d represents the time delay of the process including zero-order hold. Equation 3 can be further expressed as the sum of the contribution of the measured and unmeasured disturbances to the actual process variance:

10.1021/ie020956d CCC: $27.50 © 2004 American Chemical Society Published on Web 12/20/2003

590

Ind. Eng. Chem. Res., Vol. 43, No. 2, 2004

yut is the contribution of unmeasured disturbances to n yi,m is the contribution of meathe process while ∑i)1 t sured disturbances to the process, and (a) eut ) f(a) 0 at + ... + fd-1 at-(d-1) a eˆ ut ) f(a) d at + ... + fd+1 at-(d+1) + ...

(5)

(bi) i) ) f(b ei,m t 0 bi,t + ... + fd-1 bi,t-(d-1)

eˆ i,m t

)

i) f(b d

bi,t-d +

(bi) fd+1

bi,t-(d+1) + ...

(6)

The output under minimum variance feedback and feedforward control can be expressed as n

u ymv t ) et +

ei,m ∑ t i)1

(7)

This is the foundation for assessment of the control-loop performance not only for LTI systems but also for LTV systems from routine operating data and will be used throughout this paper. In eq 4, eˆ ut represents the inflation to eut due to nonoptimal feedback control to the unmeasured disturbances and eˆ i,m is the inflation t due to nonoptimal feedforward/feedback control to ei,m t of the measured disturbances. For a more detailed discussion, readers are referred to work by Huang and Shah.10 3. Properties of LTV Feedforward and Feedback Control Loops Because of the fact that most processes have a certain degree of time-varying behavior, it is necessary to develop performance assessment methods for timevariant control loops, particularly for those processes that have a “strong” time-varying property. The definition of the strength of the LTV property will be discussed shortly. For “strong” time-variant systems, care has to be taken with the noncommutativity associated with manipulations of LTV transfer functions. It is shown1,3 that the multiplication or division of two LTV polynomials, u(q-1,t) and v(q-1,t), is noncommutative; that is, u(q-1,t) v(q-1,t) * v(q-1,t) u(q-1,t). For example, depending on the multiplication sequence, multiplication of a time-variant term, say θ(t) q-1, by a backshift operator q-1 will produce two different results; i.e., q-1 × θ(t) q-1 ) θ(t-1) q-2, but θ(t) q-1 × q-1 ) θ(t) q-2. This is important and has to be taken into consideration in calculating the minimum-variance term when, for example, an LTV ARMA model is transferred to an LTV MA model. If the time shift of the coefficient is ignored in the multiplication as for the case of LTI systems, i.e., q-1 × θ(t) q-1 ) θ(t) q-2, then the multiplication is called pointwise multiplication and is commutative. It has been shown that pointwise multiplication may yield theoretically erroneous results especially when the parameter change in the plant or disturbance dynamics is relatively fast.3 Therefore, normal multiplication, which is noncommutative, is recommended when handling LTV transfer functions and should therefore be used in the analysis of variance

Figure 1. Schematic of a time-variant multiple-input singleoutput process under feedforward plus feedback control.

for LTV feedforward and feedback control systems whenever it is possible. To summarize, the basic property of the backshift operator q-1 when it is multiplied by an LTV transfer function, λ(q-1,t), is given as follows:

q-nλ(q-1,t) ) λ(q-1,t-n) q-n qnλ(q-1,t) ) λ(q-1,t+n) qn

(8)

A natural question is, when is this time-varying property “strong” enough so that the conventional pointwise multiplication would break down? This question is dependent on how the calculation is performed. For example, consider that a model has two parameters: one is time-varying and the other is not. If the time-varying parameter does not (or weakly) participate(s) in the calculation, then the time-varying property of the model does not have any (or much) effect. Therefore, this question should be addressed together with the objective of the computation, and in our case, this objective is the minimum-variance term of feedforward and feedback control. Consequently, the “strength” of the time-varying property is determined by how the minimum-variance term is affected by it. We will elaborate the solution to this problem in section 6 after we have discussed how the minimum-variance term is calculated. Now consider the LTV process (see also Figure 1):

˜ (q-1,t) ut + Na(q-1,t) at + Nb(q-1,t) Dt yt ) q-dT

(9)

yt is the measured process output, ut is the manipulated variable, and T ˜ ) K-1(q-1,t) J(q-1,t) is the delay-free plant transfer function represented by coprime fractions. Na ) M-1(q-1,t) L(q-1,t) represents the time-variant transfer function of the unmeasured disturbances represented by coprime fractions, while at is the “driving force” for realization of the unmeasured disturbances. Nb ) X-1(q-1,t) P(q-1,t) is the transfer function from the measured disturbance, Dt, to the output variable, yt, represented by coprime fractions. Note that T ˜ , Na, and Nb given above represent an LTV Box-Jenkins model,13 which is the most general linear dynamic model. For the LTV process in eq 9, the closed-loop response under LTV feedforward and feedback control can be represented by an LTV ARMAX model stated in the following lemma. Lemma 1. The closed-loop transfer function of an LTV Box-Jenkins model shown in Figure 1 under feedback and feedforward control can always be represented by an LTV ARMAX model.

Acl(q-1,t) yt ) Bcl(q-1,t) Dt + Ccl(q-1,t) at

(10)

Ind. Eng. Chem. Res., Vol. 43, No. 2, 2004 591

Equation 10 has the form of an LTV ARMAX model, where Acl(q-1,t), Bcl(q-1,t), and Ccl(q-1,t) are LTV polynomials in the backshift operator q-1. Proof. This follows a proof similar to that of Huang.3 Similarly, a general form of the closed-loop response to more than one measured disturbance can be derived as

Acl(q-1,t) yt )

n

-1 -1 B(i) ∑ cl (q ,t) Di,t + Ccl(q ,t) at i)1

(11)

where

Acl(q-1,t) ) 1 + R1(t) q-1 + R2(t) q-2 + ... + (i) -1 -td(i) -td(i) Rk(t) q-k + ... B(i) + b(i) + cl (q ,t) ) b1 (t) q 2 (t) q -td(i) + k - 1 + ... Ccl(q-1,t) ) 1 + ... + b(i) k (t) q

∂h(i)(q-1,t) Di,t ) η j (i)(q-1,t) bi,t

where bi,t is the driving force for the measured disturbance Di,t. 4. Performance Assessment of LTV Feedforward and Feedback Control Systems 4.1. Analysis of Minimum Variance Contributed by the Unmeasured Disturbance. For the analysis of the LTV feedforward and feedback control scheme, we adopt a second-order ARMAX(2,[nb1,nb2,...,nbn],2, (2) (n) [t(1) d ,td ,...,td ]) model, where nb1 ) nb2 ) ... ) nbn ) 2. This is a representative model applicable to most physical processes and is considered here for the sake of a simpler presentation. This time-series modeling will give a model of the following form:

1 + c1(t) q-1 + c2(t) q-2 + ... + ck(t) q-k + ... (12) and n is the number of measured disturbance variables, yt is the output variable, Di,t is the ith measured disturbance, and, t(i) d is the time delay it takes for a change in Di,t to affect the output variable yt. at is whitenoise sequence with zero mean and variance σ2a(t). Note that eq 11 is the key model for the subsequent performance assessment. It can be estimated from representative routine closed-loop operating data using any time-series analysis package with the measured disturbances Di,t as the inputs and yt as the output. However, because of its time-varying property, a recursive algorithm with the forgetting factor should be used. The forgetting factor should be chosen between 0 and 1, with 0 indicating that only the most current data point is used to estimate the model and 1 indicating that all data are used. The selection of the forgetting factor depends on the rate of time-varying of the parameters and is typically a tradeoff between the variance error and the bias error of the estimations. If the parameter uncertainty is relatively small (or there is small unexplained process noise), then a small forgetting factor would be preferred. A detailed discussion on the selection of the forgetting factor and statistical properties of the recursive algorithm can be found from work by Ljung and Soderstrom.14 From eq 11, the closed-loop response of the process output to the measured and unmeasured disturbances can further be written as

n

yt )

(m) yt,i + y(u) ∑ t i)1

(14)

(m) yt,i represents the closed-loop response to the ith measured disturbance while y(u) t represents the closedloop response to the unmeasured disturbance, and

-1 (m) (i) -1 ) A-1 yt,i cl (q ,t) Bcl (q ,t) Di,t

(15)

-1 -1 -1 y(u) t ) Acl (q ,t) Ccl(q ,t) at

(16)

Time-series analysis may be applied to each measured disturbance, Di,t, to obtain an LTV ARMA model

(17)

A ˆ cl(q-1,t) yt )

n

-1 B ˆ (i) ˆ cl(q-1,t) at ∑ cl (q ,t) Di,t + C i)1

(18)

For the sake of easiness of presentation, we will drop all circumflexes in eq 18 in the sequel. The model can then be transferred to the form of eq 13. The effect of measured and unmeasured disturbances can be calculated from eqs 15 and 16. In this subsection, we shall discuss the calculation of the minimum-variance term contributed by the unmeasured disturbance. For the unmeasured disturbance, we have

Acl(q-1,t) ) 1 + R1(t) q-1 + R2(t) q-2 Ccl(q-1,t) ) 1 + c1(t) q-1 + c2(t) q-2

(19)

Equation 16 can be written in the impulse response form -1 + f2(t) q-2 + f3(t) q-3 + ...)at y(u) t ) (f0(t) + f1(t) q (20)

Substituting eq 20 into eq 16 yields

Acl(q-1,t) [f0(t) + f1(t) q-1 + f2(t) q-2 + f3(t) q-3 + ...]at ) Ccl(q-1,t) at (21) The time-variant impulse response coefficients are estimated by equating the coefficients on the right- and left-hand sides of eq 21:

{

{

f0(t) ) 1 f1(t) ) c1(t) - R1(t) f0(t-1) f2(t) ) c2(t) - R1(t) f1(t-1) - R2(t) f0(t-2) fk(t) ) -R1(t) fk-1(t-1) - R2(t) fk-2(t-2)

k>2 (22)

and from eq 22, it follows that

f0(t) ) 1 f1(t) ) c1(t) - R1(t) f2(t) ) c2(t) - R1(t) c1(t-1) + R1(t) R1(t-1) - R2(t) f3(t) ) -R1(t) c2(t-1) + R1(t) R2(t-1) + R1(t) R1(t-1) c1(t-2) - R1(t) R1(t-1) R1(t-2) R2(t) c1(t-2) + R2(t) R1(t-2) l (23)

592

Ind. Eng. Chem. Res., Vol. 43, No. 2, 2004

The corresponding minimum variance for the LTV process with a time delay of 3 (as an example) can be calculated:

I is an n × n identity matrix. ∂h1 is an n × n matrix that consists of the autoregressive coefficients of the ARMA model in eq 29, which corresponds to q-1, while ∂h2 is an n × n matrix that consists of the autoregressive coefficients corresponding to q-2. Similarly, η j 1 is an n × n matrix that consists of the moving-average coefficients of the ARMA model in eq 29, which corresponds to q-1, while η j 2 is an n × n matrix that consists of the movingaverage coefficients corresponding to q-2. Substituting eqs 31 and 32 into eq 29 yields

σ2mv(t) ) [1 + f21(t) + f22(t)]σ2a(t)

[I + ∂h1(t) q-1 + ∂h2(t) q-2][F(D) 0 (t) Bt +

Therefore, the process response to the unmeasured disturbance can be expressed as -1 y(u) + [c2(t) - R1(t) c1(t-1) + t ) {1 + [c1(t) - R1(t)]q

R1(t) R1(t-1) - R2(t)]q-2 + ...}a (24)

-1 (D) F(D) + ... + Fd-1 (t) q-d+1 + ...]Bt ) 1 (t) q

) {1 + [c1(t) - R1(t)]2 + [c2(t) - R1(t) c1(t-1) +

j 2(t) q-2]Bt (33) [I + η j 1(t) q-1 + η

R1(t) R1(t-1) - R2(t)]2}σ2a(t) (25) However, direct long division gives -1 + [c2(t) - R1(t) c1(t) + y(u) t ) {1 + [c1(t) - R1(t)]q

R1(t)2 - R2(t)]q-2 + [-R1(t) c2(t) + 2R1(t) R2(t) + R1(t)2 c1(t) + R1(t)3 - R2(t) c1(t)]q-3 + ...}at (26) and the corresponding LTI minimum variance is given by

σ2mv(t) ) {1 + [c1(t) - R1(t)]2 + [c2(t) - R1(t) c1(t) + 2

}σ2a(t)

2

R1(t) - R2(t)]

(27)

which is not correct. 4.2. Analysis of Minimum Variance Contributed by the Measured Disturbance. We consider a general case that accounts for cross correlation among the measured disturbances, Di,t. The approach adopted in this section is to treat the measured disturbances as a single identity (a vector form), and the final variance contribution can be distinguished between measured and unmeasured disturbances but cannot be traced to individual measured disturbances because of the correlation of the measured disturbances.9 Let Dt be the vector of n measured disturbance variables, which can be represented as

Dt ) [D1,t, D2,t, ..., Dn,t]

T

(28)

{

(D) (D) The coefficients F(D) 0 (t), F1 (t), ..., Fk (t) are calculated from eq 33 and are given by

F(D) 0 (t) ) I F(D) j 1(t) - ∂h1(t) F(D) 1 (t) ) η 0 (t-1)

F(D) j 2(t) - ∂h1(t) F(D) h2(t) F(D) 2 (t) ) η 1 (t-1) - ∂ 0 (t-2)

h1(t) F(D) h2(t) F(D) F(D) 3 (t) ) -∂ 2 (t-1) - ∂ 1 (t-2) l (D) (D) h1(t) Fk-1 (t-1) - ∂h2(t) Fk-2 (t-2) F(D) k (t) ) -∂

{

k g3 (34)

It follows from eq 34 that

F(D) 0 (t) ) I F(D) j 1(t) - ∂h1(t) 1 (t) ) η

j 2(t) - ∂h1(t) η j 1(t-1) + ∂h1(t) ∂h1(t-1) - ∂h2(t) F(D) 2 (t) ) η

F(D) h1(t) η j 2(t-1) + ∂h1(t) ∂h1(t-1) η j 1(t-2) + 3 (t) ) -∂ ∂h1(t) ∂h2(t-1) - ∂h1(t) ∂h1(t-1) ∂h1(t-2) j 1(t-2) + ∂h2(t) ∂h1(t-2) ∂h2(t) η l (35)

Equation 11 can be expressed as -1 (2) -1 Acl(q-1,t) yt ) B(1) cl (q ,t) D1,t + Bcl (q ,t) D2,t + ... + -1 -1 B(n) cl (q ,t) Dn,t + Ccl(q ,t) at (36)

and Dt can be modeled by

∂ˆ (q-1,t) Dt ) ηˆ (q-1,t) Bt

(29)

Bt is the vector of driving force bi,t for the measured disturbances, which can be expressed as

Equation 36 can be rearranged to give

[]

Dt can be expressed in impulse response form as

-1 (2) -1 (n) -1 Acl(q-1,t) yt ) [B(1) cl (q ,t), Bcl (q ,t), ..., Bcl (q ,t)] D1,t D2,t + Ccl(q-1,t) at (37) l Dn,t

(D) Dt ) F(D) 0 (t) Bt + F1 (t) Bt-1 + ... +

Therefore

Bt ) [b1,t, b2,t, ..., bn,t]T

(D) Fd-1 (t)

(30)

Bt-d+1 + ... (31)

Consider eq 29 to be a vector ARMA(2,2) process, i.e.

{

∂ˆ (q-1,t) ) I + ∂h1(t) q-1 + ∂h2(t) q-2 ηˆ (q-1,t) ) I + η j 1(t) q-1 + η j 2(t) q-2

(32)

Ind. Eng. Chem. Res., Vol. 43, No. 2, 2004 593

Equation 38 can be expressed as

yt ) y(m) + y(u) t t

(39)

where -1 (1) -1 (2) -1 (n) -1 ) A-1 y(m) t cl (q ,t) [Bcl (q ,t), Bcl (q ,t), ..., Bcl (q ,t)]

[D1,t, D2,t, ..., Dn,t]T (40)

[

Figure 2. Jacketted stirred-tank heater.

that is

if θ1 > 0, then the first θ1 coefficient(s) of the moving average model is (are) 1 × n zero row vector(s). That is,

F(m) 0 (t) ) [0, 0, ...] l (t) ) [0, 0, ...] Fθ(m) 1-1

Substituting eq 12 into eq 41 yields -td(1) [1 + a1(t) q-1 + a2(t) q-2]y(m) ) [lb(1) + b(1) t 1 (t) q 2

(t) q-(td

(1)+1)

Let

-td -(td l ... l b(n) + b(n) 1 (t) q 2 (t) q

(n)+1)

(n)

l]Dt (42)

{

(i) θ ) maximum(t(i) d ) - minimum(td ) + 2 and θ1 ) minimum(t(i) d ) (for i ) 1, 2, ... n)

(43)

Then, eq 42 can be rearranged and expressed as

(D) (t) ) B(cl) Fθ(m) 1 (t) F0 (t-θ1) 1

(t) ) -a1(t) Fθ(m) (t-1) - a2(t) Fθ(m) (t-2) + Fθ(m) 1+1 1 1-1 (D) (cl) (D) B(cl) 2 (t) F0 (t-θ1-1) + B1 (t) F1 (t-θ1) (m) (m) F(m) k (t) ) -a1(t) Fk-1(t-1) - a2(t) Fk-2(t-2) + (cl) (D) Bθ-j (t)Fk-θ (t-θ1-θ+j+1) 1-θ+j+1 for j ) 0, 1, 2, ..., θ - 1, k > θ1

The minimum feedforward variance can be calculated from eq 47: T

-θ1 [1 + a1(t) q-1 + a2(t) q-2]y(m) ) [B(cl) + t 1 (t) q -(θ1+1) -(θ1+θ-1) B(cl) + ... + B(cl) ]Dt (44) 2 (t) q θ (t) q

where θ represents the total number of terms on the right-hand side of eq 44. B(cl) 1 (t) is a 1 × n row vector that consists of coefficients on the right-hand side (rhs) of eq 42, which corresponds to q-θ1. Similarly, B(cl) 2 (t) is a 1 × n row vector that consists of coefficients on the rhs of eq 42, which corresponds to q-(θ1+1) until B(cl) θ (t), which is also a 1 × n row vector that consists of coefficients on the rhs of eq 42 corresponding to q-(θ1+θ-1). y(m) can be expressed in impulse response form as t

y(m) t

)

F(m) 0 (t)

Bt +

F(m) 1 (t)

Bt-1 + ... + (m) Fd-1 (t) Bt-d+1 + ... (45)

Substituting eqs 31 and 45 into eq 44 yields (m) [1 + a1(t) q-1 + a2(t) q-2][F(m) 0 (t) Bt + F1 (t) Bt-1 + -θ1 (m) ... + Fd-1 (t) Bt-d+1 + ...] ) [B(cl) + 1 (t) q -(θ1+1) -(θ1+θ-1) + ... + B(cl) ][F(D) B(cl) 2 (t) q θ (t) q 0 (t) +

F(D) 1 (t)

q

-1

+

F(D) 2 (t)

-2

q

+ ...]Bt (46)

The impulse response coefficients (which are 1 × n row vectors) are obtained by equating the terms on the leftand right-hand sides. However, it should be noted that

(47)

T

(m) (m) (m) (m) σ2(m) mv (t) ) F0 (t) ΣBF0 (t) + F1 (t) ΣBF1 (t) + ... + T

(m) (m) (t) ΣB Fd-1 (t) (48) Fd-1

where

[]

b1,t b2,t ΣB ) cov l bn,t

(49)

and d is the time delay of the process including the delay due to a zero-order-hold device. 5. Simulation Let us consider a stirred-tank heater shown in Figure 2. When the work done by the impeller is neglected, energy balance around the tank is used to obtain the model given by

dT F Q ) (Ti - T) + dt V VFCp

(50)

where the rate of heat transfer from the jacket to the tank, Q, is given by

Q ) UA(Tj - T)

(51)

T is the tank temperature, F the volumetric flow rate, F the density, Cp the heat capacity, U the overall heattransfer coefficient, and A the area for heat transfer. The subscripts i, j, and ji denote inlet, jacket, and jacket inlet, respectively.

594

Ind. Eng. Chem. Res., Vol. 43, No. 2, 2004

In linearization of the nonlinear model in eq 50, it is assumed that the tank outlet flow rate, F, and the tank temperature, T, are the manipulated and controlled variables, respectively. The overall heat-transfer coefficient U is time-varying. In this example, the tank inlet temperature, Ti, is considered as the disturbance affecting the system, which is measured and included as the measured disturbance. The steady state is obtained by solving the dynamic equation for dT/dt ) 0. The steady-state values of the system variables and some parameters for the process are given in Bequette.15 3

3

Fs ) 1.0 ft /min, FCp ) 61.3 Btu/°F‚ft , V ) 10 ft 3, Tis ) 50 °F, Ts ) 125 °F Fjs ) 1.5 ft 3/min, FjCpj ) 61.3 Btu/°F‚ft 3, Vj ) 1 ft 3, Tjis ) 200 °F, Tjs ) 150 °F Let the temperature measurement have a time delay of 4 and a sampling time of 1. Consider that UA is timevarying and is given by UA ) 183.9[1 + 0.5 sin(t/x)] for illustration purposes, where x/2π is the oscillation period. Let us consider that the ambient temperature, Ta, is also measured. The time delay it takes for this measured disturbance to affect the tank temperature is assumed to be three samples. The measured ambient temperature is considered to be correlated with the measured tank inlet temperature, Ti. It is also assumed that a random-walk disturbance is the unmeasured disturbance added to the system and the white noise sequence, at, is its driving force. Through linearization and discretization, T(t) can be finally expressed as

Hence, the final process model can be expressed in a standard form as

(

)

q

... q-1

(

)(

{

To illustrate the effect of time-varying parameters on the performance assessment of LTV processes, this example considers that the disturbance has three different time-variant dynamics in ascending order of increasing parameter-varying rate, from relatively slow parameter change to relatively fast parameter change as given by

case 1: υ(t) ) (0.1 - 0.067e-0.15 sin (t/10))/[0.4 + 0.15 sin(t/10)] δ(t) ) 0.67e-0.15 sin(t/10)

(0.1 - 0.067e

)/[0.4 + 0.15 sin(t/x)] Ti q 0.1 (t) + at (52) 1 - q-1

-[0.15 sin(t/x)] -1

1 - 0.67e

Further assume that the measured disturbance, Ti, can be expressed as

Ti(t) )

1 b1,t 1 - 0.9q-1

(53)

where b1,t is white noise representing the driving force for realization of the measured disturbance, Ti. The ambient temperature is given by

Ta(t) )

1 b2,t 1 - 0.9q-1

(54)

where b2,t is white noise which is correlated to b1,t by

b2,t ) 0.5b1,t + 0.5e e in eq 55 represents white noise.

(55)

(58)

case 2: υ(t) ) (0.1 - 0.067e-0.15 sin(t))/[0.4 + 0.15 sin(t)] δ(t) ) 0.67e-0.15 sin(t)

-0.15 sin(t/x)

)

φ(t) ) (5.03e-0.15 sin(t/x) - 7.5)/[0.4 + 0.15 sin(t/x)] υ(t) ) (0.1 - 0.067e-0.15 sin(t/x))/[0.4 + 0.15 sin(t/x)] δ(t) ) 0.67e-0.15 sin(t/x) (57)

-0.15 sin(t/x)

- 7.5)/[0.4 + 0.15 sin(t/x)] F(t) + 1 - 0.67e-[0.15 sin(t/x)]q-1 0.02 q-3 Ta(t) + 1 - 0.8q-1

)

where the time-variant process and disturbance dynamics are given by

T(t) ) -4(5.03e

(

φ(t) υ(t) yt ) q-4 u + q-1 -1 t 1 - δ(t) q 1 - δ(t) q-1 1 0.02 1 b + q-3 b2,t + -1 1,t -1 1 - 0.9q 1 - 0.8q 1 - 0.9q-1 0.1 at (56) ... 1 - q-1

(59)

case 3: υ(t) ) (0.1 - 0.067e-0.15 sin(2t))/[0.4 + 0.15 sin(2t)] δ(t) ) 0.67e-0.15 sin(2t)

(60)

The feedforward and feedback minimum-variance control benchmark is calculated using the developed LTV algorithm for correlated disturbances and compared with the conventional minimum-ariance control benchmark that is calculated using pointwise multiplication for the three different time-variant dynamics given in eqs 58-60. The simulation results for the correlated disturbances in Figure 3 show a comparison of the difference between the developed algorithm and the conventional method. It can be seen that the difference between the two methods increases as the parameter-varying rate increases from the top to the bottom subplots. This result shows that it is important to consider noncommutativity associated with the manipulation of LTV operators in calculating the LTV feedforward and feedback minimumvariance control benchmark to evaluate the performance of LTV feedforward and feedback control systems. 6. Assessment of the Time-Varying Effect Industrial data can have a variety of different timevarying properties, which may or may not affect the conventional performance calculation significantly. It is

Ind. Eng. Chem. Res., Vol. 43, No. 2, 2004 595

Figure 3. Comparison of the time-variant feedforward/feedback minimum-variance terms using normal multiplication and pointwise multiplication for correlated disturbances.

often found that the LTI assessment method may suffice.16 In this section, we will address the question on when the conventional pointwise multiplication/ division will break down for the calculation of the minimum feedforward and feedback term using routine

Figure 4. Plot of the operating data for the adaptive controller.

Figure 5. Plot of the measured disturbance.

operating data. Given a time series, the following are steps to take in order to assess how much the timevarying effect is. (i) Fit the time series by a recursive second-order ARMAX model. A second-order ARMAX model would be sufficient for such a preliminary analysis purpose. (ii) If the fitted ARMAX model has coefficients that are all sufficiently close to their previous ones, the conventional pointwise multiplication should suffice. The closeness of the coefficients can be tested by checking whether there is an overlap of their confidence intervals. The confidence intervals are provided by most of the time-series analysis packages. (iii) Otherwise, one would carry on the analysis using both pointwise multiplication and normal multiplication as discussed in section 4. Note that the calculation according to the pointwise multiplication is relatively minor and the calculation according to the normal multiplication for second-order ARMAX model has a standard procedure and has been discussed in the previous sections. (iv) The relative difference (RD) between the feedforward and feedback minimum-variance term using normal multiplication (noncommutative) and pointwise multiplication (commutative) is defined as follows:

RD )

|NM result - PM result| |PM result|

596

Ind. Eng. Chem. Res., Vol. 43, No. 2, 2004

Figure 6. RD between LTV feedforward/feedback minimum-variance terms using normal multiplication and point multiplication.

where NM stands for normal multiplication and PM stands for pointwise multiplication. One can therefore, according to the RD, determine whether the normal multiplication should be adopted for the remaining performance monitoring. To illustrate the methodology just discussed, an industrial data set has been taken over a 5-h period with a sampling interval of 1 min. The time delay of the process is approximately 2 min including the delay due to the zero-order-hold device. The measured disturbance, F1, is considered in this analysis with the time delay of td ) 1. The operating data are shown in Figure 4, while Figure 5 is the plot of the measured disturbance. The data set used in this example clearly reveals that there is significant time-varying behavior occurring over a certain limited time period (approximately 8% time) within the given data, and the confidence intervals calculated from the recursive algorithm confirm this observation. Figure 6 shows the RD between the minimum-variance terms using NM and PM. One can see that the RD between the LTV performance method and the traditional recursive method can be as high as 50%, indicating breakdown of the PM. The result shows that it is important to take noncommutativity into account in the manipulation of LTV transfer functions for this loop.

(10) Huang, B.; Shah, S. L. Control Loop Performance Assessment: Theory and Applications; Springer-Verlag: London, 1999.

7. Conclusion

(11) Huang, B.; Shah, S. L.; Miller, R. Feedforward plus Feedback Controller Performance Assessment of MIMO systems. IEEE Trans. Control Syst. Technol. 2000, 8 (3), 580-587.

The LTV feedforward and feedback minimum-variance control benchmark for time-variant processes has been derived in this paper. This result is a significant extension of the recent results of Huang3 from LTV feedback control loops to LTV feedforward and feedback control systems. The proposed method provides a way to calculate the minimum-variance control benchmark from routine operating data and has been illustrated through a simulated example and application to a set of industrial data. Acknowledgment The financial support from NSERC and Syncrude Canada Ltd. for this work is greatly acknowledged. B.H. also acknowledges the support from Alexander von Humboldt Foundation during the final preparation of this manuscript.

Literature Cited (1) Li, Z.; Evans, R. J. Minimum-Variance Control of Linear Time-Varying Systems. Automatica 1997, 33 (8), 1531-1537. (2) Huang, B. Performance Assessment of Processes with Abrupt Changes of Disturbances. Can. J. Chem. Eng. 1999, 77 (5), 1044-1054. (3) Huang, B. Minimum Variance Control and Performance Assessment of Time-Variant Processes. J. Process Control 2002, 12, 707-719. (4) Harris, T. Assessment of Closed Loop Performance. Can. J. Chem. Eng. 1989, 67, 856-861. (5) Astrom, K. J. Introduction to Stochastic Control Theory; Academic Press: New York, 1970. (6) Box, G. E. P.; Jenkins, G. M. Time Series Analysis: Forecasting and Control; Holden-Day: 1976. (7) Sternad, M.; Soderstrom, T. LQG-optimal Feedforward Regulators. Automatica 1988, 24 (4), 557-561. (8) Stanfelj, N.; Marlin, T.; MacGregor, J. Monitoring and Diagnosing Process Control Performance: The Single Loop Case. Ind. Eng. Chem. Res. 1993, 32 (2), 301-314. (9) Desborough, L.; Harris, T. Performance Assessment Measures for Univariate Feedforward/Feedback Control. Can. J. Chem. Eng. 1993, 71 (4), 605-616.

(12) Kadali, R.; Huang, B.; Tamayo, C. A Case Study on Performance Analysis and Trouble Shooting of an Industrial Model Predictive Control System. Proceedings of the American Control Conference, San Diego, CA, 1999; pp 642-646. (13) Ljung, L. System IdentificationsTheory For the User; PTR Prentice Hall: Englewood Cliffs, NJ, 1999. (14) Ljung, L.; Soderstrom, T. Theory and Practice of Recursive Identification; MIT Press: London, England, 1983. (15) Bequette, B. W. Process Dynamics, Modeling, Analysis and Simulation; Prentice Hall: Englewood Cliffs, NJ, 1998. (16) Desborough, L. D.; Miller, R. M. Increasing Customer Value of Industrial Control Performance MonitoringsHoneywell’s Experience. Proceedings of CPC VI, Tucson, AZ, 2001.

Received for review December 3, 2002 Revised manuscript received May 8, 2003 Accepted November 11, 2003 IE020956D