I n d . Eng. Chem. Res. 1987,26, 555-564 Prausnitz, J.; Anderson, T.; Grens, E.; Eckert, C.; Hsieh, S.; O’Connell, J. Computer Calculationsfor Multicomponent Vapor-Liquid and Liquid-Liquid Equilibria; Prentice Hall: Englewood Cliffs, NJ, 1980. Renon, H.; Prausnitz, J. M. AZChE J. 1968, 14, 135. Rosenbaum, J. B.; McKinney, W. A,; Beard, H. R.; Crocker, L.; Nissen, W. I. Bureau of Mines R.I.7774, Salt Lake City Metallurgy
555
Research Center, Salt Lake City, UT, 1973. Sano, H.; Nakomoto, Y. Nippon Kaguku Zasshi 1968,89, 362. Schiff, H. Ann. 1866, 240, 125. Wilson, G. M. J . Am. Chem. SOC.1964,86, 127.
Received for review November 27, 1985 Accepted September 12, 1986
Evaluation of Schemes for Detecting and Identifying Gross Errors in Process Data Joseph Rosenberg, Richard S. H. Mah,* and Corneliu Iordache Department of Chemical Engineering, Northwestern University, Evanston, Illinois 60201
Hitherto the detection of gross errors has been focused on the u’se of three separate statistical tests. These tests sometimes give conflicting predictions, and none of them takes into consideration bound violations. In this paper we propose two new composite tests, DMT and EMT, which make use of more than one statistical test to reduce mispredictions and to account for bound violation. These tests were evaluated along with global and measurement tests (GT and MT). Although the composite tests are not restricted to linear constraints and a single gross error, the evaluation was carried out for mass flow networks with at most one gross error present in the data. The effects of various factors on the performance of the tests are summarized in the form of basic rules and guidelines. Guidelines for choosing an appropriate detection scheme are also developed. The composite tests are found to have superior overall performance. Before process data are reconciled to enhance accuracy and consistency, it is important to detect and eliminate gross errors contained in them. Hitherto the principal thrust of research has been in the development of statistical tests, namely, the global test (Almasy and Sztano, 1975), the constraint test (Mah et al., 1976), and the measurement test (Mah and Tamhane, 1982). As shown in a recent review by Tamhane and Mah (19851, these tests sometimes give conflicting predictions. Furthermore, if the measured values are upper and lower bounded and if the bounds are reliable, the reconciled and estimated values should not violate these bounds. Bound violation is an important signal pointing to the possible presence of one or more gross errors. The above-mentioned tests do not make use of this important information inherent in the bounds of measured values. In our attempt to incorporate these signals, we propose two composite tests which make repeated use of one or more statistical tests to enhance the power or reduce the probability of mispredictions (type I errors). As a comparison, we evaluated their performances against the measurement test, which has been previously evaluated by Iordache et al. (19851, and the global test. Seven flow networks have been chosen as benchmark problems in our evaluation. For simplicity, in the following treatment we shall assume that all process variables are mass flow rates, which are directly measured, and that the random measurement errors have independent normal distributions with zero mean and equal standard deviations. Under our assumptions, the basic model of process measurements in the presence of gross errors is y=x+w+g
(1)
where y is an (s X 1)vector of process measurements and x is an (s x 1) vector of true values of the measured
* Author t o whom correspondence should be addressed.
quantities. The vectors w and g represent the random error and the gross error contributions, respectively. The random vector w has an (s X s) diagonal covariance matrix Q,assumed known or estimated. In general, not all process variables are measured. Let u be an ( m X 1)vector containing unmeasured variables or measurements deleted from the set of measured variables. Then, the solution to the reconciliation problem is the set of estimates of vectors x and u satisfying the least-squares estimation problem min (y - x)Q-l(y - x) (2) x,u
subject to the material conservation constraints A~X + A ~ =u 0 (3) Assuming that none of the elements in vectors x and u are precisely known, eq 3 represents a system of homogeneous linear equations, where A, and A2 are (n X s) and (n X m) matrices of coefficients in the balance equations, respectively. The general strategy to solve the reconciliation problem with unmeasured parameters and linear constraints is to decompose the problem into two parts. First, reconciliation of a modified network involving only measured quantities is used to estimate the true values of those measured variables. Second, some of the unmeasured variables may be estimated from the reconciliation estimates and the original process constraints. This second estimation step is known as coaptation. The unmeasured variables which can be uniquely estimated are referred to as observable (Stanley and Mah, 1981). As shown by Tamhane and Mah (19851, a unique solution to the coaptation problem exists if and only if A2 is of full column rank, i.e., rank(A2) = q = m. To decouple the reconciliation and coaptation problems, the general approach is to construct a projection matrix, P, whose rows span the n-q dimensional null space of A2 (Crowe et al., 1983). Then PA2 = 0 (4)
0888-588518712626-05%$01.50/0 0 1987 American Chemical Society
556
Ind. Eng. Chem. Res., Vol. 26, No. 3, 1987
Premultiplication of eq 3 by P yields PA,x = B x = 0
(5)
The projection matrix is not unique. Diffsrent methods of construction have been reviewed by Rosenberg (1985). Equation 5 represents the projected process constraints which are subsequently used in the reconciliation step. The general reconciliation solution, giving the estimate of vector x and the adjusted vector 9 ,was obtained by Mah and Tamhane (1982). For our model, we have
3=f
= (I - QB~(BQB~)-'B)Y
(6)
Basic Relations There are two important quantities used in all statistical tests for gross error detection in the literature to date. The first quantity is known as the vector of balance residuals and is defined by r = By
(7)
In the absence of gross errors, the vector r has a normal distribution with zero mean and a covariance matrix given by
H = BQBT
(8)
The balance residual vector r and its covariance matrix H have been used to construct the global test of Almasy and Sztano (1975) and the constraint test of Mah et al. (1976). The global test statistic is defined as 7
= rTH-lr
(9)
and follows as x2 distribution with degrees of freedom equal to the rank of matrix B. A serial elimination process must be used to discard measurements in order to locate the sources of possible gross errors. In constraint tests, the normally distributed balance residuals are standardized to form standard normal test statistics for each constraint in the matrix B,
vj = (Hjj)-l/zrj for 1 I j I n - q
(10)
The test statistics qj can be used in hypothesis testing to infer which rows in the constraint matrix B are violated due to gross error corruption. Mah et al. (1976) proposed a graph-theoretical search algorithm using constraint (nodal) test statistics. Another important quantity is the vector of measurement adjustments which is the difference between the measured parameters and their adjusted (or reconciled) values, i.e.,
a=y-f
(11)
A transformed form of this vector d = &-'a (12) has been used by Mah and Tamhane (1982) to define the measurement test for any nonsingular matrix Q. The vector d is also normally distributed, with zero vector mean and a covariance matrix,
W = BT(BQBT)-'B (13) By standardizing the transformed measurement adjustments, one gets the measurement test statistics, z, =
(Wii)-l/zdj for I I i I s
(14)
which are standard normal variables. Hypothesis testing is used to infer locations of gross errors. If lzil > h, a critical
value for the test, a gross error may be present in stream i. Mah and Tamhane (1982) suggested that k should be equal to zPl2(the upper 012 point of the standard normal distribution), where
p
= 1 - (1 - 4 1 1 5
(15)
and a is the level of significance. No search algorithm is needed to conjunction with the measurement test. Performance of this test has been evaluated by Iordache et al. (1985). The measurement test statistic has been used to construct two new detection schemes, which are presented in this paper. Our goal is to reduce the frequency of mispredictions by supplementing the basic statistical test with additional decision criteria, for instance, through the use of bound violation of variables in data reconciliation or coaptation. Two different strategies based upon this idea have been implemented in the dynamic measurement test (DMT) and the extended measurement test (EMT). One more detection scheme, the double-pass search algorithm (DP), which is based upon the constraint test but which eliminates the detection inconsistencies due to error cancellations, was also evaluated (Rosenberg, 1985). But it will not be discussed further in this paper, because it is found to have less power than the other schemes. Note that DMT and EMT are not restricted to mass flow networks with directly measured variables and independent random measurements errors. These assumptions were made in our evaluation only to simplify the simulation. If nonlinear constraints are involved, they must be linearized before applying the detection schemes. Only measurement biases (no leaks) are allowed as gross errors in our schemes.
Serial Elimination Procedure Different gross error identification schemes based on the serial elimination of streams have been suggested by Ripps (1965),Nogita (1972),and Romagnoli and Stephanopoulos (1981). Although the notion of serial elimination is simple, the exact details of the elimination process have not been explicitly defined. This is done in the present algorithm, GT. The serial elimination process is coupled with the global test statistic in order to locate the possible gross errors. The first step in the serial elimination procedure is the construction of a candidate set of streams. In the GT algorithm, the global test is applied first to the set of all process measurements. If an error is indicated by the test statistic, all process streams are included in a candidate set. The measurements are deleted sequentially from the network in groups of size t , t = 1 , 2 , ..., tmW After each deletion, the process constraints are projected to reflect the reduced set of measurements, and the global test is again applied. The serial elimination of measurements proceeds until among all stream combinations at some level t we can find a set which can be declared suspect or until t exceeds tau. The parameter t,, is the maximum number of measurements that can be deleted from the original network to leave the degree of freedom of the projected constraints equal to 1 or, equivalently, to leave the rank of the projected matrix B equal to 1. Note that the rank of B may be reduced to 1 with fewer than t,, deletions in some instances. If the constraints are imposed by flow networks, t,, = (number of streams) - (minimum degree of the process nodes). Three stream sets are used in the serial elimination. Let S be the set of all measured stream flow rates, L be a temporary set which contains the measurements being deleted in a particular step, and C be the current set of
Ind. Eng. Chem. Res., Vol. 26, No. 3, 1987 557 measurements suspected of containing gross errors. The complete algorithm is as follows: (1)Determine tmw Set t = 0 and C = 0. (2) Determine the degrees of freedom, V, for the network containing measurements in set S, which, if all streams are measured, is equal to the number of physical nodes for the process graph. Calculate the global test statistic, 7. If 7 < x,,,2, declare no errors and stop. Otherwise go to step 3. (3) Set t = t + 1. If t > t,,,, declare all process measurements in S suspect and stop. Otherwise set P- = 1- a. For each possible combination L of t measurements, do the following: (a) Delete the set o f t measurements, L, from the network. (b) Attempt to obtain the projected process constraints, B, to reflect the reduced set of measurements, S - L. If one or more columns of the new matrix B are zero, discard this set L and choose the next set of t measurements in step 3a. 0therwise proceed to step 3c. (c) Determine the degrees of freedom, v, for the projected constraints, Le., the rank of B. Calculate the global test statistic, 7,for the projected constraints. (d) If p(x,2 < 7)< Pmin, replace C by L and reset Pmin = p(x; < 4. Otherwise choose the next set of t measurements in step 3a. (4) If Pmin < 1- a , declare all measurements in C in suspect and stop. Otherwise go to step 3. Step 3b accounts for the fact that the intended deletion of a particular set of measurements may give rise to one or more columns of zeros in the projected constraint matrix B. This fact has already been revealed by Tamhane and Mah (1985). In terms of graph operations, the aggregation of nodes to delete some measurements will also delete other measurements. The measurements whose corresponding columns in B are identically zero will have no effect upon the balance residuals, r , or the global test statistic. Therefore, a gross error in such streams cannot be detected. For this reason, constraint matrices with zero columns are discarded in step 3b. Furthermore, for the deletion of two sets of t measurements, the degrees of freedom for the resulting networks may not be the same. Step 3c determines the degrees of freedom (i.e., the rank of projected matrix B) for each acceptable B matrix. It is possible that at a certain value oft, there exists two or more combinations of measurements whose deletion enables the resulting networks to pass the global test. A decision criteria to select the final set to be declared in error is needed. One possibility is to choose the set L whose deletion produces the minimum global test statistic. But, as shown above, the degrees of freedom of the projected matrices B may not be the same, and therefore a comparison among the corresponding global test statistics may not be adequate. This particular algorithm chooses that combination which yields the lowest cumulative probability p(xy2 < 7),Pmin, at the particular value of t. This is a conservative criterion, because this choice of Pmh maximizes the probability that the remaining measurements S - C will be free of gross errors. Extended Measurement Test A major deficiency in most detection schemes (including the GT algorithm above) is that the reconciliation and coaptation estimates are not checked for violation of their respective upper and lower bounds. For instance, the measurement test may predict a set of streams as being in gross error. Upon the deletion of these streams from the network, the reconciliation and coaptation may yield negative or unreasonable values. This problem was raised by Crowe et al. (1983) and further explored in a recent paper by Tamhane and Mah (1985). A solution to this more complex problem is suggested in the following two
algorithms, which are modified serial elimination detection schemes, this time coupled with the measurement test statistic. We shall describe first the EMT, where the serial elimination procedure is very similar to that previously described for the GT. The difference between the two is the way of constructing the initial and updated candidate set. In EMT this is done by the measurement test and the bound violation checking. Let C’be the set of initial candidates and D be the set of unmeasured variables. Keeping the rest of the notation the same as before, the algorithm is as follows: (1) Determine tmu. Set t = 0, C’ = 0 and C = 0. (2) Perform the measurement test on each measurement in the set S. Place measurements which fail the test in candidate set C! If C’ is not empty, proceed to step 3. Otherwise determine reconciliation estimates of set S and coaptation estimates of set D. If some of the estimates violate their respective upper or lower bounds, place all measured variables in C’, i.e., C’ = S, and go to step 3. Otherwise predict no gross errors and stop. (3) Set t = t + 1. If t > t, declare all measurements suspect and stop, i.e., C’ = S. Otherwise set Pmin = 1. For each possible combination L of t measurements in C’, do the following: (a) Delete set L of t measurements from the measured set S and add it to set D. (b) Attempt to obtain the projected constraint matrix B to reflect the reduced set of measurements S - L. If one or more columns of matrix B are zero, proceed to step 3h. Otherwise proceed to step 3c. (c) Attempt to perform the measurement test on measurements in S - L. If any of the measurements fails the test, proceed to step 3h. Otherwise go to step 3d. (d) Obtain reconciliation estimates for set S - L. If some of these estimates violate their respective upper or lower bounds, proceed to step 3h. Otherwise continue to step 3e. (e) Attempt to obtain coaptation estimates for the set D + L. If some of these measurements are unobservable, proceed to step 3h. Otherwise continue to step 3f. (f) If any coaptation estimates of step 3e violate upper or lower bounds, proceed to step 3h. Otherwise alculate the objective function f = pix; < aTQ-la}.(g) Iff < Pmh, replace Pmin by f and C by L. (h) Clear set L , choose the next set of t measurements, if possible, and go back to step 3a. If no other set L is possible, go to step 4. (4) If Pmin < 1, declare all measurements in C in gross error and stop. Otherwise go back to step 3. It is assumed that all originally unmeasured variables are observable. Step 3f uses the lowest p(x,2 < aTQ-la) as a tie breaker to choose among possible subsets L. This criteria is the same as that used in the GT, because the quantity a TQ-lais identical with the global test statistic, eq 9. Although the global test statistic is used, passing the test is not a requirement of the EMT. From the equations developed by Iordache et al. (19851, one can easily infer that for streams whose corresponding columns in B are zero, the measurement test statistic of eq 14 is undefined. Therefore, a gross error in such streams cannot be detected. For this reason such constraint matrices are discarded in step 3b. Dynamic Measurement Test An important feature of the EMT is that the candidate set is fixed at an initial step. If in searching for the minimal subset C it becomes apparent that some measurements not present in C should be included, there is no mechanism for augmenting the candidate set. An alternative scheme, DMT, uses a different approach to constructing the candidate set. It attempts to enlarge the candidate set a t each step until all suspected measurements are included.
558 Ind. Eng. Chem. Res., Vol. 26, No. 3, 1987
With the same notation as in the GT and EMT, the DMT algorithm is as follows: (1)Determine tmm. Set t = 0 and clear sets L and C. (2) Add set L to set C. Clear set L. Remove the measurements in C from the measured set S. (3) Set t = t 1. Attempt to perform the MT on the remaining measurements S - C. If t > t, or the projected constraint matrix has one or more columns equal to zero, declare all measurements suspect and stop. Otherwise proceed to step 4. (4) Place the measurements having the largest magnitude MT statistic, z-, from step 3 in set L. If zmax> zBl2,go to step 2. Otherwise proceed to step 5. (5) Reconcile the measurements present in S - C. If any estimates of these measurements violate their upper or lower bounds, go to step 2. Otherwise proceed to step 6. (6) Attempt coaptation of the unmeasured parameters in set C + D. If any of these parameters are unobservable, declare all measurements suspect and stop. If all values are observable but some violate their bounds, go to step 2. Otherwise go to step 7. (7) If the deletion of measurements in set C passes steps 2-6, declare measurements in C to be suspect and stop. In the DMT algorithm the candidate set C is always augmented with the measurement which has the maximum absolute measurement test statistic. If more than one stream tie for the maximum test statistic, then all streams associated with that maximum test statistic are placed in the list of suspect streams.
+
Performance Criteria In our notation, the overall null hypothesis, Ho, is the hypothesis that no gross error is present anywhere in the process data. The alternative hypothesis, HI, is the converse of Ho: at least one gross error is present. In gross error detection, a type I error is defined as the prediction of the presence of errors while Ho is true. A type I1 error is the prediction of no error given that HIis actually true. The power of test is the probability of correctly predicting the presence of gross errors. The gross error detection schemes presented in this work indicate not only the presence of errors but also the location of these errors. For this reason, the definition of error types and of the power above should be enlarged. In evaluating the performance of the measurement test, Iordache et al. (1985) proposed two different definitions of the power of that test. However, both definitions are specific to the measurement test and, therefore, inappropriate for other gross error detection schemes. A new framework is needed for our work. We shall assume that at most one gross error is present in the measurements at any time. Let us consider two categories of events: G,, which denotes the prediction of a gross error in stream i (for 1 5 i 5 s) or the prediction of no error (for i = 0), and Lj, which denotes the occurrence of a gross error in stream j (for 1 5 j 5 s) or the event that no gross error is present in the data (for j = 0). From these events we can define a probability P(i,j) as
For the purpose of illustration let us construct a table of such probabilities whose columns correspond to the location of a gross error and whose rows correspond to predictions made by detection scheme:
0
actual location of gross error 1 2 ...
S
P r
0
P(0,O)
P(0,l)
P(0,2)
...
P(0,s)
e d
1
P(1,O)
P(1,l)
P(1,2)
...
P(1,s)
i
2
P(2,O)
P(2,l)
P(2,2)
...
P(2,s)
.
3
P(3,O)
P(3,l)
P(3,2)
...
P(3,s)
0
.
.
n s
.
. P(s,O)
C
t
.
s
.
...
...
. .
...
P(s,s)
... P(s,l)
P(s,2)
.
For instance P(3,2)is the probability of predicting stream 3 in the error when an error is present in stream 2 only. Note that the probabilities in a column do not correspond to mutually exclusive events. A detection scheme may predict multiple gross errors even when at most one gross error is present. Therefore, the sum of a particular column in the table above may be greater than unity. The power of the scheme on a particular measurement shall be defined as the probability of predicting that measurement in error given a gross error in that measurement only. Thus, the powers are simply P G j ) , for l Ij’ Is, i.e., the diagonal elements of the table above, excluding the entry P(0,O). The elements P ( i j ) , with i # j and i # 0, represent probabilities of type I error. Under H,, the event of committing at least one type I error is simply the complement of predicting no errors when no errors are present anywhere. Thus, pitype I errorlH,) = 1 - P(0,O)
(17)
Moreover, the average number of type I errors committed under Ho is the sum of all P(i,O)’s,for any i # 0, in the first column of the table above. Furthermore, the probability of committing a type I1 (TII) error with respect to measurement j given only one gross error present, p{TII),is the quantity P(O,j),for any j # 0. A particular detection scheme may have high powers at the cost of making many mispredictions. The number of incorrect predictions is not reflected directly in any of the definitions above. Therefore, an additional quantity, the selectivity of the detection scheme on a measurement, is defined. The selectivity of the detection scheme for measurement i is the expected number of times a gross error is correctly predicted in stream i, normalized by the expected total number of predictions made given an error in measurement i only. The detection scheme that makes a lot of mispredictions given the presence of one gross error will have low selectivities. The three performance criteria selected for the evaluation of various gross error detection schemes are (a) the power, (b) the selectivity, and (c) the probability of type I1 error of the detection scheme for an average number of type I errors of 0.1. Since they are, for the most part, difficult to evaluate, their values are estimated by using computer simulation.
Simulation Procedure Each simulation run consists of performing NT realizations. Each realization consists of two basic steps. The first step is to generate an appropriate random error vector, w . These errors are then added to the true values, x,to simulate a set of measured data, in the absence of any gross error. The appropriate detection scheme is applied to this
Ind. Eng. Chem. Res., Vol. 26, No. 3, 1987 559
Figure 1. Network 1.
Figure 2. Network 2.
4
Figure 4. Network 4.
-L@L 4 < 1
(
-
/
4
L
Figure 5. Network 5. Figure 3. Network 3.
generated set of data, and the predictions are recorded. The second step is to serially place a gross error in each successive measurement in turn and apply the detection scheme, again recording the predictions. Therefore, in a given realization, the detection scheme is applied (s + 1) times. The vector w remains the same in a given realization; it changes only between successive realizations. This method of performing realizations was also used by Iordache et al. (1985) for the same purpose of yielding more precise comparisons between powers associated with different measurements in the same simulation. In essence, we adopted their simulation technique. In addition to the vector of true values of the flow rates, x,the input of each simulation program also contains the constraint matrix Al, the standard deviations of the measurement errors, and the sign and magnitude of the gross errors. All computer programs require CY, a specified level of significance, for each hypothesis test performed. Note that the detection schemes based upon the measurement test (Le., EMT and DMT) use 0 as the level of significance as defined in eq 15. In every scheme the input variable a is chosen such that the average number of type I errors committed per simulation (under Ho) is approximately 0.1. This is done by trial and error. The output from e5ch of the simulation programs includes (a) estimates, P ( i j ) ,of the probabilities P ( i j )defined in eq 16 and (b) a comparison of the powers associated with different _measurementsin the same simulation. In a simulation run,P ( i j )is determined by a ratio of N ( i j ) , the number of occurrences of an event (GJL,),to NT,the total number of realizations performed,
P(i,.j) = N(i,j)/NT
(18)
In all but one network considered, the number of realizations was set to 1OOOO. Therefore, as shown by Iordache et al. (1985), in the corresponding simulation results, there is approximately 95% confidence that the absolute error E in approximating P(i,j) by P ( i j )is less than 0.01. For network 7, to reduce the computational expenses, NT was lowered to 2500, which corresponds to an error E < 0.02 in approximately 95% of the estimations. The comparison of powers of different measurements within the same simulation has been performed by the McNemar’s x2 test, again used by Iordache et al. (1985). In our tables of computed results, the outcome of McNemar’s test is not shown, but in general, it is significant wherever the powers of two different streams differ by more than the estimation error of powers. Many simulation runs were performed for the comparison of the above-mentioned detection schemes. The run conditions for a selection of these simulation runs are presented in Table IX. They cover seven networks shown in Figures 1-7. For easy reference, each simulation run
7
I
I
Figure 6. Network 6.
I
I
1
I
13
Figure 7. Network 7.
is labeled with two numbers. The first one indicates the network being considered. The second one refers to the method of gross error identification. The GT, MT (pure measurement test), EMT, and DMT are numbered 1-4, respectively. Whenever appropriate, a third label indicates a particular input data set, as illustrated in Table 1X. For instance, runs 1.3.i and 1.3.j, both for EMT on network 1, differ in the true values of x. The computed results, containing the estimates for the power, the selectivity, and p(type I1 error) for every position of the gross error and for all simulation runs associated with a particular network are summarized in Tables I-VII. The tables of results and their respective networks are labeled in a similar fashion. For example, the results of runs for network 5 are summarized in Table V. Discussion of Simulation Results Various factors influence the performance of the gross error detection schemes. They may be classified in the following categories: (1)magnitudes of the measurements and their physical bounds; (2) magnitude of the gross error and the ratio g/u; and (3) constraints, network size and configuration, and stream position, which are summarized by the constraint matrices Al and A2. The impact of these factors on the performance of various detection schemes is shown below using the simulation results presented in Tables I-VII. We shall refer mainly to GT, MT, and EMT schemes, since the performance behavior of the DMT is very similar to that of EMT. A comparison of the two schemes is illustrated in
560 Ind. Eng. Chem. Res., Vol. 26, No. 3, 1987 Table I. Simulation Results for the Runs Associated with Network 1 perform. measurements run det sch ratio g / u criteria 1 2 3 1.1 GT 3 power 0.591 0.598 0.595 kelectivity 0.789 0.796 0.777 0.303 0.299 0.292 p{TIIl power 1.3.a EMT 3 0.590 0.593 0.580 selectivity 0.824 0.819 0.813 0.305 0.298 0.305 PIT111 1.3.b EMT 3 0.599 0.581 0.593 power selectivity 0.829 0.820 0.819 0.297 0.312 0.293 PIT111 3 power 0.642 0.648 0.643 1.3.c EMT selectivity 0.807 0.819 0.807 0.312 0.310 0.313 PIT111 1.3.d EMT 3 0.673 0.660 0.650 power selectivity 0.739 0.753 0.758 0.317 0.329 0.340 PIT111 3 0.596 0.588 0.587 power 1.3.e EMT selectivity 0.823 0.821 0.819 0.298 0.301 0.303 PITW 3 power 0.592 0.585 0.592 1.3.f EMT selectivity 0.820 0.820 0.828 0.297 0.302 0.303 PIT111 power 3 0.646 0.645 0.643 1.3.g EMT selectivity 0.807 0.799 0.813 0.310 0.310 0.313 P r w 3 0.660 0.657 0.660 power 1.3.h EMT selectivity 0.741 0.733 0.741 0.310 0.330 0.329 pITII1 3 0.421 0.416 0.409 power 1.3.i EMT selectivity 0.696 0.680 0.686 0.516 0.509 0.520 PWIJ 3 0.579 0.576 0.572 power 1.3.j EMT selectivity 0.825 0.837 0.833 0.317 0.324 0.325 p(TII1 Table 11. Simulation Results for the Runs Associated with Network 2 measurements perform. run det sch ratiog/a criteria 1 2 3 4 2.1 GT 3 power 0.610 0.602 0.610 0.604 selectivity 0.812 0.815 0.820 0.811 p(TII} 0.302 0.314 0.310 0.312 2.3.a EMT 3 power 0.548 0.549 0.551 0.541 selectivity 0.751 0.755 0.754 0.758 plTII) 0.397 0.390 0.394 0.403
the runs associated with network 4. 1. Magnitudes of Measurements and Bounds. The sensitivity of DMT and EMT to bound violation is a function of the magnitude of measured values and their standard deviations, the process constraints, and the significance level of the statistical test. In runs 1.3.a-1.3.d, Q = I, while in runs 1.3.e-1.3.h Q = 0.11. Four different bounds were chosen in both cases as follows: (i) x f loa; (ii) x f 50; (iii) x f 20.; and (iv) x f 1.250. From runs 1.3.a-1.3.h, it was found that if the bounds are within x f 2a, their effect on the performance of EMT and DMT is significant. The increase in power is accompanied by a less pronounced decrease in selectivity as the interval between bounds becomes narrower. A change in the philosophies of EMT and DMT schemes can be invoked to explain this effect. For tight bound intervals, the bound violations occur frequently, which results in placing very often all measurements in the candidate set and subsequently in declaring all of them in gross error (in steps 3 of both DMT and EMT). This explains the relative increase in power with the decrease in selectivity. If the magnitudes of the true values of the measured variables are close to only one of the bounds, both power and selectivity of EMT and DMT schemes deteriorate. In
that case, the frequency of the bound violations leading to declaring all measurements in error is reduced, and hence there is also a decrease in the power. Consider, for instance, the pairs of runs 1.3.i and 1.3.j and 2.3.a and 2.3.b. The true values are far away enough from their common upper bound. But in runs 1.3.i and 2.3.a, x = 1for all flow rates and the true values are very close to the lower bound (which was set to zero), while in runs 1.3.j and 2.3.b the values of x were placed further from the lower bound (x = 5 for all flow rates). In runs 1.3.i and 2.3.a, the power and the selectivity are lower, while the probability of the type I1 error is higher. As shown by Tamhane and Mah (1985), for constraints leading to proportional columns in B, even when both bounds lie outside the range x f 2a, the performance of the detection scheme is seriously affected. The relationship between the statistical part and the bound violation checking part of the EMT and DMT becomes more visible by taking different input values for the level of significance a. This is illustrated by runs 3.3.a-3.3.c for the EMT scheme. In run 3.3.a, a = 0.1 and the statistical part is significant. Then, the performance measures for all streams are comparable. On the other extreme, if the input level (Y is very low, for instance, a = as in run 3.3.c, the statistical part of the detection algorithm becomes inactive and the performance of the scheme depends mainly on the bound violation part. In the latter case, the performance of EMT for streams 1and 5 becomes much worse than that for streams 2-4. A gross error present in a relatively large flow rate is more difficult to detect by bound violation only than is the same error associated with a small flow rate. For this reason, although the bound violation checking is important in obtaining reasonable reconciled and estimated values, the input level should be adjusted such that the statistical test is also active. The situation in run 3.3.b is an intermediate one. Although the results are not given here, the same is true for DMT also. 2. Magnitude of the Gross Error. This parameter has a remarkable effect upon the performance of all detection schemes. A better measure of its effect is the ratio g / u . Every detection scheme is highly sensitive to this ratio, as shown by runs 3.l.a-3.3.a, 4.1-4.4, and 5.1-5.3. As g / a increases, the powers increase and the probabilities of committing type I1 errors decrease. Selectivities also increase for all schemes except the MT, where a maximum value at some value g/a is observed. As suggested by Iordache et al. (19851, for large g/a, the smearing effect of the gross error onto adjoining measurements becomes significant and in such a case the MT is more susceptible to make mispredictions, which results in lowering the selectivity. But in GT, DMT, and EMT, the serial elimination and additional testing reduce the frequency of type I errors, increasing the selectivity. 3. Constraints and Network Structure. If two or more columns of the constraint matrix A, are proportional, then the powers, the selectivities, and the type I1 error probabilities for the corresponding measurements are equal. This effect has been already shown by Iordache et al. (1985) for MT, and it is valid for EMT and DMT also. If two streams in a network are parallel, then their associated columns in matrix AI are identical. This situation is found in networks 5 and 6. Measurements 3 and 4 in runs 5.1-5.3 exhibit identical performances and so do also measurements 5 and 6 in runs 6.2 and 6.3. For the EMT, the selectivities of streams 3 and 4 in run 5.3 remain low and fairly constant for the three values of g / a simulated, while the selectivities of streams 1 and 2
Ind. Eng. Chem. Res., Vol. 26, No. 3, 1987 561 Table 111. Simulation Results for the Runs Associated with Network 3 measurements run
3.1.a
det sch GT
ratio .g/o 2 3
perform. criteria power selectivity PIT111 power selectivity
4 3.2.a
MT
2 3 4
3.3.a
EMT
2 3 4
3.3.b
EMT
4
3.3s
EMT
4
P{TIIl
power selectivity
PlTIIl
power selectivity
PlTIIl
power selectivity
PPI11
power selectivity PIT111 power selectivity
PWIl
power selectivity
PlTIIt
power selectivity
PlTIIl
power selectivity
PWIl
power selectivity PlTIIl
1 0.234 0.575 0.616 0.511 0.725 0.346 0.767 0.801 0.129 0.232 0.541 0.689 0.512 0.585 0.425 0.796 0.567 0.177 0.220 0.599 0.655 0.497 0.741 0.374 0.775 0.819 0.137 0.283 0.456 0.690 0.057 0.348 0.932
2 0.242 0.584 0.610 0.516 0.725 0.340 0.779 0.801 0.120 0.236 0.548 0.686 0.522 0.590 0.417 0.807 0.563 0.164 0.225 0.598 0.648 0.488 0.737 0.388 0.760 0.814 0.148 0.528 0.543 0.360 0.500 0.540 0.476
3 0.171 0.490 0.668 0.392 0.640 0.427 0.640 0.732 0.209 0.185 0.489 0.729 0.435 0.559 0.488 0.699 0.537 0.250 0.172 0.552 0.693 0.385 0.693 0.462 0.647 0.776 0.221 0.566 0.712 0.397 0.458 0.487 0.482
4 0.232 0.571 0.619 0.509 0.719 0.345 0.774 0.805 0.129 0.228 0.542 0.695 0.532 0.592 0.405 0.804 0.564 0.169 0.226 0.598 0.647 0.495 0.720 0.375 0.763 0.784 0.142 0.534 0.529 0.349 0.447 0.513 0.475
5 0.246 0.590 0.608 0.513 0.721 0.341 0.774 0.810 0.127 0.230 0.539 0.693 0.528 0.587 0.410 0.803 0.564 0.166 0.220 0.587 0.652 0.489 0.724 0.385 0.758 0.792 0.150 0.276 0.467 0.699 0.051 0.336 0.938
Table IV. Simulation Results for the Runs Associated with Network 4 run 4.1
det sch GT
4.3
EMT
ratio g / o 2
perform. criteria Dower selectivity PlTIIl power selectivity PWl power selectivity
PlTIIl
power selectivity
PIT111
power selectivity
PWIl
4.4
DMT
power selectivity PlTIIl power selectivity PIT111 power selectivity
PlTIIl
power selectivity
PlTIIl
1 0.176 0.540 0.692 0.423 0.707 0.438 0.708 0.809 0.129 0.189 0.572 0.701 0.438 0.689 0.436 0.709 0.772 0.186 0.183 0.573 0.701 0.442 0.725 0.439 0.714 0.812 0.192
increase as g / u increases. The DMT scheme exhibits the same kind of result. Due to the structure of the graph, stream 3 or 4 cannot be deleted without a loss of observability. Stream 3 or 4 can be predicted in error only if all process streams are declared in error in step 6 of EMT or step 3 of DMT. Thus, the selectivities of streams 3 and 4 must be less than or equal to (which is 11s). The same rule can be observed at streams 5 and 6 of run 6.3.
2 0.181 0.540 0.685 0.430 0.703 0.429 0.700 0.802 0.192 0.184 0.560 0.701 0.446 0.693 0.433 0.714 0.760 0.189 0.185 0.555 0.702 0.435 0.706 0.443 0.714 0.776 0.192
3 0.183 0.553 0.688 0.426 0.706 0.435 0.701 0.809 0.195 0.175 0.531 0.711 0.433 0.653 0.440 0.720 0.727 0.182 0.171 0.522 0.709 0.427 0.671 0.453 0.717 0.749 0.197
measurements 4 0.120 0.444 0.743 0.298 0.610 0.536 0.548 0.736 0.302 0.133 0.531 0.748 0.331 0.634 0.523 0.575 0.731 0.288 0.124 0.475 0.753 0.322 0.645 0.531 0.573 0.752 0.293
5 0.151 0.499 0.711 0.384 0.680 0.469 0.655 0.785 0.220 0.159 0.517 0.722 0.399 0.638 0.462 0.674 0.717 0.211 0.160 0.516 0.715 0.402 0.666 0.464 0.673 0.746 0.222
6 0.122 0.450 0.740 0.301 0.626 0.541 0.549 0.737 0.301 0.128 0.486 0.754 0.324 0.614 0.525 0.571 0.718 0.293 0.133 0.480 0.741 0.311 0.614 0.542 0.573 0.736 0.294
7 0.183 0.554 0.688 0.429 0.704 0.430 0.703 0.810 0.195 0.178 0.545 0.710 0.439 0.666 0.438 0.711 0.730 0.186 0.175 0.544 0.707 0.427 0.674 0.450 0.707 0.755 0.196
Although not shown here, it is obvious that the rule is true also for the GT. Additional aspects of the impact of the network structure on the performance parameters are discussed next. 4. Position of the Measurement Containing a Gross Error. Our observations show that the performance of the gross error detection schemes on a stream in a process digraph is generally influenced by the degrees of its ad-
562 Ind. Eng. Chem. Res., Vol. 26, No. 3, 1987
the other adjacent nodes. This can be seen in all runs reported in this paper. (ii) The performance of the detection schemes is increased if nodes of degree 2 form a sequence. In general, the longer the sequence, the better the performance of the scheme for the streams incident to those nodes. A comparison of runs 1.1-1.3.a and 2.1-2.3 illustrates this effect. In network 7, streams 14-16 belong to a sequence of two degree 2 nodes, while streams 1 and 2 or 6 and 7 are each incident to single degree 2 nodes. Streams 14-16 enjoy slightly higher powers than do streams 1 and 2 or 6 and 7 . Selectivities of streams 1 and 2 are, however, markedly higher than those of streams 14-16. It is difficult to get any consistent rule for the streams not adjacent to degree 2 nodes or which are not parallel. In general, the lower the sum of the degrees of the adjacent nodes for a stream containing a gross error, the higher the power and the selectivity for that stream. However, as found also by Iordache et al. (1985) for the measurement test, complex networks do not obey this rule entirely. It should be noted that even when for some network the position rules hold for the MT, they may be disturbed in EMT and DMT by the bound violation checking. The relationship between the statistical test and the bound violation testing determines the behavior of EMT and DMT for different streams. For instance, in run 3.3.a where the statistical part is significant, the EMT behaves similarly to the MT, while in run 3.3.c where the statistical part is relaxed, the bound violation checking changed the position rule drastically. Comparison of Different Detection Schemes. Table VI11 summarizes the average performances of different detection schemes for some selected runs. Although variations of the performance measures among the network streams are in many cases significant and should be taken into account, we constructed this table for a more convenient comparison and interpretation of the general trends in the performances of the detection schemes dis-
Table V. Simulation Results for the Runs Associated with Network 5 perform. measurements run det sch ratio g / u criteria 1 2 3 4 5.2 MT 2 power 0.176 0.172 0.372 0.372 selectivity 0.433 0.434 0.418 0.418 p(TII} 0.727 0.736 0.583 0.583 3 power 0.426 0.441 0.615 0.615 selectivity 0.521 0.522 0.420 0.420 p(TII1 0.486 0.469 0.349 0.349 4 power 0.730 0.733 0.827 0.827 selectivity 0.526 0.527 0.420 0.420 p(TI1) 0.216 0.212 0.156 0.156 5.3 EMT 2 power 0.216 0.217 0.122 0.122 selectivity 0.626 0.638 0.222 0.222 p(TI1) 0.735 0.729 0.817 0.817 3 power 0.464 0.466 0.286 0.286 selectivity 0.759 0.773 0.233 0.233 p(TI1) 0.462 0.464 0.629 0.629 4 power 0.730 0.736 0.519 0.519 selectivity 0.850 0.853 0.240 0.240 p(TI1) 0.201 0.199 0.391 0.391
jacent nodes. The degree of a node is defined as the total number of streams incident with that node. Iordache et al. (1985) stated some rules concerning the effect of the degrees of the adjacent nodes on the power of MT, which are generally valid also for the detection schemes in the present paper. The only exception is the DP algorithm, which is not analyzed here. Let d.. be the degree of the node for which j is an inflow and do, be the degree of the node for which stream j is an outflow. The most notable effects are as follows: (i) The power associated with a stream adjacent to a node of degree 2 is higher than that of any other stream with the same d, + do; for which d , # 2 and doj # 2. This is shown by streams 1 and 2 and stream 5 in runs 4.1 and 4.2. Furthermore, two streams adjacent to a common node of degree 2 have statistically the same power, selectivity, and type I1 error probability regardless of the degree of
Table VI. Simulation Results for the Runs Associated with Nc!twork 6 measurements run 6.2
det sch MT
ratio g / u 3
6.3
EMT
3
perform. criteria power selectivity PlTW power selectivity PWI
1 0.524 0.647 0.426 0.447 0.730 0.499
2 0.481 0.609 0.460 0.383 0.706 0.536
3 0.474 0.605 0.466 0.389 0.707 0.530
4 0.523 0.650 0.424 0.437 0.727 0.509
5
0.316 0.419 0.627 0.236 0.137 0.698
6 0.316 0.419 0.627 0.236 0.137 0.698
7 0.404 0.576 0.522 0.319 0.664 0.600
Table VII. Simulation Results for the Run Associated with Network 7 measurements run 7.3
run 7.3
det sch EMT
det sch EMT
ratio g / u 4
ratio g / u 4
perform. criteria power selectivity DITIII perform. criteria power selectivity pl'W
1
0.66 0.82 0.28 9 0.56 0.50 0.36
2 0.66 0.83 0.29 10 0.60 0.50 0.32
3 0.56 0.48 0.37 11
0.60 0.51 0.33
4 5 0.53 0.51 0.51 0.48 0.38 0.41 measurements
6 0.67 0.51 0.28
12 0.61 0.53 0.32
14 0.71 0.54 0.25
13 0.52 0.49 0.41
Table VIII. Comparison of the Average Performances of Different Detection Schemes detection scheme for run GT EMT 1.1 3.1.a 4.1 MT, 3.2.a 1.3.a 1.3.d 3.3.a power 0.595 0.488 0.384 0.506 0.588 0.471 0.661 0.787 0.706 selectivity 0.677 0.583 0.819 0.723 0.750 0.360 0.468 0.397 0.429 0.303 0.329 plTI1) 0.298
7 0.67 0.51 0.28
8 0.58 0.53 0.35 ~~
4.3 0.401 0.655 0.465
15 0.71 0.51 0.25
DMT, 4.4 0.395 0.671 0.475
16 0.72 0.52 0.24
Ind. Eng. Chem. Res., Vol. 26, No. 3, 1987 563 Table IX. List of Runs key features of data run 1.1 1.3.a 1.3.b 1.3.c 1.3.d 1.3.e 1.3.f 1.3.g 1.3.h
det sch
GT EMT EMT EMT EMT EMT EMT EMT EMT EMT EMT GT EMT EMT GT MT EMT EMT EMT GT EMT DMT MT EMT MT EMT EMT
1.3.i 1.3.j 2.1 2.3.a 2.3.b 3.1.a 3.2.a 3.3.a 3.3.b 3.3.c 4.1 4.3 4.4 5.2 5.3 6.2 6.3 7.3 acu
= 0.10,
b a
= 0.01.
ccu
network 1 1 1 1 1 1 1 1 1 1 1 2 2 2 3 3 3 3 3 4 4 4 5 5 6 6 7
Q I I I I I 0.11 0.11 0.11 0.11 I I
I
I I I I I I I 1
I I I I I I 0.11
X
111 20 20 20 20 20 20 20 20 20 20 20 20 20 20 20 20 20 20 20 20 20 20 20 20 111 555 1111 1111 5555 21112 21112 21112" 21112b 21112c 5 15 155 1055 5 15 155 105 5 5 15 155 1055 101055 101055 10 15 15 105 55 10 15 15 10555 15 1525 105 10 105 5 5 1055 10 10 10
Lb
Ub
10 15 18 18.75 16.838 18.419 19.368 19.605 0 0
30 25 22 21.25 23.162 21.581 20.632 20.395 100 100
0 0
100 100
0 0 0
100
0
0
100 100
0
100
0 0
100
100 100
100
computed results in table
I I I I I I I I I I I I1 I1 I1 I11 I11 I11 I11 I11
IV IV IV V V VI VI VI1
= 0.00001.
cussed in this work. The values in Table VI11 represent averages of the computer results over all network streams and for a medium value of g / a taken in our runs g/u = 3). Our findings show that in general MT has the highest power and the lowest selectivity in comparison with the other detection schemes. This is illustrated by runs 3.l.a, 3.2.a, and 3.3.a. GT exhibits good overall performance. As shown by runs 3.1.a and 3.2.a, it has higher selectivity and lower type I1 error probability than MT, while the power is insignificantly lower. EMT and DMT exhibit performance comparable to GT for the network with loose upper and lower bounds (runs 1.1, 1.3.a, 3.l.a, 3.3.a, 4.1, 4.3, and 4.4). For tight bounds, however, the selectivity of EMT is reduced, but its power is increased. This can be seen by comparing runs 1.1and 1.3.d. Although not shown in Table VIII, the same is true for DMT as well. The EMT vs. DMT performance is very similar, as observed, for instance, from runs 4.3 and 4.4. The major difference between the two detection schemes is the time required for simulation. On the average, EMT takes roughly 20% more time in simulation than DMT for the same run. DMT constructs its candidate set faster than does the EMT. EMT uses a serial elimination step which is slow, whereas DMT does not. The simulation times can be used as a guide to infer relative speeds of the different detection schemes, because the time used for performing the tests and identifying gross errors is predominant. MT, GT, DMT, and EMT take increasingly more time for a given network. For instance, in the runs associated with network 4, for g / a = 4, the corresponding times per simulation on a CDC Cyber 1701730 range from 0.01 s for MT to 0.3 s for EMT. The relative speed of GT, DMT, and EMT with respect to that of MT is 11.8, 23.5, and 30, respectively.
Conclusions On the basis of the simulation results, we can infer the
following conclusions: (1)The performance of any gross error detection scheme depends on the magnitude of the ratio g/u, the network configuration and process constraints, and the position of the gross error. The rules associated with the effects of these factors on the performance criteria are generally valid for any detection scheme considered in this paper, when the statistical part of the detection scheme plays an active role. They may be significantly modified, when bound violation is the predominant mechanism in DMT and EMT. This is especially true of the position rules. (2) None of the detection schemes discussed in this paper is appropriate for networks with parallel streams or constraint matrices with proportional columns. Should this situation occur, additional constraints to eliminate this degeneracy are recommended. (3) For the same average number of type I errors, GT, DMT, and EMT have superior selectivity in comparison with that of MT, with a less significant decrease in the power. DMT and EMT have the useful feature of checking the reconciliation and coaptation estimates for possible bound violation. For most chemical processes, this aspect is very important. Thus, although more expensive, the EMT or DMT is a likely choice for the majority of chemical processes. The performance of these two schemes is quite similar in all networks studied in this work. The choice between them favors DMT which requires less computing time.
Acknowledgment This work is supported by National Science Foundation Grant CBT-8519182 and by the Du Pont Education Foundation.
Nomenclature AI, A2 = matrices of coefficients for the measured and the unmeasured variables, respectively, in balance eq 3 B = (s x 1) vector of reconciliation adjustments, eq 11
564
Ind. Eng. Chem. Res. 1987,26, 564-573
B = (n - q X s) projected constraint matrix equal to PA, C = current set of measurements suspected of containing gross errors C' = set of initial candidate measurements in EMT D = set of unmeasured variables d = (s X 1)vector of transformed measurement adjustments, defined by eq 12 E = error in the estimation of the power by an average value g = (s x 1) vector of gross errors g = magnitude of the gross error H = ( n - 4 x n - q ) covariance matrix of balance residuals r Ho = the null hypothesis used in hypothesis testing H I = the alternative hypothesis in hypothesis testing k = critical value for the measurement test statistic L = temporary set of measurements in serial elimination Lb = lower bound for the measurements in Table IX N ( i j ) = frequency of the event (GilLj) NT = total number of realizations in a simulation n = rank of matrix A, P = ( n - q X n ) projection matrix, defined in eq 5 p(i,j)= probability defined by eq 16 P(i,j) = estimated value of P ( i j ) P- = minimum cumulative probability used in GT and EMT p i ) = probability Q = (s X s) covariance matrix of the measurement errors q = rank of matrix A, r = (n - q x 1)vector of balance residuals, defined by eq 7 S = set of all measured variables s = number of measured process variables t = number of deleted measurements in the serial elimination t , = maximum number of deleted measurements allowable in a given network Ub = upper bound for the measurements in Table IX W = (s x s) covariance matrix of the transformed measurement adjustments, defined by eq 13 w = (s x 1)vector of random measurement errors x = general notation for the true values of measurements
(s X 1) vector of true values of measured variables 2 = (s x 1) vector of estimates of true values x y = (s x 1) vector of measured or simulated variables 9 = (s x 1) vector of adjusted values of measured variables
x=
z = measurement test statistic z, = maximum absolute value among all measurement test
statistics Greek Symbols = input level of significance for the statistical tests @ = modified level of significance for the multiple tests q = constraint test statistic u = degrees of freedom for a x 2 random variable cy
u
= standard deviation for the measurement error
= global test statistic x,,,' = critical value for the T
x2 global test statistic
Subscripts i = for process streams j = for process nodes or the location of gross error
Literature Cited Almasy, G. A.; Sztano, T. Probl. Control Inf. Theory 1975, 4, 57. Crowe, C. M.; Garcia Campos, Y.; Hrymak, A. AlChE J . 1983, 29, 881. Iordache, C.; Mah, R. S. H.; Tamhane, A. C. AIChE J . 1985,31,1187. Mah, R. S. H.; Stanley, G . M.; Downing, D. M. Ind. Eng. Chem. Process Des. Dev. 1976, 15, 175. Mah, R. S. H.; Tamhane, A. C. AIChE J. 1982,28, 828. Nogita, S. Ind. Eng. Chem. Process Des. Dev. 1972, 11, 197. Ripps, D. L. Chem. Eng. Prog. Symp. Ser. 1965,55, 8. Romagnoli, J. A.; Stephanopoulos, G. Chem. Eng. Sci. 1981,36,1849. Rosenberg, J. M.S. Thesis, Northwestern University, Evanston, IL, 1985. Stanley, G. M.; Mah, R. S. H. Chem. Eng. Sci. 1981, 36, 259. Tamhane, A. C.; Mah, R. S. H. Technometrics 1985, 27, 409.
Received for review January 8, 1986 Accepted August 29, 1986
Two-Liquid-Phase Extractive Distillation for Aromatics Recovery Fu-Ming Lee* and Daniel M. Coombs Phillips Research Center, Phillips Petroleum Company, Bartlesville, Oklahoma 74004
Two-liquid-phase extractive distillation using sulfolane as the solvent has been shown to be an effective process for recovering high-purity aromatics with high yield from a variety of hydrocarbon feedstocks. The process variables such as the water content in the sulfolane solvent, the solventto-hydrocarbon feed ratio, and t h e solvent-to-hydrocarbon reflux ratio were extensively studied in a 4-in.-diameter extractive distillation column with 74 sieve trays. The laboratory equilibrium data indicated that, under a single liquid phase, water improves the selectivity of the sulfolane solvent. However, the reverse was observed under two liquid phases. Water substantially reduces the boiling point of sulfolane. The vapor pressure of the sulfolane solvent was determined as a function of temperature and water content. The laboratory data also showed that two-liquid-phase extractive distillation is, in fact, a combination of liquid-liquid extraction and simple distillation of the boiling raffinate phase. Aromatics have become a major product stream in mordern petroleum and petrochemical operation. Currently, aromatics are much in demand as blending stocks for unleaded gasolines. High-purity aromatics are required as raw materials for plastics, synthetic fibers, and many other petrochemical products. In commercial operation, aromatics are mostly recovered by liquid-liquid extractions using sulfolane as the extractive solvent (Deal et al., 1959; DeGraff and Perga, 1969; Thompson, 1973). Extractive distillation has been a widely practiced unit operation for many years. The basic principle is to add
a low-volatile solvent to the distillation column to increase the relative volatility between close boiling components. The characteristics, design, and operation of an extractive distillation column (EDC) have been thoroughly discussed in the literature (Atkins and Boyer, 1949; Chambers, 1951; Hackmuth, 1952; Butler and Bichard, 1963). Due to the insolubility of sulfolane with nonaromatic hydrocarbons, sulfolane was previously considered entirely unsuitable for extractive distillation. Nevertheless, a novel two-liquidphase extractive distillation process for aromatics recovery using sulfolane was discovered by Cines (1977). Because
0888-588518712626-0564$01.50/0 0 1987 American Chemical Society
