Value-Optimal Sensor Network Design for Steady-State and Closed

Subscriber access provided by UNIV OF ESSEX

Article

Value-Optimal Sensor Network Design for Steady-State and Closedloop Systems Using the Generalized Benders Decomposition Jin Zhang, and Donald J. Chmielewski Ind. Eng. Chem. Res., Just Accepted Manuscript • DOI: 10.1021/acs.iecr.7b01865 • Publication Date (Web): 19 Sep 2017 Downloaded from http://pubs.acs.org on October 3, 2017

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Industrial & Engineering Chemistry Research is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Value-Optimal Sensor Network Design for Steady-State and Closed-loop Systems Using the Generalized Benders Decomposition Jin Zhang, Donald J. Chmielewski∗ Department of Chemical & Biological Engineering, Illinois Institute of Technology, Chicago, IL 60616

Abstract The problem of value-optimal sensor network design for linear systems has been shown to be of the nonconvex mixed integer programming class. While the branch and bound search procedure can be used to obtain a global solution, such a method is limited to fairly small systems. The bottleneck is that during each iteration of the branch and bound search, a fairly slow Semi-Definite Programming (SDP) problem must be solved to its global optimum. In this paper, it is demonstrated that an equivalent reformulation of the nonconvex mixed integer programming problem and subsequent application of the Generalized Benders Decomposition (GBD) algorithm will result in massive reductions in computational effort. While the proposed algorithm has to solve multiple mixed integer linear programs, this increase in computational effort is significantly outweighed by a reduction in the number of SDP problems that must be solved. Keywords: Sensor Placement; Optimization; Generalized Benders Decomposition; Linear Matrix Inequalities

1. Introduction In the seminal work of Bagajewicz, [1], the notion of cost optimal Sensor Network Design (SND) was introduced (please see [2–5] for extensive discussions and literature reviews of sensor network design problem). The cost optimal framework was simply stated as: ∗

Corresponding author. phone: 312-567-3537, fax: 312-567-8874, email: [email protected]

Preprint submitted to Elsevier

July 21, 2017

ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 29

minimize the capital cost of sensors subject to constraints on the performance of the resulting estimator. One class of performance criteria was a set of bounds limiting the allowable estimation error variance. The intuitive trade-off was that more/higher quality sensors would reduce error variances, but do so at a higher capital cost. The solution method advocated in Bagajewicz [1], was a tree-search procedure. In Chmielewski et. al. [6], the cost optimal SND problem was reformulated using the method of SemiDefinite Programing (SDP) to arrive at a Mixed Integer Convex Program (MICP) that could be solved using the branch and bound algorithm. In addition, the work extended the steady-state framework of [1] to the case of dynamic open-loop processes. To address closed-loop processes, Chmielewski & Peng, [5], and Peng & Chmielewski, [7], created a similar formulation, but included the feedback element as a design variable. This pair of papers also extended the design question to include actuator selection, so that both hardware elements (sensors and actuators) could be selected simultaneously. In a second seminal paper from the Bagajewicz group, [8], the notion of value (or profit) based SND was introduced. In that work, it was postulated that an improved sensor network would enhance the economics of the process by impacting its operating cost. Thus, the valueoptimal SND problem is similar to the cost-optimal version, but augments the capital cost objective function with an operating cost term to arrive at a Net Present Value (NPV) based objective function. In Nguyen and Bagajewicz, [9], a cut-set based tree search algorithm is advocated to solve the value-optimal SND problem. In Peng & Chmielewski, [10], the valueoptimal SND formulation is extended to closed-loop dynamic systems by combining SND notions with the minimally backed-off operating point notions of Economic Linear Optimal Control (ELOC). For additional discussion of ELOC, see Peng et.al., [11], and Omell and Chmielewski, [12]. Since the ELOC approach is also grounded in SDP, the resulting MICP could also be solved using branch and bound. In 2015, a new solution method for the ELOC problem was developed, [13]. In this work, the previous branch and bound procedure was replaced with an application of the

2 ACS Paragon Plus Environment

Page 3 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Generalized Benders Decomposition (GBD), and resulted in orders of magnitude reductions in computational effort. In 2016, a similar, though distinct, application of the GBD to the cost optimal SND problem (including the cost optimal SND involving residual precisions) also yielded orders of magnitude reductions in computational effort, [14]. The current paper aims to apply the GBD algorithm to the value-optimal SND problem, which is essentially a combination of the ELOC and cost optimal SND problems. The remainder of this section will review both the cost-optimal and value-optimal SND problems, using a steady-state perspective. In addition, the essential elements of the GBD algorithm will be reviewed. In Section II, a re-statement of the original value-optimal SND problem will be shown to be amiable to the GBD algorithm. However, it will also be shown that the resulting relaxed master problem will contain non-convex constraints beyond the expected integer variables. Thus, two theorems will be provided that will ultimately yield a relaxed master problem in the form of a Mixed Integer Linear Program (MILP) from which solutions can be obtained with extreme efficiency. Two examples are then provided to illustrate the computational advantage of the GBD approach. The final section will investigate application of the GBD to the closed-loop version of the value-optimal SND problem.

1.1. Review of Cost-Optimal and Value-Optimal SND Consider a process flow network with material balances ss = M sp , where sp is a vector of primary flow variables and ss are secondary variables (i.e., ss are completely determined by sp and the material balances). Then, the vector of all material flows can be written as 



 sp  s=  = Csp ss

(1)

where C = [I M T ]T and the dimension of s is ns × 1. Now let δ be the set of measurements of all flow variables, each corrupted by independent measurement noise: δ = s + v. Each 3 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 29

element of v is zero mean and has a variance σv2i . Since each noise term is independent of the others, the covariance matrix of v is of a diagonal form. Given the measurement vector, δ, the optimal estimate of s is ]

