Diagnosis of Poor Control Loop Performance Due to Model−Plant

Mar 26, 2010 - Gp/Gm, called the plant model ratio (PMR) in the frequency domain is introduced ... A method to estimate the PMR from closed-loop opera...
0 downloads 0 Views 4MB Size
4210

Ind. Eng. Chem. Res. 2010, 49, 4210–4229

Diagnosis of Poor Control Loop Performance Due to Model-Plant Mismatch S. Selvanathan‡ and A. K. Tangirala*,† Department of Chemical Engineering, IIT Madras, Chennai, India, and Department of Engineering Cybernetics, NTNU, Norway

Model-plant mismatch is a crucial factor that influences the performance of a closed-loop system. A new approach is proposed for the diagnosis of poor control loop performance due to model-plant mismatch (MPM) in the internal model control framework. In particular, the objective is to identify the mismatch in specific components of a transfer function model, namely, the gain, time-constant, and delay from closed-loop operating data. A new quantity Gp/Gm, called the plant model ratio (PMR) in the frequency domain is introduced as a measure of model-plant mismatch. It is shown that there exists a unique signature in PMR for each combination of mismatch in model parameters, which is the key step in the proposed method. A method to estimate the PMR from closed-loop operating data is provided. Theoretical and practical aspects of the mapping between the type of MPM and the proposed PMR are presented. Simulation studies demonstrate the effectiveness of the proposed method. 1. Introduction Monitoring and assessment of controller performance has evoked considerable interest among academicians and practitioners in the area of process control and monitoring for two decades now.1-4 The incentive for finding solutions to associated problems is immense since poor performance affects product quality, plant economy, and safety. In addition, there is a constant drive to improve the performance of existing control schemes and to optimize the overall plant performance. However, one can note that controllers often fail to operate according to their design specifications and, in many cases, they even increase the process variability, as was reported by Ender.5 Most modern industrial plants have hundreds or even thousands of automatic control loops. These loops can be simple proportionalintegral-derivative (PID) or more sophisticated model-based linear and nonlinear control loops. It has been reported that as many as 60% of all industrial controllers have performance problems.5 Having an automated means of detecting when a loop is not performing well and then diagnosing the root cause plays a vital role in addressing the problems mentioned at the outset of this paper. A variety of techniques for closed-loop identification have been developed by many researchers. These techniques can be categorized into three groups: the direct approach, the indirect approach, and the joint input-output approach.6,7 Controllers are implemented with various design objectives, including set point tracking, disturbance rejection, constraint handling, and surge attenuation.8 Bialkowski9 and Kozub and Garcia10 pointed out in their work that the major causes of poor control loop performance are (i) improper controller tuning, (ii) poor hardware (sensors, actuators) maintenance, (iii) valve stiction, (iv) model-plant mismatch (MPM), and (v) stochastic disturbances. Ha¨gglund11 developed a procedure for the automatic detection of sluggish control loops obtained from conservatively tuned controllers. A common source of oscillation is a limit cycle caused by a control valve with a deadband or excessive static friction. A process variable oscillating for that reason can readily propagate the oscillation to other variables and disturb other control loops, hence causing a plant-wide disturbance. Several authors have * To whom correspondence should be addressed: E-mail: arunkt@ iitm.ac.in. † IIT Madras. ‡ NTNU.

addressed the detection of oscillatory measurements in process data.4,12-15 There is a wealth of information in the literature on methods to diagnose the poor control loop performance due to valve stiction, improper tuning of controllers, and sensor faults. Relatively, there are very few methods for detecting and diagnosing MPM as the root cause of poor control loop performance. Even more, the existing methods do not comprehensively address the associated issues. In plants running on advanced control schemes, MPM becomes a governing factor on the performance of the control loops. Models and their uncertainties expressed in the frequency domain have regained attention in the past few years. This is mainly due to the development of robust control theory and a clear interpretation of most identification algorithms in the frequency domain.6 A multivariable chi-squared test on the output and prediction errors has been suggested to help diagnosis MPM in Kesavan and Lee.16 In Kendra and Cinar17 and Kendra et al.18 the multivariable complementary sensitivity functions are identified in the frequency domain, and are compared to their design specifications. Any significant differences would imply the presence of MPM. A prediction index, which compares the minimum achievable prediction error variance to the actual variance, has been proposed in Kozub.19 A benchmark specific to model predictive controllers was presented in Patwardhan and Shah,20 which compares the achieved and designed objective functions. The benchmark can be significantly affected by changes in process or unmeasured disturbance characteristics. Hence, it does not distinguish between process mismatch and disturbance model mismatch. Since it is common practice to use a simple step disturbance model in the controller design, this method can always indicate mismatch even if the process models are perfect (unless the disturbance is actually a step disturbance). The above methods indicate the presence of MPM, but do not attempt to specify which subset of models is mismatched and requires reidentification. In Kammer et al.,21 the control loop is opened so that a disturbance model may be identified. The loop is reclosed, and a time series model is fit between the white noise driving the disturbance and the output error. If the time series model for the output error is equal to the disturbance model then no MPM is present. Any difference between the two models would result from MPM or changing disturbance characteristics. This method has the advantage that it does not require any external excitation, but it is sensitive to

10.1021/ie900769v  2010 American Chemical Society Published on Web 03/26/2010

Ind. Eng. Chem. Res., Vol. 49, No. 9, 2010

