Distributed Model Predictive Control of Nonlinear Systems Based on

Aug 15, 2016 - Here, a nonlinear plant is considered, which is operated by a decentralized control system. The existing system ignores the interaction...
0 downloads 0 Views 767KB Size
Subscriber access provided by Northern Illinois University

Article

Distributed Model Predictive Control of Nonlinear Systems Based on Price-Driven Coordination Bardia Hassanzadeh, Hallas Pakravesh, Jinfeng Liu, and J. Fraser Forbes Ind. Eng. Chem. Res., Just Accepted Manuscript • DOI: 10.1021/acs.iecr.6b01862 • Publication Date (Web): 15 Aug 2016 Downloaded from http://pubs.acs.org on August 20, 2016

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 43

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

Distributed Model Predictive Control of Nonlinear Systems Based on Price-Driven Coordination Bardia Hassanzadeh,∗ Hallas Pakravesh, Jinfeng Liu, and J. Fraser Forbes Department of Chemical and Materials Engineering, University of Alberta, Edmonton, Alberta, Canada, T6G 1H9 E-mail: [email protected]

Phone: +1(780)964-5141. Fax: +1(780)492-2881

Abstract Here a nonlinear plant is considered, which is operated by a decentralized control system. The existing system ignores the interactions between sub-systems, which often results in uncaptured plant-wide performance. The focus of this paper is on the design of a distributed model predictive control (DMPC) network using successively linearized internal models. In this method, all existing interactions between the sub-systems should be considered in order to enhance the performance of the current decentralized DMPC scheme. A coordination layer is added to the existing network, while minor modifications are applied to the local MPC controllers, to achieve the performance and stability of a hypothetical centralized MPC for the entire plant. In this work, an interior-point algorithm is proposed to coordinate a DMPC network via the price-driven coordination approach. In addition, the convergence of the algorithm is shown, and the necessary conditions to ensure the closed-loop stability of the system are provided for the situation when the algorithm is terminated prematurely prior to convergence. ∗

To whom correspondence should be addressed

1 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

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

1

Introduction

Interconnected processes appear in a lot of chemical and petrochemical plants; thus, efficient networked control methodologies must be used, in order to increase their efficiency and profit. Traditionally, decentralized and centralized schemes are the two primary control frameworks that are used by the industries. A decentralized scheme is easy to implement, but it may result in the degradation of plant-wide control performance or even the loss of closed-loop stability because interactions between sub-systems are completely neglected (see, for example, 1–3 and references therein). A centralized scheme can obtain the best plant-wide control performance; however, it may become too complex to be implemented as the size of the control problem grows. Moreover, a centralized scheme for a large-scale process may require substantial maintenance effort and is not favorable from the fault tolerance point of view. 4,5 Therefore, in recent years, many methods have been developed for distributed model predictive control (DMPC) to efficiently control large-scale systems. In a DMPC scheme, distributed local model predictive controllers (MPC) communicate and exchange information with each other or with a coordinator to coordinate their control actions. 6–10 As a result, DMPC is able to achieve an improved plant-wide control performance while preserving the flexibility of a decentralized scheme such as low maintenance costs, unique fault tolerance properties, and resilience in design and operation. One classification divides DMPC schemes into two categories: non-coordinated DMPC and coordinated DMPC. In the context of non-coordinated distributed MPC, local controllers include interaction models and exchange information to achieve an improved control performance. Within this scope, a classification can be made based on communication design of network: fullyconnected networks where all controllers exchange information, and partially-connected networks where only a limited number of controllers send and receive information. Some important recent studies in this area includes cooperative DMPC of linear systems, 7,11 Lyapunovbased sequential and iterative DMPC of nonlinear systems, 12,13 robust DMPC of linear 2 ACS Paragon Plus Environment

Page 2 of 43

Page 3 of 43

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

systems, 14 DMPC of linear systems coupled through the inputs 15 and dissipativity-based DMPC of linear systems. 16 In contrast, inside a Coordinated DMPC (CDMPC) scheme, distributed MPCs communicate with a coordinator to achieve improved performance. Different methods have been developed for CDMPC including the price-driven approach (i.e., 17–20 ), the primal decomposition approach, 21,22 and the prediction-driven approach (i.e., 23 ). The main differences between these coordination methods appear in the ways that the interactions between subsystems are addressed. Pertaining to nonlinear systems, there exists a limited amount of research dedicated to the design of coordinated distributed control for nonlinear systems. The coordination of non-linear networks governed by open-loop stable DAE systems are also investigated by Scheu et al. 24 and Martí et al.; 20 however, the closed-loop stability of coordinated control networks has not yet been addressed. Another challenge is a situation during which the coordination scheme terminates prior to convergence. This could occur when there is not enough time to complete the required calculations. In this work, the design of CDMPC for nonlinear systems based on the price-driven coordination approach is considered. In the price-driven coordination, the coordinator calculates “prices” for local MPCs based on the exchanged information in order to coordinate the plant-wide control actions. It has been proven that the price-driven CDMPC can converge to the corresponding optimal plant-wide performance, when the coordinator and distributed controllers have enough time to obtain control actions at each sampling time. 17–19 In the proposed approach, an interior-point based model predictive controller is designed for each sub-system and a price-driven scheme coordinates the control actions of distributed subsystem controllers to achieve the centralized control performance and the stability of the entire system. The proposed scheme assumes that a hypothetical centralized MPC which successively linearizes a nonlinear system model, at every sampling time, can be designed to stabilize the closed-loop nonlinear system. A key idea behind this scheme relates to the hierarchical dual decomposition of a hypothetical centralized MPC controller with a minor

3 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

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

modification applied to the existing network of decentralized MPC controllers. The modification to the decentralized MPC controllers is equivalent to add interaction equations to the objective function of local controllers by assigning a price vector as a penalization coefficient. In fact, a coordinator in the upper-level optimization problem receives information from the distributed sub-system MPCs which belong to the lower-level optimization problem. During coordination cycles between the coordinator and the plant-wide control system, the price value is sent to the distributed controllers to coordinate their control actions in an iterative fashion and to achieve the control performance and stability of a centralized MPC. The necessary conditions that ensure the stability of the proposed CDMPC are derived. Note the closed-loop stability is considered when a system is stopped prematurely prior to convergence. To summarize, the main contributions of this paper are as follows: • This paper proposes an on-line interior-point CDMPC algorithm for the open-loop stable nonlinear systems via the successive linearization of a nonlinear model. Contrary to the CDMPC approaches that utilize active-set optimization algorithms, 17,19 the interior-point CDMPC obviates the need for identifying the correct active-set when the algorithm is stopped prematurely prior to convergence. • A Lyapunov-based stability criterion is derived for the proposed CDMPC approach under the premise that the equivalent centralized MPC would be able to stabilize the closed-loop system, specifically , when the algorithm is terminated prior to convergence. • The convergence of the proposed CDMPC algorithm is demonstrated via choosing a proper merit function. The control performance of the proposed CDMPC algorithm is illustrated using a benzene alkylation process.

4 ACS Paragon Plus Environment

Page 4 of 43

Page 5 of 43

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

1.1

System Description

In this paper, a class of nonlinear systems composed of m interconnected open-loop stable sub-systems is considered. Each of the sub-systems can be described by the following discretetime model:

xi (k + 1) = Fi (x(k), u(k))

(1)

where: i = 1, . . . , m, xi (t) ∈ Rnxi denotes the vector of state variables of sub-system i, x = [xT1 · · · xTi · · · xTm ]T ∈ Rnx denotes the state of the entire nonlinear system, and u = [uT1 · · · uTi · · · uTm ]T ∈ Rnu is the vector of control inputs with ui denoting the vector of control inputs associated with sub-system i. The states of the m sub-systems, xi (i = 1, . . . , m), are assumed to be sampled synchronously at time instants tk = kT , where: k = 0, 1, . . . and T is the sampling period. Note that in the remainder of this paper, k is used to denote tk in the discrete-time model. Furthermore, it is assumed that ui ∈ Ci ⊂ Rnui with Ci being a nonempty convex set defined as follows:

ineq Ci , {ui ∈ Rnu,i | Gineq i,u ui ≤ gi,u }

(2)

ineq for i = 1, . . . , m. Here, Gineq i,u ui ≤ gi,u is a linear inequality. Note that, no bounds are con-

