Two-layer Hierarchical Predictive Control via Negotiation of Active Constraints1

Preprints, 5th IFAC Conference on Nonlinear Model Predictive Preprints, IFAC Conference on Nonlinear Model Predictive Control 5th Preprints, 5th IFAC Conference on Nonlinear Model Predictive Control September 17-20, 2015. Seville, SpainAvailable online at www.sciencedirect.com Control September 17-20, 2015. Seville, Spain September 17-20, 2015. Seville, Spain

ScienceDirect

IFAC-PapersOnLine 48-23 (2015) 404–409

Two-layer Hierarchical Predictive Control Two-layer Hierarchical Predictive Control Two-layer Hierarchical Predictive Control 1 via Negotiation of Active Constraints 1 via Negotiation of Active Constraints via Negotiation of Active Constraints 1 Seyed-Amirreza Shahidi ∗∗ Radoslav Paulen ∗∗ Seyed-Amirreza Shahidi Paulen ∗ Radoslav ∗ Seyed-Amirreza ShahidiEngell Radoslav Paulen ∗ Sebastian ∗ Sebastian Engell ∗ Sebastian Engell ∗ ∗ ∗

Process Dynamics and Operations Group, Department of Biochemical Process Dynamics and Group, Department of Process DynamicsEngineering, and Operations Operations Group, Universit¨ Department of Biochemical Biochemical and Chemical Technische at Dortmund, and Chemical Engineering, Technische Universit¨ a t Dortmund, andEmil-Figge-Str. Chemical Engineering, Universit¨ at Dortmund, 70, 44221Technische Dortmund, Germany (e-mail: Emil-Figge-Str. 70, Germany Emil-Figge-Str. 70, 44221 44221 Dortmund, Dortmund,{radoslav.paulen, Germany (e-mail: (e-mail: [email protected], [email protected], {radoslav.paulen, [email protected], {radoslav.paulen, sebastian.engell}@bci.tu-dortmund.de). sebastian.engell}@bci.tu-dortmund.de). sebastian.engell}@bci.tu-dortmund.de).