changing disturbances, and the control loop must be initially opened to identify the disturbance model. Certain subsets of the inputs can be held constant to determine which specific inputs are mismatched.21 A detection method based on a twomodel divergence algorithm was proposed for the validation of dynamic models under closed-loop conditions by Huang.22 Jiang et al.23 proposed a new scheme to detect and isolate MPM for multivariate dynamic systems. In their work, the MPM problem is formulated in the state-space domain, as is widely done in the design and implementation of model-predictive controllers (MPC). The specific issue addressed therein was to identify which among the state-space matrices had to be re-estimated in order to account for significant plant deviations from its nominal state. Three MPM detection indices (MDIs) were proposed to detect the MPM for that purpose. A shortcoming of their work is that changes in state-space matrices cannot be directly translated to changes in gain, time-constant, or delay of the process. In addition, delay mismatches become difficult to detect with a state-space representation since such mismatches cause either an increase or decrease in the order of the system depending on an increase or decrease in delay of the process. The input-output representation, on the other hand, provides a suitable framework for directly attributing poor loop performance to changes in process characteristics such as gain, time-constant, and delay. To the best knowledge of the authors, there is no method in the literature that can diagnose and quantify the effects of gain, timeconstant, and delay mismatches on control loop performance from the closed-loop routine operating data. It may be noted that the effect of MPM on controller design and control loop performance is an evolved area of study and implementation under robust control. However, there is hardly any evidence of an analysis of the direct impact of process characteristics on loop performance in model-based control loops. In this article, the diagnosis of poor control loop performance due to model-plant mismatch in control loops is addressed. The emphasis of this work is on providing measures to detect changes in process characteristics that lead to performance degradation. The problem setting is in the IMC framework. The diagnosis methodology is developed in the frequency-domain. A new measure, plant-model ratio (PMR), Gp(jω)/Gm(jω) is introduced with the purpose of quantifying MPM. Conventionally MPM has been quantified as ∆G ) Gp - Gm where Gp is process transfer function and Gm is model transfer function. While this definition of MPM (and its norms) is useful for robust control design, its utility in the problem of interest in this work is limited. On the other hand, the PMR allows one to obtain a unique mapping between its signatures and changes in process characteristics. On the basis of this fact, it is proposed to detect signatures in PMR which can be then uniquely identified with changes in gain, time-constant, and/or time-delay. The ratio measure can be estimated from closed-loop operating data. The estimation method uses the cross-spectral analysis of set-point, plant output, and model output with the assumption that the set-point contains at least a pulse change. The foregoing assumption is satisfied to a large extent in industrial loops given that the set-points experience changes due to changing market demands and productivity constraints. The emphasis of this work is on single-input, single-output (SISO) systems. The rest of the paper is organized as follows. Section 2 begins with the problem formulation followed by a presentation of the concept of the proposed method to diagnose the model-plant mismatch. Section 3 discusses the procedure for estimating the PMR which is a key step in the proposed method. Results from simulation studies to demonstrate the proof-of-concept are

4211

Figure 1. Schematic representation of internal model control.

discussed in section 4. A few important remarks and recommendations on extensions of this work are given in section 5. Concluding remarks appear in section 6. 2. Problem Statement and Proposed Methodology The task of diagnosing poorly performing controllers is a challenging one. The major causes of poor control loop performance in model-based control loops are (i) improper controller tuning, (ii) actuator/process nonlinearities, (iii) stochastic disturbances, (iv) sensor faults, and (v) MPM. Keeping inline with the focus of this work, the last of the aforementioned factors is taken up for study. The question of interest is that given a poor performance of a (model-based) control loop and that it is solely due to MPM, then which among the process characteristics, namely, gain, time-constant, and delay have undergone a significant change? Although the probe begins with the assumption of a unstructured mismatch, the method developed here can be applied to situations with structure mismatch as well. In the latter case, the question will then translate to the following: Given poor performance of a (model-based) control loop and that it is solely due to MPM, then which parts of the model, namely, gain, time-constant(s), and delay require an update? This task is of immense value in industry and arises in the exercise of updating the model. The problem setting is in the IMC framework. The idea is that the results here can be extended to the more commonly applied MPC-based control schemes. The internal model control (IMC) framework provides an elegant way of reparametrizing the conventional feedback controller.24 An increasing number of practitioners have started using the IMC tuning rules for designing PID controllers. A primary benefit of using IMC-type schemes in the context of performance assessment is that in situations where the performance degradation is due to changing plant conditions, the presence of the model online enables us to obtain an estimate of such deviations from routine operating data. The key point is that the model predictions are also available along with the plant output, from where one can obtain an estimate of the MPM. It is not to be forgotten that the difference between plant output and the model prediction also contains the effect of disturbances. Figure 1 shows a schematic of IMC configuration where Gp(z-1), Gm(z-1), Q(z-1), and Gd(z-1) denote process, model, controller, and disturbance transfer functions, respectively. The observed output is denoted by yp[k], while the prediction of the model is denoted by ym[k] respectively. The problem statement and methodology is presented for discrete-time systems. However, these results are equally valid for continuous time IMC systems as well with some appropriate minor changes. It is reiterated here that any performance degradation is attributed

4212

Ind. Eng. Chem. Res., Vol. 49, No. 9, 2010

to MPM, that is, there is no loss in performance due to valve stiction, stochastic disturbances, improper tuning of controllers, and sensor faults. Traditionally, MPM has been characterized by ∆G ) Gp - Gm, which can never be zero because of the lack of accurate process knowledge (modeling uncertainties). The uncertainties can be either in the structure, whence then ∆G represents a structural mismatch, or in the parameters, namely, gain, time-constant(s), and delay, in which case ∆G is said to represent unstructured mismatch. We consider primarily the latter situation in the treatment of the problem under study. To begin with, replacing the conventional definition of mismatch with a measure, known as the plant-model ratio (PMR), based on the frequency response functions (FRF) of the plant and model, is proposed: GPMR(e ) )

Gp(ejω)



(1)

Gm(ejω)

The quantity GPMR(ejω) is clearly complex-valued. The measure defined in eq 1 has a few advantages over the traditional ∆G(ejω). The key advantage is understood by first writing a polar representation for GPMR, j

GPMR(ejω) )

|Gp(ejω)|ej∠Gp(e )e-jDpω jω

j m(ejω) -jDmω j∠G

|Gm(e )|e

e



∆P(ω)

M(ω)|ω)0 ) 1 M(ω)|ω)0 * M(ω)|ω*0