sidered over state variables of the sub-systems. This is mainly due to the methodology used in this work, which is based on successive linearization of the nonlinear system. Considering bounds over states could limit the applicability of this study, especially when the CDMPC algorithm is stopped prematurely prior to convergence. In this situation, the feasibility of the linearized problem does not necessarily reflect the feasibility of the nonlinear system. Therefore, in the remainder of this paper, no inequality constraints are considered for the state variables. Here, it is assumed that the sub-systems are interconnected through states and manip-

5 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

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 43

ulated variables, i.e. the manipulated input variable of one sub-system might affect other sub-systems as an internal variable in the models. The dynamics of the entire system can be described as follows:

x(k + 1) = F (x(k), u(k)))

(3)

where: F = [F1T · · · FiT · · · FmT ]T , and F is a twice differentiable Lipschitz vector function. Using a first-order Taylor expansion of (3) around the current measured state (xk ) and control input (uk ) yields: xˆ(k + 1) = Ψ(k)ˆ x(k) + Γ(k)ˆ u(k)

(4)

where xˆ and uˆ are deviation variables. Additionally, the values of the matrices Ψ, Γ are determined at each sampling time as: ∂F |x=xk ,u=uk ∂x ∂F Γ(k) := |x=xk ,u=uk ∂u

Ψ(k) :=

(5a) (5b)

In this work, it is assumed that all elements of Jacobian matrices Ψ(k) and Γ(k) are uniformly bounded. At each control interval, the prediction model is constructed from Ψ(k) and Γ(k), as a linear time-invariant state-space model along the prediction horizon. Thus, the control actions are calculated using the successively linearized model (5), and these actions are applied to the nonlinear plant (3).

1.2

Structure of the Centralized MPC

In this section, a monolithic (centralized) MPC is designed based on the successively linearized model (4) that satisfies the input constraints, for all x inside a compact set S containing the origin. The stable steady-state operating point of the system is transformed to the origin in order to have a unified definition of S. Note that, in the successive linearization 6 ACS Paragon Plus Environment

Page 7 of 43

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

approach, the predicted trajectories inside the controllers are constructed based on timeinvariant linear models at the beginning of each sampling period. This could limit the scope of this study to nonlinear systems with low to moderate degrees of nonlinearity; so that, the dynamic behavior of the system can be captured, with acceptable accuracy, within the chosen prediction horizon. Remark 1. 25 A generic type finite-time MPC, based on the successively linearized model (5) at time instant k, can be formulated as:  1 1 T T min J = X(k) QX(k) + U (k) RU (k) + xˆ(k + Hp |k)T P xˆ(k + Hp |k) X,U 2 T s.t. xˆ(k + l + 1|k) = Ψ(k)ˆ x(k + l|k) + Γ(k)ˆ u(k + l|k),

(6a) (6b)

uˆ(k + l|k) = uˆ(k + Hu |k), for Hu ≤ l ≤ Hp − 1

(6c)

uˆ(k + l|k) ∈ C

(6d)

for l = 0, · · · , Hp ; where: Hp is the prediction horizon and Hu is the control horizon; and, xˆ and uˆ are states and manipulated input variables inside the controller, respectively; and S C = i=1,··· ,m Ci . Additionally, X(k) = [ˆ x(k + 1|k)T , · · · , xˆ(k + Hp |k)T ]T is the vector of the predicted state trajectory; U (k) = [ˆ u(k|k)T , · · · , uˆ(k + Hp − 1|k)T ]T is the vector of the calculated manipulated variable moves; Q is a positive definite block-diagonal weighting matrix for the states (i.e., Q = diag{Qii }); R is a positive definite block-diagonal weighting matrix for the manipulated variables of the overall system (i.e., R = diag{Rii }); and P is a positive definite block-diagonal weighting matrix for the terminal cost of the overall system (i.e., P = diag{Pii }). In other words, the optimization problem (6) is meant to ensure: 25 xˆ(k + Hp |k) ∈ {x|xT P x ≤ δ}

for some δ > 0; where: 0 ∈ {x|xT P x ≤ δ}. 7 ACS Paragon Plus Environment

(7)

Industrial & Engineering Chemistry Research

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 43

In the remainder of this paper, it is assumed that the discrete-time system (3) is openloop stable; thus, the generic MPC formulation (6) can be simplified to (8) with a large prediction horizon Hp in order to ensure the closed-loop stability of (3):

min J = X,U

 1 X(k)T QX(k) + U (k)T RU (k) 2

(8a)

s.t. xˆ(k + l + 1|k) = Ψ(k)ˆ x(k + l|k) + Γ(k)ˆ u(k + l|k),

(8b)

uˆ(k + l|k) = uˆ(k + Hu |k), for Hu ≤ l ≤ Hp − 1

(8c)

uˆ(k + l|k) ∈ C

(8d)

In this work, an interior-point approach is used to design the centralized MPC controller. Using the interior-point method, unlike the active-set approach, 26 it is not required to identify the correct set of active constraints. The finite-time centralized MPC formulation based on the successively linearized model (5), at time instant k, is:

min JC = X,U

 1 X(k)T QX(k) + U (k)T RU (k) + µΩU (U (k)) 2

s.t. xˆ(k + l + 1|k) = Ψ(k)ˆ x(k + l|k) + Γ(k)ˆ u(k + l|k),

(9a) (9b)

uˆ(k + l|k) = uˆ(k + Hu |k), for Hu ≤ l ≤ Hp − 1

(9c)

for l = 0, · · · , Hp . In (9), µ is the barrier weighting parameter and ΩU is a generic barrier function, which is defined as follows:

ΩU (U (k)) = −

m nuX i ×Hp X i=1

ineq (j) − Gineq ω gi,U i,U (j, :) × Ui



(9d)

j=1

where: ω is a generic indicator function 27 that it’s domain satisfies the set of inequality constraints defined in (2) strictly. The choice of the indicator function ω must guarantee 8 ACS Paragon Plus Environment

Page 9 of 43

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

that the barrier function Ω(y) remains convex, increasing in ‘y’, and differentiable. 27 A key assumption of this work is that a centralized MPC (9) can be designed for the closed-loop system (3) which is locally exponentially stable around the origin. Then, the following assumption is an immediate conclusion from the converse Lyapunov theorem for discrete-time systems 28 : Assumption 1. [Lyapunov stability certificates 28 ] Given T , a family of continuous Lyapunov functions VT , and some positive real numbers (c1 , c2 , c3 , c4 ), there exists some strictly positive real numbers (c5 ≥ 1, σ, σ0 ≤ σ/c5 ) such that for x(k) ∈ Bσ and x(0) ∈ Bσ0 the following holds:

c1 ||x||2 ≤ VT (x) ≤ c2 ||x||2

(10a)

VT (F (x, u)) − VT (x) ≤ −T c3 ||x||2

(10b)

where: for a strictly positive constant σ, set Bσ is defined as: Bσ = {x ∈ Rn |0 ≤ ||x|| ≤ σ}. In addition, for all xa , xb ∈ Bσ0 : VT (xa ) − VT (xb ) ≤ c4 ||xa − xb ||(||xa || + ||xb ||)

(10c)

where: VT (x) is a Lyapunov function for the closed-loop system under the centralized MPC scheme (8). Then, the closed-loop system (3) under the control of the centralized MPC (9) is locally exponentially stable around the origin with a continuous Lyapunov function VT (x). In the following section, the idea of price-driven CDMPC approach for linear systems is extended to the control of nonlinear systems based on successive linearization of (3) at each sampling time. The centralized MPC will provide performance benchmarks for evaluation of the CDMPC for nonlinear systems.

9 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

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

2

Page 10 of 43

The Proposed Coordination Algorithm

In this section, a price-driven CDMPC algorithm is proposed for nonlinear systems, based on successive linearization. The main idea is to form a multilevel network consisting of a coordination level interacting with an existing network of MPC controllers, in order to achieve a higher performance compared to the fully decentralized scheme. In this method, distributed controllers are modified to account for interactions inside their objective functions using a price vector. The information flow between the coordinator and the local controllers provides each MPC with the optimal price vector and as a result the optimum plant-wide trajectory. This structure is depicted in Figure 1, in which: p denotes the price vector; Θi denotes the local interaction matrix of sub-system i, defined in (16); and Zi is the optimal trajectory of the local MPCs with respect to (w.r.t.) the current price value, defined in (21d).

Coordinator 1   p

Z

1

 i   p

Z

i

 m   p

CDMPC Controller 1

CDMPC Controller i

CDMPC Controller m

sub-system 1

sub-system i

sub-system m

Zm

Interactions

Interactions

Process

Figure 1: Architecture and information flow of the proposed CDMPC

2.1

Price-Driven Coordination Formulation