[ sˆ =

−1 T −1 C(C T Σ−1 v C) C Σv

δ

(2)

−2 −2 −2 −2 2 Notice that Σ−1 v = diag{σv1 σv2 . . . σvn }. Thus, if σvi = 0 (or σvi → ∞), then the estimator

will interpret as stream i is unmeasured, even though δi remains in the formulation and the dimension of the matrices is unchanged. As such, the inverse variance will be used to indicate the presence or absence of a sensor at stream i : σv−2 = αi σ ¯v−2 , where σ ¯v2i is the i i known variance of the sensor and αi is the zero-one attendance variable. The performance of the sensor network will be indicated by the estimation error: q = s−ˆ s. While this error is unknownable, one can calculate the variance of each element of q

ζl = ρl Σq ρTl ,

l = 1 . . . nq

(3)

where ρl is the lth row of the nq identity matrix and −1 T Σq = C(C T Σ−1 v C) C

(4)

Now define the performance criteria as ζl ≤ σ ¯l2 , l = 1 . . . nq (nq ≤ ns ) where σ ¯l is the maximum allowable standard deviation of the estimation error. Then, one can employ the Schur Complement Theorem to arrive at the following set of convex Linear Matrix Inequality (LMI) constraints 

  

ζl

ρl C

C T ρTl C T Σ−1 v C

 ¯l2 ,  ≽ 0, ζl ≤ σ

l = 1 . . . nq

(5)

where the ≽ sign in (5) indicates that the matrix must be positive semi-definite. In addition, 4 ACS Paragon Plus Environment

Page 5 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


one finds that { Σ−1 v

= diag

} α1 σ ¯v−2 1

α2 σ ¯v−2 2

. . . αns σ ¯v−2 ns

=

ns ∑

αi Θi

(6)

i=1

where Θi = σ ¯v−2 ρTi ρi . If the cost of a sensor at location i is ci , then the SDP version of the i (s)

cost optimal SND problem for steady-state systems is stated as

min

ns ∑

αi ∈{0,1},ζl

(s)

ci αi

(7)

i=1

s.t. ζl ≤ σ ¯l2 , l = 1 . . . nq   ζl ρl C  (n )  ≽ 0 l = 1...n   q ∑s T T T αi Θi C C ρl C

(8) (9)

i=1

The branch and bound method can be used to solve problem (7). However, it is noted that enforcement of the constraints of (9) requires the use of a SDP solver to solve the subproblem of each iteration of the branch and bound search. While SDP solvers will provide global solutions, they tend to converge slowly and require large amounts of memory. Thus, a reduction in the number of calls to the SDP solver is expected to provide a computational advantage. In Bagajewicz et.al., [8], a method for assessing Downside Expected Financial Loss (DEFL) incurred by production loss was developed. Due to the inaccuracy of the estimate of a product stream flow rate, there is a finite probability that the estimate is above the target but in fact the real flow is below it. In such a situation it is assumed that the operator will not make any correction to production. As a result, production output will be less than the target and a financial loss will be incurred. The financial loss was determined to be DEFL= 0.19947Ks T σl , where Ks is the cost of the product, T is the time window of the analysis and σl is the standard deviation of stream l. If one of the performance criteria is DEFL, then one must introduce a new set of decision variables, σl , and required each to be



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 29

such that ζl ≤ σl2 . As an extension of the cost optimal SND problem, the value-optimal SND problem is stated as

min

ns ∑

αi ∈{0,1} i=1 σl ,ζl

(s) ci αi

+

nq ∑

0.19947Ks T σl

(10)

l=1

s.t. ζl ≤ σl2 , l = 1 . . . nq   ζl ρC  )  ≽ 0 l = 1...n (n l   s q ∑ αi Θi C C T ρTl C T

(11) (12)

i=1

It should be noted that constraint (11) serves to connect the LMI constraints with the objective function. Since objective function will always strive to make the σl s small, we can expect the constraints of (11) to always be active at the solution. It should also be highlighted that it may seem more natural to assume the activity of this constraint and replace the ζl s of (12) with σl2 . However, this would result in the constraints of (12) becoming non-convex since they would no longer be LMIs. Similarly, one could replace the σl s of the objective function √ with ζl s. While this non-convexity will be about the same as the original formulation, it will be incongruities with the GBD approach to be discussed next.

1.2. Review of the Generalized Benders Decomposition Consider the optimization problem

min F (x, y) s.t. G (x, y) ≤ 0

x∈X,y∈Y

(13)

∧ { where y is a vector of complicating variables in that for a given y¯ ∈ V = y | G(x, y) ≤ } 0 for some x ∈ X , the following primal problem is easily solved

min F (x, y¯) s.t. G (x, y¯) ≤ 0 x∈X


(P)

Page 7 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


If X is a convex set and both F (x, y) and G (x, y) are convex in x, then the following relaxed master problem will converge to the global solution of (13), see Geoffrion [15] for details.

min y0

(14)

y∈Y,y0

( ) L∗ y, µ(k1 ) ≤ y0 k1 = 1 . . . K1 ( ) L∗ y, λ(k2 ) ≤ 0 k2 = 1 . . . K2

(15) (16)

where L∗ and L∗ are defined as { } L∗ (y, µ) = min F (x, y) + µT G (x, y)

(17)

{ } L∗ (y, λ) = min λT G (x, y)

(18)

x∈X

x∈X

and µ(k1 ) and λ(k2 ) are selected by the following procedure: 1. Select a point y¯ ∈ V ∩ Y . Solve the primal problem (P) and obtain its solution, x∗ , along with the multipliers associated with the constraints G(x, y¯) ≤ 0, µ∗ . Set K1 = 1, ( ) K2 = 0, µ(1) = µ∗ , the upper bound U BD = F (x∗ , y¯) and determine L∗ y, µ(1) . 2. Solve the relaxed master problem (14) and obtain its global solution, (y0∗ , y ∗ ). Set the lower bound LBD = y0∗ . If U BD ≤ LBD + ϵ, the problem has converged. 3. Set y¯ = y ∗ . If the primal problem (P) is feasible go to step 3a; if infeasible go to 3b. 3a. Solve the primal problem and obtain the solution and multipliers (x∗ , µ∗ ). Set ( ) K1 = K1 + 1, µ(K1 ) = µ∗ and determine L∗ y, µ(K1 ) . If F (x∗ , y¯) < U BD set U BD = F (x∗ , y¯). Return to step 2. 3b. Solve the following problem (denoted Problem Q): { } { T } max min λ G (x, y¯) λ∈Λ

x∈X

(19)

where Λ = {λ ≥ 0 | λT 1 = 1} and 1 = [1 1 . . . 1]T . Obtain its solution (x∗ , λ∗ ). ) ( Set K2 = K2 + 1, λ(K2 ) = λ∗ , determine L∗ y, λ(K2 ) and return to step 2. 7 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 29

As shown in references [13], [16] and [17] that the solution to Problem Q can be found from the following problem − as its solution along with its associated multipliers (x∗ , λ∗ ). Specifically, λ∗ are the optimal dual variables associated with the constraints G(x, y¯) − π1 ≤ 0. min {π} s.t. G (x, y¯) − π1 ≤ 0

(Q)

x∈X,π

If the solution to (Q) is such that π ∗ ≤ 0, then the primal problem would have been feasible. If F and G are linearly separable (i.e., F (x, y) = F1 (x)+F2 (y) and G(x, y) = G1 (x)+G2 (y)), then both L∗ and L∗ can be calculated as explicit functions of y, [15].

2. Application of GBD to Value-optimal SND for Steady-state Systems To employ the GBD approach, the attendance variable αi must be split into two parts (cc)

- one part for the cost of the component, αi , and a second for its impact on system (sp)

performance, αi

, [14]. Using an extension of Theorem 1 of [14], the equivalent form of

problem (10) is obtained as

min (cc)

αi

ns ∑

(sp)

,αi ζl ,σl

+

i=1 (cc)

s.t. αi

(sp)

αi

(s) (cc) ci αi

nq ∑

0.19947Ks T σl

(20)

l=1

∈ {0, 1} (cc)

− αi

(21)

≤ 0, i = 1 . . . ns

ζl ≤ σl2 , l = 1 . . . nq   ζl ρl C  )  ≽ 0, (n   ∑s (sp) T T T C ρl C αi Θi C

(22) (23) l = 1 . . . nq

(24)

i=1 (sp)

0 ≤ αi

≤ 1, i = 1 . . . ns


(25)

Page 9 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


It is clearly observed that constraints (24) and (25) are convex. Thus, the non-complicating (sp)

variables, x, are selected as αi

(sp)

and ζl , and the set X is defined as {αi

, ζl | (24) and (25) are

(cc)

satisfied}. Then, the complicating variables, y, should be selected as αi

and σl , and the set

(cc)

Y is defined as {αi , σl | (21) is satisfied}. The connecting constraints G(x, y) ≤ 0 consist of constraints (22) and (23). Clearly, F (x, y) and G(x, y) and X are convex with respect to x and thus a globally optimal solution will be guaranteed. In addition, F (x, y) and G(x, y) are linearly separable, indicating that L∗ and L∗ can be determined as explicit functions of y. Given these selections, it is noted that the objective is a function of only the complicating variables y (i.e. F (x, y) = F2 (y)). Hence, the primal problem is just a feasibility problem: ns ∑

min (sp)

αi

,ζl

(s) (cc) ci α ¯i

+

i=1 (sp)

s.t. αi

nq ∑

0.19947Ks T σ ¯l

(P1)

l=1 (cc)

−α ¯i

≤ 0, i = 1 . . . ns

ζl − σ ¯l2 ≤ 0, l = 1 . . . nq (24), (25)

The first step of the GBD algorithm requires y¯ ∈ V ∩ Y . An obvious choice is to select (cc)

αi

= 1 for all i, and σl to be a sufficiently large number for all l. If the resulting primal

problem is infeasible, then terminate, as the original optimization problem is infeasible. If the primal is feasible then the optimal multipliers, µ∗ , will be zero. In this case, UBD=F2 (¯ y) and L∗ (y, µ(1) ) = F2 (y). Then, in the first execution of step 2, problem (14) with (K1 = 1 and K2 = 0), gives the solution y ∗ = 0 and LBD=0. In the first execution of step 3 with y¯ = 0, the primal is obviously infeasible. Thus, the following version of problem Q should



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 29

be solved to obtain the solution x∗ and the associated multiplier λ∗ .

min {π} s.t.

(sp)

αi

(sp)

αi

(Q1)

,ζl ,π (cc)

−α ¯i

− π ≤ 0, i = 1 . . . ns

(26)

ζl − σ ¯l2 − π ≤ 0, l = 1 . . . nq

(27)

(24), (25)

Then, returning to step 2 to solve (14), with K2 incremented by 1, one will need to calculate L∗ (y, λ(1) ) where λ(1) are the optimal multipliers of (Q1) associated with (26) and (27). If λ is defined as λ = [γ1 γ2 . . . γns ω1 ω2 . . . ωnq ]T , then the solution and optimal multipliers ∑ (k ) (cc) (sp,k ) (k ) (k ) (k ) (x(k2 ) , λ(k2 ) ) are denoted as (αi 2 , ζl 2 , γi 2 , ωl 2 ) and L∗ (y, λ(k2 ) ) = βk2 − i γi 2 αi − ∑ (k2 ) 2 ∑ (k2 ) (sp,k2 ) ∑ (k2 ) (k2 ) ω σ where β αi + l ωl ζl . Therefore, the relaxed master problem = k 2 l l l i γi is found to be: {∑ ns

min (cc)

αi

∈{0,1},δl

βk2 −

ns ∑

(s) (cc) ci αi

(k ) (cc) γi 2 αi

i=1

+

i=1

nq ∑

} 0.19947Ks T σl

s.t.

(RM1)

l=1

−

nq ∑

(k2 ) 2 σl

ωl

≤ 0, k2 = 1 . . . K2

(28)

l=1

If at any iteration the solution to (RM1) is such that the primal problem, (P1), is feasible, then the solution to the primal will have an objective value equal to the previous relaxed master. In this case, the algorithm terminates since UBD=LBD. To summarize, the GBD algorithm implements an iteration between two subproblems: the relaxed master problem (RM1) and problem Q (Q1). Problem Q is a convex semi-definite program without integer variables, and requires the use of a SDP solver. The SDP solver (k2 )

used must be capable of generating the optimal dual variables (γi

(k2 )

, ωl

) associated with

constraints G(x, y) − π1 ≤ 0. Since constraint (28) is nonconvex (due to the term σl2 ), the relaxed master is a Mixed


Page 11 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Integer Nonlinear Program (MINLP). As such, the branch and bound method must be used to solve (RM1), which will greatly increase computational effort. The following subsection will illustrate how to remove this non-convexity from the relaxed master to arrive at a MILP.

2.1. Convexification of the relaxed master problem Consider the following restatement of the value-optimal optimal SND problem min F (x, y) s.t. G′ (x, y) ≤ 0

(29)

x∈X,y∈Y

where G′ (x, y) consists of αi

(sp)

(cc)

− αi

and

√ ζ l − σl .

Theorem 1 Problem (29) is equivalent to problem (20), and no duality gap exists when implementing GBD. √ Proof. For all ζ ≥ 0 and σ ≥ 0, ζ ≤ σ if and only if ζ ≤ σ 2 . Moreover, for all σ ¯ ≥ 0 and { } { } α ¯ (cc) , XG′ = ζ ≥ 0, 0 ≤ α(sp) ≤ 1 | G′ (x, y¯) ≤ 0 = ζ ≥ 0, 0 ≤ α(sp) ≤ 1 | G(x, y¯) ≤ 0 is a convex region. F (x, y) is also convex in x. Thus no duality gap exits when implementing

the GBD approach.

In this case, problem (Q) will need to be solved with G (x, y¯) − π1 ≤ 0 replaced by √ (sp) (cc) G′ (x, y¯) − π1 ≤ 0, or αi − α ¯ i − π1 ≤ 0 and ζl − σ ¯l − π1 ≤ 0. Clearly, this second group of constraints is non-convex. To address this issue, consider problem (19) using G′ (x, y¯): { max γi ≥0 ωl ≥0

s.t.

{ min x∈X

(sp)

γi (αi

(cc)

−α ¯i ) +

i

∑ i

∑

γi +

∑

}} √ ωl ( ζ l − σ ¯l )

(Q2)

l

∑

ωl = 1

l

The following theorem provides a method of determining the solution to this non-convex problem via the solution to a convex problem. (sp)∗

Theorem 2 The solution to problem (Q2), (αi

, γi∗ , ζl∗ , ωl∗ ), is given by αi


(sp)∗

′(sp)∗

= αi

, γi∗ =


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

γi′∗ /a, ζl∗ = ζl′∗ and ωl∗ = ωl′∗

√

ζl′∗ /a where a =

∑

′∗ i γi +

∑

′∗ l ωl

√

Page 12 of 29

′(sp)∗

ζl′∗ and (αi

, γi′∗ , ζl′∗ ,

ωl′∗ ) is the solution to: { max ′

{

∑

min

x′ ∈X

γi ≥0 ωl′ ≥0

∑

s.t.

′(sp) γi′ (αi

−

(cc) α ¯i )

+

i

γi′ +

∑

√ ¯l ) ωl′ (ζl′ − ζl′ σ

}} (Q3)

l

∑

i

ωl′ = 1

l

ˆ ∗ /a where λ ˆ ∗ is Proof. For any positive constant a, the solution to problem (Q2) is λ∗ = λ the solution to { max

x∈X

∑ {

= max

{

∑

= max ′

∑

{ min

∑ i

(sp) γî (αi

−

(cc) α ¯i )

−

(cc) α ¯i )

ω ˆl = a

∑

∑

(Q4)

}} ∑ ω √ ˆl √ (ζl − ζl σ + ¯l ) ζl l

l

i

γi′ +

}} √ ω ˆ l ( ζl − σ ¯l )