∆P(ω)|ω)0, π ) 0 |∆P(ω)|max > 0

The discrete domain transfer functions of FOPTD plant and FOPTD model are expressed as Gp(z-1) )

Gm(z-1) )

Kp(1 - e-Ts/τp) -1 -Ts/τp

1-z e

Km(1 - e-Ts/τm) -1 -Ts/τm

1-z e

z-Dp )

-

|Gm(e )|e

e

) A(ω)ejB(ω) (3)

It is clear from the above equation that the quantities A(ω) and B(ω) are both affected by mismatches in delay and timeconstant(s). Moreover, the manner in which they affect the amplitude and phase of ∆G(ejω) is not straightforward. Compare this situation with the one arising out of the proposed measure. It is clear that GPMR(ω) is a better representation of the mismatch than ∆G(ω) in the context of diagnosis of loop performance. In the presentation to follow, we show that a unique mapping exists between changes in process characteristics and the amplitude and phase of GPMR(ω). The theoretical results are presented for first-order plus time delay (FOPTD) processes and models for the sake of simplicity and brevity. The ideas, nevertheless, as shown later, turn out to be valid for all higherorder processes in a straightforward way so long as there is no mismatch in the structure.

1 - apz-1

z-Dp

Km(1 - am) -1

1 - amz

(4)

z-Dm (5)

where Ts is sampling time (Ts ) 1 throughout), Kp, τp, and Dp are process gain, process time-constant, and process time delay (in number of samples) and corresponding model parameters are Km, τm, and Dm. 2.1. Gain Mismatch. When the loss in performance is only due to gain mismatch, that is, Kp * Km and τp ) τm, Dp ) Dm, eq 1 simplifies to GPMR(ω) )

) M(ω)ej(∆P(ω)) (2)

j m(ejω) -jDmω j∠G

Kp(1 - ap)

z-Dm )

e

j p(ejω) -jDpω j∠G

∆G(e ) ) |Gp(e )|e jω

M(ω)



j p and G j m represent the delay-free parts of the plant where G and model transfer functions, respectively. Equation 2 is central to the proposed method from where conditions to diagnose the type(s) of mismatch are formulated in the later part of this section. From eq 2, M(ω) contains the information on mismatch in the magnitude of the FRFs of the plant and model, while ∆P(ω) contains the mismatch in the phase of the delay-free parts and the delays. Equivalently, in the absence of any mismatch between the plant and the model, GPMR(ω) ) 1, implying M(ω) ) 1 and ∆P(ω) ) 0. An immediate utility of the above representation is that delay mismatch, denoted by ∆D ) Dp - Dm only influences ∆P(ω), but not M(ω). Similarly, any mismatch in the gain is contained in M(ω) only, while timeconstant mismatches affect both quantities. Thus, the delay mismatch problem is decoupled from the time-constant mismatch. On the other hand, writing a polar form of representation for the traditional definition of mismatch yields the following jω

Table 1. Theoretical Conditions Derived for the Diagnosis of Time-Constant Mismatch

Kp Km

(6)

From eq 6, we can formulate the following condition for diagnosing the gain mismatch M(ω) * 1∀ω

(7)

∆P(ω) ) 0∀ω

(8)

The presence of gain mismatch in the control loop is diagnosed if the estimates of the magnitude ratio and phase difference satisfy eqs 7-8. The right-hand side constants of these equations are constants (independent of frequency) for gain mismatch. It may be noted that this result is true regardless of the order of the process/model and in the absence of any structural mismatch. Thus, this is a general result for all situations with unstructured uncertainties. Also, to determine the nature of gain mismatch, one further notes that

2.2. Time-Constant Mismatch. In the presence of timeconstant mismatch, that is, τp * τm and Kp ) Km, Dp ) Dm, eq 1 reduces to GPMR(ω) )

(1 - ap)(1 - e-jωam) (1 - am)(1 - e-jωap)

) M(ω)ej∆P(ω) (10)

The expressions for M(ω) and ∆P(ω) then follow: (1 - ap) M(ω) ) (1 - am) ∆P(ω) ) tan-1

[



1 - 2am cos ω + am2 1 - 2ap cos ω + ap2

]

[

am sin ω ap sin ω - tan-1 1 - am cos ω 1 - ap cos ω

(11)

]

(12)

The expressions in eqs 11 and 12 at a first glance appear complicated; however, they lead to a set of simple necessary and sufficient conditions listed in Table 1 for the time-constant

Ind. Eng. Chem. Res., Vol. 49, No. 9, 2010 Table 2. Theoretical Conditions Derived for the Diagnosis of Gain-Time-Constant Mismatch

4213

Table 4. Theoretical Conditions Derived for the Diagnosis of Gain, Time-Constant, and Delay Mismatches

M(ω)

∆P(ω)

M(ω)

∆P(ω)

M(ω)|ω)0 * 1 M(ω)|ω)0 * M(ω)ω)0

∆P(ω)|ω)0, π ) 0 |∆P(ω)|max > 0

M(ω)|ω)0 * 1 M(ω)|ω)0 * M(ω)ω*0

∆P(ω)| ) f(τ, ω) + Rω, f(τ, ω) f 0 for ω f π

and time-constant). Similar situations arise in other combinations of parametric mismatches as seen in what follows.

Table 3. Theoretical Conditions Derived for the Diagnosis of Time-Constant-Delay Mismatch M(ω)

∆P(ω)

M(ω)|ω)0 ) 1 M(ω)|ω)0 * M(ω)|ω*0

∆P(ω) ) f(τ, ω) + Rω, f(τ, ω) f 0 for ω f π

mismatch to be the cause of overall MPM. The presence of time-constant mismatch is diagnosed from either M(ω) or ∆P(ω) as provided in Table 1. Once the source of mismatch has been identified as due to mismatch in time-constants, one can then determine the nature of deviation of τp from τm using the condition (can be derived from eq 11)

