9th International Symposium on Advanced Control of Chemical Processes 9th Symposium on June 7-10, 2015. Whistler, British Columbia,Control Canadaof 9th International International Symposium on Advanced Advanced Control of Chemical Chemical Processes Processes 9th International Symposium on Advanced Control of Chemical Processes June Canada Available online at www.sciencedirect.com June 7-10, 7-10, 2015. 2015. Whistler, Whistler, British British Columbia, Columbia, Canada June 7-10, 2015. Whistler, British Columbia, Canada
ScienceDirect IFAC-PapersOnLine 48-8 (2015) 026–031
Distributed Economic Model Predictive Distributed Economic Model Predictive Distributed Economic Model Distributed EconomicReactor: Model Predictive Predictive Control of a Catalytic Evaluation Control of a Catalytic Reactor: Evaluation Control of a Catalytic Reactor: Evaluation Control of a Catalytic Reactor: Evaluation of Sequential and Iterative Architectures of Sequential and Iterative Architectures of Architectures of Sequential Sequential and and Iterative Iterative Architectures ∗ ∗
Timothy L. Anderson ∗ Matthew Ellis ∗ ∗ Matthew ∗ Timothy L. Anderson Ellis ∗,∗∗,1 Timothy L. Ellis ∗ Matthew Panagiotis D. Christofides ∗,∗∗,1 Timothy L. Anderson Anderson Matthew Ellis ∗ ∗,∗∗,1 Panagiotis D. Christofides Panagiotis Panagiotis D. D. Christofides Christofides ∗,∗∗,1 ∗ ∗ Department of Chemical and Biomolecular Engineering, University of ∗ Department of Chemical and Biomolecular Engineering, University of of Biomolecular Engineering, University ∗ Department California, Los and Angeles, CA, 90095-1592, USA. Department of Chemical Chemical and Biomolecular Engineering, University of of California, Los Angeles, CA, 90095-1592, USA. ∗∗ California, Los Angeles, CA, 90095-1592, USA. Department of Electrical Engineering, University of California, Los ∗∗ California, Los Angeles, CA, 90095-1592, USA. ∗∗ Department of Electrical Engineering, University of California, Los University ∗∗ Department of Electrical Angeles, Engineering, CA, 90095-1592, USA. of Department of Electrical Engineering, University of California, California, Los Los Angeles, CA, 90095-1592, USA. Angeles, CA, 90095-1592, USA. Angeles, CA, 90095-1592, USA. Abstract: The development and application of distributed economic model predictive control Abstract: The development application of model predictive control Abstract: The and application of distributed economic model predictive (DEMPC) methodologies to aand catalytic reactor isdistributed considered.economic Two DEMPC methodologies are Abstract: The development development and application of is distributed economic model methodologies predictive control control (DEMPC) methodologies to a catalytic reactor considered. Two DEMPC are (DEMPC) to a catalytic reactor is considered. Two DEMPC methodologies are designed formethodologies sequential and iterative implementation, respectively. The DEMPC architectures (DEMPC) methodologies to a catalytic reactor is considered. Two DEMPC methodologies are designed for sequential and iterative implementation, respectively. DEMPC architectures designed for and implementation, respectively. The DEMPC architectures are evaluated on the basis of the closed-loop performance and The on-line computation time designed for sequential sequential and iterative iterative implementation, respectively. The DEMPC architectures are evaluated on the basis of the closed-loop performance and on-line computation time are evaluated on the basis of the closed-loop performance and on-line computation time requirements compared to a centralized EMPC approach. For the catalytic reactor considered, are evaluatedcompared on the basis of the closed-loop performance andcatalytic on-line reactor computation time requirements to a centralized approach. For the considered, requirements compared to EMPC approach. For catalytic considered, DEMPC proves to be a viable option as EMPC it is able to give similar closed-loop performance while requirements compared to a a centralized centralized EMPC approach. For the theclosed-loop catalytic reactor reactor considered, DEMPC proves to be a viable option as it is able to give similar performance while DEMPC proves to a it similar while reducing the on-line computation timeas requirements relative to aclosed-loop centralizedperformance EMPC strategy. DEMPC proves to be becomputation a viable viable option option asrequirements it is is able able to to give give similar performance while reducing the on-line time relative to aaclosed-loop centralized EMPC strategy. reducing the on-line computation time requirements relative to centralized EMPC strategy. reducing the on-line computation time requirements relative to a centralized EMPC strategy. © 2015, IFAC (International Federation of Automatic Control) Hosting by Elsevier Ltd. All rights reserved. Keywords: Economic model predictive control, distributed model predictive control, process Keywords: Economic model Keywords: Economic model predictive predictive control, control, distributed distributed model model predictive predictive control, control, process process control, nonlinear systems Keywords: Economic model predictive control, distributed model predictive control, process control, nonlinear systems control, nonlinear systems control, nonlinear systems 1. INTRODUCTION formulations, implementation strategies, and theoretical 1. INTRODUCTION formulations, implementation strategies, and theoretical 1. INTRODUCTION formulations, implementation strategies, and theoretical results have been developed within the context standard 1. INTRODUCTION formulations, implementation strategies, and of theoretical results have been developed within the context of standard results have been developed within the context of DMPC (e.g, Liu et al. (2009, 2010); see, also, the Operating chemical processes in an economically-optimal tracking results have been(e.g, developed within the 2010); contextsee, of standard standard tracking DMPC Liu et al. (2009, also, and the Operating chemical processes in an economically-optimal tracking DMPC (e.g, Liu et al. (2009, 2010); see, also, the reviews of Scattolini (2009); Christofides et al. (2013) Operating chemical processes in an economically-optimal manner while maintaining closed-loop stability and satistrackingofDMPC (e.g,(2009); Liu et al. (2009, 2010); see, also, and the Operating chemical processes in an economically-optimal reviews Scattolini Christofides et al. (2013) manner while maintaining closed-loop stability and satisreviews of Scattolini (2009); Christofides et al. (2013) and the references therein). manner while maintaining closed-loop stability and satisfying thewhile process constraints is an important issue within of Scattolini (2009); Christofides et al. (2013) and manner maintaining closed-loop stability andwithin satis- reviews the references therein). fying the process constraints is an important issue references therein). fying the process process constraints constraints isaccomplish an important important issue within the chemical control. Tois thisissue objective, the referencesprocess therein). fying the an within To integrate (dynamic) optimization and control, chemical process control. To accomplish this objective, chemical process control. To accomplish this objective, model predictive control (MPC) has proven to be an To integrate process (dynamic) control, chemical process control. To accomplish this objective, To integrate process (dynamic) optimization and control, economic MPC (EMPC), whichoptimization optimizes a and general cost model predictive control (MPC) has proven to be an To integrate process (dynamic) optimization and control, model predictive control (MPC) has proven to be an attractive way in many industrial applications (e.g., Qin economic MPC (EMPC), which optimizes a general cost model predictive control (MPC) has proven to be an economic MPC (EMPC), which optimizes a general cost function that represents the process economics instead of attractive way in industrial applications (e.g., Qin MPC (EMPC),the which optimizes a general cost attractive way(2003)). in many many industrial applications (e.g., that Qin economic and Badgwell MPC is a control methodology function that represents process economics instead of attractive way in many industrial applications (e.g., Qin that represents the process economics instead of afunction quadratic cost function, has been proposed as a conand Badgwell (2003)). MPC is a control methodology that thatcost represents thehas process economics instead of and Badgwell (2003)). MPC MPCcriterion is aa control control methodology that accounts for performance by optimizing a that cost afunction quadratic function, been proposed as a conand Badgwell (2003)). is methodology a quadratic cost function, has been proposed as a control methodology that mayhas help to enable future manaccounts for performance criterion by optimizing a cost a quadratic cost function, been proposed as a conaccounts for performance criterion by optimizing a cost function over a finite-time prediction horizon subject to trol methodology that may help to enable future manaccounts for performance criterion by horizon optimizing a cost trol methodology that may manufacturing tasks like processfuture operations function aa finite-time subject to methodology thatdemand-driven may help help to to enable enable future manfunction over finite-time prediction horizon subject to trol a processover model (to predictprediction the future behavior of the ufacturing tasks like demand-driven process operations function over a finite-time prediction horizon subject to ufacturing tasks like demand-driven process operations (e.g., Huang et al. (2011); Angeli et al. (2012); Heidarineaaprocess), process model (to predict the future behavior of the ufacturing tasks like demand-driven process operations process model (to predict the future behavior of the process constraints, andfuture stability constraints. Huang et al. see, (2011); et al. (2012); a process process model (to predict the behavior of the (e.g., (e.g., Huang et (2011); Angeli et Heidarinejad al. (2012); also,Angeli the review Ellis et Heidarineal. (2014) process), stability (e.g.,et Huang et al. al. see, (2011); Angeli et al. al. (2012); (2012); Heidarineprocess), process constraints, and stability constraints. Traditionally, the constraints, cost functionand used within constraints. MPC is a jad et al. (2012); also, the review Ellis et al. (2014) process), process constraints, and stability constraints. jad et al. (2012); see, also, the review Ellis et al. and etthe references therein). Recently, significant effort Traditionally, the cost function used within MPC is a jad al. (2012); see, also, the review Ellis et al. (2014) (2014) Traditionally, the cost function used within MPC is a quadratic cost that is positive definite with respect to an and the references therein). Recently, significant effort Traditionally, the cost function used with within MPCtoisana within and the therein). Recently, effort thereferences control community has focusedsignificant on (centralized) quadratic cost that is positive definite respect and the references therein). Recently, significant effort quadratic cost that is positive definite with respect to an operating steady-state of a process. within the control community has focused on (centralized) quadratic steady-state cost that is positive definite with respect to an EMPC. within the control community has focused on (centralized) Since EMPC may use a general (nonlinear) ecooperating of a process. within the control community has focused on (centralized) operating steady-state of a process. Since EMPC may use aa general (nonlinear) ecooperating of a process.in a receding horizon EMPC. EMPC. Since EMPC may use general (nonlinear) economic cost function and may dictate a time-varying operGiven thatsteady-state MPC is implemented EMPC.cost Since EMPCand may usedictate a general (nonlinear)opereconomic function may a time-varying Given that MPC is implemented in a receding horizon nomic cost function and may dictate a time-varying operation strategy, the on-line computation required to solve Given that MPC is implemented in a receding horizon fashion (i.e., an optimization problem is solved on-line nomic cost function and may dictate a time-varying operGiven that MPC is implemented in a is receding horizon ation strategy, the on-line computation required to solve fashion (i.e., an optimization problem solved on-line ation strategy, the computation required solve may be significant for large-scale process fashion (i.e., an optimization optimization problem is solvedactions), on-line EMPC at each (i.e., sampling time to compute the is control ation strategy, the on-line on-lineespecially computation required to to solve fashion an problem solved on-line EMPC may be significant especially for large-scale process at each sampling time to compute the control actions), EMPC may be significant especially for large-scale process networks. Thus, distributed EMPC (DEMPC) may be at each sampling time to compute the control actions), significant computation delay may result when computing EMPC may be significant especially for large-scale process at each sampling time to compute the control actions), networks. Thus, distributed EMPC (DEMPC) may be significant computation delay may result when computing networks. Thus, distributed EMPC (DEMPC) may one choice to significantly reduce the on-line computasignificant computation delay may result when computing control actions for process systems of high dimension (i.e., networks. Thus, distributedreduce EMPCthe(DEMPC) may be be significant computation delay may result when computing one choice to significantly on-line computacontrol actions for process systems of high dimension (i.e., one choice to significantly reduce the on-line computationalchoice burden. date, onlyreduce a limited work on control actions forinputs) processwhich systems ofaffect high closed-loop dimension (i.e., (i.e., many states andfor mayof sta- one to To significantly theamount on-lineof computacontrol actions process systems high dimension tional burden. To date, only a limited amount of work on many states and inputs) which may affect closed-loop stational burden. To only aa limited amount of on for linear systems al. (2012); M¨ uller many states and inputs) inputs) In which may affect affect closed-loop sta- DEMPC bility and performance. the context of control of largetional burden. To date, date, only[Driessen limitedet amount of work work on many states and which may closed-loop staDEMPC for linear systems [Driessen et al. (2012); M¨ u ller bility and performance. In the context of control of largeDEMPC for linear systems [Driessen et al. (2012); M¨ u ller and Allg¨ o wer (2014)] and for nonlinear systems [Chen bility and performance. In the context of control of largescale nonlinear chemical process networks, an alternative DEMPC for linear systems [Driessen et al. (2012); M¨ u ller bility and performance. In the context of control of largeand Allg¨ werLee (2014)] and for for nonlinear systems [Chen scale process alternative Allg¨ ooower (2014)] and nonlinear systems [Chen et (2012); and Angeli (2012)] has been completed. scale nonlinear chemical process networks, anarchitecture alternative and is to nonlinear employ a chemical distributed MPCnetworks, (DMPC)an andal. Allg¨ werLee (2014)] and for nonlinear systems [Chen scale nonlinear chemical process networks, an alternative et al. (2012); and Angeli (2012)] has been completed. is to employ a distributed MPC (DMPC) architecture et al. (2012); Lee and Angeli (2012)] has been completed. While these works on distributed EMPC (DEMPC) have is to employ a distributed MPC (DMPC) architecture (e.g., Christofides et al. (2013)). DMPC hasarchitecture the ability et al. (2012); Lee and Angeli (2012)] has been completed. is to Christofides employ a distributed MPC DMPC (DMPC) While these works on distributed EMPC (DEMPC) have (e.g., et al. (2013)). has the ability While these works on distributed EMPC (DEMPC) have shown some promising results on DEMPC, more work in (e.g., Christofides et al. (2013)). DMPC has the ability to control large-scale multiple-input multiple-output with While these works on distributed EMPC (DEMPC) have (e.g., Christofides et multiple-input al. (2013)). DMPC has the ability shown some promising results on DEMPC, more work in to control large-scale multiple-output with shown some promising results on DEMPC, more work this direction is in order. to control large-scale multiple-input multiple-output with input and state constraints while remaining computationsome promising results on DEMPC, more work in in to control large-scale multiple-input multiple-output with shown this direction is in order. input and state constraints while remaining computationthis direction is in order. input and state state constraints while remaining remaining computationally feasible to constraints be implemented on-line through a dis- this direction is in order. input and while computationIn the present work, sequential and iterative distributed ally feasible to on-line aa disally feasible to be be implemented implemented on-line through through dis- In tributed implementation of the computations. Numerous the present work, distributed ally feasible to be implemented on-line through a disIn the work, sequential and iterative distributed EMPC strategies are sequential developed and and iterative applied to a benchtributed implementation of the computations. Numerous In the present present work, sequential and iterative distributed tributed implementation of the computations. Numerous EMPC strategies are developed and applied to a benchbench1 Corresponding tributed implementation of the computations. Numerous EMPC strategies are developed and applied to a mark catalytic reactor where time-varying operation of author: Panagiotis D. Christofides. Tel.:+1 310 794 EMPC strategies are developed and applied to a bench1 Corresponding author: Panagiotis D. Christofides. Tel.:+1 310 794 mark catalytic reactor where time-varying operation of 1 mark catalytic reactor where time-varying operation of the reactor gives greater yield of the product compared 1015; fax: +1 310 206 4107; e-mail:
[email protected]. Financial Corresponding author: Panagiotis D. Christofides. Tel.:+1 310 794 1 Corresponding author: Panagiotis D. Christofides. Tel.:+1 310 794 mark catalytic reactor where time-varying operation of the reactor gives greater yield of the product compared 1015; fax: +1 310 206 4107; e-mail:
[email protected]. Financial the reactor gives greater yield of the product compared support from Science Foundation and the Department 1015; fax: +1the 310National 206 4107; e-mail:
[email protected]. Financial to steady-state operation. A description of the DEMPC the reactor gives greater yield of the product compared 1015; fax: +1the 310National 206 4107; e-mail:
[email protected]. Financial support from Science Foundation and the Department to steady-state operation. A description of the DEMPC of Energy is gratefully acknowledged. support from the National Science Foundation and the Department to support from the National Science Foundation and the Department to steady-state steady-state operation. operation. A A description description of of the the DEMPC DEMPC
of of Energy Energy is is gratefully gratefully acknowledged. acknowledged. of Energy is gratefully acknowledged. Copyright © 2015, 2015 IFAC 26 Hosting by Elsevier Ltd. All rights reserved. 2405-8963 © IFAC (International Federation of Automatic Control) Copyright 2015 26 Copyright ©under 2015 IFAC IFAC 26 Control. Peer review© of International Federation of Automatic Copyright © 2015 responsibility IFAC 26 10.1016/j.ifacol.2015.08.152
IFAC ADCHEM 2015 June 7-10, 2015. Whistler, BC, Canada Timothy L. Anderson et al. / IFAC-PapersOnLine 48-8 (2015) 026–031
implementation strategies is provided. Several closed-loop simulations are performed to evaluate the approaches. Two key performance metrics are considered in the evaluation: the closed-loop economic performance under the various DEMPC strategies and the on-line computation time required to solve the EMPC optimization problems.
The stage cost le (x, u1 , . . . , um ) of the EMPC is one of the design/tuning elements of the EMPC. It is chosen to reflect the process economics and need not be a quadratic stage cost like that typically used with standard tracking MPC. The stage cost of the EMPC is referred to as the economic cost function. The computed input profile optimizes the economic cost (2a) over the prediction horizon while accounting for the following constraints. The constraint (2b) is the dynamic model of the process initialized with a state measurement (2c) received at sampling instance tk . The nominal dynamic model predicts the future behavior of the process under any input trajectories u1 (τ ), . . . , um (τ ) for τ ∈ [0, N ∆) and allows for the EMPC to compute the optimal input trajectories. The optimal input trajectories are denoted: u∗1 (τ |tk ), . . . , u∗m (τ |tk ) for τ ∈ [0, N ∆). The bounds on the inputs are given by the constraints of (2d). Lastly, the constraint (2e) represents economics-based constraints which are typically integral constraints.
1.1 Class of Nonlinear Systems The class of nonlinear systems considered are described by the following system of first-order ordinary differential equations: (1) x(t) ˙ = f (x(t), u1 (t), . . . , um (t), w(t)) nui nx for where x(t) ∈ R denotes the state vector, ui (t) ∈ R i = 1, . . . , m denotes the ith manipulated (control) input vector, w(t) ∈ Rnw denotes the disturbance vector. The (full) input vector has been divided into m input vectors given that m distributed controllers will be designed to control each of the m input vectors. The input vectors are bounded in a convex set denoted as Ui := {ui ∈ Rnui | uij,min ≤ uij ≤ uij,max , j = 1, . . . , nui } for i = 1, . . . , m where uij,min and uij,max denote the minimum and maximum bound on the jth element of the ith input vector, respectively. Additionally, the disturbance vector is assumed to be bounded: w(t) ∈ W := {w ∈ Rl | |w| ≤ θ} where θ > 0 bounds the norm of the disturbance vector. The vector field of the system of Eq. 1 is assumed to be a locally Lipschitz vector function of its arguments, and the origin of the unforced system is the equilibrium point of Eq. 1 (i.e., f (0, 0, . . . , 0, 0) = 0). A state measurement of the system of Eq. 1 is assumed to be available synchronously at sampling instances denoted as tk := t0 + k∆ where t0 is the initial time, k ∈ I≥0 and ∆ > 0 is the sampling period.
EMPC, like standard tracking MPC, is implemented in a receding horizon fashion. At a sampling instance tk , the controller receives the current state measurement x(tk ), computes the optimal input trajectories u∗1 (τ |tk ), . . . , um (τ |tk ) for τ ∈ [0, N ∆) (which corresponds from tk to tk+N ), and implements the control action computed for the first sampling period in the prediction horizon on the process: ui (t) = u∗i (0|tk ) for t ∈ [tk , tk+1 ). The process is repeated at the next sampling time by rolling the horizon one sampling period. 2. CATALYTIC REACTOR EXAMPLE A catalytic reactor example is considered to evaluate various DEMPC implementation strategies. A non-isothermal continuous stirred tank reactor (CSTR) where ethylene is catalytically converted to ethylene oxide is considered. Besides the oxidation reaction, two combustion reactions occur that consume ethylene and ethylene oxide. The three reactions are given by: 1 r1 C2 H4 + O2 → C2 H4 O (R1) 2 r2 C2 H4 + 3O2 → 2CO2 + 2H2 O (R2) 5 r3 C2 H4 O + O2 → 2CO2 + 2H2 O (R3) 2 The reactor has a cooling jacket to remove the heat generated by the three exothermic reactions. The catalytic reactor has three manipulated inputs: the volumetric flow rate of the reactor feed, the ethylene concentration in the reactor feed, and the coolant jacket temperature.
1.2 Economic Model Predictive Control In a centralized approach, one can design an EMPC system that computes control actions for all m input vectors. EMPC, implemented in a centralized approach for the system of Eq. 1, is formulated as follows: N∆ le (˜ x(τ ), u1 (τ ), . . . , um (τ )) dτ (2a) maximize u1 ,...,um ∈S(∆)
subject to
0
x ˜˙ (τ ) = f (˜ x(τ ), u1 (τ ), . . . , um (τ ), 0) x ˜(0) = x(tk ) ui (τ ) ∈ Ui , ∀ τ ∈ [0, N ∆) N∆ g(˜ x(τ ), u(τ )) dτ ≤ 0
27
(2b) (2c) (2d)
The gaseous mixture contained in the reactor is assumed to be an ideal gas. By employing other standard modeling assumptions, a dynamic model can be developed for the catalytic reactor, and the resulting dynamic model has four states: the reactor gas mixture density, the reactor ethylene concentration, the reactor ethylene oxide concentration, and the reactor temperature. The states are denoted as x1 , x2 , x3 , and x4 , respectively, and the inputs are denoted u1 , u2 , and u3 , respectively (dimensionless variable form is used for all variables). The complete ¨ ul¸sen et al. (1992) which uses model can be found in Ozg¨ the nonlinear Arrhenius reaction rate laws of Alfani and Carberry (1970). The admissible input values are given by the following sets:
(2e)
0
where i = 1, . . . , m and the notation S(∆) denotes the family of piecewise constant functions with period ∆, x ˜(τ ) denotes the predicted state trajectory under the piecewise constant input profiles, u1 (τ ), . . . , um (τ ), which are the decision variable of the dynamic optimization problem, ∆ is the sampling period of the EMPC, and N is a positive integer that denotes the prediction horizon (i.e., number of sampling periods in the prediction horizon). To distinguish between real-time and the prediction time of the EMPC, t denotes the real (continuous) time, tk denotes the discrete sampling instances where state feedback is obtained, and τ ∈ [0, N ∆) denotes the predicted time in the controller. 27
IFAC ADCHEM 2015 28 June 7-10, 2015. Whistler, BC, Canada Timothy L. Anderson et al. / IFAC-PapersOnLine 48-8 (2015) 026–031
u1 ∈ U1 := [0.0704, 0.7042] , u2 ∈ U2 := [0.2465, 2.4648] , u3 ∈ U3 := [0.6, 1.1] .
x1
1.1 1 0.9 0
50
100
150
200
250
300
350
400
450
50
100
150
200
250
300
350
400
450
50
100
150
200
250
300
350
400
450
50
100
150
200
250
300
350
400
450
2
x2
The economic performance is characterized by the average yield of ethylene oxide which is given by: tf u1 (t)x4 (t)x3 (t) dt (3) Y (tf ) = 0 tf u1 (t)u2 (t) dt 0 where tf is the length of operation of the catalytic reactor (t = 0 is taken to be the initial time of operation). The average molar flow rate of ethylene that may be fed to the reactor is fixed owing to practical considerations: 1 tf u1 (t)u2 (t) dt = 0.175 . (4) tf 0 The economic cost used within EMPC is: le (x, u) = u1 x4 x3 (5) which only consists of the numerator of the yield (3) given that the denominator is fixed by the constraint (4). The steady-state xTs = [0.998 0.424 0.032 1.002] corresponding to a steady-state input of uTs = [0.35 0.5 1.0] is within the range of operation of interest, satisfies constraint (4), and is open-loop asymptotically stable (i.e., stability is not an issue within the range of operation and the objective of applying EMPC is used to optimize the average yield). In all the simulations below, the process was initialized at the initial condition: xT0 = [0.95 0.50 0.21 1.00]
1 0 0
x3
0.4 0.2 0 0
x4
1 0.5 0
Time
Fig. 1. State trajectories under C-EMPC.
u1
1
0.5
0 0
50
100
150
200
250
300
350
400
450
50
100
150
200
250
300
350
400
450
50
100
150
200
250
300
350
400
450
u2
3 2 1 0 0
u3
1.2 1 0.8 0.6 0
Time
3. EVALUATION OF DEMPC METHODS TO THE CATALYTIC REACTOR
Fig. 2. Input trajectories computed by the C-EMPC. horizon approach was used within EMPC: at the beginning of the jth operating window, the prediction horizon was set to Nk := tp /∆ and the horizon was decreased by one at every subsequent sampling time (Nk = Nk−1 − 1 at the sampling instance tk ). At the beginning of the (j + 1)th operating window, the prediction horizon was set to tp /∆.
Several implementation strategies (centralized and distributed) are applied to the catalytic reactor. Studying the benefits of applying EMPC to the catalytic reactor has already been considered. In Ellis and Christofides (2014), improved average yield of ethylene oxide resulted by applying EMPC to the process compared to operating the reactor at steady-state as well as operating the reactor with an open-loop optimal periodic switching of the inputs ¨ ul¸sen et al. (1992). The closedu1 and u2 considered in Ozg¨ loop simulations below were programmed using C++ on a desktop computer with an Ubuntu Linux operating system R and an Intel CoreTM i7 3.4 GHz processor. To recursively solve the catalytic reactor dynamic model, the explicit Euler method was used. A step size of 0.00001 was used to simulate the closed-loop dynamics of the reactor, while a step size of 0.005 was used to solve the model within the EMPC problem; both led to stable numerical integration.
3.1 Centralized EMPC For this computational study, a centralized EMPC strategy, which is denoted C-EMPC, was considered to compare the two distributed implementation strategies and an EMPC of the form (6) was formulated for the catalytic reactor. The C-EMPC formulation is given by: Nk ∆ maximize u1 (τ )˜ x4 (τ )˜ x3 (τ ) dτ (6a) u1 ,u2 ,u3 ∈S(∆)
0
subject to x ˜˙ (τ ) = f (˜ x(τ ), u1 (τ ), u2 (τ ), u3 (τ )) (6b) ui (τ ) ∈ Ui ∀ τ ∈ [0, Nk ∆) (6c) Nk ∆ 1 u1 (τ )u2 (τ ) dτ tp 0 1 tk u∗ (t)u∗2 (t) dt (6d) = 0.175 − tp t0 +jtp 1 where i = 1, 2, 3 and u∗1 and u∗2 denote the optimal control actions applied to the reactor from the beginning of the current operating window to current sampling time, tk .
Regarding the implementation details of the EMPC systems below, a sampling period of ∆ = 1.0 was used. The optimization problems were solved using the interior point solver Ipopt (W¨ achter and Biegler (2006)). To account for real-time computation considerations, the solver was forced to terminate after 100 iterations and/or after 100 seconds of computation time. The tolerance of the solver was set to 10−5 . To satisfy the constraint on the amount of ethylene that may be fed to the reactor, this constraint was enforced over operating windows of length tp = 47, that is the average molar flow rate of ethylene must be equal to 0.175 at the end of each operating window. A shrinking
Figs. 1-2 depict the closed-loop state and input trajectories under the C-EMPC scheme over ten operating windows. Similar to the results of Ellis and Christofides (2014), the 28
IFAC ADCHEM 2015 June 7-10, 2015. Whistler, BC, Canada Timothy L. Anderson et al. / IFAC-PapersOnLine 48-8 (2015) 026–031
1.1
u∗1 (0|tk ), u∗2 (0|tk )
x1
EMPC-1
1 0.9 0
u∗1 (τ |tk ), u∗2 (τ |tk )
29
50
100
150
200
250
300
350
400
450
50
100
150
200
250
300
350
400
450
50
100
150
200
250
300
350
400
450
50
100
150
200
250
300
350
400
450
Process 2
EMPC-2
x2
u∗3 (0|tk )
1 0 0
x(tk )
x3
0.4
Fig. 3. A block diagram of the S-DEMPC 1-2 scheme.
0.2 0 0
1
x4
C-EMPC distributes the ethylene in a non-uniform fashion with respect to time to optimize the yield of ethylene oxide. The average yield of ethylene oxide of the reactor under the C-EMPC is 10.22%. On the other hand, the average yield of ethylene oxide of the reactor over the same length of operation under constant steady-state input values is 6.38%, and the average yield under EMPC is 60% better than that of steady-state operation.
0.5 0
Time
Fig. 4. State trajectories under the S-DEMPC 1-2.
u1
1
3.2 Sequential DEMPC
0.5
0 0
A sequential implementation strategy computes the control actions for the process by computing a series of distributed controllers in succession. The first controller computes an input trajectory for the first input trajectory (i.e., u1 of system (1)). The input trajectory u1 (t) is sent to the next controller to solve for the input trajectory u2 (t). The input trajectory u3 (t) is computed by the third controller after the input trajectories u1 (t) and u2 (t) are received from the previous controllers. The process is repeated until control actions for all m input vectors have been computed.
50
100
150
200
250
300
350
400
450
50
100
150
200
250
300
350
400
450
50
100
150
200
250
300
350
400
450
u2
3 2 1 0 0
u3
1.2 1 0.8 0.6 0
Time
Fig. 5. Input trajectories computed by the S-DEMPC 1-2.
For this process example, which has three inputs, a reasonable choice of input grouping can be made as a consequence of the integral input constraint (4) (i.e., u1 and u2 should be computed by the same EMPC, while it is worth investigating if the input u3 can be placed on another EMPC system). This input pairing will be used in all of the DEMPC schemes below and the resulting EMPC system that computes control actions for u1 and u2 is denoted as EMPC-1, and the other EMPC that computes control actions for u3 is denoted as EMPC-2. The formulations of each EMPC system follows from the centralized EMPC formulation (6) and are omitted due to space constraints.
variable is u3 only, u1 (τ ) and u2 (τ ) are set to the values computed by EMPC-1, and there is no integral input constraint (6d)). A block diagram of the resulting control architecture is given in Fig. 3. Figs. 4-5 show the closed-loop state and input trajectories under the S-DEMPC 1-2, respectively, and the trajectories are similar to those under the C-EMPC (Figs. 1-2). For this closed-loop simulation, the average yield was 10.20% (recall, the average yield under the C-EMPC was 10.22%). The difference between the average yield under the CEMPC and under the S-DEMPC 1-2 is small and likely numerically insignificant given the solver parameters used. Some differences in the state trajectories are observed from Fig. 1 and Fig. 4 (e.g., x1 (t) and x4 (t)). It is important to note that given the nonlinear nature of the process considered, there is no guarantee, in general, that the centralized EMPC and sequential EMPC scheme will lead to the same optimal input trajectories.
Sequential DEMPC 1-2 The first configuration considered, which is referred to as the sequential DEMPC 1-2 (abbreviated as S-DEMPC 1-2), first computes EMPC-1 for the optimal input trajectories u∗1 (τ |tk ) and u∗2 (τ |tk ) for τ ∈ [0, Nk ∆). Then, EMPC-2 computes the input trajectory u∗3 (τ |tk ) after receiving u∗1 (τ |tk ) and u∗2 (τ |tk ) from EMPC-1. Since the input trajectory u3 (τ ) has not been determined when EMPC-1 is computed, it is set to be the resulting input trajectory under a PI controller implemented in a sample-and-hold fashion over the prediction horizon (other methods for the assumed profile of u3 (t) within EMPC-1 could be considered). The input constraints are accounted for in the computed PI input trajectory (e.g., if the PI controller computed a control action greater than the upper bound on u3 , it was set to u3,max ). The optimal input trajectories u∗1 (τ |tk ) and u∗2 (τ |tk ) are used in EMPC-2 to predict the behavior of the reactor over the prediction horizon (i.e., the optimization problem of EMPC-2 is similar to (6) except the decision
Sequential DEMPC 2-1 Another sequential implementation of EMPC-1 and EMPC-2 may be considered by reversing the execution of EMPC-1 and EMPC-2. In this case, EMPC-2 computes its optimal input trajectory u∗3 (τ |tk ) first. The sequential DEMPC approach is referred to as sequential DEMPC 2-1 (S-DEMPC 2-1). To solve EMPC-2, the trajectories u1 (τ ) and u2 (τ ) are set to the input trajectories resulting from two PI controllers implemented in sample-and-hold fashion. While the bounds on admissible input values are accounted for in the PI input trajectories, the input average constraint of Eq. 4 is 29
IFAC ADCHEM 2015 30 June 7-10, 2015. Whistler, BC, Canada Timothy L. Anderson et al. / IFAC-PapersOnLine 48-8 (2015) 026–031
1.1
x1
x1
1.1 1 0.9 0
50
100
150
200
250
300
350
400
1 0.9 0
450
x2
x2
1 50
100
150
200
250
300
350
400
x3
x3
250
300
350
400
450
50
100
150
200
250
300
350
400
450
100
150
200
250
300
350
400
50
100
150
200
250
300
350
400
450
50
100
150
200
250
300
350
400
450
0.2 0 0
450
1
x4
50
1
x4
200
0.4
0.2
0.5 0
150
1 0 0
450
0.4
0 0
100
2
2
0 0
50
50
100
150
200
250
300
350
400
0.5 0
450
Time
Time
Fig. 9. State trajectories under the I-DEMPC (1 iteration).
Fig. 6. State trajectories under the S-DEMPC 2-1.
u1
1
0.5
u2
0 0
50
100
150
200
250
300
350
400
0 0
450
3
3
2
2
1 0 0
50
100
150
200
250
300
350
400
100
150
200
250
300
350
400
450
50
100
150
200
250
300
350
400
450
50
100
150
200
250
300
350
400
450
1.2
u3
1 0.8 0.6 0
50
1 0 0
450
1.2
u3
0.5
u2
u1
1
1 0.8 0.6
50
100
150
200
250
300
350
400
450
0
Time
Time
Fig. 7. Input trajectories of computed by the S-DEMPC 2-1. EMPC-1 u∗,c−1 (τ |tk ), 1 (τ |tk ) u∗,c−1 2
Fig. 10. Input trajectories of computed by the I-DEMPC (1 iteration).
∗,f u∗,f 1 (0|tk ), u2 (0|tk )
(τ |tk ) u∗,c−1 3
solution (even after infinite iterations). Moreover, there is no guarantee that the input solution computed at each iteration improves upon the closed-loop performance over the previous iteration. An iterative DEMPC (I-DEMPC) scheme is designed for the catalytic reactor and a block diagram of the I-DEMPC control architecture is given in Fig. 8. The computed input trajectories at each iteration of the I-DEMPC is denoted as u∗,c i (τ |tk ), i = 1, 2, 3 where c is the iteration number. At the first iteration, the input trajectory u3 in EMPC-1 is initialized with the sample-and-hold input trajectory computed from the same PI controller used in the S-DEMPC 1-2 scheme, and similarly, the input trajectories u2 and u3 for EMPC-2 are computed from the PI controllers of the S-DEMPC 2-1 scheme. The control action applied to the reactor is denoted as u∗,f i (tk |tk ) for i = 1, 2, 3 where f is the number of iterations of the iterative DEMPC scheme (f is a design parameter of the scheme). When f = 1, the I-EMPC scheme is decentralized in the sense that there is no communication between EMPC-1 and EMPC-2.
Process u∗,f 3 (0|tk )
EMPC-2 x(tk )
Fig. 8. A block diagram of the I-DEMPC scheme. not accounted for in the PI input trajectories. The block diagram describing this DEMPC architecture is similar to that of Fig. 3 with communication between EMPC1 and EMPC-2 in the opposite direction. Figs. 6-7 are the closed-loop state and input trajectories under the SDEMPC 2-1 approach. Compared to the other trajectories more noticeable differences are observed. 3.3 Iterative DEMPC Instead of sequential computation of the distributed EMPC schemes, parallel computation may be employed. Given the control actions are computed without the knowledge of the control actions computed by the other distributed EMPC schemes, an iterative approach may be used to (ideally) compute control actions closer to the centralized solution. It is important to note that given the nonlinearity and non-convexity of the optimization problems, it is difficult, in general, to guarantee that an iterative DEMPC strategy will converge to the centralized
For this example, no closed-loop performance benefit was observed after iterating more than once through the IDEMPC scheme. In fact, using the previous iterate solution to compute the next iterative gave worse closed-loop performance than applying the first computed iteration to the process. One method considered to compensate for this problem was to use the best computed input solution over all iterations to compute the next iteration. However, minimal closed-loop performance benefit was 30
IFAC ADCHEM 2015 June 7-10, 2015. Whistler, BC, Canada Timothy L. Anderson et al. / IFAC-PapersOnLine 48-8 (2015) 026–031
Table 1. The average yield and computation time under the EMPC strategies. Strategy Sequential DEMPC 1-2 Sequential DEMPC 2-1 Iterative DEMPC (f = 1) Centralized EMPC
Yield (%)
Comp. time (s)
10.20 9.92 10.05 10.22
1.039 2.969 0.832 4.244
31
Angeli, D., Amrit, R., and Rawlings, J.B. (2012). On average performance and stability of economic model predictive control. IEEE Transactions on Automatic Control, 57, 1615–1626. Chen, X., Heidarinejad, M., Liu, J., and Christofides, P.D. (2012). Distributed economic MPC: Application to a nonlinear chemical process network. Journal of Process Control, 22, 689–699. 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 & Chemical Engineering, 51, 21–41. Driessen, P.A.A., Hermans, R.M., and van den Bosch, P.P.J. (2012). Distributed economic model predictive control of networks in competitive environments. In Proceedings of the IEEE 51st Annual Conference on Decision and Control, 266–271. Maui, HI. Ellis, M. and Christofides, P.D. (2014). Optimal timevarying operation of nonlinear process systems with economic model predictive control. Industrial & Engineering Chemistry Research, 53, 4991–5001. Ellis, M., Durand, H., and Christofides, P.D. (2014). A tutorial review of economic model predictive control methods. Journal of Process Control, 24, 1156–1178. Heidarinejad, M., Liu, J., and Christofides, P.D. (2012). Economic model predictive control of nonlinear process systems using Lyapunov techniques. AIChE Journal, 58, 855–870. Huang, R., Harinath, E., and Biegler, L.T. (2011). Lyapunov stability of economically oriented NMPC for cyclic processes. Journal of Process Control, 21, 501– 509. Lee, J. and Angeli, D. (2012). Distributed cooperative nonlinear economic MPC. In Proceedings of the 20th International Symposium on Mathematical Theory of Networks and Systems. Melbourne, Australia. Liu, J., Chen, X., Mu˜ noz de la Pe˜ na, D., and Christofides, P.D. (2010). Sequential and iterative architectures for distributed model predictive control of nonlinear process systems. AIChE Journal, 56, 2137–2149. Liu, J., Mu˜ noz de la Pe˜ na, D., and Christofides, P.D. (2009). Distributed model predictive control of nonlinear process systems. AIChE Journal, 55, 1171–1184. M¨ uller, M.A. and Allg¨ower, F. (2014). Distributed economic MPC: a framework for cooperative control problems. In Proceedings of the 19th World Congress of the International Federation of Automatic Control, 1029– 1034. Cape Town, South Africa. ¨ ul¸sen, F., Adomaitis, R.A., and C Ozg¨ ¸ inar, A. (1992). A numerical method for determining optimal parameter values in forced periodic operation. Chemical Engineering Science, 47, 605–613. Qin, S.J. and Badgwell, T.A. (2003). A survey of industrial model predictive control technology. Control Engineering Practice, 11, 733–764. Scattolini, R. (2009). Architectures for distributed and hierarchical model predictive control – A review. Journal of Process Control, 19, 723–731. W¨achter, A. and Biegler, L.T. (2006). On the implementation of an interior-point filter line-search algorithm for large-scale nonlinear programming. Mathematical Programming, 106, 25–57.
observed with this method as well. Thus, f = 1 for this case given that using more than one iteration did not improve the closed-loop performance. The resulting closedloop trajectories are given in Figs. 9-10. The trajectories have similar characteristics as the centralized case. 3.4 Evaluation of DEMPC Approaches The average yield and average computation time required to solve the optimization problem at each sampling time over the entire simulation were considered for all the cases. The sequential DEMPC computation time is computed as the sum of the computation time of EMPC-1 and EMPC-2 at each sampling time because the sequential DEMPC schemes are computed sequentially. The iterative DEMPC is computed as the maximum computation time of any one EMPC at each sampling time (recall only one iteration was used). The average yield and average computation time for all the cases is given in Table 1. The centralized EMPC, sequential DEMPC 1-2, and iterative DEMPC schemes all gave similar closed-loop performance. The sequential DEMPC 1-2 and iterative DEMPC result in approximately a 70% reduction in computation time over the centralized EMPC. The sequential DEMPC 21 scheme not only had the worst performance of all the strategies considered (albeit still better than steady-state operation), but also, required a comparable amount of time to solve the optimization problems as the centralized case, thereby implying a strong dependence of closed-loop performance on controller calculation sequence. DEMPC was able to yield comparable closed-loop performance while substantially reducing the on-line computation time. This demonstrates that a distributed implementation may allow EMPC to be used on processes where centralized control is not feasible due to the solve time. This example illustrates another key point within the context of DEMPC. Specifically, the inclusion of integral constraint in EMPC may be an important consideration for input selection in DEMPC. From the sequential DEMPC results, the computed u3 profile is impacted by the assumed input profiles u1 and u2 (Fig. 7), while u1 and u2 are not affected as much by the assumed profile u3 (Fig. 5) compared to the centralized EMPC case (Fig. 2). This behavior may be due to the enforcement of the integral input constraint, and for this example, there may only be one method to distribute a fixed amount of ethylene to the reactor that maximizes the yield that is independent of u3 . REFERENCES Alfani, F. and Carberry, J.J. (1970). An exploratory kinetic study of ethylene oxidation over an unmoderated supported silver catalyst. La Chimica e L’Industria, 52, 1192–1196. 31