ω ˆl = a

∑

x′ ∈X

λ ≥0

+

∑ l

i

γî +

{

−

(cc) α ¯i )

l

x∈X

i

s.t.

∑

min

ˆ λ≥0

(sp) γî (αi

i

γî +

i

s.t.

∑

min

ˆ λ≥0

s.t.

{

(sp) γi′ (αi

+

∑

ωl′ (ζl′

}} √ − ζl′ σ ¯l )

l

√ ωl′ ζl′ = a

l

∑ √ ∑ Since a can be selected as any positive constant, the constraint i γi′ + l ωl′ ζl′ = a may be ∑ √ ∑ removed as long as i γi′ + l ωl′ ζl′ is guaranteed to be positive. Since ζl′ > 0 for all l, this ∑ ∑ guarantee is achieved by enforcing i γi′ + l ωl′ = 1. Thus, given the solution, (α′(sp)∗ , ζ ′∗ , λ′∗ ) of (Q3), the solution to (Q2) is αi √ ω ˆ l∗ /a = ωl′∗ ζl′∗ /a.

(sp)∗

′(sp)∗

= αi

, γi∗ = γî∗ /a = γi′∗ /a, ζl∗ = ζl′∗ and ωl∗ =

Given Theorem 2, if the solution to problem (Q3) is (α′(sp,k2 ) , ζ ′(k2 ) , γ ′(k2 ) , ω ′(k2 ) ), then


Page 13 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


the appropriate relaxed master problem is easily obtained as {n s ∑

min (cc)

αi

∈{0,1},σl

s.t. βk2 −

∑

(s) (cc) ci α i

} 0.19947Ks T σl

(RM2)

l

(k2 )

γi

(cc)

αi

i

∑

+

i

nq ∑

−

∑

(k2 )

ωl

σl ≤ 0, k2 = 1 . . . K2

l

′(k2 ) ′(sp,k2 ) ∑ γ αi + l i i

′(k ) ′(k ) ωl 2 ζl 2 )/a,