2.3. Delay Mismatch. If the only source of mismatch is the mismatch in delay, that is, Dp * Dm and τp ) τm, Kp ) Km, eq 1 simplifies to GPMR(ω) ) e-(Dp-Dm)jω

M(ω) ) 1∀ω

(14)

∆P(ω) ) Rω∀ω

(15)

where R ) ∆D ) Dm - Dp. Equations 14 and 15 show that delay mismatches in the model-based control loops do not cause deviations in the magnitudes of the FRF of the process and model but introduce a phase difference in the FRFs that varies linearly with respect to the frequency. Since Dm and Dp represent the sample delays, the minimum value of |R| is unity for any mismatch in delays to manifest in the sampled data. Therefore, it is appropriate to impose |R|min ) 1 in addition to the conditions specified in eqs 14 and 15. We now turn to establishing conditions for various other combinations. 2.4. Gain and Time-Constant Mismatch. In this case, τp * τm, Kp * Km, and Dp ) Dm. Equation 1 then reduces to Kp (1 - ap)(1 - e-jωam) GPMR(ω) ) Km (1 - a )(1 - e-jωa ) m p

(16)

M(ω) and ∆P(ω) are now written as follows, Kp (1 - ap) Km (1 - am)

∆P(ω) ) tan-1

[



1 - 2 cos ωam + am2 1 - 2 cos ωap + ap2

]

GPMR(ω) )

Kp -(Dp-Dm)jω e Km

(19)

From eq 19 we can formulate the following condition for diagnosing the gain-and-delay mismatch M(ω) * 1∀ω ∆P(ω) ) Rω∀ω

(20)

where R ) ∆D ) Dm - Dp. To determine the nature of gain and delay mismatch, one further notes that

(13)

leading to the following condition for the diagnosis of delay mismatch:

M(ω) )

2.5. Gain and Delay Mismatch. This situation arises when Kp * Km, Dp * Dm, and τp ) τm, reducing eq 1 to

[

am sin ω ap sin ω - tan-1 1 - am cos ω 1 - ap cos ω

(17)

]

2.6. Time-constant and delay mismatch. With the above assumption, τp * τm, Dp * Dm, and Kp ) Km; then eq 1 takes the form GPMR(ω) )

(1 - ap)(1 - e-jωam) (1 - am)(1 - e-jωap)

(21)

M(ω) and ∆P(ω) are now written as follows: (1 - ap) M(ω) ) (1 - am)

[

∆P(ω) ) tan-1



1 - 2 cos ωam + am2 1 - 2 cos ωap + ap2

] [

am sin ω 1 - am cos ω tan-1

(22)

]

ap sin ω + ω(∆D) (23) 1 - ap cos ω

From eqs 22 and 23, one arrives at the diagnostics for this case listed in Table 3. 2.7. Gain, Time-Constant, and Delay Mismatch. When there exists a mismatch in gain, time-constant, and delay that is, Kp * Km, τp * τm, and Dp * Dm, eq 1 simplifies to GPMR(ω) )

Kp (1 - ap)(1 - e-jωam) -j∆Dω e Km (1 - a )(1 - e-jωa ) m p

(24)

M(ω) and ∆P(ω) are now written as follows:

(18)

From eqs 17 and 18, we can formulate conditions for diagnosing the presence of gain and time-constant mismatch (Table 2). The conditions in this case (Table 2) are a combination of the conditions for mismatches in individual parameters (gain

e-j∆Dω

M(ω) )

Kp (1 - ap) Km (1 - am)



1 - 2 cos ωam + am2 1 - 2 cos ωap + ap2

(25)

From eqs 25 and 26, the conditions listed in Table 4 are

4214

Ind. Eng. Chem. Res., Vol. 49, No. 9, 2010

[

∆P(ω) ) tan-1

] [

am sin ω 1 - am cos ω tan-1

]

ap sin ω + ω(∆D) (26) 1 - ap cos ω

formulated for diagnosing the presence of gain, time-constant, and delay mismatch together. The signatures of M(ω) and ∆P(ω) are unique for every type of mismatch. For an easy understanding of the proposed

approach, pictorial representations of the theoretical conditions for the diagnosis of different combinations of mismatches are given in Figure 2 and Figure 3. It turns out that the diagnostic conditions listed above can be largely simplified to a condensed stepwise procedure, which examines three critical aspects, namely, (i) the value of M(ω)|ω)0 which provides an estimate of Kp/Km, (ii) the flatness of M(ω) which is indicative of the presence of ∆τ, and (iii) the linearity of ∆P(ω) which is indicative of a delay mismatch. The procedure is presented in

Figure 2. Conditions to diagnose the type of mismatch among gain, time-constant, and delay.

Figure 3. Conditions to the diagnosis of combinations of mismatches.

Ind. Eng. Chem. Res., Vol. 49, No. 9, 2010

Kp * Km else Kp ) Km if flat, τp ) τm, else τp * τm

conventional definition of MPM ∆G ) Gp - Gm would not lead to a formulation of such simple signatures as has been obtained with the proposed PMR. In the following section, a method to estimate PMR from routine operating data is given.

linear: Dp * Dm, else Dp ) Dm

3. Estimation of PMR

Table 5. Theoretical Conditions for the Diagnosis of MPM in Model-Based Control Loops assessment procedure step 1 step 2 step 3

M(ω)|ω)0 * 1 test of flatness (zero slope) of M(ω) linearity check of ∆P(ω)

4215

diagnosis of MPM

Table 5. The condensed procedure can be applied to higherorder processes/models as well with the only difference being that the second aspect, that is, the flatness of M(ω) will be indicative of a mismatch in one or more time-constants of the process. The reader may note with interest that the use of the

In this section, we provide a method to estimate PMR from routine operating data in model-based control loops. For this purpose, consider a linear time invariant (LTI) SISO feedback control system operating under an IMC scheme. The set point, process input, process output, model output, and disturbance