Abstract: We study a constraint negotiation approach for coordination of distributed model We a negotiation for coordination of model Abstract: Abstract: We study study The a constraint constraint negotiation approach forrepresent coordination of distributed distributed model predictive controllers. constraint negotiationapproach techniques feasible-path coordination predictive controllers. The constraint negotiation techniques represent feasible-path coordination predictive controllers. Theof constraint represent feasible-path coordination mechanisms having a goal achievingnegotiation an optimaltechniques distribution of shared resources, that represent mechanisms having a achieving an distribution shared represent mechanisms having among a goal goal of ofthe achieving an optimal optimal distribution of shared resources, resources, that represent physical couplings, subsystems. We present recent of developments of thethat coordination physical couplings, among the subsystems. We present recent developments of the coordination physical couplings, the subsystems. We the present recentalgorithm developments of the coordination scheme that enhanceamong its convergence. We apply presented on a benchmark problem scheme that convergence. apply presented algorithm on scheme that enhance enhance its convergence. We apply the the of presented algorithm on aa benchmark benchmark problem of hierarchical MPC its to illustrate the We effectiveness the proposed coordination scheme.problem of of hierarchical hierarchical MPC MPC to to illustrate illustrate the the effectiveness effectiveness of of the the proposed proposed coordination coordination scheme. scheme. © 2015, IFAC (International Federation of Automatic Control) Hosting by Elsevier Ltd. All rights reserved. Keywords: Hierarchical control, Decentralized systems, Predictive control. Keywords: Keywords: Hierarchical Hierarchical control, control, Decentralized Decentralized systems, systems, Predictive Predictive control. control. 1. INTRODUCTION cation, i.e. chemical industry and electrical power sys1. cation, i.e. industry electrical power sys1. INTRODUCTION INTRODUCTION cation, i.e. chemical chemical industry and electrical power systems (Christofides et al., 2013; and Halvgaard et al., 2012). tems (Christofides et al., 2013; Halvgaard et al., 2012). (Christofides et al., 2013; Halvgaard et al., 2012). The technological growth motivates the development of tems We can distinguish two types of methods in this framework The growth motivates the of The technological growthprocessing motivates plants. the development development of We distinguish types of in more technological and more complex This usually We can distinguish two typestheory of methods methods in this this framework framework thatcan originate fromtwo market in economics. more and more complex processing plants. This usually more more complex processingsystems plants.ofThis usually that resultsand in complex interconnected interacting originate from market theory in economics. that originate from market theory in economics. results in complex interconnected systems of interacting (1) Price-based approaches based on Lagrangian relaxresults complex 2009). interconnected systems of difficult interacting entities in (Scattolini, Such systems can be (or (1) approaches based entities (Scattolini, Such can (or (1) Price-based Price-based approaches based on on Lagrangian Lagrangian relaxrelaxation of the primal feasibility. entities (Scattolini, 2009). Such systems caninbe bea difficult difficult (or even impossible) to 2009). manage andsystems to control centralized ation of the primal feasibility. even impossible) to manage and to control in a centralized ation of the primal feasibility. (2) Constraint-negotiation, that can be regarded as a even impossible) to manage and to control in a centralized fashion so distributed and hierarchical solutions are often (2) that can as fashion and hierarchical (2) Constraint-negotiation, Constraint-negotiation, that can be be regarded regarded as aa subset of auction-based, approaches that seek fashion so distributed andmanagers. hierarchical solutions solutions are are often often favored so by distributed engineers and subset of auction-based, approaches that seek a favored by engineers and managers. subset of auction-based, approaches that a market-clearing by iteratively reallocating the seek shared favored by engineers and managers. market-clearing by iteratively reallocating the shared In case of strong coupling among the constituting systems, market-clearing by iteratively reallocating the shared resources among the subsystems. In case of strong coupling among the constituting systems, In of strong the constituting systems, resources e.g.case in the form coupling of sharedamong resources, various hierarchical resources among among the the subsystems. subsystems. e.g. in of hierarchical e.g. in the the form formmethods of shared shared resources, various hierarchical decomposition canresources, be appliedvarious in order to divide In price-based coordination, the goal is to to determine the In price-based coordination, the is determine the decomposition methods can be applied in order to divide price-based coordination, the goal goal is to to to tothat determine the equilibrium prices of all shared resources will result decomposition methods can be applied in order to divide the original (centralized) problem in subproblems whose In equilibrium prices of all shared resources that will result the original (centralized) problem in subproblems whose prices (no of all shared that will in market clearing agent can resources further increase its result profit the original (centralized) in subproblems solutions depend on someproblem (coordination) variables,whose such equilibrium in clearing agent further increase solutions on some variables, market clearing (no (no agent can further increase its its profit profit at market the equilibrium price). To can solve the optimization probsolutions depend onavailable some (coordination) (coordination) variables, such as prices, depend that are as an information to such each in at the equilibrium price). To solve the optimization probas prices, that are available as an information to each To solvemethods, the optimization problemthe byequilibrium price-basedprice). coordination a Lagrangian as prices, that are available as upper an information to layer each at subsystem and computed in the coordination lem by price-based coordination methods, aa Lagrangian subsystem and computed in the upper coordination layer lem by price-based coordination methods, Lagrangian relaxation of the coupling constraints is performed. The subsystem and computed in the upper coordination layer in order to satisfy the global system constraints and to the constraints in order to the global and relaxation ofmultipliers, the coupling coupling constraints istheperformed. performed. The Lagrangianof that represent is prices for The the in orderdesired to satisfy satisfy the performance global system system constraints and to to relaxation achieve optimal of constraints the overall system. multipliers, that the prices for achieve Lagrangian multipliers, that represent the local prices for the the different shared resources, are represent added to the economic achieve desired desired optimal optimal performance performance of of the the overall overall system. system. Lagrangian A characteristic feature of plant complexes in the process- different resources, are to local economic different shared resources, are added added to the the for localproduction economic objectiveshared functions in order to account A characteristic feature of plant complexes in the processA plant are complexes thenetworks process- objective functions order account production ingcharacteristic industry is feature that theof units coupledin by objective functions in inand order to values account for production and/or consumption theirto arefor iterated until ing industry is the units coupled by ing industrye.g. is that that theelectricity, units are are cooling coupledwater by networks networks consumption and their values are iterated until of utilities, steam, or raw and/or and/or consumption and their values are iterated until the price converges to the equilibrium price. Some recent of utilities, steam, electricity, cooling water of utilities,Ine.g. e.g. steam, electricity, cooling water oractraw raw materials. these complex networks, the units canor as the price converges to the equilibrium price. Some recent the price to the of equilibrium Some recent studies of converges the application price-basedprice. coordination can materials. In these complex networks, the units can act as materials. these complexofnetworks, units can act as of application price-based can consumers In and producers resources.the The utilization of studies studies of the the application of(2007); price-based coordination can be found in Cheng et al.of Mart´coordination ı et al. (2013a); consumers and producers of resources. The utilization of consumers and producers The utilization of be found in Cheng et al. (2007); Mart´ ı et al. (2013a); these shared resources has of to resources. be coordinated between the be found inetCheng et al. (2007); Mart´ı et al. (2013a); Stojanovski al. (2015). these shared resources has to be coordinated between the these sharedstructures resourcesinhas to be coordinated between the Stojanovski distributed order to fulfill the global demandStojanovski et et al. al. (2015). (2015). distributed structures order to the demanddistributed structures in order to fulfill fulfill the global global demand- The price-based techniques possess an inherent downside; supply constraints for in each shared resource network. The price-based techniques possess an supply constraints for each shared resource network. The possess an inherent inherent downside; sinceprice-based they relax techniques the constraints of the problem,downside; they are supply constraints for each shared resource network. they relax the constraints of the problem, they are The distributed optimization methods are especially suit- since since they relax the constraints of the problem, they are of infeasible-path type. This means that the only feasible The distributed optimization methods are suitThe distributed optimization methods are especially especially suit- of able for the application in model predictive control (MPC) infeasible-path type. This means that the only feasible of infeasible-path type. This the only feasible solution is the optimal one means that isthat achieved upon the able for the application in model predictive control (MPC) able for the application model predictive control (MPC) of interconnected units.inFor a variety of challenging con- solution is optimal one is upon the solution is the the optimal one that that is achieved achievedThis upon the convergence of the price-based coordination. might of units. For variety of conof interconnected units. For a a been variety of challenging challenging con- convergence of the price-based coordination. This might trolinterconnected tasks, MPC has already adopted by industry, convergence of the price-based coordination. This might be regarded as an unwanted feature, especially when an trol tasks, has adopted by trol tasks,itMPC MPC has already already been adopted by industry, industry, however, is mostly used tobeen control subsystems with be regarded an unwanted especially when be regardedinas asthe an real-world unwanted feature, feature, especially when an an application (industrial) environment is however, it is mostly used to control subsystems with however, it is mostly usedhierarchical to controlMPC subsystems weak couplings. Therefore schemeswith are application (industrial) environment is application in the real-world (industrial) environment is in question. in Forthe thisreal-world reason, we study constraint negotiaweak couplings. Therefore hierarchical MPC schemes are weak hierarchical schemes are in question. For this reason, we study constraint negotiaof thecouplings. interest Therefore in the main fields of MPC the MPC appliin question. For this reason, we study constraint negotiation algorithms, in this work, that represent feasible-path of the interest in the main fields of the MPC appliof the interest in the main fields of the MPC appli- tion algorithms, in work, represent feasible-path tion algorithms, in this this work, that representwas feasible-path solution approaches. Such typethat of algorithm presented 1 The authors acknowledge the support of European Commission solution approaches. Such type of algorithm was presented 1 The authors acknowledge the support of European Commission solution approaches. Such type of algorithm was presented in Riggs and Bitmead (2010) whose algorithm is adapted 1 under agreement number (DYMASOS). The grant authors acknowledge the611281 support of European Commission in Riggs and Bitmead (2010) whose algorithm is under grant agreement number 611281 (DYMASOS). in Riggs and Bitmead (2010) whose algorithm is adapted adapted under grant agreement number 611281 (DYMASOS). 2405-8963 © 2015, IFAC (International Federation of Automatic Control) Copyright 2015 IFAC 405 Hosting by Elsevier Ltd. All rights reserved. Copyright 2015 responsibility IFAC 405Control. Peer review© of International Federation of Automatic Copyright ©under 2015 IFAC 405 10.1016/j.ifacol.2015.11.312

2015 IFAC NMPC September 17-20, 2015. Seville, SpainSeyed-Amirreza Shahidi et al. / IFAC-PapersOnLine 48-23 (2015) 404–409

and extended in this work. The presented algorithm uses a principle of multi-agent negotiation over the set of active coupling constraints. The presented algorithm is applied on a benchmark problem and proves to achieve a promising coordination performance. 2. PROBLEM DEFINITION

Our goal is to solve the problem (1) while performing limited communication among the subsystems. Referring to Eq. (2c) as a global or joint constraint, it is evident that each subproblem should have an access to some (partial) information on it. The termination of the coordination procedure is reached when the optimality conditions (4a) are satisfied i.e. when the Lagrange multipliers associated with joint constraints become equal λ1 = λ2 = · · · = λn .

We study the problem of distributed MPC for n systems where the centralized control problem reads as min

X,U

s.t.

Np n  

J˜i (xi (tj ), ui (tj ))

(1a)

i=1 j=1

f˜ i (xi (tj ), ui (tj )) ≤ 0, n  g˜ i (xi (tj ), ui (tj )) ≤ 0,

(1b)

hi

(1d)

(1c)

i=1



  dxi  , xi (tj ), ui (tj ) = 0, dt tj

with i ∈ {1, . . . , n}, j ∈ {1, . . . , Np } where X and U are vectors of state and input variables and Np stands for the ˜i length of a prediction horizon. We assume that J˜i , f˜ i , g and hi are differentiable sufficiently many times.

3. PROPOSED METHODOLOGY In this section, we present the constraint negotiation procedure that is adapted from Riggs and Bitmead (2010). The procedure represents a coordination algorithm that mimics a negotiation between multiple agents. It is assumed that each system knows its own dynamics, cost function and local constraints and actively contributes to the balance (fulfillment) of the global constraints. 3.1 Coordinate Descent Search f,T T Consider a vector U f := (uf,T which satisfy 1 , . . . , un ) Eqs. (2b) and (2c).

The following steps are performed: (1) Starting with S1 , we compute ucd 1 such that

The optimization problem can be condensed into the problem of finding a vector U := (uT1 , . . . , uTn )T such that min U

s.t.

n 

ucd 1 := arg min J1 (u1 ) g 1 (u1 ) ≤ −

(2a) f 1 (u1 ) ≤ 0. (2b) (2c)

(2) For i ∈ {2, . . . , n − 1} compute ucd i

i=1

where ν ≥ 0 and λ ≥ 0 represent the vectors of Lagrange multipliers associated with local and global constraints respectively.

The optimization problem (2) can be decomposed into n local optimization problems. Let us assume that the local constraints are satisfied with strict inequalities at the optimum (i.e. optimal νij = 0, ∀i, j ∈ {1, . . . , n}) and that only the global constraints determine (constrain) the n global optimum (i.e. i=1 g i (ui ) = 0 at the optimum). We will refer to each subproblem as a local problem of subsystem Si (i ∈ {1, . . . , n}). The optimality conditions for the subproblem Si can be written as 0 =∇ui Ji (ui ) + ∇Tui g i (ui )λi , ∀i ∈ {1, . . . , n}, (4a) n  gij (ui ), ∀j ∈ {1, . . . , n}. (4b) 0 =λij

g i (ufi ),

(5b) (5c)

ucd i

such that

:= arg min Ji (ui )

s.t.

(6a)

ui

i=1

with possibly nonlinear functions Ji , f i and g i (for any i ∈ {1, . . . , n}). This formulation, equivalent to (1), is employed further in this paper to simplify the notation. The optimality conditions of the problem (2) may be formulated using the Lagrange function n  L(U , ν, λ) := Ji (ui ) + ν Ti f i (ui ) + λT g i (ui ), (3)

n  i=2

i=1

i ∈ {1, . . . , n}, f i (ui ) ≤ 0, n  g i (ui ) ≤ 0.

(5a)

u1

s.t. Ji (ui ),

405

g i (ui ) ≤ −

n 

g j (ufj ) −

(3) Compute

ucd n such that ucd n := arg min un s.t.

g j (ucd j ), (6b)

j=1

j=i+1

f i (ui ) ≤ 0.

i−1 

Jn (un )

g n (un ) ≤ −

n−1 

(6c) (7a)

g j (ucd j ),

(7b)

j=1

f n (un ) ≤ 0.

(7c)

It is evident that the obtained U cd = (ucd,T , . . . , ucd,T )T is n 1 both locally and globally feasible, i.e. it satisfies Eqs. (2c) and (2b). Besides, the following holds n n   Ji (ufi ). (8) Ji (ucd ) ≤ i i=1

i=1

If the Lagrange multipliers, found from (4a) by evaluating it at U cd , are equal to each other, then the global optimum was reached; otherwise, we proceed with the constraint negotiation procedure. 3.2 Constraint Negotiation Procedure

i=1

Here λi refers to the vector of local Lagrange multipliers of the ith subsystem.

406

The procedure of constraint negotiation is to distribute the relevant information among the subsystems which

2015 IFAC NMPC 406 September 17-20, 2015. Seville, SpainSeyed-Amirreza Shahidi et al. / IFAC-PapersOnLine 48-23 (2015) 404–409

leads to the optimal satisfaction of the global (coupling) constraints (2c). The skeleton of the presented algorithm is taken from Riggs and Bitmead (2010) and is further adapted and extended in this work. The methodology can be shown to be similar to sensitivity-based coordination by Scheu and Marquardt (2011). Also a similar proof, as in Scheu and Marquardt (2011), can be presented that shows the convergence of the scheme to the global optimum under the assumptions of convexity of functions Ji , f i , g i and linear independence constraint qualification. The procedure can be conceived as a circular information flow, in which subsystem Si exchanges information with subsystem Si+1 and so forth. The basic idea of the procedure of constraint negotiation is to pick a set of active constraints and perturb a member of this set in the direction that leads to the best decrease of the global cost. Let us consider the case when U cd satisfies (2c). Then a subsystem i ∈ {1, . . . , n} is picked to perform a perturbation on ui such that ui = ucd i + pi . The influence of this perturbation on the other subsystems is then examined in a sequential manner where the subsystems j (∀j �= i) are communicated the perturbation on ui as well as the response to this perturbation by the rest of the systems. We assume that the resulting vector U cn , after this perturbation, does not disturb the set of active constraints. We further require U cn to satisfy local and global constraints and to attain a decrease of the global objective function. This vector is then passed to the coordinate descent algorithm to evaluate whether the global optimality conditions, particularly the equality of (local) Lagrange multipliers, are satisfied. In fact the perturbations can be iterated with different i ∈ {1, . . . , n} and the system Si∗ whose perturbation (negotiation proposal) results in the biggest decrease of the overall objective function can be selected as the driver of the negotiation. The proposed distribution of shared resources, i.e. the one which results in a biggest decrease of the global cost, is accepted and implemented. If the Lagrange multipliers for the new point do not match, λ1 �= λ2 �= · · · = � λn , the new point can be used for further improvement, i.e. the negotiation procedure is repeated, so the iterative procedure is established and continued until the point such that the resulting solution from constraint negotiation results in equality of Lagrange multipliers. To perform the constraint negotiation procedure, the following steps are executed: (1) Compute the Lagrange multipliers corresponding to global constraints to each subsystem. (2) Pick a subsystem Si , i ∈ {1, . . . , n}. (3) Solve the following optimization problem T � � � 1 min  λTj ∇g i (ui )�ucd + ∇Ji (ui )�ucd  pi pi i i n−1 

j�=i

� 1 + pTi ∇2ui Ji (ui )�ucd pi i 2 cd s.t. f i (ui + pi ) ≤ 0.

(9a)

(9b)

cd (4) Compute ∆Ji = Ji (ucd i + diag(α)pi ) − Ji (ui ) for cd some 0 < αi ≤ 1 and also g i (ui + diag(α)pj ).

407

(5) Pick a subsystem Sj , j �= i. Initialize the index sets Ip := j and Ie := {1, . . . , n} \ {i, j}. (6) Solve the following optimization problem min Jj (ucd j + pj ) pj

s.t 0 ≥f j (ucd j + pj ), 0

≥g i (ucd i

(10)

+ diag(α)pi ) +

�

g k (ucd k )

k∈Ie

+

�

g k (ucd k

+ pk ).

k∈Ip cd and compute ∆Jj = Jj (ucd j + pj ) − Jj (uj ). (7) Repeat step 6 with j taken as an element of Ie (when j is picked perform Ie := Ie ⊖ j and Ip := Ip ⊕ j) until Ie = ∅. �n (8) (Optional) Go to step 2 to find i∗ such that i=1 ∆Ji is minimal. (9) Input U cn to the coordinate descent algorithm.

Hence that step 3 performs an SQP-like calculation with the second-order approximation of the objective function and a first-order optimality inclusion via a term that consists of constraints gradients and an average of the Lagrange multipliers from other subsystems. This optimization results in the enforcement of a satisfaction of global optimality conditions exactly when the objective function is quadratic, the joint constraints are linear and Lagrange multipliers are equal among the subsystems. A vector α is regarded in this scheme as a proposal execution rate or a negotiation parameter. In the hierarchical control notation, α stands for a degree of freedom of the negotiation coordinator and represents, in fact, a step size of the negotiation similarly to the step size applied by the coordinator in price-based schemes. Riggs and Bitmead (2010) propose to adjust this parameter to counteract the possibly �n non-decreasing value of the objective function (i.e. i=1 ∆Ji ≥ 0) or to suppress a possible violation of local or global feasibility. In this study we propose a novel approach for the calculation of α. Calculation of a Negotiation Coordination Parameter. The original constraint negotiation algorithm (Riggs and Bitmead, 2010) involves a multi-step process of evaluation and adaptation of the parameter α to arrive at an acceptable performance in terms of feasibility and decrease of the objective value. The choice of the parameter α might not be straightforward. Riggs and Bitmead (2010) propose a use of scalar negotiation parameter. However, for multivariate decision vector ui , it might turn out that, in order to achieve a faster convergence of the coordination procedure, taking the negotiation parameter as a vector, i.e. allowing more degrees of freedom to negotiation parameter, would be necessary. This would apply in cases when the initial point from which the negotiation starts is very far from a global optimum or when a substantial diversity appears in the magnitudes of the constraints, decisions and local objective functions that is not properly regularized. In practice one often encounters these situations. We suggest to formulate the choice of α as an optimization problem. As mentioned above, the coordination procedure is expected to enforce the optimality conditions (4a). The objective function for the identification of the best value

2015 IFAC NMPC September 17-20, 2015. Seville, SpainSeyed-Amirreza Shahidi et al. / IFAC-PapersOnLine 48-23 (2015) 404–409

407

of the negotiation parameter α, with respect to the rate of negotiation convergence, can be described by two goals: • Global (primal) optimality. To this end we consider the following optimization problem � cd ∗ (u ∗ + diag(α)pi∗ ) + min J Jj (ucd i i j + pj ), α,p j

j�=i∗

1≤j≤n, j�=i∗

(11a)

s.t. 0 0 0

≥f i∗ (ucd i∗ + diag(α)pi∗ ), ≥f j (ucd j + pj ), 1 ≤ j ≤ n, ≥g i∗ (ucd ∗ + diag(α)pi∗ ) �i + g j (ucd j + pj ),

(11b) j �= i∗ ,

(11c) (11d)

j�=i∗ ∗

where the index i can be selected using the procedure described in Section 3.2 with some a priori chosen value of αi ∈ (0, 1). We remark that if the functions J (·), f (·), or g(·) are nonlinear, an approximation of the constraint or objective functions can take place in order to simplify the problem (11) and to express it in the form of sequential linear (or quadratic) programming. Hence in this approach the subsystems communicate the sensitivities of the constraints and objectives w.r.t. shared resources. Note that similar schemes based on sharing the derivative information have been applied in price-based coordination in order to enhance its convergence (Mart´ı et al., 2013b). • Optimal resource allocation (market clearing). The aim here is to enforce the equality of the local Lagrange multipliers of the subsystems. Practically, we can formulate an optimization � � problem similar T � � to (11) with the objective �λ − 1 nλ � that minimizes 1 the difference of Lagrange multipliers w.r.t. the mean value of the local Lagrange multipliers among the subsystems. It is clear that the second proposal involves a bigger number of degrees of freedom and constraints. On the other hand, if the values of the underlying objectives are not available to the coordinator or the problem is highly nonlinear, an approximation is required (as discussed above), the proposed approach will provide a reasonable trade-off between computational complexity and convergence of constraint negotiation. 4. CASE STUDY - FOUR-TANK BENCHMARK We study a four-tank benchmark (Johansson, 2000) whose scheme is shown in Figure 1. The control task is to track the time-varying references in order to drive the system in the different steady-states. The cross-sectional areas of all tanks are S = 0.06 m2 . The considered values of the discharge areas of tank outlets are a = [1.31 × 10−4 , 1.51 × 10−4 , 9.27 × 10−5 , 8.82 × 10−5 ] m2 and the distributions ratios of three-way valves are [γ1 , γ2 ] = [0.3, 0.4]. Maximum allowed heights of the tanks are hmax = [1.36, 1.36, 1.3, 1.3] m for each tank. Pumps of the plant can provide maximum flow rates of [q1max , q2max ] = [3.26, 4] m3/h. It is assumed that the liquid levels in tanks are initially set to be at h0 = [0.65, 0.66, 0.65, 0.66] m that stands for a steady state together which the flow rates q 0 = [1.63, 2.00] m3/h. 408

Fig. 1. Schematic diagram of a quadruple-tank system. 4.1 Dynamic Model The dynamics of the plant can be described as a following set of differential equations a1 � dh1 a3 � γ1 =− (12a) 2gh1 + 2gh3 + q1 , dt S S S a2 � a4 � γ2 dh2 =− (12b) 2gh2 + 2gh4 + q2 , dt S S S a3 � dh3 (1 − γ2 ) =− q2 , (12c) 2gh3 + dt S S a4 � dh4 (1 − γ1 ) =− q1 . 2gh4 + (12d) dt S S

As discussed in Alvarado et al. (2011), this model describes the behavior of the plant reasonably well when the heights in tanks are over 0.2 m. Therefore, the physical lower bounds of the tank heights have been set to 0.2 m. A continuous linear model can be formulated via a linearization of (12) around the initial steady state (h0 , q 0 ) as dx = Ac x + Bc u, dt y = Cc x,

(13a) (13b)

where x = (x1 , x2 , x3 , x4 ) indicate the states variables and u = (u1 , u2 ) stand for control inputs, i.e. deviations of the tank levels and flow rates from the steady state (h0 , q 0 ) of the system. The linearization approach is selected to fulfill the assumptions of the coordination scheme. The matrices Ac , Bc and Cc are defined as   −T1 0 T3 0  0 −T2 0 T4  Ac :=  , (14a) 0 0 −T3 0  0 0 0 −T4   1 − γ1 γ1 0 0  S  BcT :=  S γ 1 − γ (14b) , Cc := I, 2 2 0 0 S S � in which Ti = ai /(S 2h0i /g) with i ∈ {1, 2, 3, 4}. The discrete-time system is formulated using zero-order hold with sampling time of 90 s.

2015 IFAC NMPC 408 September 17-20, 2015. Seville, SpainSeyed-Amirreza Shahidi et al. / IFAC-PapersOnLine 48-23 (2015) 404–409

h1 [m]

2 1

h2 [m]

0 0 2

q1 [m3 /h]

2000

3000

4000

5000

6000

7000

8000

9000

1000

2000

3000

4000

5000

6000

7000

8000

9000

1000

2000

3000

4000

5000

6000

7000

8000

9000

1000

2000

3000

4000

5000

6000

7000

8000

9000

1 0 0 4 2 0 0 4

q2 [m3 /h]

1000

2 0 0

time [s] Fig. 2. Comparison of the simulation results obtained for centralized (blue solid line), decentralized (green dotted dashed line) and hierarchical (black dashed line) control schemes. The references are plotted using red solid line. 4.2 Control Objective As in Alvarado et al. (2011), a set of references for heights in the first and the second tank is considered. The overall duration of the control task is 9000 s which is divided into three consecutive intervals of 3000 s with set-points of set hset 1 = {0.3, 0.5, 0.9} m and h2 = {0.3, 0.75, 0.75} m. The set steady-state flow-rates q1 and q2set that correspond to these set-points are calculated according to the nonlinear model of the system. The centralized model predictive controller minimizes the cost J :=

N p−1 j=1

A distributed MPC scheme is realized such that the system is decomposed into two subsystems, one concerned with the heights in the first and the third tank and the other one concerned with the heights in the second and the fourth tank. Therefore, the states in these subsystems are: (16) xs1 := (x1 , x3 )T , xs2 := (x2 , x4 )T . Upon observing (13), it is clear that this decomposition leaves no state couplings to be present. On the other hand, both subsystems are given the freedom to adjust (optimize) the flow-rates q1 and q2 which gives rise to the notion of shared resources, i.e. input couplings, as the global (overall) system optimum requires qi1 (j) = qi2 (j), i = 1, 2, j = 1, . . . , Np − 1, (17) where the superscript assigns the decision to the respective subsystem.

 2  2 set h1 (j) − hset    1 (j) Qx + h2 (j) − h2 (j) Qx

 2  2 + q1 (j) − q1set (j)Qu + q2 (j) − q2set (j)Qu , (15)

where N p = 10 is the same as control horizon. The weighting matrices with appropriate dimensions are Qx = I and Qu = 0.01 I. Note that the stability of MPC is out of the scope of this work, thus, we do not impose terminalcost term into the objective which generally results in a non-decomposable objective. 409

The cost function for subsystem k can be stated as Jk :=

N p−1 j=1

 2  2 set hk (j) − hset    k (j) Qx + qk (j) − qk (j) Qu .

(18)

The overall measure of optimality, Joverall , is given as a weighted sum (the same weighting as in (15) is used) of the tracking errors during the whole control task.

2015 IFAC NMPC September 17-20, 2015. Seville, SpainSeyed-Amirreza Shahidi et al. / IFAC-PapersOnLine 48-23 (2015) 404–409

for coordination of the negotiation process adopting, for example, the techniques of real-time optimization.

Table 1. Evaluation of the performance index for different control schemes. Control scheme

Joverall

Centralized Decentralized (system 1) Decentralized (system 2) Hierarchical (nDOF = 1)

0.5527 1.9240 1.7881 0.5579

409

REFERENCES

4.3 Results We use MATLAB and YALMIP toolbox (L¨ ofberg, 2004) to perform simulations of the centralized, decentralized and the proposed hierarchical control schemes. The costs that the different control schemes reach are shown in Table 1. We can clearly see that the decentralized solution, when the system 1 (or system 2) optimizes its own objective function results in a big deterioration of the global cost. This behavior evidently justifies the usage of a hierarchical control scheme. We implement the constraint negotiation algorithm with different degrees of freedom assigned to coordinating layer. This is reflected in the number of elements, nDOF , of the vector α. Hence, if nDOF = Np × nu , where nu is the number of control inputs, the centralized solution is replicated as the coordinator possesses the same degrees of freedom as the local subsystems. In the performed simulation studies, we found out that using nDOF = 1 is sufficient in this case study and that the constraint negotiation scheme in one iteration converges to a point where the equality of Lagrange multipliers is satisfied with a reasonable threshold. This scheme, that performs one iteration of constraint negotiation and that uses nDOF = 1, achieves 0.9% loss of optimality compared to the centralized solution. Figure 2 shows and compares the results obtained for the centralized, decentralized (minimization of the objective (18) for k = 1) and the hierarchical control schemes. We can observe only negligible differences between centralized and hierarchical control schemes. 5. CONCLUSIONS In this study we presented an adaptation of the constraint negotiation algorithm taken from the literature to the management of physically coupled systems of systems that uses a principle of multi-agent negotiation over the set of active coupling constraints. We also presented some extensions of the algorithm that enhance its convergence rate. The presented algorithm was applied on a benchmark problem and proved to achieve a promising coordination performance. Moreover, the presented algorithm stands for a feasiblepath alternative to price-based (Lagrangian relaxation) methods that reach a primal feasibility only upon the convergence, thus it is well-suited for the robust implementation in conservative industrial environment. In our future work, we will further test the algorithm on more complex case studies taken from the operational practice of the processing industry. An interesting area for further research is a usage of iterative building of the surrogate subsystems’ models 410

Alvarado, I., Limon, D., de la Pe˜ na, D.M., Maestre, J., Ridao, M., Scheu, H., Marquardt, W., Negenborn, R., Schutter, B.D., Valencia, F., and Espinosa, J. (2011). A comparative analysis of distributed MPC techniques applied to the hd-mpc four-tank benchmark. Journal of Process Control, 21(5), 800 – 815. doi: http://dx.doi.org/10.1016/j.jprocont.2011.03.003. Cheng, R., Forbes, J., and Yip, W. (2007). Price-driven coordination method for solving plant-wide mpc problems. Journal of Process Control, 17(5), 429–438. Christofides, P.D., Scattolini, R., Mu˜ noz de la Pe˜ na, D., and Liu, J. (2013). Distributed model predictive control: A tutorial review and future research directions. Computers and Chemical Engineering, 51, 21–41. doi: 10.1016/j.compchemeng.2012.05.011. Halvgaard, R., Poulsen, N.K., Madsen, H., and Jorgensen, J.B. (2012). Economic Model Predictive Control for building climate control in a Smart Grid. In 2012 IEEE PES Innovative Smart Grid Technologies (ISGT), 1–6. IEEE. doi:10.1109/ISGT.2012.6175631. Johansson, K.H. (2000). The quadruple-tank process: a multivariable laboratory process with an adjustable zero. IEEE Transactions on Control Systems Technology, 8(3), 456–465. doi:10.1109/87.845876. L¨ofberg, J. (2004). Yalmip : A toolbox for modeling and optimization in MATLAB. In Proceedings of the CACSD Conference. Taipei, Taiwan. URL http://users.isy.liu.se/johanl/yalmip. Mart´ı, R., Sarabia, D., Navia, D., and de Prada, C. (2013a). A method to coordinate decentralized NMPC controllers in oxygen distribution networks. Computers & Chemical Engineering, 59, 122–137. doi: 10.1016/j.compchemeng.2013.05.023. Mart´ı, R., Sarabia, D., Navia, D., and de Prada, C. (2013b). Coordination of Distributed Model Predictive Controllers Using Price-Driven Coordination and Sensitivity Analysis. In H. Michael (ed.), 10th IFAC International Symposium on Dynamics and Control of Process Systems, 215–220. doi:10.3182/20131218-3-IN2045.00035. Riggs, D.J. and Bitmead, R.R. (2010). Negotiation of coupled constraints in coordinated vehicles. In 49th IEEE Conference on Decision and Control (CDC), 2010, 479– 484. doi:10.1109/CDC.2010.5717221. Scattolini, R. (2009). Architectures for distributed and hierarchical Model Predictive Control A review. Journal of Process Control, 19(5), 723–731. doi: 10.1016/j.jprocont.2009.02.003. Scheu, H. and Marquardt, W. (2011). Sensitivity-based coordination in distributed model predictive control. Journal of Process Control, 21(5), 715 – 728. doi: http://dx.doi.org/10.1016/j.jprocont.2011.01.013. Stojanovski, G., Maxeiner, L.S., Kr¨amer, S., and Engell, S. (2015). Real-time shared resource allocation by price coordination in an integrated petrochemical site. In 2015 European Control Conference (ECC). IEEE.