1694
I n d . Eng. C h e m . Res. 1992,31,
Walker, D.; Koros, W. J. Gas Separation Membrane Materials Selection Criteria: Weakly and Strongly Interacting Feed Component Situations. Polym. J. 1991,23,481. Wankat, P. C. Rate-Controlled Separations; Elsevier Applied Science: New York, 1990. Yang, R. T.Gas Separation by Adsorption Processes; Butterworth Boston, MA, 1987.
1694-1704
Yen, L.; McKetta, J. J. A Thermodynamic Correlation of Nonpolar Gas Solubilities in Polar, Nonassociated Liquids. NChE J. 1962, 8 (51,501. Receiued for review October 25, 1991 Revised manuscript received March 17, 1992 Accepted April 5, 1992
Probabilistic Approach to Robust Process Control Charles D. Schaper,* Dale E. Seborg, and Duncan A. Mellichamp Department of Chemical and Nuclear Engineering, University of California, Santa Barbara, California 93106
A probabilistic approach to robust process control is developed. First, a statistical measure of a controller’s ability to reject disturbances is introduced. Next, a new robust control framework of characterizing model uncertainty descriptions by probability distributions is developed. The statistical measure of disturbance rejection is then incorporated within the framework. In the proposed probabilistic approach, process knowledge can be incorporated in the design procedure and controller performance can be analyzed by probability measures. Several simulation examples demonstrate the advantages of the new approach.
Introduction An important objective in designing a process control system is robustness to modelling error. Previous approaches to robust process control design have generally used bounds around the parameters or frequency response of a nominal plant model to describe model uncertainty. The control system is then designed to minimize the effects of a worst-case situation. Current design approaches for robustness are described by Morari and Zafhiou (1989). Process control applications of these design techniques include those of Agamennoni et al. (1988)and Skogestad et al. (1988).Advantages of existing design techniques for robustness include the following: (1)closed-loopstability is guaranteed over the entire range of model uncertainty (robust stability); (2) an upper bound on a given performance measure is guaranteed (robust performance). Because the controllers are generally designed for worst-case situations that may have a low probability of occurring, the resulting robust controllers may be very conservative for more typical operating conditions that have a much higher probability of occurring. In this paper, a new approach to robust process control design is developed in which model uncertainty is characterized by probability distributions. This approach allows closed-loop performance tradeoffs to be analyzed as a function of the likelihood of controller performance; that is, performance can be characterized by a probability measure for all situations between nominal and worst-case conditions. The result is a more complete analysis strategy that can result in better controller design. In the subsequent development, a general linear representation of the plant description is used in which modeling error is described by probability distributions. Modeling error due to both parameter uncertainty and the linear approximation of a nonlinear plant can be included within this probabilistic framework. It should be noted that the error resulting from the approximation of a nonlinear system by a linear model may be greater than any model parameter uncertainty. For example, this situation
* Present address: Department of Electrical Engineering, Stanford University, Stanford, CA 94305.
could occur when a fundamental physical model of the process does not exist or is too complex for controller design, and consequently, an empirical linear model (e.g. a transfer function model) is developed from experimental data. In this instance, the parameters of the linear approximation can be represented by probability distributions. Although we describe some methods and examples of approximating this type of modeling error, it is not the intent of this paper to provide a well-formulated description of how to identify model uncertainty descriptions. However, we note that probabilistic descriptions of modeling error can be developed from a wide variety of sources, including statistical information on phenomenological model parameters, empirical model parameters, or frequency response (Cloud and Kouvaritakis, 1987;Correa, 1989;Goodwin and Salgado, 1989;Stengel and Ryan, 1989). Also, process knowledge is usually available in the form of engineering heuristics and information about the range of operating conditions. The probabilistic model description is sufficiently general to capture such prior process knowledge and incorporate it within the design procedure. In addition to the development of a general probabilistic framework, a statistical measure of closed-loop disturbance rejection capabilities is introduced for process control applications. A disturbance rejection measure is generally more appropriate for process control applications because the set-point remains constant for long periods of time. In the development of this measure, it is important to note that performance specifications for outputs or inputs can be formulated in terms of statistical moments. For example, a typical product specification is expressed in terms of a mean and standard deviation (also referred to as root mean square). Well-known control design strategies have been developed to minimize statistical moments of the outputa and inputs. These controller design strategies include qlassical methods such as minimum variance control (Astrijm, 1970;Box and Jenkins, 1976;Kucera, 1979),in addition to current methods such as robust linear quadratic Gaussian (LQG) control strategiea (Stengel, 1986; Bernstein and Haddad, 1990) and constrained minimum variance control (Makila et al., 1984;Hotz and Skelton, 0 1992 American Chemical Society
Ind. Eng. Chem.Res.,Vol. 31, No. 7,1992 1695 r
1986; Hsieh et al., 1989). In addition, statistical quality control applications are often formulated in terms of the statisticalmomenta of the inputa and outputs (MacGregor, 1988). Since statistical moments already are used for design specification,there is strong motivation to develop an analysis methodology for process control that is based on the statistical momenta of the process signals. In this paper, the covariance matrix of the inputs, states, and outputs, is used to characterize the disturbance rejection properties of a process control system. Results are presented which allow the computation of bounds on variances of the important process signals for closed-loop conditions. These bounds are then analyzed to determine whether the disturbance rejection characteristics of the controller are adequate. The paper is organized by first introducing the statistical measure of a controller’s ability to reject disturbances. The probabilistic framework for robust process control is then developed, and the statistical measure is incorporated within the framework. The results are then generalized in the frequency domain.
Statistical Measure of Performance The process is assumed to be described by a linear state-space representation in which model uncertainty can be present. This description is given by S
= (A + AA)x + (B + AB)u + Dlw y = (C + AC)X
+D ~ w
(3)
u=&+&
(4)
where f is the state vector for the controller, S E B*,the initial condition is f ( 0 ) = 0, and the set point is assumed to be zero. This linear description of the control law is UBed to facilitate computational solutions to the analysis and synthesis problems (see Bernstein and Haddad (1990)); Proportional control action results if B # 0 and A ,B, C are zero. Integral control !ction is specified by setting appropriate eigenvalues of A to zero. The controller matrices are represented in partitioned form as (5)
By combining eqs 1-4 and performing some algebra, one can express the closed-loop system relations as the augmented state-space model
where
(6)
J
r
1
(11)
We now consider w(t) to be a zero mean, white noise, stochastic signal with covariance matrix R. This representation corresponds to a signal that is uncorrelated with respect to time. If this signal is passed through an appropriate filter, a wide variety of zero mean disturbances can be effectively modeled. The covariance matrix for i is defined by
Q = lim ~ [ i ( tt )( t ) T ] t--
(2)
i=Af+&
i = (A + A A ) +~ (B + AB)W
(7) L
(1)
where u E B5 is the control input vector, w E Bb is the disturbance vector (which might include elements accounting for process and measurement noise), x E 934 is the state vector, y E W . y is the observation vector. The matrix perturbations are considered to belong to the set U where (AA,AB,AC)E U. For clarity of presentation, the time variable t will often be omitted. A general class of linear feedback controllers can be described by the following state-space model
.
(12)
where C is the expectation operator which is taken with respect to the disturbance w(t). In the following development, we will not consider the more general case of &w(t)# 0. Thus, the disturbance is characterized solely by ita range of possible values, statistically quantified by R. With this representation, the analysis (and design) of controllers is then based on the system’s long term statistical behavior. The covariance matrix can be represented in a partitioned form
where the submatrix Qll =
lim &[x(t)x(t)T] t--
(14)
is the covariance matrix for x(t). Also, Qlz = QTl with definition and dimension that follow accordingly. Provided that A + AA is asymptotically stable for all AA, AB, AC) E U, then the closed-loop covariance matrix is a positive semidefinite solution to a Lyapunov equation (Bryson and Ho, 1975)
6
0=
(A + d)Q + Q(A +
+
(B + Afi)R(8 +
(16)
We wish to develop the following bound Q on the closedloop covariance matrix
Q
IQ, V(AA, AB,AC) E U
(16)
implying that a positive semidefinite matrix must result after subtracting Q from Q. To determine this bound, we first partition Q as
1696 Ind. Eng. Chem. Res., Vol. 31, No. 7, 1992
1
Qzi
Qzz
J
with dimensions that correspond with those of Q. The following properties of eq 16 are useful and can be obtained from basic properties of positive semidefinite qatrices. Property 1. Each principal submatrix of Q, after subtraction from the corresponding submatrix of Q, must yield a nonnegative determinant. Property 2. Each diagonal element of Q must be greater than or equal to the corresponding element of Q. Consequently, the variance of each state variable of x(t), determined by the corresponding diagonal element of Qll, is bounded by the corresponding diagonal element of Qll. The input energy expended is also of concern in practical control applications. One suitable measure of this property, the covariance of the input during closed-loop operation, is given by
Q,
=
lim e[u(t) u(t)*] t+-
(18)
By substituting eq 2 into eq 4, for the case where AB and AC equal zero, and then substituting the resulting relationship into eq 18,we arrive at a bound on the covariance matrix of the input energy that is given by
Q, =
cQz2eT + cQZIC%T + DCQlzeAT+ DCQllC’%T
+ DDzRDmT (19)
Similarly, we can determine a bound on the covariance of the observation signal defined by
Q, * lim e[y(t) yWTI t--=
(20)
For the case where AB and AC equal zero, this bound is found by substituting eq 2 into eq 20 to arrive at
Q, = CQ1lCT + DzRDT
(21)
Using property 1,it is straightforward to show that eqs 19 and 21 provide a bound on the covariance of the input and the covariance of the output, respectively. In addition, we have the following. Property 3. The covariance bounds determined from eqs 15,19, and 21 are valid for any zero mean disturbance with covariance matrix bounded by R. Consequently, the bounds hold for a class of disturbances (rather than a single one) with the characteristics of zero mean, stable, and covariance matrix bounded by R. One approach that can be used to compute the bound on the covariance matrix is_to ap ly a grid over the region defined by A + AA and B + A and solve a Lyapunov equation at each node of the grid. A necessary condition for stability is that all of these Lyapunov solutions are positive definite. A bound on the variance of each state variable can then be determined by finding the maximum corresponding diagonal element over the grid. Bounds on the input energy and output covariance matrices can be determined using a similar approach that involvea a search over partitioned matrices. The disadvantage of this approach is that it can only approximate the bound since the actual bound may lie between the nodes of the grid. Practical consequences can be reduced by using a fine grid, but the computational requirements will increase. Another approach for obtaining an upper bound on the closed-loop state covariance matrix is to replace the uncertainty terms in eq 15 by a bounding function with the result that the network of Lyapunov equations then reduces to a single modified Lyapunov or Riccati equation. The disadvantage
of this approach is that the bounding functionsare actually approximations; thus substantial loss in accuracy can result. These issues are discussed in the following sections. Example 1. The use of the covariance matrix to characterize the disturbance rejection capabilities of a linear feedback controller is demonstrated using a nonlinear model of a chemical reactor. Material and energy balances for a continuous stirred tank reactor with a first-order, irreversible, exothermic reaction are given by
where the notation used is standard. Consider the problem of controllingthe outlet concentration CAby manipulating the coolant temperature T,. V, QF, and all model parameters are assumed to be constant. These equations can be expressed in dimensionless quantities as i1 = -xl + &(I i z= -xz
- xl)exz/(1+4y)
+ BhDu(l - xl)ex2/(1+xz/y) - @x2 + j3u + d
(24) (25)
Feed temperature TFis considered to be a disturbance and modeled as a firsborder process driven by zero mean white noise with unit intensity T d d = -d + K ~ w (26) To perform the transformation from this nonlinear structure to a linear structure, bounds must first be imposed on the state variables. An operating window in the xl-xz state plane is defined by these bounds with the nominal point denoted by (x10,x2,)and a window width denoted by (xlr,x2,). The nonlinear state equations are transformed to linear state relationships by first expressing eq 24 and eq 25 in the form, x1 = x < + xloand x 2 = x i + xzo, where the prime denotes a deviation variable. The nonlinearity of the reactor model is dependent on the single term (1 - x’,, - x’~)e(xh+x’z)/(l+(xh+1’2)/y) (27) Consequently, only this term needs to be bounded in order to define an uncertainty description. With some algebraic manipulation, this nonlinear term can be shown to be bounded by
The function 6 quantifies the approximation error which is bounded between -1 and +l. Furthermore, the error due to the approximation of the nonlinear term is expressed through C,,the center of the nonlinear term, and R,, the range around the center:
L+ + Lc, = 2
Ind. Eng. Chem. Res., Vol. 31, No. 7,1992 1697 1.0 I
.
I.
.
0.15 I
i
0.1
0.5
d
xi
-0.5
-1.5
I
OPen-lWP
0
10
I
0.0s
0
-1.0
I
1
'1
I
0
-0.M
I . . . . 20 . , 30. . 40 . I 50 10
0
dmc Figure 1. Disturbance signal for the example.
20
30
40
so
dmc
Figure 2. Dimensionless concentration for open-loop conditions and closed-loop conditions.
The resulting state equations expressed in deviation form are then given by
0.8
open-Imp
y' = Cr'
where x' = [ x ' ~
and -0.8
I
10
0
30
20
40
dmc
1
so
Figure 3. Dimensionless temperature for open-loop conditions and closed-loop conditions. 3
(36)
2
PIT
(37)
c = p 01
(38)
B=(O
1
u
o -1
The relationships for A and AR are given by N,,= -1 - Da(eXd(l++zo/Y))
-2 -3
Nlz = DaC, N,,= -B,J)a(e-%/(1+Xzo/Y))
Numerical parameters chosen for this example are Bh = 1, Da = 0.072,Kd = 3, xlo = 0.3,xz0 = 1.96,xl, = 0.3,xZr = 1.96,0 = 0.3,y = 20, Td = 1.5. Alternative uncertainty descriptions for this procese are given in Doyle and Morari (1989). A controller is then designed to place the closed-loop poles for the nominal model (6 = 0) at -10. (The real parts for the open-loop poles range from -1 to -1.4 for this particularly benign system.) The parametrization of the resulting controller, represented in partitioned form, is given by
0'
]
.
10
20
30
40
I
50
dnr Figure 4. Dimensionless input response for closed-loopconditions.
Nzz = -1 - 0 + B,J)aC, R11 0 R12= DaR, R11= 0 Rz2 = BJlaR,
-20 0 K=[l 0 1.4~10~1.2~10~ -810
0
(39)
It can be shown that this controller contains integral action and that robust stability is satisfied.
The performance of the controller is evaluated with the random variable w described as white noise with unit intensity. The disturbance d is expressed by eq 26. The time response for a single realization is shown in Figure 1. The open-loop and closed-loop responses for this disturbance are compared in Figures 2 and 3, where the initial steady-state conditions are xlo = 0.3 and xzo = 1.96. The corresponding input for closed-loop control is presented in Figure 4. The controller is seen to provide adequate disturbance rejection. The statistical techniques described in the previous section are then used to analyze the performance properties of this controller in terms of the variance of the input and state and output variables, as well as to check for robust stability. Becaw only a single term represents the model uncertainty, a grid approach is used to compute the covariance bounds. The uncertainty parameter 6 is discretized over the region [-1,+1], and a Lyapunov equation is solved at each of the 20 nodes. For open-loop conditions, the bounds on the variance of x , and x z are shown in Figure 5 as a function of 6. The Lyapunov relationship is then used to show that the controller achieves robust stability during closed-loop operation. The bounds on the variances of xl and x z under
1698 Ind. Eng. Chem. Res., Vol. 31, No. 7, 1992
1.2
0.8 lY6
,
i
t
i
5 -0.6
-1
-0.2
0.2
6
Figure 5. Bounds on the variances of the temperature and concentration dimensionless variables for open-loop Conditions. 0.4
0.3
.I
0.2
0.1
0
-1
0.2
-0.2
-0.6
0.6
1
s Figure 6. Bounds on the variances of the temperature and concentration (XlO) dimensionless variables for closed-loopconditions. 35
30
,I
23 20 1s 10
Using this strategy, upper bounds on the variance of the inputs, states, and outputs can be determined. These guarantees are useful to the control design engineer since they provide conditions on time-domain, observable characteristics of the process. Since the output specifications are often expressed in terms of statistical moments, this methodology for process control can be used directly to determine whether the controller performs adequately when subjected to disturbances of a given type class. Furthermore, the guarantees on the variances can be employed on line to evaluate the functionality of the controller and to detect for possible faults. Consider that traditional fault detection strategies often use variances of the inputs, states, and outputs to assess whether a fault has occurred (Himmelblau, 1978; Isermann, 1984). These variances are generally computed on line and compared to a priori limits. A key parameter in the success of these fault detection strategies is therefore the choice of a priori limits on the variances. If the limits are too tight, numerous false alarms may be detected. Conversely, if the limits are too conservative, fault detection may be significantly delayed. Despite their importance, these variance limits are generally determined on an ad hoc basis. The approach presented above can be used to compute a realistic bound on the covariance matrix (and therefore the variances) and thus provide a viable method of initializing the limits employed in fault detection strategies.
-1
-0.2
-0.6
s
0.2
0.6
1
Figure 7. Bound on the variance of the input energy.
open-loop and closed-loop conditions are presented in Figures 5 and 6, respectively, as a function of 6 using the grid computational approach. The improvement in variances due to closed-loop control is evident by comparing these figures. The bound on the dimensionless input energy under closed-loop control is shown in Figure 7. The asymmetric behavior that results from use of the grid approach demonstrates that certain values of model uncertainty improve control system performance, as would be expected. Note the low value for the concentration variance bound for values of 6 near -1. As a consequence, the input energy required for this situation is low. The global bounds on the variances are found easily from Figures 6 and 7; they are 0.36 for rl, 0.02 for rz,and 31 for u. These values can be used to assess whether the controller achieves adequate disturbance rejection. The variance bounds are valid for any zero mean disturbance having a variance Qd satisfying
Kd2/2rd (40) This value was determined by the Lyapunov relationship for white noise signals passing through linear filters and eq 26. Qd
5
Probabilistic Representation of Uncertainty We now consider the matrix elements of the real-valued perturbations (AA,AB,AC) to be functions of elements belonging to a p-dimensional random vector, 8. This vector is modeled as a random process b = g(8,e) (41) where e ( t ) is a stochastic process with a particular underlying probability distribution. The initial conditions for the parameter vector are represented by a probability density function f(&t=O). We consider e(t)such that 8 is a stationary process whose statistics are completely determined by f@;t=O). Realizations for this vector of uncertain parameters can thus be represented by a set 0 C wp.
The proposed formulation of robust performance based on probabilistic uncertainty ensures that all performance objectives are guaranteed to a certain probability level of modeling error. This definition is used to develop probability distributions for the performance bounds expressed by the covariance matrices presented above. We now express the covariance matrix for the closedloop state vector as
Q
=
lim e [ t ( t )%(t)Tle] t-c-
(42)
where the expectation is evaluated with respect to w(t)for a given realization of the uncertain parameter vector, 8 . We assume that a bound on the closed-loop covariance matrix can be found such that (43) Pr (QIQ) 1 8 where the distribution function is taken with respect to the random variable 8. The probability distribution function of Q is denoted by F(Q). Consider that the probability that the event Q IQ consists of all outcomes such that Q = g[8] IQ. Hence, (44) F(Q) = Pr (Q IQ) = Pr @[e] IQ) For a specific Q, all realizations that satisfy g[8] I Q form
Ind. Eng. Chem. Res., Vol. 31, No. 7, 1992 1699 a set in RP denoted by bg. We can then determine F(Q) = Pr (6 E bg) (45) From this relationship,the probability level associated with eq 43 can be djtermined. Similar relationships can be developed for Q, and 8,. A method to approximate the probability distribution on the bound of the covariance matrix is now developed. We extend the concept of a bounding function (Petersen and Hollot, 1986) to the probabilistic case. Theorem 1. Define a bounding function Q E P,the set of symmetric matrices, for all positive definite Q E Pfl (where P denotes positive definite matrix) such that Pr ( d Q + QAdT IQ) = 0 (46) Now suppose there exists Q E Pa for a given controller that satisfies o = AQ QAT+ BRBT+Q (47) Then (A + d , B ) is st+bilizable for all (AA,AB,AC)E U if and only if A AA is asymptotically stable for all (AA,AB,AC) E U. It then follows that Pr (Q I Q) 1 (48) Proof. Asymptotic stability follows from Bernstein and Haddad (1990). Since asymptotic stability is achieved for all (AA,AB,AC) E U, subtraction of eq 47 from eq 15 yields the following form
with the Cholesky decomposition of V-' = UD1/2D1/WT. Proof. Given that
0 5 [BiK1/21- Ai/K1/2]Q[BiK1/2- A ~ / K ' / ~ (57) ]~ the proof follows by algebraic manipulation after comto yield puting V-' = UD1/2D1/WT
aTa
*
I (58) where B = D1/W4. The scalar K can be used to reduce the conservativeness of the bounding function. The closed-loop system then has the form P
(59)
+
+
(49) Hence the inequality in eq 48 follows from eq 46. A linear bounding function can be developed for uncertainty of the form
P
AB =
C&Pi
is1
(51)
P
AG = CdiCCi i=l
(52)
with the structure of the model uncertainty given through the matrices (A*i,B*i,C*i).The vector 4 = [41...4p,1T.is a real parameter random vector with a particular probabllrty distribution describing the uncertainty. We consider the case where the uncertain terms are confiied to an ellipsoid given by
f$Tv-'c$ I $I (53) where V E P P defines the general orientation of an ellipsoid bounded by the real-valued scalar $. A factorization of V-' is denoted by LTDL where L is a lower diagonal matrix of dimension p X p and D is a diagonal matrix of dimension p X p . Theorem 2. A bounding function for the uncertainty structures of eq 50 to eq 52 is given by P
Q = $KQ+ KCA,QAT is1
J
L
for i = 1 ... p . Now consider the case where a normal distribution with zero mean is specified for the uncertain elements of 4. In this instance, considering V to be the covariance matrix, we have where the scalar $I now denotes the 1- u probability level of the x2 probability distribution with d degrees of freedom. A simple example involving only a single uncertain term demonstratas the results of this section. For multiple uncertain terms, the computations required to map the statistical description of model uncertainty to the closedloop state vector are more complex. However, the example clearly illustrates the benefit of analyzing system performance in terms of probability. (A more realistic example with two sources of uncertainty is considered in the following section.) Example 2. Consider the following simple cascaded linear system
(i:) (: $1;) =
+u):(
+
(:).
(62)
where u is an uncertain parameter with a mean of unity and a standard deviation of 0.2. A normal distribution for u is assumed. The disturbance d is modeled as
d=-d+w (64) where w is a zero mean, white noise stochastic process. The controller is to be designed by placing the closed-loop poles so that constraints on the variances of the states and the input energy are satisfied. A two-state controller is employed which can be expressed in partitioned form by
(54)
where K is a positive real scalar and A satisfies
where the parameter X denotes the position of the closed-loop poles under nominal conditions. It is noted that this controller contains integral action.
1700 Ind. Eng. Chem. Res., Vol. 31, No. 7, 1992 1
0.9
0.9 0.8 0.7 0.6 0.5 0.4
0.3
O.q//.
.
.
0.01
I
.
.
.
,
0.03 0.04 9 Figure 8. Probability distribution of the variance for xl. O O
0.02
0.2 0.1
j
0
0.05
9
Figure 10. Probability distribution of the variance for u.
1
0.9 0.8
0.9 0.8
0.7 0.6
0.6
0.7
03
0.5 0.4
0.4
::;I/ 0.3
.
.
,
.
0
0.2
0.3
0.4
O.S
0.6
0.7
1
""II 0.3
O
0.8
9 Figure 9. Probability distribution of the variance for x2.
A statistical analysis of the system variances was conducted for X = 0.1 and 0.05. In Figure 8,the probability distribution (PD) expressed by the bound on the variance of state variable x1 is shown. Note that the PD is evaluated with respect to the model uncertainty term u while the variance is evaluated with respect to the PD for the noise w(t). The sensitivity of the variance bounds can be analyzed through these results. It is seen that the performance as measured by the variance bound on x1 degrades relatively rapidly for low likelihood (or worst-case) situations. In this context, worst case is taken to be results that are outaide of the 0.95 probability level. In addition, the difference in the variance bound at the 0.95 probability level is seen to be a factor of 3 greater for X = 0.1 than for X = 0.05. The PD expressed by the bound on the variance of state variable x 2 is shown in Figure 9. Note that the bound on the variance for x2 is smaller for X = 0.1 than for X = 0.05; while the opposite is true for xl. The PD expressed by the bound on the variance of the input u shown in Figure 10 exhibits a shape similar to the variance for the controlled variable xl. The linear bounding function is then applied with the results for the variance bound on x1 shown in Figure 11. These results are obtained by determining at each node the optimal value of K. Hence, this plot presents the best result that can be achieved when the linear bounding function is used,since, in a normal situation, a single value of K is used to determine the variance bound. It is seen by comparing Figure 11to Figures 8 and 9 that a fairly accurate bound is obtained for high likelihood cases. However, for lower likelihood situations (greater than 0.6), the h e a r bounding function fails to yield adequate results for this particular example.
Frequency Domain Formulation In this section, the probabilistic approach to robust control is extended to frequency domain performance measures. Single-input, single-output transfer functions, denoted by G(s;8) are considered in which the uncertain
.
.
O
0.1
.
.
.
. 0.3
0.2 9
.
1 0.4
Figure 11. Probability distribution of the variance for x1 when bounding functions are used.
parameters are taken to be the p-dimensional vector, 8. The probability that the parameter vector belongs to a certain set, 8 C JP, is given by Pr
(e E e) = Jf(8) 0
d8
(66)
where f ( 8 )is the probability density function (PDF) for
e.
The standard approach for describing model uncertainty in the frequency domain is to place bounds on each element of 8, a representation of uncertainty that corresponds to a hyperrectangle in 9. In the proposed probabilistic approach, this type of representation is equivalent to assuming a uniform PDF for each element of the p-dimensional vector 8, + abi, where 8, is a nominal value for 8 and Ibi, where bi bounds the uncertain element. One alternative probabilistic description is given by choosing a Gaussian PDF for 8. This latter representation of model uncertainty requires only the mean and the variance, V, for 8 to completely specify all other statistical moments. More complex probability descriptions of model uncertainty could also be used. For example, many parameters in process control problems cannot have negative values; an appropriate probability distribution can be defined for this and other situations. Probability density functions of model uncertainty in the frequency domain can be generated from the PDF for 8. In addition, a model uncertainty set R&b;8,u) at each frequency can be constructed from the mapping G(jw;B) such that Pr [G(jw;8) E RG(jw;8,v)Vo]= Y (67) The union of these seta over frequency is denoted by URc(jw;8,u). This set can be graphically depicted in a Nyquist diagram as the union of uncertainty regions at each frequency which contains G(jw;8) with probability Y for 8 E 8. In general, one can either specify Y or 8 so as to generate the set Rc(jw;f3,u). If 8 C YP is specified, then the
Ind. Eng. Chem. Res., Vol. 31, No. 7, 1992 1701 probability Y is determined from eq 66. This approach can be used to generate a model uncertainty set for robust stability analysis. The second approach for generating the set &(jw;e,v) uses the PDF f(e) to determine a parameter region 8 C RP such that Pr (e E e) = Y (68) The mapping G(jw;O) is then evaluated for all 8 E 8 to construct URG(jw;8,v). This approach can be used to construct a model uncertainty set for robust performance analpis. It is noted that numerous parameter uncertainty sets 8 satisfy eq 68 and thus there is more than one U&(jc@,v) that satisfies eq 67. For the general situation, generation of the set 8 to satisfy eq 67 is dependent on the PDF for 8. Care should be taken when developing 8, especially for a Gaussian PDF. This uncertainty description is unbounded and thus can admit values of B that have no practical importance. In this situation, a more meaningful PDF can be achieved by truncating the Gaussian PDF at a reasonable value while still retaining the same mean and covariance and imposing a suitable normalization of the resulting PDF. After 8 has been determined from either of these two approaches, the function G(jw;B) maps the set 8 into a complex plane “template” at each frequency (69) &(jw;8,v) = (RGGw;e)IeEe) e
c
For a noncompact 8 , model parameters outside the set 8 may be mapped by &(jw;f?) into URG(jw;e,v). The set URG(jw;e,v) will then enclose G(jw;B) with probability greater than Y. Thus for the general case Pr [G(jo;e)E & ( j w ; 8 , v ) V w ] 2 V (70) In addition, the model uncertainty set URG(jw;e,V)can exhibit an intricate shape in a Nyquist diagram and be difficult to express quantitatively. This uncertainty description can be approximated at each frequency as a disk centered at G(jo;B,) for a nominal value 0, (set equal to 60) with the radius given by L,(w;v)= SUP IP(jw) - G(jw;B,)I (71)
Po’,)
where PO’,) E RG(jw;B,v). The multiplicative representation of uncertainty is then given by
The model uncertainty disk at frequency w generated by L,(w;v) and G(jw;8,) will be denoted by DG(~w;u).
This proposed probabilistic approach to robust control uses an explicit separation of robust stability and robust performance. Roughly speaking, robust stability is ensured for one set of perturbations to the nominal model, while performance objectives are maintained for a second, inclusive set of perturbations to the nominal model. Probability levels are used as a guideline to construct these sets. Relationships for robust stability with probabilistic uncertainty descriptions expressed in the frequency domain can easily be developed from deterministic uncertainty descriptions. For example, the condition for robust stability with probability a,that is, VG(jw;B) E UD&;a), can be expressed in terms of an =-norm relation (Morari and Zdiriou, 1989) IIT(jo;e,) ~,(j~;~)ll,= sup ITGw;e,)L , ( ~ ~ ; < ~ )1I (73) w
where the T(jw;B,) is the complementary sensitivity
function at nominal conditions. This relation provides a sufficient condition for robust stability with a probability of at least a. Now a value of a less than 1implies that at most a 1- a probability of closed-loop instability exists. In some situations, a low probability of instability may be tolerable in order to achieve significantly higher performance from the control system. To guard against the improbable unstable situation, a fault detection and action strategy can be employed at the supervisory control level. This approach is very different from the overly conservative methodology that is currently employed. Robust performance has traditionally meant that certain performance specifications are met for the entire model uncertainty set (Morari and Zafiriou, 1989). In a probabilistic context, robust performance implies that specifications are maintained with probability 1or VG(jw;B) E U&(jw;e,v) with = 1. Although any performance criterion can be incorporated within the proposed approach, the m-norm with a disk representation of uncertainty will be considered here. It is easy to show that the closed-loop system achieves robust performance with probability 8, or VG(jw;B)E UDG(jw;@), if the nominal closed-loop system is stable and SUP IT(jw;e,)IL,(w;8) + Is(jw;e,) w(jw)l < 1 (74) LJ
w
where S(jo;e,) is the sensitivity function evaluated at nominal conditions and wGw) is a bound on the maximum peak of the sensitivity function. Some properties of the probabilistic robustness measures will now be summarized. Property 4. Robust performance with an H, objective is equivalent to simultaneously satisfying eq 74 with 8 = 1 and robust stability with probability 1. Property 5. For an H, performance objective and RG(jO$,@), robust performance with prob&(j@,a) ability /3 implies robust stability with probability a if a =
8.
Property 6. For an H, performance objective and a specified probability level a of robust stability, there exists a limiting value (assuming that PI < a) such that V8 < &,robust performance with probability 8 is satisfied only if robust stability with probability a is sacrificed. Hence a reasonably conservative approach to probabilistic controller design would be to guarantee robust stability with probability 1 (a = 1) and then use the probabilistic robust performance condition, expressed through 8, as a design variable. The effect of 8 on controller performance can then be analyzed. The controller is then designed using the uncertainty bounds as the principal analysis and turning parameters. The major tradeoffs in the design of such a probabilistic-based robust controller involve the closed-loop performance under high likelihood (near nominal) conditions against performance under low likelihood (near worst-case) conditions. In all cases, stability is guaranteed. Example 3. A simple example is presented to demonstrate the use of a probabilistic approach in the frequency domain. A linear model is used to exemplify the techniques described earlier. These techniques can also be used for nonlinear problems if a probabilistic model identification method (Goodwin and Salgado, 1989) is employed. Consider a classical linear CSTR problem with a firstorder reaction carried out isothermally where the reactor outlet concentration, C A ( t )is controlled by a pure inlet stream with concentration C~U.The flow rate of the pure manipulated stream F&) is much smaller than that of the primary inlet stream F with concentration CAI. The vol-
1702 Ind. Eng. Chem. Res., Vol. 31, No. 7,1992
:'i
1.2 1.0
0 H
2
,
3
9
2.8 10-~01h) Figure 12. Probabilistic uncertainty bounds. 1.4
2.1
3.5
0, 2
4.2
-
FX
Table I. Numerical Values of the Parameters Darameter value Darameter value V 2.8 x 103 L u2k 2.5 X h-2 G 5.7 x 102 L F,, 3.7 X lo3 L/h CF 2.8 X lo3 LJh Fmi. 2.0 x lo3 L/h Ck 2 h-' k,, 2.6 h-' u2F 5.2 X 10' L2/h2 kmin 1.4 h-'
ume of the reactor, V, is held constant by an overflow line. Thus, the component balance can be expressed as v dCA = CAiF + hfA(t) - CA(t)F- kCA(t)V (75) dt where k is the reaction rate coefficient and hfA(t) = CA,$&). A variable time delay is added to the measurement of cA(t) equal to $F'where $ is a constant. Uncertainties are specified for k and F. Consequently, uncertainty is present in the gain, time constant and time delay of the system. A Gaussian PDF is specified for both k and F. In order to avoid physically unrealistic values of k and F for this unbounded PDF, upper and lower bounds are placed on both: k E [kmin,kmax]and F E [Fmi,,Fm,]. The numerical parameters chosen for this example are presented in Table I. It is desired that the closed-loop system achieve robust stability with probability j3 = 1;hence closed-loop stability is required for all model parameters within the rectangular region shown in Figure 12. Also the controller must achieve the performance objective of no offset for step changes in the set point. The function l/w(jw) used to bound the maximum peak of the sensitivity function is chosen to be 2.5 for all frequencies. The performance of the controller as a function of model uncertainty can now be analyzed. First, ellipsoids for a x2&2)distribution are used to generate a joint probabilistic uncertainty region for k and F. These ellipsoids are obtained from the relationship (K - 212 ( F - 2.8 x 103)~ 57 (76) 0.025 5.2 x 104 where y is obtained from the x2@(2)distribution as a function of the probability level 8. The ellipsoid that corresponds to B = 0.7 is also shown in Figure 12. The parametric uncertainty regions are then transformed to the frequency domain. The corresponding uncertainty regions for j3 = 0.7 and j3 = 1are shown in Figure 13. The frequency domain uncertainty regions are approximated in the form of disks using eq 72. The disk radii for 6 = 0.7 and B = 1 (or a! = 1)are shown in Figure 14 as functions of frequency. It is noted that there is a great loss of accuracy (or increase in conservativeness) in approximating the uncertainty regions by disks. Consequently, the use of probability concepts to achieve better analysis methods loses some importance. However, the
-
+
-2
4
-3
-1
3
1
~
7
u
e~ ~( jxi IO-'
Figure 13. Probabilistic uncertainty regions in the frequency domain for @ = 0.7 and /3 = l along with the nominal model (solid l i e ) .
10-2
10''
100
10 fW-7
102 (M)
103
Figure 14. Disk radii of model uncertainty.
loss in accuracy at low likelihood regions is seen to be much greater than that at high likelihood regions. (Methods to minimize this problem include graphical approaches or the time-domain approaches described earlier.) The analysis concepts are used with a pole placement controller designed for the nominal conditions. For this example, the controller is given by 8.4 X 1O3(O.33s 1) CAS) = (77) As + 1 which is chosen to have a Smith-predictor type structure. The closed-looppole at nominal conditions is denoted by
+
A.
Parameter X was then tuned such that both the robust stability and robust performance criteria are satisfied for their respective uncertainty sets. Furthermore, at least one inequality for describing robust stability with probability a or robust performance with probability j3 reduces to an equality. A search showed that X must at least be to guarantee robust stability with probability 1. The value of A that converts the probabilistic-based robust performance inequality into an equality is shown in Figure 15 as a function of 0. The limiting value j31 is approximately 0.6. The characteristics of the controllers corresponding to 0 = 1 and 0 = 0.7 can be compared by analyzing their responses to a step change in the set point. For the nom-
Ind. Eng. Chem. Res., Vol. 31, No. 7, 1992 1703 1.2
14
-
2
.
0.7
1
0.9
0.8
1
‘“t
0.8 0.6
0.4
0.08 ‘0.6
1
--
1
O*’ 0 0
1
10
20
0.9 0.8
r-
1
I
Q
0.5 0.4
0.3
.
0.2 0.1
-
“ 0
60
1.7 1.6
l.4 8:: :u
50
40
Figure 18. Comparison of responses to a set-point change at conditions within the 0.95 confidence interval.
Figure 16. Filter parameter as a function of @.
A
30
1 (h)
P 1.1
1
n
: 4
1.5
3
1.4
1.4
1.3
1.3
2
-
8 10
20
30
40
50
60
=:2 +
2 9
1.2
I(h)
Figure 16. Comparison of responses to a set-point change at nominal conditions. 1.1 I
1
1.1 1.0
1.0
1.2
1.4
1.6
1.8
2.0
J(0.2)
J(G;