Figure 4. Comparison of true diagnostic measures with the averages estimated from 1000 realizations.

4216

Ind. Eng. Chem. Res., Vol. 49, No. 9, 2010

are denoted by r[k], u[k], yp[k], ym[k], and d[k], respectively. The actuator and sensor dynamics are neglected for the remainder of the presentation. The expressions of process and model outputs can be written for the IMC structure as follows: yp[k] ) Gp(q)u[k] + d[k] ) Gp(q)u[k] + Gd(q)a[k] (27) ym[k] ) Gm(q)u[k] where

ˆ (ω)|ω)0. Figure 5. Analysis of the distribution of the estimate M

(28)



Gp(q) )

∑ g [k]q

-k

p

k)0

and ∞

Gd(q) )

∑ g [k]q

-k

d

k)0

In a noise-free environment (that is, when d[k] ) 0), the natural estimate of PMR is

Ind. Eng. Chem. Res., Vol. 49, No. 9, 2010

Yp(ω) Gp(ω) ) Gm(ω) Ym(ω) where Yp(ω) and Ym(ω) are the Fourier transforms of yp and ym, respectively. In practice, however, due to the presence of disturbances, this will lead to a poor estimate of the PMR. To overcome the effect of disturbances, a smoothed estimate is desired. This issue is reminiscent of the estimation of the transfer function using the ETFE versus the smoothed estimate involving the cross-spectra.6 On the basis of such ideas, eqs 27 and 28 are multiplied by r[k - l] followed by taking expectation on both sides. Subsequently, taking Fourier transforms on both sides of the equation leads to equations involving cross-spectra

ˆ (ω)|ω)π. Figure 6. Analysis of the distribution of the estimate M

4217

Φypr(ω) ) Gp(e )Φur(ω) + Gd(e )Φar(ω)

(29)

Φymr(ω) ) Gm(ejω)Φur(ω)

(30)





Noting the fact that the set-point and the disturbance are ideally uncorrelated, we arrive at the expression for the estimate of PMR as ΦYpR(ω) ΦYmR(ω)

)

ˆ p(ejω) G ˆ m(ejω) G

ˆ PMR(ω) )G

(31)

where the capped variable is introduced as a natural consequence

4218

Ind. Eng. Chem. Res., Vol. 49, No. 9, 2010

Table 6. Empirical (Data-Based) Conditions for the Diagnosis of MPM in Model-Based Control Loops

step 1 step 2 step 3

assessment procedure

diagnosis of MPM

M(ω)|ω)0 g 1.05 or M(ω)|ω)0 e 0.95 test of flatness (zero slope) of M(ω): M(ω) ) Rcτω + β linearity of ∆P(ω): ∆P(ω) ) RcDω

Kp * Km, else Kp ) Km if flat, i.e., |Rcτ| e 0.001, τp * τm, else τp ) τm if |RcD | g 0.9 (linear): Dp * Dm, else Dp ) Dm

in estimation problems owing to the finite sample size and the fact that Φar(ω) is numerically nonzero (but a small value) in practice. The estimation method provided above may be misinterpreted as a method for identifying the process model from closedloop operating data since Gm(ejω) is a known quantity. However, the estimate is only used to identify signatures relevant to the detection of MPM, but not to identify the process parameters. The key idea in this work is to propose a method to detect changes (and the direction of such changes) in process parameters based on routine operating data without having to explicitly estimate the process parameters. Table 5 is reflective of this fact. The flatness of M(ω) and linearity of ∆P(ω) are the key signatures that are sought, rather than the estimates of τp and/ or Dp from operating data. The signature-based approach of the proposed method demands only a minimal excitation, that is, a pulse excitation in the set-point, whereas an explicit estimation of process parameters imposes stronger excitation requirements either in the set-point or in a dither signal as is well-known in closed-loop identification. It may be noted, however, that gain can be estimated merely from step response data. Consequently, one of the conditions listed in Table 5 makes use of the gain estimate by examining the value of M(ω)|ω)0. A few remarks on the quality of the estimate are in place here. 3.1. Properties of the Estimate. (1) The set-point is assumed to contain some amount of excitation to provide a decent estimate of the PMR. It is assumed that a pulse-type variation with a minimum of two pulses is contained in the set-point. Such an assumption is not unreasonable since in several processes, the set-point changes are not merely a step-type but rather follow a pulse-type change. (2) Statistical properties of the estimate, particularly the bias and variance of the estimate, have to be assessed. The distribu-

Figure 7. Set-point profile considered for the simulation.

tion of the estimate is necessary to place confidence intervals on the estimate. A rigorous theoretical analysis is required to arrive at answers to the aforementioned issues and is reserved for future work. In this article, Monte Carlo simulations are used to check the biasedness and to determine an empirical distribution of the estimate obtained for an FOPTD process/ model. Figure 4 panels a and b show plots of the realization averages of M(ω) and ∆P(ω) versus the corresponding true values when the MPM is due to simultaneous mismatches in gain, time-constant, and delay parameters. The disturbance corrupting the output is assumed to be a zero-mean white-noise sequence with variance such that the signal-to-noise ratio (SNR) level is set to 10. The figures are generated from 1024 samples resulting from simulations for 1000 different realizations of disturbance. The plots clearly indicate that the proposed method given in eq 31 yields an unbiased estimate of the PMR. (3) To determine the distribution of the resulting estimate, Monte Carlo simulations were performed under similar conditions as those to check for the unbiasedness of the estimate. ˆ (ω)|ω)0 obtained for the Figure 5a shows the distribution of M 1000 different realizations of the disturbance. An empirical distribution is fit using the dfittool function of MATLAB. The best fit is achieved by a log-logistic distribution as supported by plots in Figure 5b-d. Here, the survival function, also known as a survivor function or reliability function, is a property of any random variable that maps a set of events, usually associated with mortality or failure of some system, onto time. It captures the probability that the system will survive beyond a specified time. Another name for the survival function is the complementary cumulative distribution function. To check for the invariance of the distribution with the frequencies, a similar ˆ (ω)|ω)π (Figure 6). analysis was carried out for the estimate M The log-logistic distribution once again turned out to be the