For each sub-system, the MPC is formulated based on the successively linearized sub-system model. Specifically, for sub-system i, i = 1, . . . , m, the prediction model used in the design 10 ACS Paragon Plus Environment

Page 11 of 43

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

of the sub-system MPC at time instant k takes the following form:

xˆi (k + l + 1|k) =Ψii xˆi (k + l|k) + Γii uˆi (k + l|k) +

X (1 − β)Ψij xˆj (k|k) + vˆi (k + l|k) (11a) j6=i

xˆi (k|k) =xi (k)

(11b)

with:

β=

   0 l = 0

(11c)

  1 l = 1, · · · , Hp − 1 where: Ψii and Γii denote coefficient matrices for the linearized system description of subsystem i, namely the ith block diagonal elements of the plant-wide matrices Ψ and Γ given in (5); and vˆi is defined as the interacting or linking variable that contains unknown interaction information between different sub-systems. Note that the term vˆi in (11) characterizes the interaction of sub-system ‘i’ with other sub-systems. The interaction error for sub-system ‘i’ is defined as:

ei (k + l|k) , vˆi (k + l|k) −

m X

 βΨij xˆj (k + l|k) + Γij uˆj (k + l|k))

(12)

j6=i

for l = 0, · · · , Hp − 1. The coordinator finds a price for the sub-system i, such that the interaction term vˆi determined by the price, ensures that ei (k + l|k) = 0 and the CDMPC approaches the performance of the centralized MPC optimal trajectory. The overall interaction error over the prediction horizon can be described as follows:

E(k|k) =

m X

Ei (k|k)

i=1

11 ACS Paragon Plus Environment

(13a)

Industrial & Engineering Chemistry Research

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 12 of 43

where:      Ei =    

ei (k|k)



  ei (k + 1|k)    ..  .   ei (k + Hp − 1|k)

(13b)

According to (12), Ei (k|k) can be written as: 



−Γ1i uˆi (k|k)   h i     − Ψ x ˆ (k + 1|k) + Γ u ˆ (k + 1|k) 1i i 1i i     .   ..     h i    − Ψ1i xˆi (k + Hp − 1|k) + Γ1i uˆi (k + Hu − 1|k)      .   . .         v ˆ i (k|k)       vˆi (k + 1|k)   Ei (k|k) =   .   . .         vˆi (k + Hp − 1|k)     . ..         −Γmi uˆi (k|k)     h i     − Ψ x ˆ (k + 1|k) + Γ u ˆ (k + 1|k) mi i mi i     ..   .    h i − Ψmi xˆi (k + Hp − 1|k) + Γmi uˆi (k + Hu − 1|k)   Xi (k)    = Θi (k)   Ui (k)    Vi (k)

12 ACS Paragon Plus Environment

(14a)

(14b)

Page 13 of 43

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

where: Xi (k) , [ˆ xi (k + 1|k)T , · · · , xˆi (k + Hp |k)T ]T is the vector containing predictions for state variables; Ui (k) , [ˆ ui (k+1|k)T , · · · , uˆi (k+Hp |k)T ]T is the vector containing predictions for manipulated input variables; Vi (k) , [ˆ vi (k + 1|k)T , · · · , vˆi (k + Hp |k)T ]T is the vector containing predictions for linking variables; and Θi (k) is the coefficient matrix for the linking constraints defined as: T

 T T T Θi (k) = θ1,i , · · · , θm,i , · · · , θi,i

where: θj,i =

     0

Hp nxi ×Hp nxi , 0Hp nxi ×Hu nui , I

    θΨj,i , θΓj,i , 0Hp nx

j ×Hp nxi



(15a)

for j = i (15b)



for j 6= i

with ‘0’ being zero matrices of appropriate dimensions, I being an identity matrix of size Hp nxi × Hp nxi , as well as θΨj,i and θΓj,i being Hp nxj × Hp nxi and Hp nxj × Hu nui matrices, respectively, and defined as: 

θΨj,i

0 ···  0  −Ψji 0 ···  = ... ...  0   0 0 −Ψji

   0 ··· 0  0 −Γji   ..  ..  0 . −Γ .  0 ji    , θ =    ..  Γj,i  .. ... ... 0  .  .     0 0 ··· 0 −Γji

(15c)

Equivalently, matrix Θi can be written in terms of an augmented matrix of coefficients for Xi , Ui , and Vi as follows: 

θΨ1,i   ..  .   Θi =  0Hp nxi ×Hp nxi  ..  .   θΨm,i

θΓ1,i .. . 0Hp nxi ×Hp nui .. . θΓm,i

0 .. .



     IHp nxi ×Hp nxi    ..  .   0

13 ACS Paragon Plus Environment

(16a)

Industrial & Engineering Chemistry Research

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

=[ΘXi , ΘUi , ΘVi ]

Page 14 of 43

(16b)

and the overall interaction error E, defined in (13), can be written as:   Xi (k) m X    E(k) = Θi (k)   Ui (k)    i=1 Vi (k)

(17)

Based on the above, the optimization problem (9) can be re-written as follows:

min JP =

X,U,V

m X

JPi

(18a)

i=1

s.t. xˆi (k + l + 1|k) = Ψii xˆi (k + l|k) + Γii uˆi (k + l|k) +

X

(1 − β)Ψij xˆj (k|k) + vˆi (k + l|k)

j6=i

for l = 0, . . . , Hp − 1   Xi (k) m X    Θi  U (k)  i =0   i=1 Vi (k)

(18b)

(18c)

uˆi (k + l|k) = uˆi (k + Hu |k) for Hu ≤ l ≤ Hp − 1

(18d)

where:  1 T T JPi = (Xi (k)) Qii (Xi (k)) + Ui (k) Rii Ui (k) + µi ΩUi (Ui (k)) 2 nui ×Hp X  ineq ΩUi (Ui (k)) = − ω gi,U (j) − Gineq i,U (j, :) × Ui

(18e) (18f)

j=1

Optimization problem (18) is the main core of the design of the CDMPC system. Note that the cost function and constraints are separable in terms of sub-systems. Relaxing

14 ACS Paragon Plus Environment

Page 15 of 43

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

constraint (18c) into JPi yields:

min JD =

X,U,V

m X

JDi

(19a)

i=1

s.t. xˆi (k + l + 1|k) = Ψii xˆi (k + l|k) + Γii uˆi (k + l|k) +

X

(1 − β)Ψij xˆj (k|k) + vˆi (k + l|k)

j6=i

for l = 0, . . . , Hp − 1

(19b)

uˆi (k + l|k) = uˆi (Hu |k) for Hu ≤ l ≤ Hp − 1

(19c)

where:

JDi

  Xi (k)     1  (Xi (k))T Qii (Xi (k)) + Ui (k)T Rii Ui (k) + µi ΩUi (Ui (k)) + pT Θi  = U (k)  i  2   Vi (k)

(19d)

In (19), ‘p’ is a price vector provided by the coordinator and is used to drive the distributed system to the plant-wide optimal operation, see Section 2.2 for more details. Optimization problem (19) is separable w.r.t. the sub-systems i = 1, . . . , m. In the proposed CDMPC, the formulation of local MPC ‘i’ is considered as:

min JDi

(20a)

Xi ,Ui ,Vi

s.t. xˆi (k + l + 1|k) = Ψii xˆi (k + l|k) + Γii uˆi (k + l|k) +

X

(1 − β)Ψij xˆj (k|k) + vˆi (k + l|k)

j6=i

for l = 0, . . . , Hp − 1

(20b)

uˆi (k + l|k) = uˆi (Hu |k) for Hu ≤ l ≤ Hp − 1

(20c)

15 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

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 43

where:

JDi

2.2

  Xi (k)    1 (Xi (k))T Qii (Xi (k)) + Ui (k)T Rii Ui (k) + µi ΩUi (Ui (k)) + pT Θi  = Ui (k)    2   Vi (k)

(20d)

Coordinator Design

In order to simplify the description of the coordinator design, the overall formulation (19) is re-written in a more compact form:

min JD =

m X

Z

J¯Di

(21a)

i=1

eq s.t. Geq i Zi (k) = gi

(21b)

where:  1 J¯Di = (Zi (k))T Υi (Zi (k)) + µi ΩZi (Zi (k)) + pT Θi Zi (k) 2

(21c)

Zi (k) represents the vector of decision variables, namely: Zi (k) = [Xi (k)T , Ui (k)T , Vi (k)T ]T

(21d)

ΩZi (Zi (k)) is the corresponding barrier function for sub-system i, defined as:

ΩZi (Zi (k)) = ΩXi (Xi (k)) + ΩUi (Ui (k)) + ΩVi (Zi (k)) = ΩUi (Ui (k))

(21e) (21f)

with: ΩXi (Xi (k)) = 0, ΩVi (Vi (k)) = 0, and Geq i defined as: eq eq eq Geq i = [Gi,Xi , Gi,Ui , Gi,Vi ]

16 ACS Paragon Plus Environment

(21g)

Page 17 of 43

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



Geq i,Xi

Geq i,Ui



0 ··· 0   Inxi   −Ψii In  · · · 0 x   i = , ..  .. ..  0  . . .     0 0 −Ψii Inxi   0 ··· 0  −Γii  ..  .  0 −Γii . . .    = . , .. ..  .. . . 0      0 ··· 0 −Γii

Geq i,Vi = − IHp nxi ×Hp nxi   P (1 − β)Ψ x ˆ (k) ij j  j6=i      0   gieq =   .   . .     0

(21h)

(21i)

(21j)

(21k)

In (21), the price vector p can be interpreted as the Lagrange multipliers associated with interaction equality constraints. This vector is determined by the coordinator at each coordination cycle, based on the information provided by the local MPC controllers. In order to find a plant-wide solution to the distributed system (21), the primal-dual optimization problem (22) is solved:

max JD

(22a)

p

min JD = Z

m X

J¯Di

(22b)

i=1

eq s.t. Geq i Zi (k) = gi

(22c)

where: Z = [X T , U T , V T ]T . Denote L the overall Lagrangian of the local sub-problems (21),

17 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

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 18 of 43

and define Li as: 1 eq Li (Zi , ϑi , µi , p) = Zi (k)T Υi Zi (k) + pT Θi Zi (k) + ϑTi Geq i Zi (k) − gi ) + µi ΩZi (Zi (k)) (23) 2 with: ϑi being the Lagrange multiplier of the equality constraints for sub-system ‘i’. In Barrier methods, it is common to express the first-order optimality conditions via the corresponding KKT conditions and solve the resulting set of equations with a Newton’s method variant. 27 Correspondingly, the modified KKT conditions of the primal problem in (22), for sub-system ‘i’, are: nzi ×Hp

∇Zi Li (Zi , ϑi , λi , p) = Υi Zi (k) +

ΘTi p

+

T Geq i ϑi

+

X

T λi (j) Gineq (j, :) =0 i

(24a)

j=1 eq R1,i (Zi , p) = Geq i Zi (k) − gi = 0   R2,i (Zi , p) = diagj=1,··· ,nz ×Hp (λi (j)) φi + µi ρ = 0

(24b) (24c)

i

where: "

 dω −1  dω −1  −1 dω φi = ,··· , ,··· , dζi (1) dζi (j) dζi (nzi × Hp ) ζi (j) = Gineq (j, :)Zi (k) − giineq (j) i

#T (24d) (24e)

ρ = [1, · · · , 1]T1×nzi Hp

(24f)

and λi (j) can be interpreted as the Lagrange multiplier associated with the ‘j’th inequality constraint of sub-system ‘i’. The use of λi (j) simplifies the notation and can always be eliminated with the relation defined in (24c). Inside local MPCs, it is ensured that the parameter µi is monotonically decreasing to zero at each iteration of the interior-point algorithm, i.e. µi is manually decreased in every iteration by the following routine: (si +1)

µi

(si )

= α i × µi 18

ACS Paragon Plus Environment

(25)

Page 19 of 43

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

where: αi ∈ [0.05, 0.1], and si is the iteration number of the interior-point algorithm in the sub-system ‘i’. 27 As µi → 0 the KKT conditions (24) converge to the optimal solution of (21), for a fixed value of the price vector ‘p’ received from the coordinator. On the other hand, the dual problem is an unconstrained optimization problem w.r.t the price vector ‘p’. Denote the optimal value of Z at the current value of p as Z ∗ (k), then according to (22), the Lagrangian of the unconstrained problem can be written as:

¯ Z ∗ (k)) = −JD (p, Z ∗ (k)) L(p,

(26a)

based on the convex nature of the lower level problem, at Zi∗ the equality constraint (24b) eq ∗ is always satisfied, i.e. Geq i Zi (k) − gi = 0 and µi = 0. Therefore, penalizing the equality

constraint with ϑi and adding −µi ΩZi (Zi∗ (k)) to the right-hand-side of (26a) would result in:

¯ Z ∗ (k)) = − JD (p, Z ∗ (p)) − ϑT Geq Z ∗ (k) − g eq ) − µi ΩZ (Z ∗ (k)) L(p, i i i i i i

(26b)

in other words, at Z ∗ (k):

¯ Z ∗ ) = −L(p, Z ∗ ) L(p,

(26c)

The first-order optimality condition for the dual problem can be written as: nu ×Hp



eq T

¯ Z ) = − ΥZ (k) + Θ p + G ∇p L(p, ∗



T

ϑ+

X j=1



T

eq T

but ΥZ (k) + Θ p + G

ϑ+

Pmu ×Hp j=1

 dφ(j) T  dZ − ΘZ ∗ = 0 λ(j) dZ dp

(27a)

T  dφ(j) = 0, since ∇Z L = 0. This is equivalent λ(j) dZ

to say:

¯ Z ∗ ) = − ΘZ ∗ = 0 ∇p L(p,

19 ACS Paragon Plus Environment

(27b)

Industrial & Engineering Chemistry Research

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 43

in other words, the KKT condition of the dual problem is:

ΘZ ∗ = 0.

(27c)

Hence, an efficient price-updating method, such as Newton’s method, can be used to achieve the plant-wide performance for the coordination problem defined in (21):

p(k + 1) = p(k) − αH−1 (k)J (k),

(28)

where J (k) is the Jacobian and H(k) is the Hessian matrix of L, and α is the step size in Newton’s method that is determined using line search method given in. 29 The Jacobian is calculated by: ∗

J (k) = ΘZ =

m X

Θi Zi∗ (k)

(29)

i

Then the Jacobian is equal to the overall interaction error vector ‘E(k)’, defined in (14a). Accordingly, the Hessian matrix H(k) can be calculated as: m

dE(k) X dJ (k) = = Θi (k)∇p Zi∗ (k) H(k) = dp(k) dp(k) i

(30)

Thus, the coordinator calculates Θi (k) at each sampling time and is provided with the sensitivity matrix of the decision variables Zi (k) with respect to changes in the price vector information. The problem of determining the sensitivities can be written as:     T ∇p Zi (k) −Θi      = 0  Γi  ∇ ϑ (k) p i         ∇p λi (k) 0

20 ACS Paragon Plus Environment

(31a)

Page 21 of 43

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

where: 

T Geq i

Υi   Γi =   

diagj=1,··· ,nz

i

Geq i

0

  dφi (λ (j)) i ×Hp dZi

0

T Gineq i



   0    diagj=1,··· ,nz ×Hp (φi (j))

(31b)

i

The system of equations (31) is derived from differentiating the optimality conditions of the primal CDMPC problem for sub-system ‘i’. The required sensitivity information can be obtained by differentiating the optimality conditions (24). Then, the following system of equations is derived: T

T

ineq ∇2Zi ,Zi Li dZi + ∇2Zi ,p Li dp + Geq dλi = 0 i dϑi + Gi

(32a)

∇Zi R1,i dZi + ∇p R1,i dp = 0

(32b)

∇Zi R2,i dZi + ∇p R2,i dp + ∇λi R2,i dλi = 0

(32c)

Therefore, the sensitivity matrix ∇p Zi (k) is obtained by solving the linear system of equations (31), assuming that the matrix Γi is full rank. The implementation of the CDMPC algorithm at sampling time ‘k can be summarized as illustrated in Algorithm 1.

As can be seen in the proposed CDMPC scheme (22), a nonnegative merit function M can be defined as:

M = L∗ − L(p)

(33)

where: L∗ is the optimum value of the dual problem objective function L. Based on (33), the value of M can be zero iff p = p∗ . Now, consider the changes of M w.r.t. the CDMPC

21 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

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 43

Algorithm 1: the proposed successive linearization CDMPC algorithm Input: coordination cycle counter s = 0, price vector p = 0, stopping criterion  > 0, Zi0 (k) ; Local MPCs: Calculate Ψ(k) and Γ(k) using (5); Coordinator: Send the value of ps to the local MPCs; Coordinator: Calculate the overall interaction error E(k) using (17); while ||E(k)|| >  do Local MPCs: Solve problem (19) using the interior-point method and update Z¯is (k); Local MPCs: Calculate sensitivity matrices J (k) and H(k) using (29) and (30); Local MPCs: Send J (k) and H(k) to the coordinator; Coordinator: Calculate the price vector ps using (28); Coordinator: Send ps to the local MPCs; Coordinator: Calculate the overall interaction error E(k) using (17) s = s + 1; end Local MPCs: Apply the receding horizon action to the plant; Output: Z¯is (k), and ps ; cycle ‘s’:  dM T dp ∇s M = dp ds  dL(p) T dp = − dp ds

(34) (35)

which reduces to:

∇s M = −J T

dp ds

(36)

Furthermore, the rate of change of ‘p’ w.r.t. the iteration number ‘s’ is directly proportional to the gradient ascent direction J in the price update scheme (28). Thus, the rate of change of the merit function M w.r.t. the iteration number ‘s’ is:

∇s M = −δJ T J , δ ≥ 0

(37)

which means ∇s M ≤ 0. Therefore, the proposed CDMPC scheme (22) converges to local

22 ACS Paragon Plus Environment

Page 23 of 43

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

optimal plant-wide solutions (18).

3

Stability Analysis

In this section, the closed-loop stability of the proposed CDMPC scheme is discussed. The analysis is based on the premise that the coordination algorithm has been stopped prematurely prior to convergence is achieved, i.e.:

0 ≤ ||E(k)|| ≤ 

(38)

where:  is the stopping criterion for the proposed CDMPC algorithm. Proposition 1. As  → 0, the optimal solution to the CDMPC problem (22) converges to the optimal solution of the corresponding centralized MPC trajectory (9). Proof. The plant-wide formulation of the centralized MPC problem (9), stated in (18), can be re-written in terms of Zi (k) as:

min JP = Z

m X

J¯Pi

eq s.t. Geq i Zi (k) = gi m X

(39a)

i=1

(39b)

Θi Zi (k) = 0

(39c)

1  J¯Pi = (Zi (k))T Υi (Zi (k)) + µi ΩZi (Zi (k)) 2

(39d)

i=1

where:

for i = 1, · · · , m. Denote LP the overall Lagrangian of the plant-wide MPC problem (39)

23 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

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 24 of 43

defined as:

LP (Zi , ϑΘ , ϑi , µi ) =

m  X 1

m  X  Zi (k)T Υi Zi (k) + ϑTΘ Θi Zi (k)

2 m  X

i=1

+

i=1

 eq ϑTi Geq Z (k) − g ) + µ Ω (Z (k)) i i Zi i i i

