8112
Ind. Eng. Chem. Res. 2007, 46, 8112-8119
Assessing Model Prediction Control (MPC) Performance. 2. Bayesian Approach for Constraint Tuning Nikhil Agarwal and Biao Huang* Department of Chemical and Materials Engineering, UniVersity of Alberta, Edmonton, AB, Canada, T6G 2G6
Edgar C. Tamayo Syncrude Canada Ltd., Fort McMurray, Alberta, Canada, T9H 3L1
Performance assessment of model predictive control (MPC) systems has been focusing on evaluation of the variability with, for example, minimum variance or LQG/MPC tradeoff curve as benchmarks. These previous studies are mainly concerned with the dynamic performance of MPC. However, the benefit of MPC is largely attributed to its capability for economic optimization. The economic performance, on the other hand, is also dependent on the variability reduction achieved through dynamic control. There is a need to assess MPC performance by considering economic performance, variability reduction, and their relationships. One of the good indications of this relation is the constraint tuning. In practical MPC applications, the constraint setups are important whenever an MPC is commissioned, and constraint tunings are not uncommon, even when the MPC is already on-line. Thus, the questions to ask are which constraints should be adjusted, and what is the benefit to do so? By investigating the relationship between variability and constraints, problems of interest are solved under the Bayesian inference framework (namely, through the Bayesian approach for decision evaluation and decision-making). The decisions that are referenced are whether to tune the constraints to achieve the optimal economic MPC performance and which constraints should be tuned. A detailed case study for a distillation column MPC application is provided to illustrate the proposed performance assessment methods. 1. Introduction Process control has helped industries to improve the efficiency of the operations and reduce the losses. A control system for any process is required to satisfy three general classes of needs:15 (i) suppress the influence of external disturbances, (ii) ensure stability of the process, and (iii) optimize the performance of the process. Satisfying the process specifications is the key operational objective for any process. After the objective is satisfied, the next operational goal is to optimize the process and to make the operations more profitable. The process operations, to a great extent, are governed by the market forces of demand and profits; therefore, the operating conditions must change as the demands and the product price patterns change. However, the changes must be made in such a manner that the economic objective function is always optimized. The model predictive control (MPC) system performs the job of economic objective function optimization, on top of the dynamic optimization. The success of the MPC system relies on the accuracy of the model of the process for good control. Using the process model, the controller predicts the behavior of the dependent variables (the output variables, or the controlled variables (CVs)) of a dynamic system, with respect to the changes in the independent process variables (the input variables, or the manipulated variables (MVs)). The current process values and the process model are used to predict the future movements of the MVs that will result in the operation of the process honoring all the constraint limits defined for the CVs and MVs. Today, many commercial controllers are available, each with its own control philosophy11,12 and economic objective function, such as a linear objective function through linear programming * To whom correspondence should be addressed. Tel.: +1-780-4929016. Fax: +1-780-492-2881. E-mail address:
[email protected].
for optimizing the process operations14 or a quadratic function, along with the linear objective function, to optimize and maintain the process variables at the desired values.5 Having an MPC controller installed for a process does not guarantee optimization of the process operations. For optimization, the controller is required to be tuned. The handles available to tune the controller are the control and prediction horizon, the step size for the moves on the MVs, etc. Tuning of the controller, with these parameters, requires a thorough understanding of the behavior of the process and the control philosophy of the MPC application being used. The controllers are tuned with these parameters at the design stage, and it is not recommended to change them on day-to-day basis. The performance of the controllers is also affected by the constraint limits given to the CVs and MVs. Providing the controller with incorrect and/or conflicting constraint limits on the process variables is one of the factors that can prevent the controller from operating the process at the optimum operating point. According to Singh and Seto,13 “It is a common phenomenon that a lot of APC benefits are lost because the operators tend to over-constrain the CVs/MVs. The reasons for over-constraining may be many, such as poor APC models/tuning, which causes the controller to swing. Hence, the operator tends to tighten the MV limits in order to stop the swing. Another reason can be the conservative nature of operator. He or she may just feel more comfortable if the APC does not have so much leeway to move around the MVs. In any case, such an action can lead to significant lost benefits.” Because of the fact that, in daily operations of the process plant, it is common to change these parameters, the controller performance is definitely affected when they are changed. Thus,
10.1021/ie0704694 CCC: $37.00 © 2007 American Chemical Society Published on Web 11/01/2007
Ind. Eng. Chem. Res., Vol. 46, No. 24, 2007 8113
The remainder of the paper is organized as follows: Bayesian inference and network are introduced in section 2. Calculation of optimal operating points and variability distributions is given in section 3. Section 4 derives the Bayesian methods for MPC constraint analysis and tuning. The creation of the Bayesian network is illustrated in section 4.2, followed by two case studies: one for a binary distillation column (in Section 5) and the other for an industrial case study (in section 6). Conclusions are given in section 7. 2. Bayesian Inference
Figure 1. Typical Bayesian network.
it becomes imperative to assess the effect of any decision made with regard to constraint limit changes. Therefore, a methodology is required to assess the change in the controller performance that is due to these changes. In Part 1 of this paper,2 the algorithms for assessing potential of MPC economic performance and probabilistic optimization have been developed. In this part, the following two problems are considered: (1) If a decision is made regarding changes of constraint, what is the expected return? This is a decision analysis problem. (2) If a specific economic profit is desired, which constraints should be changed? This is a decision-making problem. For Bayesian analysis purposes, the following information is required: the routine operating data, process steady-state gains, the profit/loss terms associated with the CVs to be within and outside the constraint limits, and other related process information, such as which CVs and MVs are allowed to change their limits and the preference to change them (prior probabilities). The algorithm then performs optimization on the various combinations of making the changes to obtain new optimum operating points. The optimization results are then used to establish the subsequent Bayesian decision analysis/making network. A Bayesian network for the MPC application, under consideration, is created with all the variables (CVs and MVs) for which the changes can be made and quality variables that reflect economic performance. The probabilities estimated through the various optimization purposes are used to form the Conditional Probability Distribution Table (CPD or CPT) for the network.4,10 The decisions regarding the changes are the evidence for the analysis. The analysis then provides the probabilities for the CVs to be in different locations, which are then used to make an assessment of the performance of the controller. The network can also be used to obtain guidelines for making the decisions, to obtain the desired performance level of the controller. The contribution of this paper can be summarized as follows: (1) A framework for Bayesian analysis of the decisions related to the constraint changes for the CVs and the MVs, using a probabilistic optimization function, is created. (2) Guidelines for making constraint changes, to achieve target controller performance, are derived according to the statistical inferences.
2.1. Bayesian Analysis. Named after Thomas Bayes, Bayesian analysis is a branch of statistical inference that can be applied for decision-making and statistical analysis, using the knowledge of prior events to predict future events. The Bayes theorem forms the backbone of Bayesian analysis. It enables calculation of the conditional probabilities for a hypothesis4,16 and is also known as the principle of inverse probability. Probability for a hypothesis A conditional on a given evidence B is the ratio of probability of the conjunction of A and B to the probability of B,4 i.e., P(A|B) )
P(A,B) P(B|A) × P(A) ) P(B) P(B|A) × P(A) + P(B|¬A) × P(¬A)
(1)
where P(A) is the prior probability of occurrence of hypothesis A (also known as prior). P(B|A) is the likelihood of obtaining evidence B, given hypothesis A is true, and P(A|B) is the posterior probability of A to be true, given the evidence B that is obtained. The term ¬A represents the case when the hypothesis A is not true. 2.2. Bayesian Networks. For a system that is comprised of more random variables than A and B, a network connecting all variables can be built. This network, which represents the relationship between the various random variables, is called a Bayesian network. A Bayesian network is defined as4 “a graphical structure that allows us to represent and reason about an uncertain domain. The nodes in the network represent a set of random variables.” A pair of nodes are connected through directed arcs that represent a relationship between the nodes. The node through which the arc originates is called the parent node and the node where it terminates is called the child node. The nodes in a Bayesian network are the variables of interest, and the link between them represents the probabilistic dependencies among the nodes.10 To specify the probability distribution of a Bayesian network, prior probabilities are to be defined for the root nodes, i.e., the nodes with no predecessor and the Conditional Probability Distribution Table (CPD or CPT) is defined for all non-root nodes, for all possible combinations of their direct predecessors.3 The CPT quantitatively represents the relationship between the parent and the child nodes. A typical Bayesian network is shown in Figure 1. Node C has two parent nodes (A and B) and one child node (D). Nodes A and B have two states (1,2), node C has three states (X,Y,Z), and node D has two states (P,Q). The tables beside nodes A and B in Figure 1 are their prior probabilities, and those tables beside nodes C and D represent their CPT. A Bayesian network cannot have directed cycles; i.e., a node cannot be reached again by following the directed arcs from the child nodes directly. Thus, Bayesian networks are also called directed acyclic graphs (DAGs). Using DAGs, the parameters can be represented as nodes or random variables and can be associated with a prior distribution. DAGs, which include decision and utility nodes, as well as chance nodes, are known
8114
Ind. Eng. Chem. Res., Vol. 46, No. 24, 2007
as influence diagrams or decision networks and be used for optimal decision-making.8 The CPD/CPTs quantify the dependencies between the nodes. The CPDs are probability distributions P(xi|Pai), where xi is the ith node and Pai represents all of its parent nodes.8,9 There are three types of nodes, two of which are used in this paper: (1) Chance Nodes: These nodes represent random variables and are associated with CPT or prior. (2) Utility Nodes or Value Nodes: These nodes represent the value of the utility function (benefit function). The parents for these nodes are the nodes whose outcome directly affects the utility. These nodes are associated with the utility table, with the value for each possible instance of its parents perhaps including an action taken.4 Thus, if xj is the set of observed variables, xi the set of variables whose values we are interested in estimating, and xk the set of variables in the system (not included in xi and xj), then inference in a Bayesian analysis means to compute
P(xi ) a|xj ) b) )
P(xi ) a,xj ) b)
)
P(xj ) b) P(xi ) a,xj ) b,xk) ∑ x k
(2)
P(xi,xj ) b,xk) ∑ x ,x i k
where P(xi ) a|xj ) b) is the probability for node xi to take value a (provided that node xj takes the value b), P(xj ) b) the probability for node xj to take value b, and P(xi ) a,xj ) b) the unconditional probability of conjunction of xi and xj. The probability distributions form the basis for the statistical analysis of the data. In this work, we are limited to using the Gaussian distribution to represent or approximate process variability; other distributions are also possible for the Bayesian network, in principle. For a random variable y, the Gaussian distribution is represented as
p(y) )
1
x2πσ2
[
exp -
(y - yj)2 2σ2
]
(3)
where y is the data with a mean yj and a standard deviation σ. Plotting p(y) against y gives the probability density function, and plotting the integral of p(y) in the range of -∞ to y gives the cumulative distribution function. The cumulative distribution function (cdf) can mathematically be represented as
cdf(y) )
[
( )]
y - yj 1 1 + erf 2 σx2
1 xπ
∫0y
q
J)
6
∑ ∑ P(yi ∈ Ωk) × F(yi ∈ Ωk) i)1 k)1
(6)
Let (yji0,ujj0) be the base-case mean operating points and let (yji,ujj) be the optimum operating points, where the base-case operating points must move by (∆yi,∆uj) to reach the optimal operating points. The equality constraints to be satisfied for the economic optimization are given by n
∆yi )
[Kij × ∆uj] ∑ j)1
(7)
yji ) yji0 + ∆yi
(8)
ujj ) ujj0 + ∆uj
(9)
where i ) 1, 2, ..., m and j ) 1, 2, ..., n Considering the acceptable limit for constraint violation for the output variables to be 5%,6,7 a set of inequalities can be defined that also must be satisfied while optimizing the objective function defined in eq 6. For the limits change, the inequalities constraints are given by eqs 10 and 11. These inequalities define the new constraint limits for the CVs and MVs17 after the constraint changes:
Lyi + 2σi0 - yholi × ryi e yji e Hyi - 2σi0 + yholi × ryi
(10)
Luj + 2Rj0 - uholj × ruj e ujj eHuj - 2Rj0 + uholj × ruj
(11)
where i ) 1, 2, ..., m and j ) 1, 2, ..., n and yholi, uholj are half of the constraint range for yi and uj; ryi and ruj are the percentage change in the constraint limits for the process variables, and σi0 and Rj0 are the base-case standard deviation and the quarter of the range for yi and uj, respectively. The subscripts “yi” and “uj”’ represent the variable for the ith output and jth input variables, respectively. Because the constraint limits for any CV or MV cannot be changed by large values, in this paper, a change of 10% of the existing limit range for the changeable variables is adopted. Thus, for the constraint limit change, the optimization problem can be defined as
max
yj1, ..., yjmuj1, ..., ujn
J
(12)
s.t. eqs 7, 8, 9, 10, 11 (4) where J is given by eq 6.
where
erf(y) )
information and notations of this section have been presented in Part 12 in detail. Consider the following simplified economic objective function specified in Part 1:2
4. Bayesian Methods for MPC Constraint Analysis and Tuning
exp(-t2) dt
(5)
3. Calculation of Optimal Operating Points and Variability Distributions In this section, we shall consider how the means of the ith output (CV), which is denoted as yji, and the jth input (MV), which is denoted as ujj, are calculated by constraint tuning. The means yji and ujj are the optimal operating point calculated from the optimization, which is discussed next. The background
4.1. Building Bayesian Inference Network. The data and the optimization results thus obtained for all the cases of the limit change can now be used to build the Bayesian network for the process. The network comprises of the N variables available for making the changes as the parent nodes, all the q quality CVs as the child nodes, and one utility node, which represents the value of the benefit function for the process. For the constraint limits change, the parent nodes have two states (“change limits” and “do not change limits”). The child nodes in the network have six states (“Zone 1”, “Zone 2”, “Zone 3”,
Ind. Eng. Chem. Res., Vol. 46, No. 24, 2007 8115 Table 1. Coefficients L and Q, and Change Allowed for CVs CV
L coefficient
Q coefficient
allowed to change
1 2
R1 R1
β1 β2
yes no
Table 2. Coefficients L and Q, and Change Allowed for MVs MV
L coefficient
Q coefficient
allowed to change
1 2
0 0
0 0
no yes
Table 3. States for Parent Nodes for Limit Change Figure 2. Bayesian network for an m × n plant.
“Zone 4”, “Zone 5”, and “Zone 6”). The notation of states has been defined in Part 1.2 All parent nodes will be assigned prior, and all child nodes will be assigned CPT. The term “prior” defines prior probabilistic information about the parent nodes. If prior is set as uniform, the Bayesian inference will be equivalent to the likelihood estimation. This can be derived as follows. To calculate the posterior of parameter A, given the observation B, we have
P(B|A) × P(A) P(A|B) ) P(B) where P(B) is a constant, given B as an observed value. If A has a uniform prior, then P(A) is also a constant. Denote (P(A)/ P(B)) ) λ. Then,
parent nodes
MV2
CV1*
states
(yes,no)
(yes,no)
Table 4. States for Child Nodes for Limit Change child nodes
CV1
states
Zone 1, Zone 2, Zone 3, Zone 4, Zone 5, Zone 6
Table 5. Prior Probability for Making Limits Changes for Parent Nodes parent node
change
do not change
MV2 CV1*
0.5 0.5
0.5 0.5
Table 6. Possible Cases for Applying Limit Changes parent node
Case 1
Case 2
Case 3
Case 4
MV2 CV1*
no no
no yes
yes no
yes yes
Table 7. Optimum Operating Point for Each Identified Case
P(A|B) ) λP(B|A) where P(B|A) is a likelihood function. Thus, the Bayesian inference for the posterior is the same as the likelihood function in this case. The selection of prior is often not obvious unless it is uniform, but it is not a problem in the context of this work and will be explained shortly. The CPT determines the probabilities for the child nodes to be in each of the six zones or the states. For the constraint change, the 2N optimization results, when combined with the base-case standard deviation and assuming a Gaussian distribution, provides probabilities for CVs to be in each of the six zones, as illustrated in Figure 7 of Part 1.2 The utility node provides profit data for q quality variables in each of the six zones. From the profit data provided in the utility node, the expected return can be calculated according to eq 6 of Part 1 (or eq 7 of Part 1 if quality variables are independent).2 The Bayesian network thus created is shown in Figure 2. Pa, Ch, and U represent the parent nodes, the child nodes, and the utility node for the network. The network can now be used for decision evaluation and decision-making purposes. 4.2. Illustration on Building a Bayesian Network. To illustrate the procedure of the proposed method, consider four variables of a 2 × 2 system (i.e., a 2 CV and 2 MV system) for which the steady-state gain matrix (K) is represented as
[
k k K ) k11 k12 21 22
]
where CV1 is a constraint variable and CV2 is a quality variable. The economic linear and quadratic coefficients for the output and the input variables considered in the application (if available), and whether they are allowed to change the limits, are listed in Tables 1 and 2, respectively. Note that, in many applications, the economic linear and quadratic coefficients may
child node
Case 1
Case 2
Case 3
Case 4
CV1 CV2
jy11 yj21
jy12 yj22
jy13 yj23
jy14 yj24
not be available, but the maximum and minimum returns may be given directly. Interpolation among the six zones will provide returns of operating in each zone. Because there are two variables available for making the changes, the Bayes network that has been created for the aforementioned system will have two parent nodes (MV2 and CV1*) and one child node (CV2). To avoid notation confusion, we have used CV1* to represent the logical CV1, meaning that it takes only two values (“yes” or “no”), while the physical CV1 can take any of the real values. Since MVs are not considered as quality variables here, we omit the superscript * for MVs. Namely all MVs in the network stand for logical variables taking values ‘Yes’ and ‘No’ only. The parent and child nodes, along with their respective states, are described in Tables 3 and 4. The prior probability for both MV2 and CV2* is assumed to be 0.5, which means that there is no preference to make or not make a change (i.e., the prior is uniform). Table 5 defines the prior probabilities for the parent nodes. For the aforementioned system, there will be 22 ) 4 cases for which the optimizations are performed with these combinations of limits changes for (MV2, CV1*): (no, no), (no, yes), (yes, no), (yes, yes), according to eq 12, as illustrated in Table 6. The optimum operating points are then obtained for each case of the change limits, and these values are listed in Table 7. In this table, yjip is the optimal operating point for the ith CV and the pth case. However, because of variability, only the average operating point can be at the optimal value. The actual data can be distributed in any of the six zones. The probability for the ith CV to be in each of the six zones is estimated, assuming the data to be in a Gaussian distribution, with the mean being the optimal operating point calculated from optimization and
8116
Ind. Eng. Chem. Res., Vol. 46, No. 24, 2007 Table 10. Posterior of Child Nodes node
Zone 1
Zone 2
Zone 3
Zone 4
Zone 5
Zone 6
CV2
γ21
γ22
γ23
γ24
γ25
γ26
Table 11. Variables Available for Change CVs
MVs
change limit? Figure 3. Bayesian network for the considered 2 × 2 plant. Table 8. Conditional Probability Table for the Child Node (CV2)
Case 1 Case 2 Case 3 Case 4
Zone 1
Zone 2
Zone 3
Zone 4
Zone 5
Zone 6
P11 P21 P31 P41
P12 P22 P32 P42
P13 P23 P33 P43
P14 P24 P34 P44
P15 P25 P35 P45
P16 P26 P36 P46
U node
Zone 1
Zone 2
Zone 3
Zone 4
Zone 5
Zone 6
F1
F2
F3
F4
F5
F6
variance being that calculated from the base operation. These probabilities form the CPT for the network (see Table 8), where Ppk is the probability for CV2 to be in the kth zone for the pth case. In this example, because the profit function is defined only for the CV2, it is the only variable that affects the value for the expected return. The return values taken by the value node in the six zones or states are specified in Table 9. The returns in this table could be used to calculate the return according whichever zone that the mean operating lies or the expected (average) return, according to the posterior calculated. The nominal expected return for the existing system is calculated using eq 7 of Part 1,2 according to the profit values of each zone given in Table 9 and the probabilities for the variables to be in six zones (Case 1 of Table 8). With CV2 as the child node, MV2 and CV1* as parent nodes, as mentioned previously in Table 5, and the CPT for CV2 (as given in Table 8), the Bayesian network for the system is created as shown in Figure 3. The network is now ready to perform decision evaluation and decision-making. 4.2.1. Decision Evaluation. For the purpose of decision evaluation, the network is provided with the decision to be taken, which is equivalent to the evidence in the network. The evidence is then utilized to estimate the posterior probabilities for the child node to have values in the six zones, which are then used to estimate the expected return. Because, in this case, the parent nodes will be assigned values according to the decisions made, the prior of the parents nodes is not needed. As an example, if the decision is made to change the limits for MV2 only, then the probability of CV2 to be in any of the six zones is estimated using Bayes’ Theorem. For illustration purposes, one calculation is shown below:
γ21 ) P(CV2 ) Zone1|MV2 ) yes) )
P(CV2 ) Zone1,MV2 ) yes) P(MV2 ) yes)
∑
)
P(CV1*,CV2 ) Zone1,MV2 ) yes)
CV1*
∑
no yes no no yes yes yes no no no
change limit? MV1 MV2 MV3 MV4
yes yes yes no
joint probability P(CV1*,CV2,MV2) can be calculated as
Table 9. Utility Table with Profit Values CV2 state
CV1 CV2 CV3 CV4 CV5 CV6 CV7 CV8 CV9 CV10
P(CV1*,CV2,MV2 ) yes)
CV1*,CV2
where, according to Figure 3, using the multiplication rule,4 the
P(CV1*,CV2,MV2) ) P(CV1*)P(MV2)P(CV2|CV1*,MV2) The probability for CV2 to be in any of the six zones is similarly estimated and the values are given in Table 10. Based on Tables 9 and 10, the value of the expected return is then estimated using eq 7 of Part 1.2 Thus, comparing to the nominal expected return, the decision to change the limits for MV2 can be evaluated. The calculation can be done using many available Bayesian inference software programs, such as the Bayesian Net Toolbox (BNT, http://bnt.sourceforge.net/). 4.2.2. Decision-Making. For decision-making purposes, the network is provided with the target value of the return. Thus, reading from the table for the utility node, the value closest to the target value is selected and the corresponding state of the child node acts as the evidence for the analysis. For decisionmaking, the prior does matter. However, when we create the Bayesian network, optimizations are performed for all possible values of the parent nodes without any preference. Thus, the prior should be selected as uniform, and, consequently, the Bayesian inference yields a maximum likelihood estimation. If the target is set to have a specific expected return, then, from Table 9, the state for CV2 to be in a specific zone (e.g., Zone 2) is the evidence for the evaluation purposes. Thus, with CV2 to be in Zone 2 as evidence, the maximum a posteriori explanation (equivalent to the maximum likelihood estimation in the case of uniform prior) for a decision on whether to change the limits for parent node CV1*, node MV2, or both can be made. 5. Case Study of a Simulated Distillation Column Consider the MPC application for the simulated binary distillation column shown in Figure 9 of Part 1.2 The MPC controller designed for the process has 10 controlled variables (CVs) and 4 manipulated variables (MVs). The list of the CVs and the profits associated with them to be under-spec (u/s), inspec (i/s), and over-spec (o/s) are listed in Table 1 of Part 1,2 and the list of MVs is given in Table 2 of Part 1.2 CV2 and CV8 are two independent quality variables; as observed from Table 1 of Part 1,2 they are to be minimized. The process variables that are available for making the constraint changes are listed in Table 11, and the optimization problem can be defined as eq 13.
Ind. Eng. Chem. Res., Vol. 46, No. 24, 2007 8117 6
max
yj1,...,yj10uj1,...,uj4
J)
∑ ∑ P(yi ∈ Ωk) × F(yi ∈ Ωk) i)2,8 k)1
(13)
s.t. eqs 7-11 Based on the information provided in Table 11, the Bayesian network for the system can be created with 7 parent nodes (MV1, MV2, MV3, CV2*, CV5*, CV6*, and CV7*), 2 child nodes (CV2 and CV8), and 1 utility node. The parent nodes and prior probabilities are listed in Table 12. Corresponding to the 7 parent nodes for the network with two limit change states, 27 ) 128 optimizations were performed. The 128 optimization results obtained for all the child nodes were used to create the CPT for the child nodes. As for the process under consideration, CV2 and CV8 are the two variables that affect the overall economic performance of the operations; the utility node for the process will have 62 ) 36 values for various combinations of the states of these CVs. The expected (average) return from the base-case process operation was 141.59 units, and this was used as the basis for evaluating the decisions made for applying limit changes to the controller. The network was then used for decision evaluation and decision-making purposes. 5.1. Decision-Making. For the process under consideration, if the target was set to achieve return of 170.00 units, then the states (or the locations) of all the key process variables (i.e., the CVs that affect the utility node (CV2 and CV8)) were determined (see Table 13). With these states as the evidence, the maximum a posteriori states (maximum likelihood estimate) for the parent nodes that were to have their limits changed were calculated. The inference results indicate that the parent nodes 2, 3, 4, 5 and 6 (i.e., MV2, MV3, CV2, CV5, and CV6) were expected to have their constraint limits changed by 10%. The results are illustrated in Figure 4 where the variables are numbered according to sequence shown in the first column of Table 12. To check the validity of the results, we actually applied the suggested constraint changes to the MPC controller and used probability optimization to determine the optimal operating points of the new MPC control. The MPC was simulated again for the changed conditions. The average return for the new basecase process operation was then calculated as 176.05 units, which is not too much different from the targeted return (170.00 units). 5.2. Decision Evaluation. For the decision to change the limits of MV2, MV3, CV2, CV5, and CV6, the maximum a posteriori estimate of the states of CV2 and CV8 were Zone 4 and Zone 3, respectively. Corresponding to this decision, the expected (average) return was estimated to be 174.17 units. A comparison of the expected return of the controller before and after the decision was made and is shown in Figure 5, where “Existing” means the expected return of the base-case operation and “Case Study” means the expected return potential when the decision to change/not change is made. Thus, it can be inferred that the decision to increase the limits for MV2, MV3, CV2, CV5, and CV6 will increase the expected return of the operations. 6. Case Study for an Industrial Process The process for the industrial case study is a coker plant, where the bitumen feed is upgraded to synthetic crude oil. The bitumen is thermally cracked in a fluidized-bed reactor. The light hydrocarbon products leave the reactor as overhead vapors through the scrubber. The scrubber overhead effluent is then routed to a fractionator, where it is separated into raw fuel gas, naphtha, light gas oil, and heavy gas oil. The naphtha is then
Figure 4. Tuning guidelines for the limits-change case.
Figure 5. Comparison of the expected return in a simulated distillation column. Table 12. Prior Probability for Making Constraint Changes for Parent Nodes change limits
do not change limits
0.5 0.5 0.5 0.5 0.5 0.5 0.5
0.5 0.5 0.5 0.5 0.5 0.5 0.5
MV1 MV2 MV3 CV2* CV5* CV6* CV7*
Table 13. States for Key Child Nodes CV
state
CV2 CV8
Zone 4 Zone 3
fractionated into various streams, such as C3 and lighter compounds, C4 compounds, and naphtha (C5 and heavier compounds). The separated C4 compounds and naphtha are recombined before they are exported to the naphtha hydrotreating units for further processing. Figure 6 gives a schematic for the process. The MPC application for the process discussed has seven CVs and seven MVs. The list of the CVs and the profits associated with them to be under-spec (u/s), in-spec (i/s), and over-spec (o/s) are provided in Table 14; the list of MVs is given in Table 15. CV2 reflects the product (naphtha) quality and is considered for economic assessment in industry. Thus, it is taken as the quality variable for this analysis. The in-spec profit associated with CV2 has a minimum of 1 unit to a maximum of 7 units. The profits for the four in-spec zones (i.e., Zone 2 to Zone 5)were estimated by interpolation in the four zones, with the profit being from 1 unit to 7 units. The process variables that are available for making the constraint changes are listed in Table 16, and the optimization problem can be defined as eq 14. 6
max
yj1,...,yj7uj1,...,uj7
J)
∑ P(y2 ∈ Ωk) × F(y2 ∈ Ωk)
(14)
k)1
s.t. eqs 7-11 Based on the information provided in Table 16, the Bayesian network for the system was created with five parent nodes
8118
Ind. Eng. Chem. Res., Vol. 46, No. 24, 2007
Figure 6. Coker plant. Table 14. Control Objective and the Assumed Profits for the CVs Profit CV
description
objective
u/s
i/s
o/s
1 2 3 4 5 6
naphtha fractionator level C4 in naphtha C3 in naphtha depropanizer flooding debutanizer pressure sponge column top temperature presaturator bottom level
constraint maximize constraint constraint constraint constraint
0 -50.00 0 0 0 0
0 1 to 7 0 0 0 0
0 -20.00 0 0 0 0
constraint
0
7
0
0 Figure 7. Comparison of the expected return in an industrial process.
Table 15. List of MVs MV
description
1 2 3 4 5 6 7
feed temperature depropanizer bottom temperature depropanizer feed flow SP lean oil flow to presaturator lean sponge oil flow to LGO stripper lean oil flow to depropanizer reboiler debutanizer reflux
Table 16. Variables Available for Constraint Change CVs
MVs
change limit? CV1 CV2 CV3 CV4 CV5 CV6 CV7
no yes no yes yes no no
change limit? MV1 MV2 MV3 MV4 MV5 MV6 MV7
no off off yes off off yes
Table 17. Prior Probability for Making Constraint Changes for Parent Nodes
MV4 MV7 CV2* CV4* CV5*
change limits
do not change limits
0.5 0.5 0.5 0.5 0.5
0.5 0.5 0.5 0.5 0.5
(MV4, MV7, CV2*, CV4*, and CV5*), one child node (CV2), and one utility node. The parent nodes and their prior probabilities are listed in Table 17. Corresponding to five parent nodes for the network with two limit change states (“yes” and “no”), 25 ) 32 optimizations were performed. The 32 optimization results obtained for all the child nodes were used to create the CPT for the child nodes. As for the process under consideration, CV2 is the variable that affects the overall economic performance of the operations; the utility node for the process will have 61 ) 6 values. The base-case nominal expected return from the process was 2.80
units, and this will be used as the basis for evaluating the decisions made for applying limits change to the controller. The network was then used for decision-making and decision evaluation purposes. 6.1. Decision-Making. For the process under consideration, the target was set to achieve a return of 5.00 units; the state (or the location) of the key process variable, CV2, then was determined to be Zone 4. With this state as the evidence, the maximum a posteriori states (the maximum likelihood estimate) for the parent nodes that were to have their limits changed were calculated. The result indicates that the parent nodes 1, 2, 3, and 4 (i.e., MV4, MV7, CV2*, and CV4*) were expected to have their constraint limits changed by 10%. 6.2. Decision Evaluation. On the other hand, for the decision to change the limits of MV4, MV7, CV2, and CV4, the maximum a posteriori estimate of the location of CV2 was Zone 4; the return to operate in Zone 4 is 5.00 units, and the expected (average) return for the mean operating point in Zone 4 was estimated to be 3.94 units. A comparison of the expected return of the controller before and after the decision was made is shown in Figure 7. Thus, it can be inferred that the decision to relax the limits for these parent-node variables will increase the expected return from the process. It is also noted that the maximum a posteriori estimate of the state for CV2 is the same as that obtained, while using the network for decision-making purposes (as noted previously). 7. Conclusion This paper has discussed the use of Bayesian analysis for performance assessment of model predictive control (MPC) controllers, using probabilistic optimization, with regard to constraints tuning. The probabilistic optimizer performs the optimization of the control operating points, taking into account the process variability. The Bayesian inference is exploited to evaluate the decisions regarding constraint changes for MPC and to obtain the guidelines for changing the constraint limits, to achieve the target return from the process. Two case studies have also been discussed, and the utility of the tool in process industry is illustrated. The results obtained from the application can be used for day-to-day monitoring of
Ind. Eng. Chem. Res., Vol. 46, No. 24, 2007 8119
the MPC controllers and for obtaining the guidelines of decisionmaking for constraint changes. The decision guidelines obtained from the algorithm can be applied to improve the overall controller performance and thus improve the returns from the process. Although we have focused on MPC economic performance assessment, in principle, the method proposed in this paper could be applied to any constrained control with an economic performance objective. Nomenclature E(R) ) expected return F(yi ∈ Ωk) ) return associated with yi in Zone k Hkyi ) high limit of the kth zone for yi J ) objective function K ) steady-state process gain matrix Kij ) steady-state gain for yi and uj Lkyi ) low limit of the kth zone for yi Luj ) low limit for uj Lyi ) low limit for yi m ) number of output variables n ) number of input variables P(¬A) ) probability for A not to be true P(A|B) ) probability for A to be true, given B P(A) ) probability for A to be true P(B|A) ) likelihood of B to be true, given A q ) number of quality variables ruj ) change in limits (as percentage of original range) for the uj high limit ryi ) change in limits (as percentage of original range) for the yi low limit Rj0 ) base-case quarter of range for uj ujj ) optimized mean operating point for uj yji ) optimized mean operating point for yi ujj0 ) base-case mean operating point for uj yji0 ) base-case mean operating point for yi uholj ) half of limits for uj yholi ) half of limits for yi uj ) jth input variable (MV) for MPC yi ) ith output variable (CV) for MPC Greek Symbols σi ) standard deviation for yi σi0 ) base-case standard deviation for yi Ωi ) yi state matrix
Acknowledgment This work is supported in part by Natural Sciences and Engineering Research Council of Canada and Syncrude Canada, Ltd. Literature Cited (1) Agarwal, N.; Huang, B.; Tamayo, E. C. Bayesian approach for constraint analysis of MPC and industrial application. Proceedings of DYCOPS, 2007. (2) Agarwal, N.; Huang, B.; Tamayo, E. C. Tamayo. Assessing Model Prediction Control (MPC) Performance. 1. Probabilistic approach for constraint analysis. Ind. Eng. Chem. Res. 2007, 46, 8101-8111. (3) Charniak, E. Bayesian Networks without tears. Artif. Intell. Mag. 1991, (Winter), 50-63. (4) Korb, K. B.; Nicholson, A. E. Bayesian Artificial Intelligence; Chapman and Hall/CRC: London, U.K., 2004. (5) Krishnan, A.; Kosanovich, K. A.; DeWitt, M. R.; Creech, M. B. Robust model predictive control of an industrial solid phase polymerizer. Proc. Am. Control Conf. 1998, 3386-3390. (6) Latour, P. H.; Sharpe, J. H.; Delaney, M. C. Estimating benefits from advanced control. ISA Trans. 1986, 25 (4), 13-21. (7) Martin, G. D.; Turpin, L. E.; Cline, R. P. Estimating control function benefits. Hydrocarbon Process. 1991, (June), 68-73. (8) Murphy, K. P. An introduction to graphical models. Technical Report, University of California, Berkeley, CA, 2001. (9) Murphy, K. P. How to use the Bayes Net Toolbox. Available via the Internet at http://www.cs.ubc.ca/∼murphyk/Software/BNT/bnt.html, 2004. (10) Pearl, J.; Russel, S. Bayesian networks. Technical Report, Computer Science Department, University of California, Berkeley, CA, 2000. (11) Qin, J.; Badgwell, T. An overview of industrial model predictive control technology. AIChE Symp. Ser. 1997, 93 (316), 232-256. (12) Qin, J.; Badgwell, T. A survey of industrial model predictive control technology. Control Eng. Pract. 2003, 11, 733-764. (13) Singh, P.; Seto, K. Analyzing APC performance. Chem. Eng. Progress 2002, 98 (8), 60-66. (14) Sorensen, R. C.; Cutler, C. R. LP integrates economics into dynamic matrix control. Hydrocarbon Process. 1998, 77 (9), 57-65. (15) Stephanopoulos, G. Chemical Process ControlsAn Introduction to Theory and Practice; Prentice-Hall of India, Pvt., Ltd.: New Delhi, 1997. (16) Tan, S. B. Introduction to bayesian methods for medical research. Ann. Acad. Med. Singapore 2001, 30, 444-446. (17) Xu, F.; Huang, B.; Akande, S. Performance assessment of model predictive control for variability and constraint tuning. Ind. Eng. Chem. Res. 2007, 46, 1208-1219.
ReceiVed for reView April 2, 2007 ReVised manuscript receiVed July 25, 2007 Accepted August 29, 2007 IE0704694