Ind. Eng. Chem. Res., Vol. 49, No. 9, 2010

4219

Table 7. The Values of Performance Index (ISE) for Different Cases of MPM no MPM

FOPTD δK

ISE

δτ

SOPTD

δD

FOPTD

SOPTD

10%

50%

10%

50%

20%

40%

δK ) 20%, δτ ) 50%, δD ) 20%

δK ) 20%, δτ2 )20%, δτ1 )20%, δD ) 33%

46.49

36.07

47.12

56.22

46.93

48.58

52.07

59.94

53.72

40.97

best choice. The difference, however, was the variances in the estimates (hence the distributions) at these two different frequencies. The variance of the estimate at ω ) π was much larger than that at ω ) 0. This should be expected since the

Figure 8. Plant and model ouput behavior for 10% and 50% gain mismatch.

estimates of the cross-spectra at high frequencies exhibit more variability. The mean values at ω ) 0 and ω ) π were 1.491 and 0.667, respectively, while the true values were 1.5 and 0.7, respectively, once again reflective of the unbiasedness of the

4220

Ind. Eng. Chem. Res., Vol. 49, No. 9, 2010

Figure 9. Plots of M(ω) and ∆P(ω) for 10% and 50% gain mismatch.

estimate. A separate analysis of ∆Pˆ(ω) at the aforementioned frequencies revealed similar characteristics with a log-logistic distribution for the estimate. (4) The Monte Carlo simulations achieved in providing insights into the distribution of the estimates. The variance of the estimates and/or the parameters of the distributions

would be useful in statistically verifying the hypothesis arising out of the diagnostic conditions given in Table 5. However, it was not possible to arrive at a closed-form expression for these quantities in terms of sample size and properties of the data. A natural alternative, therefore, was to propose thresholds for the critical values of M(ω)|ω)0,

Ind. Eng. Chem. Res., Vol. 49, No. 9, 2010

4221

Figure 10. Plant and model ouput behavior for 10% and 50% time-constant mismatch.

flatness (slope) of M(ω), and the linearity (slope) of ∆P(ω). These thresholds are obtained by a rigorous Monte Carlo simulation for each diagnostic situation mentioned in Table 5. The resulting values and the corresponding diagnostic checks are listed in Table 6. 4. Simulations A control system consisting of a process characterized by the transfer function Gp(z-1) ) (Kp(1 - ap)/(1 - z-1ap))z-Dp and model Gm(z-1) ) (Km(1 - am)/(1 - z-1am))z-Dm is simulated with an IMC scheme, with pulse type set-point

changes as shown in Figure 7. These pulse-type set point changes are not uncommon in several processes and can be therefore considered as a part of the routine operating conditions. The tuning parameter, λ and the model parameters Km, τm, and Dm are fixed at 2, 1, 10, and 5, respectively. The closedloop system is simulated under four different situations of MPM, namely, gain mismatch, time-constant mismatch, delay mismatch, and a mismatch in all three parameters. The mismatches are expressed in terms of percentages for simplicity and are defined as follows:

4222

Ind. Eng. Chem. Res., Vol. 49, No. 9, 2010

δK )

Kp - Km 100 Km

(32)

δτ )

τp - τm 100 τm

(33)

δD )

Dp - Dm 100 Dm

(34)

To determine the degradation in performance, the ISE is used as a measure of the performance of the closed-loop system. The

ISE value is calculated for the zero MPM case and compared with the value corresponding to each case of mismatch. The extent of degradation, evidently, depends on the magnitude of MPM. For the cases discussed as follows, the performance index values are tabulated in Table 7. 4.1. Gain Mismatch. Two different levels of gain mismatch at 10% and 50% are considered. The cross spectra between process output yp[k] and set point r[k] and model output ym[k] and r[k] are estimated. The magnitude ratio M(ω) and phase difference ∆P(ω) are computed from the estimate of PMR as given in eq 31. Figure 9a,b shows the estimates of these

Figure 11. Plots of M(ω) and ∆P(ω) for 10% and 50% time-constant mismatch.

Ind. Eng. Chem. Res., Vol. 49, No. 9, 2010

4223

Figure 12. Plant and model ouput behavior for 20% and 40% delay mismatch.

measures for the two different levels of gain mismatch. Though visual inspection of Figure 8a,b clearly indicates the presence of gain mismatch, the nonexistence of time-constant mismatch and delay mismatch cannot be ruled out. Proceeding with the first step outlined in Table 6, the value of M(ω)|ω)0 is 1.1009 when δK ) 10% and 1.5001 when δK ) 50%, clearly indicating the presence of gain mismatch. From the values |R1| ) 0.0009 < Rcτ, |R2| ) 0.001 < RcD for case 1a and |R1| ) 0.0009 < Rcτ, |R2| ) 0.0002 < RcD for case 1b, it can be concluded that there exists only gain mismatch. Further, the extent of mismatch can be calculated from these values by calculating the deviation from the no MPM case, that is, M(ω)|ω)0 ) 1, which yield 10.09%

and 50.01%, respectively. These values are, as expected, in agreement with the true values. 4.2. Time-Constant. Two different levels of time mismatch at 10% (τp ) 11) and 50% (τp ) 15) are considered. Visual inspection of Figures 10a,b does not reveal any information of MPM. The magnitude ratio M(ω) and phase difference ∆P(ω) are computed from the estimate of PMR as given in eq 31. Both these quantities are plotted in Figure 11a,b. The stepwise procedure outlined in Table 6 is followed. The values of M(ω)|ω)0 ) 0.997, |R1| ) 0.018 > Rcτ, and |R2| ) 0.007 < RcD when δτ ) 10% and M(ω)|ω)0 ) 0.9983, |R1| ) 0.023 < Rcτ, and |R2| ) 0.033 < RcD