(40)

i=1

where: ϑΘ is the Lagrange multiplier of the interaction constraints of all the sub-systems. Correspondingly, the modified KKT conditions 27 of the plant-wide problem in (39) can be written as: m X

nzi ×Hp

Υi Zi (k) + ΘTi ϑΘ +

T Geq i ϑi

+

X

! λi (j)

T (j, :) Gineq i

eq Geq i Zi (k) − gi = 0   diagj=1,··· ,nz ×Hp (λi (j)) φi + µi ρ = 0 i

m X

=0

(41a)

j=1

i=1

Θi Zi (k) = 0

(41b) (41c) (41d)

i=1

for i = 1, · · · , m. Comparing the KKT conditions of the plant-wide MPC problem stated in (41) with the KKT conditions of the proposed CDMPC scheme stated in (24) shows that as P ∗ E(k) = m i=1 Θi Zi (k) → 0, the optimal price vector ‘p ’ of the CDMPC scheme converges to the optimal Lagrange multiplier ϑΘ of the plant-wide MPC problem. Furthermore, when the stopping criterion of the proposed CDMPC is chosen to be zero, i.e.  = 0, the L2 -norm P of the overall interaction E(k) becomes zero, i.e. ||E(k)|| = || m i=1 Θi Zi (k)|| = 0, based on (38). Therefore, as  → 0, the optimal solution to the CDMPC problem (22) converges to the optimal solution of the corresponding centralized MPC trajectory (9). In the following proposition, it is shown that if the proposed CDMPC algorithm is stopped prior to convergence the state variables and control actions of the closed-loop system remain bounded. Proposition 2. Consider the CDMPC scheme (22) is stopped prematurely at time instant 24 ACS Paragon Plus Environment

Page 25 of 43

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

‘k’ prior to convergence, i.e. ||E(k)|| ≤ . Then, the decision variables Xi (k), Ui (k), and Vi (k) remain bounded, for i = 1, · · · , m. Proof. Consider the KKT conditions (24), the following set of equations can be obtained for sub-system ‘i’: T

s Qii Xis + ΘTi,Xi ps + Geq i,Xi ϑi = 0

(42a)

nui ×Hp

Rii Uis

+

ΘTi ps

+

T s Geq i,Ui ϑi

T λi (j) Gineq =0 i,Ui (j, :)

(42b)

psi − ϑsi = 0

(42c)

eq eq s s Vis − Geq i,Xi Xi − Gi,Ui Ui + gi = 0

(42d)

X

+

j=1

T

T

where: psi is the i-th portion of the price vector ‘p’ and ‘s’ denotes the iteration number of the proposed CDMPC algorithm. Consequently, the overall set of equations for sub-systems i = 1, · · · , m can be written as: T

s QX s + ΘTX ps + Geq X ϑ = 0

(43a)

nu ×Hp s

T s

RU + Θ p +

T s Geq U ϑ

T λ(j) Gineq =0 U (j)

(43b)

ps − ϑs = 0

(43c)

eq s s eq V s − Geq =0 X X − GU U + g

(43d)

+

X j=1

T

T

From (43), the price vector ‘p’ and the Lagrange multiplier associated with the equality constraints ‘ϑ’ are equal. Therefore, ‘ϑ’ can be eliminated and the set of equations (42) can be rewritten as: T

eq T s T Xis = −Q−1 ii (Θi,Xi + Gi,X Θi,Ii )p

Uis

=

−Rii−1



(44a) nui ×Hp

(ΘTi,Ui

+

T T s Geq i,Ui Θi,Ii )p

+

X

T  λi (j) Gineq (j, :) i,Ui

j=1

25 ACS Paragon Plus Environment

(44b)

Industrial & Engineering Chemistry Research

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

T

Page 26 of 43

T

eq eq s s Vis = Geq i,Xi Xi + Gi,Ui Ui − gi .

(44c)

For a fixed value of the price vector ‘p’, the local MPC problem defined in (19) is convex and has a unique solution. Furthermore, according to the rule stated in (25), µi → 0 as the interior-point algorithm converges to the optimal solution.Then, λi → 0 as µi → 0 based on the KKT conditions (24), for i = 1, · · · , m. Therefore, at the optimal point of the local MPCs (19), the derivation (44) can be further simplified to: T

eq T T s Xis = −Q−1 ii (Θi,Xi + Gi,X Θi,Ii )p  s T T Uis = −Rii−1 ΘTi,Ui + Geq i,Ui Θi,Ii p

(45a) (45b)

T

T

eq eq s s Vis = Geq i,Xi Xi + Gi,Ui Ui − gi .

(45c)

P m s Since ||∆E (k)|| = i=1 Θi (k)Zi (k) ≤ , the following result is obtained from (44): s

m X s s s T ||∆E (k)|| = Θi (k)[Xi (k), Ui (k), Vi (k)] ≤  s

i=1

= Π × ps + Ξ ≤ 

(46a)

where:

Π=

m X

"   eq T T eq T T −1 T T Θi,Xi − Q−1 ii (Θi,Xi + Gi,X Θi,Ii ) , Θi,Ui − Rii (Θi,Ui + Gi,Ui Θi,Ii ) ,

i=1

# T Geq i,Xi

Ξ=

m h X



T Q−1 ii (Θi,Xi

0, 0, −gieq

iT

+

eq T T Gi,X Θi,Ii )



+

T Geq i,Ui



Rii−1 (ΘTi,Ui

+

T T Geq i,Ui Θi,Ii )

.



(46b) (46c)

i=1

Now, using the triangular inequality ||a|| − ||b|| ≤ ||a + b||, (46a) yields:

||Π|| × ||ps || − ||Ξ|| ≤ ||Π × ps + Ξ|| ≤ 

26 ACS Paragon Plus Environment

(47a)

Page 27 of 43

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

Rearranging (47a), an upper-bound can be derived for the norm of the price vector at the coordination cycle ‘s’:

||ps || ≤ dps ()

(47b)

where:

dps () =

 + ||Ξ|| ||Π||

(47c)