(k ) γi 2

′(k ) γi 2 /a,

(k ) ωl 2

′(k ) ωl 2

√

′(k )

= ζl 2 /a, and where βk2 = ( = √ ∑ ′(k2 ) ∑ ′(k2 ) ′(k ) a = i γ i + l ωl ζl 2 . Furthermore, the solution to problem (Q3) can be obtained from

min {π}

′(sp)

αi

,ζl′ ,π

′(sp)

(cc)

−α ¯ i − π ≤ 0, i = 1 . . . ns √ ζl′ − ζl′ σ ¯l − π ≤ 0, l = 1 . . . nq

s.t. αi

(24), (25) =

min

′(sp)

αi

,ζl′ ,τl >0,π ′(sp)

{π}

(Q5)

(cc)

− π ≤ 0, i = 1 . . . ns   ′  ζl τl  ζl′ − τl σ ¯l − π ≤ 0,   ≽ 0, l = 1 . . . nq τl 1

s.t. αi

−α ¯i

(24), (25)

Clearly, this problem is convex and a global solution can be readily obtained. Summarizing the GBD algorithm, it implements an iteration between the convexified relaxed master problem and problem (Q5) (see Figure 1). The convexified relaxed master (RM2) is a mixed integer linear program (MILP) without any matrix inequalities. Thus, it can be solved using a variety of MILP solvers. It is also noted that the convexified relaxed master at each iteration is the same as the previous except that a new linear constraint is



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 29

Covexified Relaxed Master (RM2) (MILP solved with GUROBI) J i( k ) , Zl( k Ek

Di( cc ) Vl

2

2)

2

Problem (Q5) (SDP solved with MOSEK) Figure 1: Illustration of GBD algorithm. added. In the examples of this paper, the MILP solver GUROBI [18] was used. Problem (Q5) is a convex semi-definite program without integer variables, and can be solved by a SDP solver. As mentioned before, the SDP solver used must be capable of generating the (k2 )

optimal dual variables (γi

(k2 )

, ωl

) associated with constraints G(x, y) − π1 ≤ 0. The SDP

solver MOSEK [19] was employed. Both solvers were used via YALMIP [20] in a MATLAB environment.

2.2. Case Studies for Steady-State Systems This subsection will demonstrate the computational advantage of the GBD approach compared to the brand and bound method. All examples are solved on 64-bit windows PC with Intel 2.5 GHz CPU and 8 GB RAM.

s2 s1

Node 1

s3 Node 2

s4

Figure 2: Process flow diagram for Example 1.

Example 1: Consider the flow diagram of Figure 2, [1]. If s2 and s3 are selected as primal


Page 15 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


variables, the matrix C of equation (1) is C = [C0T C0T C0T ]T where 

T

 1 1 0 0  C0 =   1 0 1 1

(30)

The nominal flow rates are assumed to be s = [sT0 sT0 sT0 ]T where s0 = [150.1, 52.3, 97.8, 97.8]T (kg/min). The sensor options are 1%, 2%, or 3% relative errors. The resulting measurement precisions and sensor costs are presented in Table S1 in the Supporting Information. The objective is to minimize the sum of the sensor cost and the DEFL associated with streams 1 and 4. The parameter Ks is assumed to be $3600/year, and T is assumed to be 5 years. Both the branch and bound and the GBD methods generate the same solution: 1% sensor placed at streams 2 and 3, the standard deviation of the estimation error on streams 1 and 4 are found to be 1.1163 and 1.0205, with an optimal objective function value of $12,672. The GBD approach achieves a 98.6% reduction in computational effort, as indicated in Table 1. Table 1: Comparison of computational efforts for Example 1. Method

Iterations

Branch and Bound GBD

612 11

6 5

1

7

16 2

9 18

3

10

4

2 3

5

1

7 22

20

19

17 6 21

4

Percent Reduction 98.6%

8

15 11

SDP Time Calls (sec) 1224 492 11 7

11

10 23

24

8 13

9 14

12

Figure 3: Process flow diagram for Example 2, [22] Example 2: This example is adapted from Example 3 of [21], which was adapted from 15 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 29

[22]. Consider the 24 stream network of Figure 3. The flow rate and the cost of sensors for each stream are listed in Table S2 in the Supporting Information. The precision of sensors are assumed to be 2.5%. The design objective is to minimize the sum of sensor cost and the DEFL associated with steams 3, 10, 16, 17, 20 and 24. The parameter Ks is assumed to be 10, and T is assumed to be 5. The branch-and-bound and GBD approach both obtain the same solution, which requires sensors at streams 4, 8, 10, 16, 20, 23 and 24, with the resulting standard deviation of estimation error on streams 3, 10, 16, 17, 20 and 24 being 2.7028, 2.5000, 2.4971, 0.1292, 0.7496 and 1.1242, respectively. The optimal objective function value is $188. Comparison of the computational efforts for the two approaches is given in Table 2 and illustrates a computational effort reduction of almost 80%. Table 2: Comparison of computational efforts for Example 2. Method

Iterations


2222 257



3. Value-Optimal SND for Dynamic Process This section illustrates application of the GBD to the value-optimal SND problem for closed-loop processes. The first subsection presents a short review of the nonconvex mixed integer programming formulation for the value-optimal SND problem, while the second describes the application of GBD to this class of problems and demonstrates the computational advantage of GBD compared to the branch and bound approach.

3.1. Review of Value-Optimal SND for Dynamic Process Consider a dynamic process model: s˙ = f (s, m, p), q = h(s, m, p), and q min ≤ q ≤ q max , where s, m and q are the state, manipulated and performance variables, each being a column 16 ACS Paragon Plus Environment

Page 17 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


vector with dimension of ns , nm and nq , respectively, and p is a white noise process with a mean p¯ and spectral density of Σp . The Optimal Steady-State Operating (OSSOP) can be obtained from the following optimization problem:

min

sss ,mss q min ≤q ss ≤q max

{g(q ss )} s.t.

0 = f (sss , mss , p¯)

(31)

q ss = h(sss , mss , p¯)

where g(·) is the operating cost function. The solution to problem (31) is denoted as sOSSOP , mOSSOP , q OSSOP . If the process is operated at the OSSOP, then the actual non-steady-state disturbances (i.e., p ̸= p¯) will make q deviate from q OSSOP and cause significant constraint violations. Thus, a Backed-off Operating Point (BOP) should be selected that is some distance from the operating constraints, but still close to the OSSOP, as indicated in Figure 4. The BOP is denoted as (sBOP , mBOP , pBOP ), and must satisfy the steady-state process model:

0 = f (sBOP , mBOP , p¯), q BOP = h(sBOP , mBOP , p¯)

(32)

q min ≤ q BOP ≤ q max

(33)

To characterize the Expected Dynamic Operation Region (EDOR), we start with a linearization of the nonlinear dynamic model s˜˙ = A˜ s + Bm ˜ + G˜ p, q˜ = Ds s˜ + Dm m ˜ where the deviation variables (˜ s, m, ˜ p˜, q˜) are defined as s˜ = s−sBOP , m ˜ = m−mBOP , p˜ = p−pBOP , and q˜ = q − q BOP . The manipulated variable is a linear feedback of the state estimate: m ˜ = Lsˆ˜, where L is to be determined and sˆ˜ is the optimal state estimate of s˜. If the estimation error is e = s˜ − sˆ˜, then its covariance Σe is calculated as: AΣe + Σe AT + GΣp GT − Σe C T Σ−1 v CΣe = 0


(34)


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Steady-State Operating Line

Page 18 of 29

Expected Dynamic Operating Regions Minimally BaFked-off Operating Point

Different Controller Tuning Values

Optimal Steady-State Operating Point

Figure 4: Illustration of the MBOP problem. where C and Σv are from the measurement equation δ˜ = C s˜ + v where v is a zero-mean ns ∑ white noise with spectral density of Σv . Note that Σ−1 αi Θi where Θi is as defined in = v i=1

(6). Then, using the covariance analysis of [7], the variance of lth output of q is denoted as ζl and is equal to the lth diagonal of Σq : ( ) ζl = ρl Σq ρTl = ρl Ds Σe DsT + (Ds + Dm L)(Σs − Σe )(Ds + Dm L)T ρTl , l = 1 . . . nq

(35)

where ρl is defined as the lth row of an nq × nq identity matrix, and Σs is the covariance matrix associated with s and is found as the positive define solution to

AΣs + Σs AT + GΣw GT + BL(Σs − Σe ) + (Σs − Σe )LT B T = 0

(36)

One can impose upper bounds to the EDOR in each direction by enforcing the constraints √ ζl ≤ σl . However, to guarantee the resulting BOP and its associated EDOR is contained √ in the region defined by q min and q max , the following constraints must be enforced: ζl ≤


Page 19 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


qlBOP − qlmin and

√ ζl ≤ qlmax − qlBOP , l = 1 . . . nq . These can be written equivalently as σl ≤ qlBOP − qlmin σl ≤ qlmax − qlBOP , l = 1 . . . nq √ ζl ≤ σl , l = 1 . . . nq

(37) (38)

The value-optimal SND for closed-loop dynamic processes is then a result of combining the steady-state perspective (BOP) and the dynamic perspective (EDOR). { min

sBOP ,mBOP ,q BOP

g(q BOP ) +

ns ∑

} (s)

ci αi

i=1

σ>0,ζ>0,αi ∈{0,1}

s.t. (32), (33), (34), (35), (36), (37), (38)

Theorem 1 of Peng and Chmielewski [10] indicates that the nonlinear equalities of (34)-(36) hold if and only if the LMIs of (40)-(42) are satisfied. As a result the SND can be states as: { min

sBOP ,mBOP ,q BOP ,αi ∈{0,1} σ>0,ζ>0,M0 ≥0,M2 ≥0,M1

g(q BOP ) +

ns ∑

} (s) ci α i

s.t.

(39)

i=1

(32), (33), (37), (38) (n )   ∑s T T αi Θi C − A M2 − M2 A M2 G   C i=1  ≽0 T −1 G M2 Σw  ζl ρl (Ds M0 + Dm M1 ) ρl Ds    (D M + D M )T ρ T M0 I s 0 m 1 l   I M2 DsT ρTl

(40)     ≽ 0, l = 1 . . . nq  

(41)

(AM0 + BM1 ) + (AM0 + BM1 )T + GΣw GT ≼ 0

(42)

ζl ≤ (qlmax − qlmin )2 /4, l = 1 . . . nq

(43)

Note the addition of (43). While (43) is redundant with (37), it will be useful in the following section. If αi∗ , M0∗ , M1∗ , M2∗ is part of the optimal solution, then a feasible controller under 19 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 29

this configuration can be calculated as L = M1∗ (M0∗ − (M2∗ )−1 )−1 . It should be noted that the formulation of the value-optimal SND problem for a discrete-time dynamic process can be obtained through a simple modification of constraints (32), (33), (40) and (42), [7].

3.2. Application of GBD Similar to the steady-state SND case, the attendance variable, αi is split into two parts: (cc)

the component associated with cost, αi , and the component impacting system performance, (sp)

αi

. Then the problem (39) can be equivalently expressed as: { g(q BOP ) +

min (cc) (sp) sBOP ,mBOP ,q BOP ,αi ,αi σ>0,ζ>0,M0 ≥0,M2 ≥0,M1

ns ∑

} (s) (cc) ci αi

(44)

i=1

s.t. (32), (33), (37) (cc)

∈ {0, 1}

(sp)

− αi

αi

αi √

(cc)

(45)

≤ 0, i = 1 . . . ns

(46)

ζl ≤ σl , l = 1 . . . nq (n )  ∑s (sp) T T αi Θi C − A M2 − M2 A M2 G   C i=1  ≽0 T −1 G M2 Σw

(47)



(48)

(41), (42), (43) (sp)

0 ≤ αi

≤ 1, i = 1 . . . ns

(49) (cc)

Based on this reformulation the complicating variables, y, should be selected as αi , sBOP , mBOP , q BOP and σl , and the set Y is the constraints (32), (37) and (45). (sp)

complicating variables, x, are selected as αi

The non-

, ζ, M0 , M1 and M2 , and the set X is constraints

(41), (42), (43), (48) and (49). Finally, the connecting constraints, G′ (x, y) ≤ 0, are (46) and (47). Again, G′ (x, y) is linearly separable, and thus L∗ and L∗ can be determined as explicit functions of y. Using Theorem 1 again, it is found that the GBD algorithm is guanranteed to find a global solution. However, similar to the results of Section 2.1, use 20 ACS Paragon Plus Environment

Page 21 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


of G′ (x, y) will result in Problem Q being non-covex. Fortunately, Theorem 2 can be re√ (k2 ) ′(k2 ) (k2 ) ′(k2 ) ′(k ) applied such that the solution to problem Q is γi = γi /a, ωl = ωl ζl 2 /a, and √ ∑ ′(k2 ) ∑ ′(k2 ) ′(k ) ′(k ) ′(k ) ′(k ) a = i γi + l ωl ζl 2 where ζl 2 , γi 2 and ωl 2 are the solution and multipliers from

min

′(sp)

,ζl′ >0 π,αi M0 ≥0,M2 ≥0,M1 ′(sp)

π

(Q6)

(cc)

− π ≤ 0, i = 1 . . . ns   ′  ζl τl  ζl′ − τl σ ¯l − π ≤ 0,   ≽ 0, l = 1 . . . nq τl 1 s.t. αi

−α ¯i

(41), (42), (43), (48), (49)

Figure 5: Furnace reactor system.


(50) (51)


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 29

3.3. Case Studies for Closed-loop Systems Example 3: Consider the furnace reactor system of Figure 5. The system matrices are 

−10

 0

0

0

0

      8000 −8000 0 0 0      A= 0 2000 −1500 0 0         0 0 0 −5000 0   0 0 0 0 −5000     0 0 0 10          −75 75000    0    0       G= 0  B= −25 0 0             5  0 −8500 8.5 × 10   0      5 0 0 −500 × 10 0 The process state variables 1, 2 and 3 correspond to T0 , TF and TR , respectively, and states 4 and 5 correspond to the O2 and CO concentrations in the furnace. The manipulated variables 1, 2 and 3 correspond to feed flow rate, Ff eed , fuel flow rate, Ff uel , and furnace vent position, Pv , respectively. The nominal value of the state and manipulated variables are snom = [300K; 375K; 500K; 4%; 100ppm] and mnom = [10000 bbl/day; 10 bbl/day; 0.1%]. The disturbance is a white noise process with a mean of p¯ = 300K and a spectral density of Σp = 100. The constraints on state variables and manipulated variables are: 370K ≤ TF ≤ 380K, 495K ≤ TR ≤ 505K, 1% ≤ CO2 ≤ 7%, 70ppm ≤ CCO ≤ 130ppm, 9900bbl/day ≤ Ff eed ≤ 10100bbl/day, 8bbl/day ≤ Ff uel ≤ 12bbl/day, 0.09% ≤ Pv ≤ 0.11%. It is assumed that two types of sensors are available for each state variable. The data of the sensors is provided in Table 3 and generates a matrix C equal to [cT1 cT1 cT2 cT2 cT3 cT3 cT4 cT4 cT5 cT5 ]T where ci = ρi .


Page 23 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Table 3: Sensor data for reactor furnace example. (s)

i

Type

Location

σ ¯v2i

ci ($)

1

T emp

c1

0.0625

300

2

T emp

c1

0.25

200

3

T emp

c2

0.0625

300

4

T emp

c2

0.25

200

5

T emp

c3

0.0625

300

6

T emp

c3

0.25

200

7

O2 Concentration

c4

0.1

1800

8

O2 Concentration

c4

0.4

1200

9

CO Concentration

c5

0.1

2400

10 CO Concentration

c5

0.4

1600

The economic objective function (to be maximized) is:

g = 10Ff eed − 30Ff uel − 0.1CCO −

ns ∑

(s)

ci αi

i=1

Both the branch and bound method and the GBD approach obtain the same optimal solution: two temperature sensors with variance of 0.0625 measuring the feed flow temperature and furnace temperature, with an optimal objective function value of $99,175. The GBD implementations achieves over 84% improvement in computation time, as indicated in Table 4. Table 4: Comparison of computational efforts for Example 3. Method

Iterations


78 14



Example 4a: Consider the jacketed CSTR system of Figure 6, adapted from Ref. [23].



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(TC) Temp. Controller

Page 24 of 29

VP Fvg

Fo, To, CAo

To Compressor Section

Pressure Controller

From Feed Section

Main Reactor

(PC)

Fc, Tc

VT

VL F4

Level Controller (VC)

Fco, Tco

F,T,CA

Coolant Reactor Coolant Pump

F2

To Stripper Section F3

Product Pump

Material Flow Information Flow

Figure 6: CSTR system of Example 4. The process is described by the following material and energy balances: F0 F C˙ A = CA0 − CA − k0 e−E/RT CA V V F F (−∆H) −E/RT U A(T − Tc ) 0 T˙ = T0 − T + k0 e CA − V V ρCp V ρCp U A(T − Tc ) Fc T˙c = (Tc0 − Tc ) + Vj Vj ρj Cpj V˙ = F0 − F P˙ =

RT (V k0 e−E/RT CA − Fvg ) 64 − V

The state, manipulated and disturbance vectors are s = [CA T Tc V P ]T , m = [F Fc Fvg ]T and p = [F0 CA0 ]T . The manipulated inputs F , Fc and Fvg have nominal values of 40 f t3 /h, 56 f t3 /h and 10.6 lb · mol/h, respectively. The disturbance inputs F0 and CA0 have mean values of 40 f t3 /h and 0.5 lb · mol/f t3 , respectively. The spectral density of F0 and CA0 equal to 0.12 and 0.012 , respectively. The operating constraints are as follows. 0.15 lb · mol/f t3 ≤ CA ≤ 0.35 lb · mol/f t3 , 595 ◦ R ≤ T ≤ 650 ◦ R, 560 ◦ R ≤ Tc ≤ 620 ◦ R, 40 f t3 ≤ V ≤ 55 f t3 , 1600 lbf /f t2 ≤ P ≤ 2600 lbf /f t2 , 0 ≤ F ≤ 80 f t3 /h, 0 ≤ Fc ≤ 200 f t3 /h, 0 ≤ Fvg ≤ 20.5 lb · mol/h. It is assumed that 1% precision sensors are available 24 ACS Paragon Plus Environment

Page 25 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


at each state, and each sensor has an annualized cost of $1000/yr. The economic objective function (to be maximized) is defined as: ns [ ] ∑ (s) g = 8760 × 0.375(CA0 − CA )F − 0.015Fc − 0.00225Fvg − ci αi i=1

Given the process constraints and objective function, the OSSOP is found to be sOSSOP = [0.15 lb·mol/f t3 622.5 ◦ R 610 ◦ R 40 f t3 2100 lbf /f t2 ], mOSSOP = [40 f t3 /hr 56.4 f t3 /hr 14 lb· mol/h], with optimal objective function value of $43, 724/yr. Linearizing the nonlinear dynamic model and objective function around the OSSOP, one can obtain the linear model x˙ = Ax + Bu + Gw, g = g(sOSSOP , mOSSOP , p¯) + dTg [xT uT wT ]T , where dTg is the partial derivative of the objective function evaluated at the OSSOP, and A, B, G and dg are given Table S3 in the Supporting Information. The branch and bound method and the GBD approach are used to solve this sensor selection problem. Both approaches yield the same optimal solution: 1% sensors measuring CA , V and P , respectively, with an optimal objective function value of $33, 369/yr. The branch and bound method required 660 iterations and a runtime of 2,125 s. Using the GBD approach, only 48 iterations and 71 s were required. Example 4b: Continuing Example 4a, remove the assumption of a linear steady-state model and use the nonlinear steady-state model and nonlinear objective function directly. To make the nonlinear relaxed master problem easier to solve, the nonlinear model is converted to a dimensionless variable form. Then, the relaxed master problem (with mixed integer nonlinear constraints and nonlinear objective function) is solved using SCIP [24]. In this case, the GBD approach required 75 iterations and 184 s. A summary of computational effort for the 3 cases of Example 4 is provided in Table 5.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 29

Table 5: Comparison of computational efforts for Example 4. Method

Iterations

Branch and Bound (Linear steady-state model) GBD (Linear steady-state model) GBD (Nonlinear steady-state model)

660 48 75

SDP Time Calls (sec) 1320 2125 48 71 75 184

Percent Reduction 96.7% 91.3%

4. Conclusion In this work, we have shown that the GBD algorithm can be a powerful tool for solving value-optimal SND problems. In addition to the improvement in computational efficiency, the proposed algorithm decomposed the problem in such a way that off-the-shelf solvers can be employed. Thus, significantly less effort will be required to implement as compared to previous procedures, which required an ad’hoc development of the branch-and-bound algorithm to be compatible with the SDP solver. It should also be noted that the fairly basic problem formulations of this article can be extended (adding gross-error detectablility or resilience performance criteria as well as adding actuator selection to the closed-loop formulation), and still be able to utilize the GBD algorithm. In citations [5–7], it is shown that these extensions can be made without changing the overall structure of the problem. One extension that will require a more careful consideration is when bias detection is involved. This problem is the subject of future research.

Associated Content The Supporting Information is available free of charge via the Internet at http://pubs.acs.org/. Precision, relative error and sensor cost for Example 1; Flow rate and flowmeter cost of each stream for Example 2; The values of matrices A, B, G and the vector dg for Example 4.


Page 27 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Acknowledgement Both authors thank the National Science Foundation (CBET-1511925) for financial support.

References 1. Bagajewicz M. J. Design and retrofit of sensor networks in process plants. AICHE J. 1997, 43 (9), 2300-2306. 2. Bagajewicz M. J. Process Plant Instrumentation: Design and Upgrade; Technomic Publishing Company: Pennsylvania, U.S.A., 2001. 3. Kubrusly C. S.; Malebranche H. Sensor and controllers location in distributed systema survey. Automatica, 1985, 21 (2), 117-128. 4. Padula S. L.; Kincaid R. K. Optimization strategies for sensor and actuator placement. NASA/TM-1999-209126, 1999. 5. Chmielewski D. J.; Peng J. K. Covariance-based hardware selection - part I: globally optimal actuator selection. IEEE Trans. Control Syst. Technol. 2006, 14 (2), 355-361. 6. Chmielewski D. J.; Palmer T.; Manousiouthakis V. On the theory of optimal sensor placement. AICHE J. 2002, 48 (5), 1001-1012. 7. Peng J. K.; Chmielewski D. J. Covariance based hardware selection-part II: equivance results for the sensor, actuator and simultaneous selection problems. IEEE Trans. Control Syst. Technol. 2006, 14 (2), 362-368. 8. Bagajewicz M. J.; Markowski M.; Budek A. Economic value of precision in the monitoring of linear systems. AIChE J. 2005, 51 (4), 1304-1309.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 28 of 29

9. Nguyen D. Q.; Bagajewicz M. J. New efficient breadth-first/level traversal tree search method for the design and upgrade of sensor networks. AIChE J. 2011 57 (5), 1302-1309. 10. Peng J. K.; Chmielewski D. J. Optimal sensor network design using the minimally backed-off operating point notion of profit. Proc. of the Am. Cont. Conf., Portland, OR, USA, June 2005, WeA07.4, 220-224 11. Peng J. K.; Manthanwar A. M.; Chmielewski D. J. On the tuning of predictive controllers: the minimum back-off operating point selection problem. Ind. Eng. Chem. Res. 2005, 44, 7814-7822. 12. Omell B. P.; Chmielewski D. J. On the tuning of predictive controllers: impact of disturbances, constraints, and feedback structure. AIChE J. 2014, 60 (10), 3473-3489. 13. Zhang J.; Omell B. P.; Chmielewski D. J. On the tuning of predictive controllers: application of generalized Benders decompostion to the ELOC problem. Comput. Chem. Eng. 2015, 82, 105-114. 14. Zhang J.; Wang X.; Chmielewski D. J. Covariance-based hardware selection part IV: solution using the generalized Benders decomposition. AICHE J. 2016, 62, 3628-3638. 15. Geoffrion A. M. Generalized benders decomposition. J. Optim. Theory Appl. 1972, 10, 237-260. 16. Bagajewicz M. J.; Manousiouthakis V. On the generalized benders decomposition. Comput. Chem. Eng. 1991, 15, 691-700. 17. Floudas C. A.; Aggarwal A.; Ciric A. R. Global optimum search for nonconvex NLP and MINLP problems. Comput. Chem. Eng. 1989, 13 (10), 1117-1132. 18. Gurobi

Optimizer

Reference

Manual.

Gurobi

http://www.gurobi.com


Optimization,

Inc.,

2015.

Page 29 of 29

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


19. MOSEK ApS. The MOSEK optimization toolbox for MATLAB manual. Version 7.1 (Revision 28), 2015. http://docs.mosek.com/7.1/toolbox/index.html 20. J. Löfberg. YALMIP : A Toolbox for Modeling and Optimization in MATLAB. Proceedings of the CACSD Conference, Taipei, Taiwan, 2004. 21. Gala M.; Bagajewicz M. J. Rigorous methodology for the design and upgrade of sensor networks using cutsets. Ind. Eng. Chem. Res. 2006, 45, 6679-6686. 22. Madron F.; Veverka V. Optimal selection of measuring points in complex plants by linear models. AICHE J. 1992, 38(2), 227-236. 23. Bhushan M.; Rengaswamy R. Design of sensor network fased on the signed directed graph of the process for efficient fault diagnosis. Ind. Eng. Chem. Res. 2000, 39, 999-1019. 24. Achterberg T. SCIP: solving constraint integer programs. Math. Prog. Comp. 2009, 1 (1), 1-41.

Process Monitoring

Value-Optimal SND Problem

Process Control

Solved by

Computational Efficiency

GBD

>>

Brand and Bound

Figure 7: Table of Contents.


Value-Optimal Sensor Network Design for Steady-State and Closed

Recommend Documents