4224

Ind. Eng. Chem. Res., Vol. 49, No. 9, 2010

Figure 13. Plots of M(ω) and ∆P(ω) for 20% and 40% delay mismatch.

when δτ ) 50% clearly indicate the presence of time-constant mismatch. Further, the sign of R1 in both the cases indicates that τp is greater than τm in both the cases. 4.3. Delay Mismatch. The system is simulated for the delay mismatch of 20% (Dp ) 6) and 40% (Dp ) 7). Visual inspection of Figure 12a,b does not reveal any information of MPM. The cross spectra, ΦYpR(ω), ΦYmR(ω), are estimated and subsequently used to compute M(ω) and ∆P(ω). Both these quantities are plotted in Figure 13a,b. The stepwise procedure outlined in Table 6 is followed. The values of M(ω)|ω)0 ) 0.9981, |R1| ) 0.0009 ( Dm, δK ) δτ ) 0 K p > K m, τ p > τ m , D p > D m K p > K m, τ p > τ m , D p > D m

5. Extensions to Other Scenarios: Remarks and Recommendations The theoretical derivations and simulation studies presented until this point have been set up in somewhat restrictive conditions. However, as we discuss in the following section, the proposed method can be suitably extended or modified to other scenarios as well. The ambiguity arising because of poor

Ind. Eng. Chem. Res., Vol. 49, No. 9, 2010

4227

Figure 16. Diagnostic plots for SOPTD plant/model with random walk-ype disturbance.

controller tuning possibly manifesting as significant MPM is also discussed. 5.1. Random Walk Disturbances. The processes affected by random walk-type disturbances are common in chemical plants. It is not straightforward to apply the proposed method on such processes because of the presence of the integrating effects (nonstationarities) in the outputs. One solution is to consider the differenced outputs so that the nonstationarity can be eliminated. However, one must difference the outputs only when the integrating effects are detected. There are several wellestablished methods to detect unit roots in time series.25-27 Once the differenced outputs are available, the diagnostic procedure is applied as usual. To illustrate this idea, the SOPTD system considered in section 4.5 is simulated with random walk disturbance. Figure 16a shows the diagnostics when the

procedure is applied to differenced outputs. As earlier, the diagnostics and bounds presented in Table 6 together pinpoint the sources of MPM. Figure 16b shows erratic behavior of PMR when the procedure is applied to the outputs as they are, that is, without differencing. 5.2. Poor Controller Tuning. Poor controller tuning is a relative term. There are two possibilities. Either the controller has been tuned poorly right from the start or the controller tunings are deemed poor given that MPM has become significant or that it is witnessing disturbances that it has not been designed for. The proposed method is able to correctly determine the absence or presence of significant MPM in all the aforementioned situations. In the first case where the controller has been tuned suboptimally, then MPM is truly not the cause and will be identified as not so by the method. The reason is that the

4228

Ind. Eng. Chem. Res., Vol. 49, No. 9, 2010

same input excites the plant and model. If the controller tuning is due to unforeseen disturbances, then MPM is not the cause. By way of estimation of PMR, the effect of disturbances on the estimate is minimized. Consequently, the technique will not show MPM as the cause of poor performance. Finally, if MPM is truly the cause of poor performance (although to an observer it could be the controller being unable to deliver), then PMR will correctly highlight MPM as significant since the same input is being fed to the plant and model. 5.3. Structural Mismatch. In a structural mismatch scenario, clearly one can speak of gain and (possibly) delay mismatch, but not of a time-constant mismatch. Examining mismatches in delay is appropriate only if the FOPTD (or the approximated model) model did not include delay as an approximation of the higher-order dynamics. However, one could ask if the approximation errors have become significantly larger than those which were acceptable at the beginning. Furthermore, the nature of PMR itself changes considerably for structural mismatches. A striking difference is that the PMR has a relative order of zero (not strictly proper) for unstructured mismatches, while it is strictly proper for structured mismatches. In light of these facts, the proposed method in its present form cannot be directly applied to highlight structural mismatches. Nevertheless, the principles presented in this work, with some modifications to their implementation, may be applied to detecting significant increase in the errors and pin-pointing the element of the model that needs attention. As an example, consider the situation of approximating a SOPTD process with an FOPTD model such that the gain and delay of the model match that of the process, with the time-constant selected using a least-squares technique. MPM exists right from beginning, however, that which is acceptable. While gain and delay mismatches can be detected using the existing technique, a modification is required to address the question if the timeconstant of the model needs to be reidentified. A PMR generated using the data at the time of model identification can be generated to serve as a template PMR for all subsequent operations of the plant. At any future point, the PMR can be projected against the template PMR, whereupon statistical techniques can be used to distinguish between process and chance variations and to pinpoint the elements requiring reidentification. A discussion and the nature of such modifications and their effectiveness warrants a separate treatment and is beyond the scope of the current work. 5.4. Multivariate Systems. The proposed method cannot be directly applied to a multivariate system where, in general, input and output relations involve interactions. Therefore, the extension of the proposed approach to the multivariate systems involves an additional step of decoupling these interactions. This can be carried out by means of the partial coherence function as has been illustrated in ref 28. The basic idea is to break up a multivariate system into a number of SISO systems characterized by conditioned inputs and conditioned outputs. Subsequently, the proposed method can be employed on each SISO system. The detection and diagnosis of significant MPM in each subsystem can be carried out in the domain of conditioned outputs and conditioned set points. A rigorous treatment of these ideas is a subject of future study. 6. Conclusions In this work, the problem of characterizing the deviations of process parameters from that of the model has been addressed. The primary application of this method is envisaged in the area of closed-loop performance assessment where degradation in