Therefore, according to (44) and (47b), it can be shown that the decision variables Xi (k), Ui (k), and Vi (k) are bounded, i.e.: ||Xis || ≤ dXis ()

(48a)

||Uis || ≤ dUis ()

(48b)

||Vis || ≤ dVis ()

(48c)

where: T

eq T T s dXis () = ||Q−1 ii (Θi,Xi + Gi,X Θi,Ii )|| × dp () T

(48d)

T s dUis () = ||Rii−1 (ΘTi,Ui + Geq i,Ui Θi,Ii )|| × dp ()

(48e)

eq eq s s dVis () = ||Geq i,Xi || × ||Xi || + ||Gi,Ui || × ||Ui || + ||gi ||

(48f)

Equivalently, predetermined upper-bounds exist on the norm of the decision variables Xi (k) and Ui (k) at a chosen stopping criterion ‘’ for the CDMPC scheme (22), i.e.:  + ||Ξ|| ||Π|| T  + ||Ξ|| T + Geq i,Ui Θi,Ii )|| × ||Π|| T

eq T T ||Xis || ≤ ||Q−1 ii (Θi,Xi + Gi,X Θi,Ii )|| ×

(48g)

||Uis || ≤ ||Rii−1 (ΘTi,Ui

(48h)

27 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

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 43

in which,  = 0 corresponds to the optimal solution of centralized MPC problem (9). From proposition 2, the state variables Xi (k) = [ˆ xi (k +1|k)T , · · · , xˆi (k +Hp |k)T ]T and the control actions Ui (k) = [ˆ ui (k|k)T , · · · , uˆi (k + Hp − 1|k)T ]T remain bounded if the CDMPC scheme (22) is stopped at ||E(k)|| ≤ . In the next proposition, it is shown that the closedloop nonlinear system xi (k+1) = Fi (x(k), u(k)) remains bounded when the stopping criterion is chosen to be  > 0, for i = 1. · · · , m. Proposition 3. Consider the CDMPC scheme (22) and the centralized MPC scheme (9) are applied to system (3) at sampling time ‘k’. Then, if the CDMPC scheme (22) is stopped prematurely at ||E(k)|| ≤ , the closed-loop trajectory remains bounded at sampling time ‘k + 1’. Proof. Since F is a Lipschitz function, it can be shown that:

F (xD (k), uD (k)) − F (xC (k), uC (k)) ≤ Lx (xD (k) − xC (k)) + Lu (uD (k) − uC (k)) (49) where: Lx and Lu are Lipschitz constants for x and u, respectively; (xC , uC ) denotes the trajectory and control action of the system under the centralized MPC scheme (9), and (xD , uD ) denotes the trajectory and control action of the system under the proposed CDMPC scheme (22). But, xD (k) = xC (k) and (uD (k) − uC (k)) = (ˆ uD (k) − uˆC (k)), which reduces (49) to:

F (xD (k), uD (k)) − F (xC (k), uC (k)) ≤ Lu uˆD (k) − uˆC (k)

(50)

Using ||a − b|| ≤ ||a|| + ||b||, (50) can be stated as:   F (xD (k), uD (k)) − F (xC (k), uC (k)) ≤ Lu uˆD (k) + uˆC (k)

28 ACS Paragon Plus Environment

(51)

Page 29 of 43

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

Based on (48), ||ˆ u(k)|| ≤ ||U s (k)|| ≤

qP m

i=1

d2Uis (); then, (51) can be further simplified to: v v ! u m u m uX uX t d2 s () + t d2 s (0)

F (xD (k), uD (k)) − F (xC (k), uC (k)) ≤ Lu

Ui

Ui

i=1

(52)

i=1

Rearranging (52), ||F (xD (k), uD (k))|| is bounded by: v v ! u m u m uX uX F (xD (k), uD (k)) ≤ F (xC (k), uC (k)) + Lu t d2 s () + t d2 s (0) Ui

Ui

i=1

(53a)

i=1

or:

xD (k + 1) ≤ xC (k + 1) + Lu

v v ! u m u m X uX u t d2 s () + t d2 s (0) Ui

i=1

Ui

(53b)

i=1

In the following theorem, it is shown that the closed-loop system (1) is locally exponentially stable around the origin under proposed CDMPC control scheme (22), if it is stopped prematurely prior to convergence. This result is under the premise that a centralized MPC control scheme (9) can be designed for system (1), which is exponentially stable around the origin with the Lyapunov function VT (x). Theorem 1. Consider the discrete-time system (3) is controlled by the CDMPC scheme (22). Assume that all the local sub-systems (1) are open-loop stable, and a centralized MPC control scheme (9) can be designed to stabilize the entire system (3) around the origin. Then, the closed-loop system under the proposed CDMPC control scheme (22) is locally exponentially stable around the origin. Proof. It is required to verify that VT (xD (k + 1)) − VT (xD (k)) ≤ 0 at any two consecutive sampling times ‘k’ and ‘k + 1’ in order to prove the CDMPC control scheme (22) is locally h i exponentially stable around the origin. Adding and subtracting VT (xC (k +1)) − VT (xC (k)) 29 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

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 30 of 43

to the left-hand side of the above inequality results in: h i VT (xD (k + 1)) − VT (xD (k)) = VT (xD (k + 1)) − VT (xD (k)) + VT (xC (k + 1)) − VT (xC (k)) h i − VT (xC (k + 1)) − VT (xC (k)) (54a)

Suppose xC (k) = xD (k), then VT (xC (k)) = VT (xD (k)). As a result, (54a) is simplified to: h i VT (xD (k + 1)) − VT (xD (k)) = VT (xC (k + 1)) − VT (xC (k)) h i + VT (xD (k + 1)) − VT (xC (k + 1))

(54b)

Based on Assumption 1, for all xC , xD ∈ Bσ0 − Bσ1 and 0 < σ1 < σ0 , the right-hand-side of (54b) can be rewritten as: 2 VT (xD (k + 1)) − VT (xD (k)) ≤ − T c3 xC (k) + c4 xD (k + 1) − xC (k + 1) ×   xD (k + 1) + xC (k + 1) (54c)

Then, using ||a|| − ||b|| ≤ ||a − b||, (54c) results in:  2 2 2  VT (xD (k + 1)) − VT (xD (k)) ≤ − T c3 xC (k) + c4 xD (k + 1) − xC (k + 1) (54d) or: 2 2 2 VT (xD (k + 1)) − VT (xD (k)) ≤ − T c3 xC (k) − c4 xC (k + 1) + c4 xD (k + 1) (54e)

which, according to Proposition 3, can be further simplified to:

30 ACS Paragon Plus Environment

Page 31 of 43

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

2 2 VT (xD (k + 1)) − VT (xD (k)) ≤ − T c3 xC (k) − c4 xC (k + 1) v v !!2 u m u m uX uX + c4 xC (k + 1) + Lu t d2Uis () + t d2Uis (0) i=1

i=1

(54f) Rearranging (54f) results in: 2 VT (xD (k + 1)) − VT (xD (k)) ≤ − T c3 xC (k) v v ! u m u m uX uX + 2c4 Lu xC (k + 1) t d2 s () + t d2 s (0) Ui

Ui

i=1

+ c4 L2u

i=1

v v !2 u m u m uX uX t d2 s () + t d2 s (0) Ui

(54g)

Ui

i=1

i=1

Then, (54g) can be restated as:

VT (xD (k + 1)) − VT (xD (k)) ≤ − T c3 σ12 + 2c4 Lu σ0

v v ! u m u m uX uX t d2 s () + t d2 s (0) Ui

Ui

i=1

+ c4 L2u

i=1

v v !2 u m u m uX uX t d2 s () + t d2 s (0) Ui

i=1

Ui

(54h)

i=1

that is equivalent to:

VT (xD (k + 1)) − VT (xD (k)) ≤ 0

31 ACS Paragon Plus Environment

(54i)

Industrial & Engineering Chemistry Research

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 32 of 43

when c3 and c4 are chosen such that: !2

! c3 ≥ c4

2σ0 Lu

qP m

i=1

d2Uis () +

qP m

i=1

+ L2u

d2Uis (0)

qP m

i=1

d2Uis () +

qP m

i=1

d2Uis (0)

T σ12 (54j)

Therefore, there exists a set of values for c3 and c4 which guarantees that VT (xD (k + 1)) ≤ VT (xD (k)). Since VT is decreasing along xD at any two consecutive sampling times ‘k’ and ‘k + 1’, the Lyapunov function of the centralized MPC scheme (9) can be considered as the Lyapunov function of the CDMPC scheme (22), and the proposed CDMPC scheme (22) is locally exponentially stable around the origin.

4

Alkylation Process Case Study

The proposed CDMPC is applied to a benzene alkylation example. The process consists of four continuous stirred tank reactors (CSTRs) and a flash tank separator, as shown in Figure 2. Pure benzene is fed via stream F1 and pure ethylene is fed via streams F2 , F4 and F6 . Two catalytic reactions take place in CSTR-1, CSTR-2 and CSTR-3. In reaction 1, benzene (A) reacts with ethylene (B) and produces the desired product ethylbenzene (C); in reaction 2, ethylbenzene reacts with ethylene to form a byproduct, 1, 3-diethylbenzene (D). In the flash tank separator, most of benzene is separated, flows overhead and recycled back to the plant. A portion of the recycle stream Fr1 is fed to CSTR-4 and the other portion Fr2 is fed back to CSTR-1 together with an additional feed stream F10 , which contains 1, 3diethylbenzene. Moreover, reaction 2 and reaction 3, in which 1, 3-diethylbenzene reacts with benzene to produce ethylbenzene, take place in CSTR-4. Finally, it is assumed that throughout the process, all the materials are in liquid phase and their molar volumes are constant. The definition of the process variables is given in Tables 1. A more detailed version of the process model is presented by Christofides et al.. 10 32 ACS Paragon Plus Environment

Page 33 of 43

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

F2,B

F4,B

F6,B

F1,A

F3

Q1

CSTR-1

F4

Q2

CSTR-2

Q3

CSTR-3 F7

Fr2 Fr

FP

Fr1

F10,D

F9

F8

CSTR-4

Q5

Q4

Separator

Figure 2: Process flow diagram of alkylation of benzene Table 1: Alkylation process parameters Process Parameter CA1 ,CB1 ,CC1 ,CD1 CA2 ,CB2 ,CC2 ,CD2 CA3 ,CB3 ,CC3 ,CD3 CA4 ,CB4 ,CC4 ,CD4 CA5 ,CB5 ,CC5 ,CD5 CAr ,CBr ,CCr ,CDr T1 ,T2 ,T3 ,T4 ,T5 F3 ,F5 ,F7 ,F8 ,F9 F1 ,F2 ,F4 ,F6 ,F10 Fr , Fr1 , Fr2 Q1 ,Q2 ,Q3 ,Q4 ,Q5 CpA ,CpA ,CpA ,CpA

Description Concentration of species in CSTR-1 Concentration of species in CSTR-2 Concentration of species in CSTR-3 Concentration of species in Separator Concentration of species in CSTR-4 Concentration of species in Fr , Fr1 , Fr2 Temperature in each vessel Effluent flow rates from each vessel Feed flow rates to each vessel Recycle flow rates External heat/coolant input to each vessel Heat capacity at liquid phase

In the design of the proposed CDMPC, the process is divided into two sub-systems. The first sub-system includes CSTR-1 and CSTR-2. The second sub-system includes CSTR-3, the separator, and CSTR-4. Consequently, the following sub-system decomposition of state and input variables is defined as:    X1 = [CA1 , CB1 , CC1 , CD1 , T1 , CA2 , CB2 , CC2 , CD2 , T2 ]   U1 = [Q1 , Q2 , F4 ]

33 ACS Paragon Plus Environment

(55a)

Industrial & Engineering Chemistry Research

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 34 of 43

and:    X2 = [CA3 , CB3 , CC3 , CD3 , T3 , CA4 , CB4 , CC4 , CD4 , T4 , CA5 , CB5 , CC5 , CD5 , T5 ]   U2 =

(55b)

[Q3 , Q4 , Q5 , F6 ]

The prediction horizon Hp and the control horizon Hu are taken to be 15. In addition, the sampling period is chosen to be T = 30[sec]. The decision variables of optimization are scaled to be between 0.1 and 1 in order to enhance numerical computation. The weighting matrices Q and R are I25 and 10−4 × I7 , respectively; here, ‘I’ represents the identity matrix of appropriate dimension. In this case study, the stopping criterion for the proposed CDMPC scheme is chosen to be  = 10−6 , i.e. the CDMPC scheme (22) is stopped prematurely when ||E(k)|| ≤ 10−6 . In addition, a log-barrier function is used in the formulation of the local MPC controllers, i.e.: nui ×Hp

ΩUi (Ui (k)) = −

X

ineq ln gi,U (j) − Gineq i,U (j, :) × Ui



(56)

j=1

Furthermore, the closed-loop performance of the proposed CDMPC scheme (22) is compared to a decentralized MPC scheme, in which the interactions between sub-systems are completely ignored. This scheme can be defined as:

min JDC = X,U

m X

JDCi

(57a)

i=1

s.t. xˆi (k + l + 1|k) = Ψii xˆi (k + l|k) + Γii uˆi (k + l|k) for l = 0, . . . , Hp − 1

(57b)

uˆi (k + l|k) = uˆi (Hu |k) for Hu ≤ l ≤ Hp − 1

34 ACS Paragon Plus Environment

(57c)

Page 35 of 43

3.5 JC (Centralized)

3

Overal Objective Function

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

JD (CDMPC) JDC (Decentralized)

2.5

2

1.5

1

0.5

0 0

500

time [s]

1000

1500

Figure 3: Comparison of the overall objective function for the centralized, coordinated and decentralized control schemes, based on scaled decision variables where:

JDCi =

 1 (Xi (k))T Qii (Xi (k)) + Ui (k)T Rii Ui (k) + ΩUi (U (k)) 2

(57d)

In Figure 3, the trends of the overall objective functions of the CDMPC design JD and the decentralized design JDC are compared to the centralized objective function JC . It should be noted that all the values in Figure 3 are calculated based on the scaled versions of the decision variables. As can be seen, CDMPC can efficiently track the optimal centralized MPC trajectory (9). The result obtained by decentralized MPC shows that the decentralized MPC is not able to achieve the performance of the centralized MPC (9), i.e. the mean square error between the decentralized MPC and the centralized MPC scaled objective functions is at least two orders of magnitude less than that of between the CDMPC and the centralized MPC, over the whole simulation time.

35 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

490

500

480 T [K]

470

2

T1 [K]

480

460 460

450 0

500

1000

1500

2000

2500

440 0

3000

500

1000

time [s]

1500

2000

2500

3000

2000

2500

3000

time [s]

500

490

490 480

480

470

5

T [K]

T3 [K]

460

470

450 440 430 0

500

1000

1500

2000

2500

460 0

3000

500

1000

time [s]

1500 time [s]

Figure 4: Temperature profiles in CSTR-1, CSTR-2, CSTR-3, and CSTR-4; solid lines (—): centralized scheme, circles (◦): CDMPC scheme, dotted lines (· · · ): steady-state values

Q1 [J/s]

3.5

x 10

6

4.1 4.7 5.3 0

Q2[J/s]

4

x 10

500

1000

1500

2000

2500

3000

500

1000

1500

2000

2500

3000

500

1000

1500 time [s]

2000

2500

3000

6

4.4 4.8 5.2 0

8.698

x 10

4

F4 [m3/s]

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 36 of 43

8.697 0

Figure 5: Control action inputs of sub-system 1: CSTR-1 and CSTR-2; solid lines (—): centralized scheme, circles (◦): CDMPC scheme, dotted lines (· · · ): upper and lower bounds

36 ACS Paragon Plus Environment

Page 37 of 43

4.1

x 10

6

6

9.9

x 10

9.7

Q4 [J/s]

Q3 [J/s]

4.5

4.9

6.5

9.3

8.9

5.3 0

500

x 10

1000

1500 time [s]

2000

2500

8.5 0

3000

6

8.698

x 10

500

1000

1500 time [s]

2000

2500

3000

500

1000

1500 time [s]

2000

2500

3000

4

6

Q5 [J/s]

F [m3/s]

6.1

5.7

5.3 0

500

1000

1500 time [s]

2000

2500

8.697 0

3000

Figure 6: Control action inputs of sub-system 2: CSTR-3, the separator, and CSTR-4; solid lines: centralized scheme, circles: CDMPC scheme, dotted lines: upper and lower bounds 4

CDMPC Coordination Cycles

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

3

2

1

0

0

500

1000

1500

time [s]

2000

2500

3000

Figure 7: Coordination cycles of the proposed CDMPC scheme (22) In addition, temperature profiles of CSTR-1, CSTR-2, CSTR-3, and CSTR-4 in Figure 4 show that the state trajectories obtained from the CDMPC scheme (22) tracks the optimal 37 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

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

plant-wide trajectory perfectly. Control action inputs of sub-systems 1 and 2 are also shown in Figures 5 and 6, respectively. These also demonstrate that the proposed CDMPC control scheme (22) is capable of achieving the performance of the centralized MPC control scheme (9). Furthermore, Figure 7 shows the of coordination cycles ‘s’ of the proposed algorithm, i.e. the number of iterations between the local MPCs and the coordinator. This shows that the proposed scheme requires two coordination cycles at each sampling time to satisfy the stopping criterion ||E(k)|| ≤ 10−6 .

5

Conclusion

In this study, a nested optimization approach for price-driven coordination of distributed MPC networks is proposed for nonlinear systems. The developed CDMPC scheme upgrades an existing decentralized MPC control network via adding a coordinator level and applying limited modifications to the local MPC controllers. The main goal is to improve the closed-loop performance of the currently installed decentralized scheme to have an optimal plant-wide performance via capturing the interactions between local processing units. The proposed approach is an on-line scheme that implements a price-driven CDMPC algorithm for open-loop stable nonlinear systems via successive linearization. Unlike the CDMPC approaches that utilize active-set optimization algorithms, 17,19 the proposed method uses an interior-point approach which obviates the need for identifying the correct active-set. In this algorithm, the distributed MPCs solve their local problems for a fixed value of the price vector and communicate the sensitivity information with the coordinator. Then, the coordinator updates the price vector via a variant of the Newton’s method and communicates the updated price with the local MPCs. This iterative procedure continues until an optimal pant-wide solution is obtained. The convergence of the proposed CDMPC algorithm is demonstrated via choosing a proper merit function. Assuming that the centralized MPC (9) can stabilize the entire system, a Lyapunov-based criterion is derived to ensure the

38 ACS Paragon Plus Environment

Page 38 of 43

Page 39 of 43

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

closed-loop stability of the proposed CDMPC, even when the algorithm is stopped prematurely prior to convergence. The closed-loop performance of the proposed CDPMC scheme is evaluated with a benzene alkylation benchmark. According to the simulation study, the proposed CDMPC control scheme (22) is able to improve the closed-loop performance of the previously installed decentralized MPC scheme (57) and tracks the performance of the centralized MPC control scheme (9).

6

Acknowledgments

The authors would like to thank the Natural Sciences and Engineering Research Council of Canada (NSERC) for financial support of this work.

References (1) Sandell, N. R.; Varajya, P.; Athans, M.; Safonov, M. C. Survey of decentralized control methods for large scale systems. IEEE Transactions on Automatic Control 1978, 23, 108–128. (2) Bakule, L. Decentralized control: an overview. Annual Reviews in Control 2008, 32, 87–98. (3) Sun, Y.; El-Farra, N. H. Quasi-decentralized model-based networked control of process systems. Computers & Chemical Engineering 2008, 32, 2016–2029. (4) Lu, J. Z. Challenging control problems and emerging technologies in enterprise optimization. Control Engineering Practice 2003, 11, 847–858. (5) Morari, M.; H. Lee, J. Model predictive control: past, present and future. Computers & Chemical Engineering 1999, 23, 667–682.

39 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

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

(6) Camponogara, E.; Jia, D.; Krogh, B. H.; Talukdar, S. Distributed model predictive control. IEEE Control Systems Magazine 2002, 22, 44–52. (7) Rawlings, J. B.; Stewart, B. T. Coordinating multiple optimization-based controllers: new opportunities and challenges. Journal of Process Control 2008, 18, 839–845. (8) Scattolini, R. Architectures for distributed and hierarchical model predictive control a review. Journal of Process Control 2009, 19, 723–731. (9) Christofides, P. D.; Scattolini, R.; De La Peña, D. M.; Liu, J. Distributed model predictive control: a tutorial review and future research directions. Computers & Chemical Engineering 2013, 51, 21–41. (10) Christofides, P. D.; Liu, J.; Muñoz De La Peña, D. Networked and distributed predictive control: methods and nonlinear process network applications; Advances in Industrial Control Series; Springer-Verlag, London, England, 2011. (11) Stewart, B. T.; Venkat, A. N.; Rawlings, J. B.; Wright, S. J.; Pannocchia, G. Cooperative distributed model predictive control. Systems & Control Letters 2010, 59, 460–469. (12) Liu, J.; Muñoz De La Peña, D.; Christofides, P. D. Distributed model predictive control of nonlinear process systems. AIChE Journal 2009, 55, 1171–1184. (13) Liu, J.; Chen, X.; PeÃśa, D. M. D. L.; Christofides, P. D. Sequential and Iterative Architectures for Distributed Model Predictive Control of Nonlinear Process Systems. Part I: Theory. American Control Conference. Baltimore, Maryland, 2010; pp 3148– 3155. (14) Al-Gherwi, W.; Budman, H.; Elkamel, A. A robust distributed model predictive control algorithm. Journal of Process Control 2011, 21, 1127–1137.

40 ACS Paragon Plus Environment

Page 40 of 43

Page 41 of 43

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

(15) Maestre, J. M.; Muñoz De La Peña, D.; Camacho, E. F. Distributed model predictive control based on a cooperative game. Optimal Control Applications and Methods 2011, 32, 153–176. (16) Tippett, M. J.; Bao, J. Distributed model predictive control based on dissipativity. AIChE Journal 2013, 59, 787–804. (17) Cheng, R.; F. Forbes, J.; Yip, W. Price-driven coordination method for solving plantwide MPC problems. Journal of Process Control 2007, 17, 429–438. (18) Cheng, R.; Forbes, J. F.; Yip, W. Dantzig-Wolfe decomposition and plant-wide MPC coordination. Computers & Chemical Engineering 2008, 32, 1507–1522. (19) Marcos, N. I.; Forbes, J. F.; Guay, M. Price-driven coordination of distributed MPC controllers for constrained dynamic systems. Industrial & Engineering Chemistry Research 2013, 52, 17451–17464. (20) Martí, R.; Sarabia, D.; Navia, D.; De Prada, C. A method to coordinate decentralized NMPC controllers in oxygen distribution networks. Computers & Chemical Engineering 2013, 59, 122–137. (21) Feb, N. J.; Dantzig, G. B.; Wolfe, P. Decomposition principle for linear programs: decomposition principle for linear programs. Rand Symposium on Mathematical Programming. 1959; pp 101–111. (22) Cohen, G. Optimization by decomposition and coordination: a unified approach. IEEE Transactions on Automatic Control 1978, 23, 222–232. (23) Marcos, N. I.; Forbes, J. F.; Guay, M. Prediction-driven coordination of distributed MPC controllers for linear unconstrained dynamic systems. International Journal of Control 2014, 87, 1496–1512.

41 ACS Paragon Plus Environment

Industrial & Engineering Chemistry Research

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

(24) Scheu, H.; Busch, J.; Marquardt, W. Nonlinear distributed dynamic optimization based on first order sensitivities. American Control Conference. Baltimore, Maryland, 2010; pp 1574–1579. (25) Sharma, R.; Nesic, D.; Manzie, C. Idle speed control using linear time-varying model predictive control and discrete time approximations. IEEE International Conference on Control Applications. Yokohama, Japan, 2010; pp 1140–1145. (26) Marcos, N. I. Coordinated-Distributed Optimal Control of Large-Scale Linear Dynamic Systems. PhD Dissertation, University of Alberta, 2011. (27) Boyd, S. P.; Vandenberghe, L. Convex Optimization, 1st ed.; Cambridge University Press, 2004. (28) Khalil, H. K. Nonlinear Systems, 2nd ed.; Prentice-Hall, 1996. (29) Cheng, R. Decomposition and Coordination of Large-scale operations optimization. PhD Dissertation, University of Alberta, 2007.

42 ACS Paragon Plus Environment

Page 42 of 43

Page 43 of 43

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

7

Table of Contents Graphic

Decentralized MPC: MPC 1

MPC i

sub-system 1

MPC m

sub-system i

sub-system m Interactions

Interactions

Process

Coordinated Distributed MPC: CDMPC Coordinator

CDMPC Controller 1

CDMPC Controller i

sub-system 1

CDMPC Controller m

sub-system i

sub-system m Interactions

Interactions

Process

Centralized MPC: Centralized MPC

sub-system 1

sub-system i

sub-system m Interactions

Interactions

Process

TOC: The abstract graphical representation of the decentralized MPC, the proposed coordinated distributed MPC, and the centralized MPC architectures

43 ACS Paragon Plus Environment