performance could occur due to model-plant mismatch. The study and results presented in this work has been set up in the IMC framework. The major contribution of this work is the formulation of the MPM characterization in the frequency domain for the identification of specific signatures to distinguish among mismatches in gain, time-constant and delay. The current work introduced a quantity Gp(ω)/Gm(ω), termed as the plantmodel ratio (PMR) as against the conventional Gp - Gm to quantify MPM. The theoretical results show that the advantage of the proposed PMR is the ease of representation in the complex frequency domain. The amplitude part of the PMR contains the effects of gain and time-constant mismatches, while the phase component contains the effects of mismatch in delay and timeconstant. The key outcome of this work is a simple three-step procedure to identify mismatches in gain, time-constant, and delay from closed-loop data through the use of PMR. A method to estimate PMR from the response due to minimal set-point changes has been presented. Statistical analysis of the estimate using Monte Carlo studies show that the amplitude and phase estimates follow a log-logistic distribution. Theoretical analysis, however, requires a more rigorous treatmentsa subject reserved for future study. Simulation studies on FOPTD and SOPTD processes demonstrate the effectiveness of the proposed methodology. Although the methodology and its illustration has been demonstrated on SISO systems with IMC schemes, the recommendations and ideas outlined in section 5 can be used to apply and/or extend to other situations, namely, (i) structural mismatch, (ii) multivariate systems, and (iii) MPC driven control loops. Acknowledgment The authors would like to thank the anonymous reviewers for their critical comments which have helped in enhancing the depth of the presentation. Literature Cited (1) Huang, B.; Shah, S.; Kwok, E. Performance assessment of multivariable processes. Automatica 1989, 33, 1175–1183. (2) Harris, T.; Boudreas, F.; Macgregor, J. Performance assessment of multivariable feedback controllers. Automatica 1996, 32, 1508–1518. (3) Huang, B.; Shah, S. Performance Assessment of Control Loops; Advances in Industrial Control; Springer: London, 1999. (4) Jelali, M. An overview of control performance assessment technology and industrial applications. Control Eng. Pract. 2005, 14, 441–466. (5) Ender, D. Process control performance: Not as good as you think. Control Eng. 1993, 40, 180–190. (6) Ljung, L. System Identification Theory for the User; Prentice Hall: Upper Saddle River, NJ, 1999. (7) Ljung, L.; Forssell, U. Closed-loop identification revisited. Automatica 1999, 35, 1215–1241. (8) Harris, T.; Seppala, C.; Desborough, L. A review of performance monitoring and assessment techniques for univariate and multivariate control systems. J. Process Control 1999, 9, 1–18. (9) Bialkowski, W. Dreams versus reality: A view from both sides of the gap. Pulp Pap. Can. 1993, 94, 19–27. (10) Kozub, D.; Garcia, C. Monitoring and Diagnosis of Automated Controllers in the Chemical Process Industries, Technical Report, AIChE Annual Meeting, 1993; Paper No. 145i. (11) Ha¨gglund, T. A control-loop performance monitor. Control Eng. Pract. 1995, 3, 1543–1551. (12) Thornhill, N.; Huang, B.; Zhang, H. Detection of multiple oscillations in control loops. J. Process Control 2003, 13, 91–100. (13) Choudhury, M.; Shah, S.; Thornhill, N. Diagnosis of poor controlloop performance using higher-order statistics. Automatica 2004, 40, 1719– 1728. (14) Tangirala, A.; Shah, S.; Thornhill, N. A new tool for plant-wide oscillation detection. Process Control 2005, 15, 931–941.

Ind. Eng. Chem. Res., Vol. 49, No. 9, 2010 (15) Tangirala, A.; Kanodia, J.; Shah, S. Non-negative matrix factorization for detection and diagnosis of plant wide oscillations. Ind. Eng. Chem. Res. 2007, 46, 801–817. (16) Kesavan, P.; Lee, J. Diagnostic tools for multivariable model based control systems. Ind. Eng. Chem. Res. 1997, 36, 2725–2738. (17) Kendra, S.; Cinar, A. Controller performance assessment by frequency domain techniques. J. Process Control 1997, 7, 181–194. (18) Kendra, S.; Basila, M.; Cinar, A. Control Syst. Mag. IEEE 1994, 14, 37–47. (19) Kozub, D. Controller Performance Monitoring and Diagnosis: Experiences and Challenges, Proceedings of the 5th International Conference on Chemical Process Control, Tahoe City, CA, 1996; pp 83-95. (20) Patwardhan, R.; Shah, S. Issues in performance diagnostics of model-based controllers. J. Process Control 2002, 12, 413–427. (21) Kammer, L.; Gorinevsky, D.; Dumont, G. Semi-intrusive multivariable model invalidation. Automatica 2003, 39, 1461–1467. (22) Huang, B. On-line closed-loop model validation and detection of abrupt parameters changes. J. Process Control 2001, 11, 699–715. (23) Jiang, H.; Li, W.; Shah, S. Detection and isolation of model-plant mismatch for multivariate dynamic systems. Fault Detect., SuperVision, Saf. Tech. Processes 2006, 12, 1396–1401.

4229

(24) Morari, M.; Zafiriou, E. Robust Process Control; Prentice Hall: Englewood Cliffs, NJ, 1989. (25) Dickey, D. A.; Fuller, W. Distribution of the estimators for autoregressive time series with a unit root. J. Am. Stat. Assoc. 1979, 74, 427–431. (26) Phillips, P.; Perron, P. Testing for a unit root in time series regression. Biometrica 1988, 75, 335–346. (27) Kwiatkowski, D.; Phillips, P.; Schmidt, P.; Shin, Y. Testing the null hypothesis of stationarity against the alternative of a unit root. J. Econometr. 1992, 54, 159–178. (28) Selvanathan, S.; Tangirala, A. Time-delay estimation in multivariate systems using Hilbert transform relation and partial coherence functions. Chem. Eng. Sci. 2010, 65, 660–674.

ReceiVed for reView May 13, 2009 ReVised manuscript receiVed February 1, 2010 Accepted February 1, 2010 IE900769V