A Reactive Batch Distillation Column Case Study

ILC are based on learning of the inputs of the process from previous batches for ...... E. J.; González, A. H. Frontiers in Advanced Control Systems;...
0 downloads 0 Views 2MB Size
Subscriber access provided by BUFFALO STATE

Process Systems Engineering

Model Learning Predictive Control for Batch Processes: A Reactive Batch Distillation Column Case Study Alejandro Marquez Ruiz, Marco Loonen, M. Bahadir Saltik, and Leyla Ozkan Ind. Eng. Chem. Res., Just Accepted Manuscript • DOI: 10.1021/acs.iecr.8b06474 • Publication Date (Web): 24 Apr 2019 Downloaded from http://pubs.acs.org on April 25, 2019

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 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 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.

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 42 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

Model Learning Predictive Control for Batch Processes: A Reactive Batch Distillation Column Case Study Alejandro Marquez-Ruiz, Marco Loonen, M. Bahadir Saltik, and Leyla Özkan∗ Department of Electrical Engineering, Eindhoven University of Technology, De Zaale, Eindhoven 5612 AJ, The Netherlands E-mail: [email protected] Abstract In this paper, we present the control of batch processes using Model Predictive Control (MPC) and iterative learning Control (ILC). Existing combinations of MPC and ILC are based on learning of the inputs of the process from previous batches for a fixed linear time-invariant model (LTI). However, batch processes are inherently time varying therefore, LTI models are limited in capturing the relevant dynamic behaviour. An attractive alternative is to use Linear Parameter Varying (LPV) models because of their ability to capture nonlinearities in the control of batch processes. Therefore, in this work we propose a novel method combining MPC and ILC based on LPV models, and we call this method Model Learning Predictive Control (ML-MPC). Basically, the idea behind the method is to update the LPV model of the MPC iteratively, by using the repetitive behavior of the batch process. To this end, three different application-dependant options to estimate the parameters and disturbances of the model are proposed and are compared in simulation on a nonlinear batch reactor. Finally, the ML-MPC with one of the estimation method is applied to an industrial Reactive Batch Distillation Column (RBD) in simulation.

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

Introduction Growing quality demands of costumers, environmental constraints, limited raw materials and increasing energy demand pushes process industry to continuously improve their process and their operations. Continuous processes are often developed for bulk production of a specific product and are therefore less dynamic in their capacity and product variations. Meanwhile, batch processes, which are defined by their start-up and well-defined end of the operation, are generally used for smaller quantities of specialized products, but often have the ability to change their recipes-setup between the batches to provide a wider product range.

For the control of batch processes, methods such as PID, run-to-run, iterative learning or model based control, can be used as described in Bonvin 1 . Model Predictive Control (MPC) is a control method that includes the input and output constraints explicitly in its formulation, therefore, has become an accepted standard within process industry for constrained multivariate control problems. 2 In a typical implementation, MPC solves an optimization problem subject to system dynamics and constraints. The solution of this optimization problem is a sequence of inputs that can steer the performance output to desired levels. Naturally, the performance of such a control strategy depends on the quality of the model that is used in predicting the process behaviour .

There are two basic methods for modeling of systems; modeling based on system identification which uses only input-output-data, or modeling based on first principles. In the latter, process models with a wide operating window can be obtained. This kind of models are generally nonlinear hence require nonlinear optimization methods and can be computationally inefficient in a model based control strategy. An attractive modeling approach to deal with nonlinear systems is Linear Parameter Varying (LPV) models. These models describe the nonlinearity by a set of linear equations by using a scheduling parameter. The use of LPV models in control synthesis has several challenges of which the identification of the model is 2

ACS Paragon Plus Environment

Page 2 of 42

Page 3 of 42 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

the most important. 3–6 In this direction, we have shown that a well structured LPV model for reacting systems can be obtained by considering the reaction contribution as a disturbance term and using the extent transformation. 7–9 The resulting model has a parameter dependent state space form with a diagonal state matrix. An important advantage of such representation is that the model parameters have physical meaning.

Another important characteristic of batch processes is the repetitive nature of the operation. This could be an advantage for control purposes. While humans have the ability to improve after repetitive mistakes, classical feedback controllers applied in the process industry do not have this capability. If a controller can store information from the control actions calculated in the previous batches and use this information to improve the future and the current performance, then this is called a learning action. In the 1980’s, the first learning applications were introduced in the field of robotics, 10 and is now widely known in the field of systems and control as the Iterative Learning Control (ILC) methods, see Wang et al. 11 for a survey. ILC has the ability to obtain asymptotic convergence over iterations of the same repetitive task while being subject to model uncertainties and iteration- invariant disturbances. The conventional applications of ILC is operated in open loop, where the error of the previous iteration is used to update the trajectory of the next iteration, hence it is unable to reject real time disturbances. Therefore, some methods have been proposed that have implemented ILC as a feed-forward action in the closed loop control problem, see Wang et al. 11 for more details.

Certainly, the combination of MPC and ILC (IL-MPC) is a good approach in the control of batch processes. This combination has a lot of advantages such as double integration mechanism (in the batch and between batches), constraints satisfaction among others. The general learning procedure with MPC is achieved by correcting the MPC input with past information. 12 The mathematical formulations for IL-MPC can differ in the incorporation of integrators and constraints; see for example. 2,13 Despite of all these advantages, the model

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

used for prediction is still a limitation, also, the most of the processes in industry frequently show nonlinear behaviour, 14 , 15 . 16 For nonlinear processes, such as the batch reactors, several studies have been conducted where a type of IL-MPC is implemented with nonlinear model. For example 16,17 propose a iterative nonlinear model predictive control (INMPC) strategy. INMPC has been applied to two different batch reactions abd convergency, stability and robustness of the closed-loop system have been analyzed and proved. 18 Another nonlinear method, called as the nonlinear model predictive iterative learning control (NMPILC), is described in 19 which uses a fuzzy model that contains local linear models. In general terms, to utilize a nonlinear model in IL-MPC is the ideal case, nevertheless, as it is usual in nonlinear MPC, there are a lot of open issues when it comes to implementation such as computational load. Therefore, to close the gap between accuracy in the prediction and implementation in real time, the combination of IL-MPC with LPV models arises as a good alternative.

In this work, we propose a method not only to improve the performance of an MPC controller by utilizing the repetitive behavior of batch processes in an iterative learning fashion but also to overcome the limitations of MPC approaches using LPV models by updating the MPC and disturbance models. Similar strategies have been proposed for classic ILC, 20,21 however in this work, all these ideas are combined in just one method called Model Learning Predictive Control (ML-MPC). This method inherits the advantages of ILC and MPC while a suitable model in structure (linear) and accuracy is used for prediction in every batch. For the estimation of the parameters and disturbance, we evaluate three different methods that cover a broad range of scenarios. The performance of each method is evaluated implementing them on a simple nonlinear batch reactor. Furthermore, one of the proposed learning methods is used for the control of an industrial Reactive Batch Distillation (RBD) column.

This paper is organized as follows: In Chapter 2 we provide the basic principles of MPC,

4

ACS Paragon Plus Environment

Page 4 of 42

Page 5 of 42 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

MPC for LPV models, ILC and the current methods of Iterative Learning Model Predictive Control (IL-MPC). Chapter 3 presents the problem statement. In Chapter 4, we describe the proposed method of model learning predictive control, including implementations on a simple example of a batch reactor. In Chapter 5 the method is applied to a more advanced process in a reactive batch distillation column. Finally, in Chapters 6, the conclusions are presented.

Preliminaries MPC for LPV systems MPC is a well known control technique that uses a model to predict the future states and outputs. It solves an optimization problem, in which, usually, a quadratic cost function is minimized. In this way, the controller can steer the process outputs towards the reference profiles. The main advantage of MPC is that it can consider explicitly constraints on inputs, states, and outputs in the optimization problem. The optimization problem is solved at −1 every time sample to determine the best sequence of inputs {δuk+j }N j=0 for a given prediction

horizon N ∈ Z>0 . The first input from this sequence is applied, and at the next time sample the same optimization problem is solved again. This is known as the receding horizon strategy. Despite these advantages of MPC, the use of this technology for the control of batch processes is limited or non existing. This is mainly due to the nonlinearities, and absence of steady states. Therefore typical linear MPC techniques could not be implemented. For batch processes, MPC based on a first-principles nonlinear model is a suitable approach. Certainly, there are different approaches to deal with these issues such as multiple models, successive linearization, or empirical nonlinear models such as artificial neural networks, however, these approaches are often subject to performance limitations 22 and computational issues. An attractive approach to deal with batch processes and nonlinear systems is to use LPV models 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 42

in the MPC formulation. With such a representation we have the potential to capture nonlinear and time varying characteristics of systems but also use linear control theory which is well developed. In the available MPC methods for LPV systems, the future behavior of the parameter (scheduling variable) is assumed unknown over the prediction horizon. One of the most used method to solve the MPC problems for this kind of uncertainty is the worst case scenario which can be described by Equation (1). 23–26

min

δuk+j

kxk+N k2P

max

+

θk+j ∈ Θ

N −1 X

kek+j k2Q + kδuk+j k2R

j=0

dk+j ∈ D

subject to: xk+j+1 = A(θk+j )xk+j + B(θk+j )uk+j + Bd dk+j yk+j = Cxk+j + Duk+j

(1)

ek+j = rk+j − yk+j uk+j = uk+j−1 + δuk+j xk+j ∈ X, uk+j ∈ U, xk+N ∈ Ω dk+j ∈ D, θk+j ∈ Θ where xk+j ∈ Rnx denotes the system state, yk+j ∈ Rny is the system output, uk+j ∈ Rnu is the control vector, dk+j ∈ Rnd are disturbances, and kxk2Q = x| Qx with Q  0. The system is subject to hard constraints on state xk+j ∈ X and input uk+j ∈ U, where X ⊂ Rnx , and U ⊂ Rnu are assumed convex and compact. Finally, Ω ⊆ X is an invariant set of the system. In this work, the sets U, X and Ω are assumed polyhedrons, i.e., described by linear inequalities of the form,

U = {uk+j |Au uk+j 6 bu } X = {xk+j |Ax xk+j 6 bx } Ω = {xk+N |AN xk+N 6 bN } 6

ACS Paragon Plus Environment

(2)

Page 7 of 42 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

The main issue with the formulation given by optimization problem (1) is the computational demand. 27 Therefore, it has led to approaches in which the problem is approximated (see 23 and 28 ). We also find multiparametric approaches described as explicit MPC for LPV systems 26 and tube based MPC (TMPC). 29 In 26 the high computational issues are reduced by pre-computing the optimal inputs as a piecewise affine function of the state. The precomputed solutions are stored in a look-up table, such that only this look-up table can be used online. 26 However, this method only works well if the size of the system is low. 30 On the other hand, TMPC is based on the assumption that, while the current value of the scheduling parameter is measured exactly, the future variations in the scheduling variable belong to a sequence of sets Θk+j (tube) that describe the expected deviations from the nominal trajectory. 25 With this, the controller is able to derive the worst-case cost for expected future changes only, resulting in a single linear program. The sequence Θk+j can be constructed in several ways, namely classical, anticipative and oracle, each corresponding to a different MPC scenario. 25 An additional interested approach called MPC for Quasi-LPV systems (QMPC) has been proposed by. 31,32 In this formulation it is assumed that the scheduling parameters θk+j and disturbances dk+j can be calculated by means of explicit functions of the states and inputs as dk+j = fd (xik+j , uik+j ) and θk+j = fθ (xik+j , uik+j ). With this assumption, the min-max optimization problem given in (1) can be avoided. The main drawback of this method is to find the functions fd (·) and fθ (·). However, for chemical process, most of the time such functions are given by constitutive equations of the physical parameter of the process, which is advantageous if a first principle model is utilized in the controller.

In general, the MPC methods for LPV systems proposed in the literature are well defined mathematically, nevertheless, real time implementations are difficult due to their complexity and high computational load. In spite of all these issues, TMPC and QMPC combined with Iterative Learning Control (ILC) have inspired some of the methods proposed in this paper.

7

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 8 of 42

Iterative Learning Control (ILC) ILC is the control method that uses data from previously performed iterations. With this data, the controller is able to learn from past "mistakes" and correct for this. Figure 1 illustrates the ILC procedure. In this Figure P is the plant, and L and Q are filters. Notice in Figure 1 how the control signal induces an error signal through the plant with disturbance, and how this error signal of iteration j is used for the determination of the control signal in iteration j + 1.

L Q ui−1 j Memory rj

ej

yj C

uij

P

Figure 1: Illustration of an ILC control algorithm.

In classical ILC, the goal is to generate the optimal feed-forward input. The first solution for this has been applied by Arimoto on robotics, 10 which is the also called as P-Type ILC, 33 first-order learning algorithm 12 or Arimoto algorithm, 34 and represented by + Lei−1 uij = ui−1 j j

(3)

∀i = 1, . . . , ∞, ∀j = 1, . . . , N where uij ∈ Rnu is the input, ei−1 = rj − yji−1 , ej ∈ Rny is the error between the output and j the reference, N is the number of samples in batch, and L is an LTI filter called learning filter. 12 The indices i ∈ Z≥0 and j ∈ Z≥0 indicate the batch index and the time sample within a batch, respectively. To ensure that the ILC algorithm is learning the system dynamics and not the disturbances, the filter L is designed to be the inverse of the process sensitivity. 35 The 8

ACS Paragon Plus Environment

Page 9 of 42 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

drawback of this function is that high-frequency terms of the input can increase continuously causing instability across batches. To overcome this problem and ensure convergence to a steady input profile, a robustness filter Q is added to this algorithm as can be seen in Figure 1 and Equation 4. This robustness filter usually is a type of low-pass filter. 35

  i−1 uij = Q ui−1 + Le j j

(4)

The ILC algorithm can now be executed as shown in Algorithm 1. Algorithm 1 Iterative Learning Control algorithm 1: procedure Learning procedure 2: Set i = 1, and an initial input sequence uij ; 3: Run batch, with input sequence uij ; 4: Save input as ui−1 and the error trajectory as ei−1 j j ; i−1 i−1 5: Filter ej using filter L, such that uILC = Lej ; + uILC ; 6: Set uij = ui−1 j 7: Filter uij using robustness filter Q, such that uij = Quij ; 8: Set i = i + 1. Return to step 3. 9: end procedure

Feedback control can also be included in the closed loop operation with ILC in the so-called Current-Iteration Iterative Learning Control. 36,37 Then the control update algorithm is given by,

i−1 uik+j = Q[ui−1 k+j + Lek+j ] + | {z } Learning control

Ceik+j | {z }

Feedback control

Here, the feedback control has the ability to react to batch specific disturbances. Usually, in the corresponding control scheme, as shown in Figure 1, only PD controllers are considered, 36 as ILC has a natural integration component from one batch to the next. ILC strategy requires repetitive processes like robotics and batch processing and if this is the case it can reject repetitive errors using past data. It adds a batch to batch integrator to the control loop and the non-causal behaviour can overcome the delay in the error suppression.

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

Page 10 of 42

Additionally, model information can be used to improve the convergence speed. However, it requires error and input filtering to prevent instabilities. 36

Iterative Learning Model Predictive Control (IL-MPC) Different methods have already been developed to combine the benefits of the MPC control with iterative learning capabilities, as seen in Lee and Lee 12 and Adam and González 13 . In general, the IL-MPC problem considers additive disturbances, and mathematically can be written as,

min

∆δuik+j

kxik+N k2P

+

N −1 X

keik+j k2Q + k∆δuik+j k2R

j=0

subject to: xik+j+1 = Axik+j + Buik+j + Bd dik+j i yk+j = Cxik+j + Duik+j i i eik+j = rk+j − yk+j

(5)

uik+j = uik+j−1 + δuik+j ∆δuik+j = δuik+j − δui−1 k+j xik+j ∈ X, uik+j ∈ U xik+N ∈ Ω, dik+j ∈ D In this formulation ∆ describes the input change rate from one batch to the next. Similar to how δ includes an integrator in the in-batch time direction, the ∆ introduces an extra integrator in the batch-to-batch direction. Another feature of this ∆ is that the speed of changes over batches can be constrained. Now, in order to calculate the optimal control solution ∆δuik+j , the problem (5) can be transformed into a quadratic programming (QP) problem.

10

ACS Paragon Plus Environment

Page 11 of 42 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

Problem statement In the previous sections, we have described MPC and ILC strategies separately and also presented information on control approaches which combine both technologies. In IL-MPC approaches, the optimal input is calculated using the predefined model as shown in Equation 5. Changes in the plant from batch to batch causing model degradation, has not been considered. Besides, MPC-LPV methods are difficult for online implementation. Motivated by these observations, we present the following problem:

min

∆δuik+j

kxik+N k2P

+

N −1 X

keik+j k2Q + k∆δuik+j k2R

j=0

subject to: i i xik+j+1 = A(θk+j )xik+j + B(θk+j )uik+j + Bd dik+j i yk+j = Cxik+j + Duik+j

(6)

i i eik+j = rk+j − yk+j

uik+j = uik+j−1 + δuik+j ∆δuik+j = δuik+j − δui−1 k+j xik+j ∈ X, uik+j ∈ U, xik+N ∈ Ω i dik+j ∈ D, θk+j ∈Θ

Our principal objective, therefore, is to design an MPC controller to achieve constrained tracking of the LPV representation of the batch process for a given reference trajectory. As i this is a complicated problem to solve due to the need of knowledge about θk+j and dik+j , we

propose an iterative method to derive these parameters from previous batches. We call this strategy model learning predictive control.

11

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

Model Learning Predictive Control (ML-MPC) In this paper, we propose a method where we update the model used for control, such that the MPC controller performance can be improved over batches. To this end, we use the available data from the previous batch to (re-)estimate the model parameters, as illustrated in Figure 2. yk

uk

i=0 k

1 θˆk+j dˆ1k+j estimation

MPC

i=1

Batch

MPC

k i=1

i θˆk+j dˆik+j estimation

i−1

i=0

i−1

Batch

Figure 2: Illustration of the proposed control structure

In this illustration, we can see two time directions, the first one being the in batch time k, and the second being the sequential batches denoted by i. After the initial batch, data is 1 used to get estimations for θk+j , d1k+j for the next batch. The MPC controller then computes

the optimal inputs for the next batch using these updated parameters. And this sequence is then further repeated as described step by step in Algorithm 2.

12

ACS Paragon Plus Environment

Page 13 of 42 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

Algorithm 2 Model Learning Predictive Control Algorithm 1: procedure Model update process 2: Define the LPV model structure; 3: Set the initial parameter θ1 4: Set the initial disturbance vector d1k + j 5: Set i = 1; 6: Run a single batch, with an MPC controller using the LPV model; 7: Save all the available data, such as the trajectories of inputs ui−1 , outputs yi−1 and error ei−1 ; 8: Set i = i + 1 9: Perform the estimation to obtain the new θi and dik ; 10: Update the MPC model with these new parameters and disturbance; 11: Return to step 6. 12: end procedure

This method is shown in the control block diagram in Figure 3.

Memory

Estimation

P

MPC

Figure 3: Model Learning Predictive Control scheme.

Controller formulation As explained before we utilize an MPC controller with an LPV model in this ML-MPC method. Additionally, we include two integrators by the use of both δ and ∆ in the MPC formulation as described by Equation (7).

13

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

min

∆δuik+j

kxik+N k2P

+

N −1 X

Page 14 of 42

keik+j k2Q + k∆δuik+j k2R

j=0

subject to: i i xik+j+1 = A(θˆk+j )xik+j + B(θˆk+j )uik+j + Bd dˆik+j i = Cxik+j + Duik+j yk+j

(7)

i i eik+j = rk+j − yk+j

uik+j = uik+j−1 + δuik+j ∆δuik+j = δuik+j − δui−1 k+j i−1 ∆δuik+j = uik+j − uik+j−1 − ui−1 k+j + uk+j−1

xik+j ∈ X, uik+j ∈ U, xik+N ∈ Ω i In this problem, we include the estimation of the trajectory of the parameter θˆk+j , and the

time-varying disturbance dˆik+j , which are determined from previous batches. In the learning method proposed in this paper, the parameters and disturbances can be estimated using multiple methods. Particularly, we propose three different methods inspired in the MPC for LPV systems approaches that interpret the available data in different ways, such methods are presented in the next section. Remark 1 Note that optimization problem (7) is equivalent to the traditional IL-MPC when the parameter θk+j is constant during the batch. This fact shows that IL-MPC is a particular case of the method proposed in this paper (ML-MPC), also the performance of the controller given by (7) is at least the same of (5)

Parameter and Disturbance Estimation Methods We have considered three different methods for estimating the model parameters and disturbances. The first one is inspired by the Quasi-LPV approach. The second one is formulating

14

ACS Paragon Plus Environment

Page 15 of 42 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

the estimation problem as an identification problem, while, the third method is inspired by adaptive control in the batch direction. All these methods are explained in the following.

Method 1 In the case of chemical processes it is common to have constitutive relations. Therefore, it makes sense to approach the estimation problems in a similar fashion as seen in. 31 In this method, we use the prior knowledge of the system (constitutive relations) to define a function of the parameters and disturbances. We can describe this method by defining the following set of equations for θk+j and dk+j : i−1 i θˆk+j = fθ (xi−1 k+j , uk+j )

(8)

i−1 dˆik+j = fd (xi−1 k+j , uk+j )

To illustrate this method, let us consider a simple batch reactor given by the following set of nonlinear equations: 2 UA ∆HV dT =− (T − Tj ) − k0 e−E/RT CA2 dt M Cp M Cp dCA = −k0 e−E/RT CA2 dt Be defining x1 = T , x2 = CA , u = Tj , θ =

UA M Cp (T )

(9)

(T ) and d = − ∆HV k e−E/RT CA2 the model M Cp (T ) 0

can be written as:

x˙ = −θ(x − u) + d

(10)

where θ and d can be calculated by means of the following functions, UA M Cp (x) ∆HV (x) −E/Rx1 2 k0 e x2 d = fd (x, u) = − M Cp (x) θ = fθ (x, u) =

Note that, the parameters and disturbances can be updated in the ith batch by means of 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 42

fθ (·) and fd (·). More detail about this example are given in the next section (Illustrative Example). Naturally, in some cases, it is not possible to have such function for both the parameters and the disturbance, however, the method can still be used if we have a function only for either the parameters or the disturbance. An advantage of Method 1 is the ability to estimate a complete trajectory for time-varying parameters. However, it requires additional measurements besides the available input-output data, namely state information. Finally, the original optimization problem can be written as,

min

∆δuik+j

kxik+N k2P

+

N −1 X

keik+j k2Q + k∆δuik+j k2R

j=0

subject to: i i xik+j+1 = A(θk+j )xik+j + B(θk+j )uik+j + Bd dik+j i yk+j = Cxik+j + Duik+j i i eik+j = rk+j − yk+j

(11)

uik+j = uik+j−1 + δuik+j ∆δuik+j = δuik+j − δui−1 k+j xik+j ∈ X, uik+j ∈ U, xik+N ∈ Ω i−1 i θk+j = fθ (xi−1 k+j , uk+j ) ∈ Θ i−1 dik+j = fd (xi−1 k+j , uk+j ) ∈ D

Method 2 In this method we use an approach similar to the data driven identification as seen in Forgione et al. 38 . We use an optimization problem to estimate the parameters and disturbance from the input and output data. The cost function of this problem is user defined, and depends on the a priori knowledge the user has about the process.This problem in general can be described by:

16

ACS Paragon Plus Environment

Page 17 of 42 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

 i θˆk+j

 dˆik+j

= arg

min

θk+j ∈θ,dk+j ∈D

i−1 J(xi−1 k+j , uk+j )

subject to: i−1 i−1 i i i xi−1 k+j+1 = A(θk+j )xk+j + B(θk+j )uk+j + Bd dk+j

(12)

i−1 i−1 yk+j = Cxi−1 k+j + Duk+j i dik+j ∈ D, θk+j ∈Θ

where J(·) is usually the typical quadratic cost function used in Least Squares problems. Note that, problem 12 can be converted into a Quadratic Programming problem. This is a significant advantage compared with the proposal of 38 since both the parameter and disturbance can be estimated by solving a single optimization problem. Despite its advantages, this method is suitable for real time applications only if the parameter θk+j does not vary in time (θk+j = θ). Otherwise, the problem (12) is a nonlinear optimization problem, which is usually more difficult to solve. Another problem with this method is that by using inputoutput data only, the optimal solution could be affected by the noise in the system resulting in undesirable behavior by the controller. On the other hand it can happen that the number of parameters and disturbances are larger than the size of the input-output data. These results in an ill posed problem. These issues are analyzed in the illustrative example presented in the next section (Illustrative Example). Finally, with Method 2, the original optimization problem does not change its formulation.

Method 3 This third method contains some similarities to the original ILC method, where the error trajectory of the batch is utilized to estimate the parameters and disturbance. Method 3 is inspired in adaptive control strategies in the presence of unknown input gain. To this end, a projection of the error as described in Hovakimyan et al. 39 is used to update the parameter 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 42

and disturbance, mathematically the parameters and disturbances can be calculated as, i−1 i−1 i θˆk+j = θˆk+j + Γθ Proj(θˆk+j , −xi−1 ei−1 k+j ¯ k+j )

(13)

ˆi−1 ei−1 ) dˆik+j = dˆi−1 k+j + Γd Proj(dk+j , ¯ k+j where, ¯ e is the difference between the output data and the predicted output with the input data of the previous batch. Furthermore, Proj(·, ·) is the projection operator that ensures boundedness of the parametric estimates by definition, and it is defined as: 39 Proj(θ, y) ,    y, if f (θ) < 0    y, if f (θ) ≥ 0 and 5 f T y ≤ 0   D E   5f  y − 5f , y f (θ), if f (θ) ≥ 0 and 5 f T y > 0 k5f k k5f k

f (θ) ,

(14)

T (θ + 1)θT θ − θmax θmax T θ θmax θmax

Additional information of the projection operator could be found in more detail in Hovakimyan et al. 39 . Similar to Method 1, the parameter and disturbance can be estimated separately, as long as they can both be derived by a good projection of the state and/or error data. Because different projections of the same data points can give information about both the parameters and disturbance, this method is suitable for time-varying parameters.

Illustrative Example: Nonlinear batch reactor We apply the estimation methods described in the previous section to a relatively simple example of a nonlinear batch reactor. This example has also been used for the evaluation of an IL-MPC method in Oh and Lee 2 . In the reactor, an exothermic reaction takes places to form the required product C by means of the reaction A + B ⇒ C. To achieve the desired product quality, a specific temperature profile should be tracked. The jacket temperature Tj is used to achieve this goal. The dynamics of the reactor is described by the set of differential 18

ACS Paragon Plus Environment

Page 19 of 42 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

equations presented in Equation (15). UA ∆HV dT =− (T − Tj ) − k0 e−E/RT CA2 dt M Cp M Cp dCA = −k0 e−E/RT CA2 dt

(15)

where T is the reactor temperature, and Tj is the temperature of the heating/cooling jacket. The parameters of the process dynamics are shown in Table 1. Table 1: Nonlinear batch reactor Parameters. Parameter Values Units UA 0.09 l/min M Cp ∆HV −1.64 Kl/mol M Cp E/R 13, 550 K k0 2.53 l/mol min CA (0) 0.9 mol/l ◦ T (0) 20 C

As we stated in a previous sections, the nonlinear reactor model equations (15) can be written as,

x˙ = −θ(x − u) + d

(16)

In this example, θ is assumed as a constant parameter and d a time-varying disturbance. Although θ is considered constant during the batch, this parameter could change over batches due to fouling, as it includes the heat transfer coefficient of the heating jacket. In comparing the three methods, we consider the same MPC settings. We set a desired temperature profile as shown in Figure 4. 2

19

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

Figure 4: Reactor temperature: reference profile. The MPC is tuned with the control horizon N = 30, the cost matrices Q = 106 , R = 106 and the inputs are in the sets u ∈ [10, 50], δu ∈ [−1, 1] and ∆u ∈ [−10, 10]. For initialization of this reactor problem we apply as the initial parameter θˆ0 = 3 and initial  > 0 disturbance vector d = 0 · · · 0 . With these parameter and disturbance values, the performance of the MPC in tracking of the reference trajectory is not satisfactory. This is not surprising due to the plant-model mismatch, as the actual value of θ is 0.09 and the disturbance trajectory is shown in Figure 5.

Figure 5: Reaction rate profile.

In the following, we study the performance of these methods in simulation on this batch reactor.

Method 1 For this specific example we do not have a function to estimate θˆ as in Equation 8. Therefore, we assume that we are able to estimate the parameter perfectly such that θˆ = θ. Applying this parameter learning with the ML-MPC algorithm results in the reference tracking improvement as shown in Figure 6 20

ACS Paragon Plus Environment

Page 20 of 42

Page 21 of 42 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

Figure 6: Input and Outputs with method 1.

To improve the reference tracking further, Method 1 can be implemented on the batch reactor to estimate the disturbance vector based on the measurable temperature T and concentration CA in the reactor. From the differential equation of the reactor, we know that the reaction rate can be described as the disturbance as shown in Equation (17). i−1 ∆HV 2 dˆik = fd (Tk , CAk ) = − k0 e−E/RTk CA i−1 , k M Cp

(17)

which results in the improvement of reference tracking, from t = 0 min till t = 6 min, as shown in Figure 7. In this figure, the reference tracking improvement is clearly seen at the beginning of the batch. This is expected because the reaction rate has the highest amplitude in the start of the batch.

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 42

Figure 7: Input and Outputs with method 1: From t = 0 min till t = 6 min .

Method 2 In this method, we use Euler approximation to represent the differential equation as a difference equation. This derivation is explained as follows, x˙ ≈

xk+1 − xk = θˆ (T j − T ) +dˆk , ∀k = 0, · · · , N | {z } ∆t | {z } αk

βk

 ˆ θ   I 0 0   ˆ   d1  .. . .  . 0  .  .    ..  I ··· I   {z } dˆN a | {z } 









 β1   α1  .   .  ..  =  ..       βN αN | {z } | b

(18)

γ

which we can calculate γ using the following least squares problem after substituting the input and output data in a and b,

minkaγ − bk, γ

22

ACS Paragon Plus Environment

(19)

Page 23 of 42 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

This optimization problem, however, has two issues. First of all, it is ill-posed as there are more parameters than the rank of the matrix a. The second problem is that any noise in the input or output affects the estimation of the disturbance as shown in Figure (8).

Figure 8: Least squares estimation With noise: T − T (0) If we evaluate the result of the dˆ for this estimation, we see that the noise in the estimated dˆ is much larger than actual d as shown in Figure 9.

Figure 9: Least squares estimation with noise: Disturbance d.

A solution to overcome such problems is utilizing a regularization method. In the case of the batch reactor, we know that the reaction rate converges to zero, and therefore the L1 or "lasso" regularization is a good option. The least squares problem then becomes:

minkaγ − bk + λkγk1 , γ

(20)

where λ is a tuning parameter with λ = 100. The solution of this problem for the first batch is shown in Figure 10. We can notice how the regularization term filters the solution from the system noise.

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

Figure 10: L1 least squares estimation: T − T (0). The result of the estimation for dˆ is shown in Figure 11. In this figure, we see that there is ˆ These fluctuations can, however, still some fluctuations in the least squares estimation of d. ˆ This is be filtered using a lowpass filter to obtain a good estimation of the trajectory d. shown by the "filtered L1 LS" result, compared with d.

Figure 11: L1 least squares estimation: Disturbance d. Now we can implement the estimated θˆ in the controller model for the next batch. The result for the updated parameter is shown in Figure 12. In this result, we see that the estimation of θˆ results in a direct improvement step towards the reference tracking performance. By implementing the estimated dˆ for the next batch we see an improvement in the reference tracking performance of the controller for the start of the batch as shown by Figure 13.

24

ACS Paragon Plus Environment

Page 24 of 42

Page 25 of 42 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

Figure 12: Input and Outputs with method 2.

Figure 13: Input and Outputs with method 2: From t = 0 min till t = 6 min .

Method 3 Method 3 is implemented using the learning function in Equation (21) based on the error. i−1 i−1 i θˆk+j = θˆk+j + Γθ Proj(θˆk+j , −xi−1 ei−1 k+j ¯ k+j )

dik

=

di−1 k

+

λei−1 k , 25

ACS Paragon Plus Environment

(21)

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

With λ = 0.25 we obtain the resulted presented in Figure 14 and Figure 15.

Figure 14: Input and Outputs with method 1.

Figure 15: Input and Outputs with method 3: From t = 0 min till t = 6 min .

An interesting observation of these results is related to the stability of this learning method. In the original ILC algorithm using Equation 3 it is known that high frequent input oscillations can reduce the controller performance over iterations, whereas in the case of this 26

ACS Paragon Plus Environment

Page 26 of 42

Page 27 of 42 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

implementation in the MPC controller the δ formulation prevents these high frequent fluctuations by the penalty on the input change δu, and therefore helps to stabilize this learning algorithm. It is important to highlight that some methods are better in the disturbance learning, and others better in parameter learning, and therefore we will compare both aspects separately for the specific batch reactor. We compare the reference tracking performance over batches by taking the kek2 as the performance measure. In Figure 16, we present the performance of each method for the learning of the parameter. In this figure, 10 batches are displayed where the parameter of the plant is changed from θ = 0.09 to θ = 0.05 between the sixth and the seventh batch. This causes an increase in the error during the seventh batch.

Figure 16: kek2 convergence for parameter learning using different methods.

From this plot we can see that Method 1 and Method 2 have an almost equal performance, meaning that the least squares optimization is well defined for this reactor problem as Method 1 is able to learn the parameter and disturbance exactly. Method 3 has a clearly slower learning behavior for this application. When we compare these results to the paper of Oh and Lee, 2 we see that their IL-MPC algorithm converges in approximately five to six batches, whereas ML-MPC converges in approximately two batches, such that the learning speed over batches is more than twice as good compared to the work of Oh and Lee 2 . Of course this also depends on the quality of the estimation method used. For the learning of the disturbance, we can observe that Method 1 gives a better performance over Method 2 as seen in Figure 17. This can be explained by the remaining mismatch in 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

disturbance estimation as shown previously in Figure 11, whereas Method 1 estimates the exact reaction rate. Method 3, in the way it is implemented, includes errors caused by the controller in the disturbance trajectory. And with this property, it is able to improve the controller performance over batches compared to the other two methods. The drawback of the method 3 is that the disturbance trajectory can not be interpreted as the reaction rate.

Figure 17: kek2 convergence for disturbance learning using different methods.

From the results we can conclude that the proposed methods offer good controller improvement over batches, and thus being able to take into account changes in the system. The results, however, also show that each estimation method have its own strengths and limitations. To implement the ML-MPC it remains important to have the good knowledge of the process in order to select the right estimation method for each application. We have also observed that for the batch reactor Method 1 would be the best option for implementation, under the assumption that both temperature and concentration can be measured. If the concentration measurement in the plant would not be possible, Method 2 offers a good alternative. In the next section, we apply the ML-MPC algorithm to a more advanced process in a reactive batch distillation column using only one of the learning methods.

Reactive Batch Distillation (RBD) Processes A reactive batch distillation column is a unit operation utilized generally for chemical processes involving equilibrium reactions. With the use of a distillation column, light products 28

ACS Paragon Plus Environment

Page 28 of 42

Page 29 of 42 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

are removed so that the reaction proceeds in the products direction. A general form of RBD is illustrated in Figure 18. The reactor of such RBD is enclosed by a jacket which can heat or cool the reactor. On top of the reactor is the distillation column with multiple trays. The vapor rising from the reactor and the liquid flow (reflux) coming from the top of the column get into contact at the tray with liquid holdup.

Figure 18: Illustration of a reactive batch distillation column

Rigorous Dynamic Model of RBD A general rigorous dynamic model for the i-th stage of the RBD is formulated and given by the following set of equations, dnli,j = Li−1 xi−1,j − Li xi,j − ξi,j + nLi ri,j dt dngi,j = Vi+1 yi+1,j − Vi yi,j + ξi,j dt NR X dnT,i hi = Li−1 hi−1 − Li hi − Vi Hi + Vi+1 Hi+1 + nLi (−∆Hi,l )ri,l + Qext,i dt l=1 ∀i = 1, . . . , NT , ∀j = 1, . . . , NC 29

ACS Paragon Plus Environment

(22)

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

where nli,j , ngi,j , ri,j ∈ R are the number of moles in the liquid and gas phase, and the reaction rate of j-th component. Vi ∈ R and Li ∈ R are the outlet gas and liquid flow, and Vi+1 ∈ R and Li−1 ∈ R are the inlet gas and liquid flow from the (i+1)-th and (i−1)-th tray. Furthermore, hi , Hi , −∆Hi,j , Qext,i ∈ R are the liquid, gas and reaction molar enthalpies, and the external heat flow respectively. Besides, nT,i is the total number of moles, and is given by nT,i = nLi + nGi , where nLi and nGi are the total liquid and gas molar holdup. The liquid mole fraction xi,j and the vapor mole fraction yi,j are given by,

nLi =

NC X

nlj,i , nGi =

j=1

xi,j

NC X

ngj,i ,

j=1

ng nl = i,j , yi,j = i,j , nLi nGi

∀i = 1, . . . , NC , ∀j = 1, . . . , NT , The model given by the set of equations (22) is subject to the following conditions presented in Table 2. Table 2: Conditions of the rigorous dynamic modeling of RBD processes Li−1 Li Vi+1 i=1 LD 0 0 i=2 0 Li + LD Vi+1 i = NT Li−1 0 0

Vi 0 0 Vi

ξi,j xi−1,j ri,j Qext,i 0 xi+1,j 0 0 ξi,j 0 0 0 ξi,j xi−1,j rNT ,j Qext,NT

LPV Representation of RBD Column The extent decomposition for RBD processes can be achieved finding a linear transformation Tli , and Tgi , similar to the one computed in Amrhein et al. 40 , such that the system (22) can be represented in terms of new independent states for the gas, and the liquid phase as follows:

30

ACS Paragon Plus Environment

Page 30 of 42

Page 31 of 42 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

Liquid phase: 



  T  xr,li  T1e ,li      x  = T T  nl  in,li   2e ,li  i     T λli T3e ,li | {z }

(23)

Te,li

Gas phase:     T xin,gi  T1e ,gli   ngi  = T2Te ,gi λgi | {z }

(24)

Te,gi

where xr,li ∈ RNR is the extent of reaction, xin,{l,g}i ∈ Rp+1 are the extent of inlet flow, xinv,li ∈ RNC −NR −p−1 , and xinv,gi ∈ RNC −p−1 are the extent of reaction and inlet flow invariants. Finally, we obtain the following LPV representation for the i-th tray,

        0 0   uli   Bd,li  0 0   xe,li  Be,li  x˙ e,li  Ae,li (θli )                  x˙  =  0 Ae,gi (θgi ) 0   xe,gi  +  0 Be,gi 0   ugi  +  Bd,gi  di  e,gi             ˙ Bd,Ti Qext,i Ti αLi α ¯ Vi α Q i Ti 0 0 θTi      nli  Cli 0 xe,li     = 0 1 Ti Ti 



(25) Details about the extent transformation for the RBD can be found in Marquez-Ruiz et al. 9 .

Case Study: Synthesis of unsaturated polyester in a RBD Consider an industrial reactive batch distillation process for the synthesis of unsaturated polyester from maleic anhydride and propylene glycol. This process involves four types of reactions. 41 First Maleic anhydride reacts with propylene glycol and produces a maleic 31

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

acid end group and a propylene glycol end group with an ester bridge. This is a fast and exothermic reaction (∆H = −40 kJ/mol). Esterification proceeds by the reaction of different acid and alcohol end groups to form new ester bridges and water, or by reaction of a glycol hydroxyl group with an acid end group to form an ester bridge and water. The double bond in maleic anhydride can be isomerized or saturated. Saturation of the double bond causes crosslinking in the polymer, and approximately 10-20 % of the double bonds are saturated in the preparation of the polyester. The three reactions, esterification, isomerization and saturation form a network of nine reactions. These reactions are summarized in Table 3. Table 3: Basic Reactions in the Polyesterification of Unsaturated Carboxylic Acids with Diols 41 Reaction 1 2 3 4 5 6 7 8 9

Stoichiometry 0 + H2 O RCOOH1D + R0 OH  RCOOR1D 0 0 RCOOH2D + R OH  RCOOR2D + H2 O RCOOHS + R0 OH  RCOORS0 + H2 O RCOOH1D  RCOOH2D 0 0 RCOOR1D  RCOOR2D RCOOH1D + 0.5R0 OH  RCOOHS RCOOH2D + 0.5R0 OH  RCOOHS 0 + 0.5R0 OH  RCOORS0 RCOOR1D 0 + 0.5R0 OH  RCOORS0 RCOOR2D

The RBD consists of 6 stages (NT = 6) with only 3 internal trays. Stage equilibrium is assumed in the process with NRTL activity coefficient model for the liquid phase; molar vapor holdup is negligible with respect to the molar liquid holdup, an isothermal total condenser is assumed, and the reaction is limited to the reboiler (Reactor). Only four species are assumed to evaporate (p = 4). In total, the chemical process has 8 species (NC = 8) with 9 reactions, but 5 independent (NR = 5), i.e., the stoichiometric matrix N T ∈ R8×9 but has rank 5. The following are mole and temperature initial conditions in the reactor: nlNT ,RCOOH1D (0) = 20 kmol, nlN

0 T ,R OH

(0) = 20 kmol, and TNT (0) = 373 K.

32

ACS Paragon Plus Environment

Page 32 of 42

Page 33 of 42 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

Control Objective In the RBD, the objective is to steer the temperature of the reactor to a reference temperature. This reference is a predefined trajectory to obtain optimal product qualities from a batch in the RBD. This reference profile is shown in Figure 19. With this temperature profile, the water in the reactor vessel is eliminated as shown by the bottom plot in Figure 19.

Figure 19: The optimal reference profiles for the temperature and the water holdup in the reactor.

We have two inputs, the heat-flow to the reactor and the reflux of distilled liquid in the top of the reactor for the control of the reactor temperature. In a recent paper by 42 , several conventional control strategies have been studied for this system. In this work, we consider a model based control strategy. The tracking the temperature reference is difficult since when the boiling point of the products in the reactor is reached, the temperature cannot be increased anymore by the heat flow input. This behaviour can be observed by the optimal profile in Figure 19 where the temperature only increases when the temperature is below the boiling point. The reflux input of the RBD has also interesting dynamic behavior that has to be taken into account. The reflux input only attenuates the temperature in the reactor when 33

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

there is already some product distilled. When this is the case, the reflux cools the reactor when the distillation tray above the reactor has a lower temperature than the reactor, while it would heat the reactor when the tray above the reactor has a higher temperature than the reactor.

Simulation results In this work, we will consider the control of the reference temperature by manipulating the heat-flow only. We assume that the optimal input trajectory of the reflux is used for the input, which can be seen in Figure 20

Figure 20: Optimal reference profile for the input. By considering the manipulation of the heat flow in the reactor jacket only, the LPV model can be reduced to: dTi = −θT i Ti + αQi Qext,i + αLi Li−1 + βi ri + κi Vi | {z } |{z} |{z} | {z } | {z } dt A(θ)

x

B(α)

u

(26)

d

First, we apply the ML-MPC with a fixed parameter θ for the batch time. In this case we learn the disturbance by using Method 1, where we assume the liquid and vapor flows and the mass in the reactor measurable. With these measurements the disturbance trajectory can be calculated. The results of learning the parameter and disturbance trajectories for the controller are shown in Figure 21. In these results, we see that due to a small cooling action around 5, 000 seconds in the first batch, the temperature in the reactor drops and it is not able to reach the reference anymore. 34

ACS Paragon Plus Environment

Page 35 of 42 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

Figure 21: Reference tracking after one batch of learning the disturbance trajectory using Method 1. At the same time, the heat flow input is high, because the boiling temperature of the reactor is decreased. This boiling point decrease has been caused by the reflux that returns water back to the reactor. In the second batch, we learn the disturbance caused by the reflux, therefore, the controller cools the reactor slightly less around 5, 000 seconds, such that the temperature in the reactor does not decrease, and the reference tracking performance is improved for the next 10, 000 seconds. While this method improves the performance of RBD operation, it does not consider changes in the plant as the parameters are not re-estimated. We, therefore, implement another controller with the time-varying parameter. We start in the initial batch with the poorly chosen parameter constant for the full batch. In the second batch the varying trajectory is estimated. The results for this application are shown in Figure 22.

35

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

Figure 22: Reference tracking after one batch of learning the parameter and disturbance trajectory using estimation Method 1.

While in the first batch the temperature is only increased by the exothermic reaction in the reactor as the heat flow input is zero, in the second batch the controller is able to track the reference better. Similar to the previous controller implementation, cooling of the reactor around 5, 000 seconds causes a temperature drop in the reactor in the second batch. And by learning this undesired behavior, the controller for the following batches then increases the temperature around 5, 000 seconds to be able to remain at the reference temperature. We have also tested the ML-MPC algorithm for 15 batches, with changes in the initial conditions of the plant at the eleventh and fourteenth batch. We obtain the convergence of the Euclidean norm of the error over batches as seen in Figure 23.

Figure 23: Convergence of the Euclidean norm over batches, with changes in the initial conditions at batch 11 and 14.

36

ACS Paragon Plus Environment

Page 36 of 42

Page 37 of 42 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

Conclusions In this paper we have proposed a Model Learning Predictive Control (ML-MPC) method, based on the repetitive behavior of the batch processes. To this end, the LPV model used in the controller is updated using information from the previous batch. Inspired by the MPC for LPV systems in literature, three different application-dependant options to estimate the parameters and disturbance of the model have been proposed and compared in simulation on a nonlinear batch reactor. The best controllers are able to converge to their best performance within approximately two batches, which is better than the IL-MPC method in Oh and Lee 2 . Finally, by applying the ML-MPC on a reactive batch distillation column, we have shown the ability to adjust the controller to complex nonlinear behaviors. The parameter and disturbance converge to their actual values in limited number of batches. In the case of changes in the initial conditions the controller is also able to adapt quickly.

References (1) Bonvin, D. Control and optimization of batch processes. IEEE control systems 2006, 26, 34–45. (2) Oh, S.-K.; Lee, J. M. Iterative learning model predictive control for constrained multivariable control of batch processes. Computers & Chemical Engineering 2016, 93, 284–292. (3) Bamieh, B.; Giarre, L. Identification of linear parameter varying models. International journal of robust and nonlinear control 2002, 12, 841–853. (4) Bruzelius, F. Linear parameter-varying systems-an approach to gain scheduling; Chalmers University of Technology, 2004.

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

(5) Casella, F.; Lovera, M. LPV/LFT modelling and identification: overview, synergies and a case study. Computer-Aided Control Systems, 2008. CACSD 2008. IEEE International Conference on. 2008; pp 852–857. (6) Tóth, R. Modeling and identification of linear parameter-varying systems; Springer, 2010; Vol. 403. (7) Marquez-Ruiz, A.; Mendez-Blanco, C.; Özkan, L. Control of homogeneous reaction systems using extent-based LPV models. IFAC-PapersOnLine 2018, 51, 548–553. (8) Marquez-Ruiz, A.; Mendez-Blanco, C.; Porru, M.; Özkan, L. Computer Aided Chemical Engineering; Elsevier, 2018; Vol. 44; pp 583–588. (9) Marquez-Ruiz, A.; Méndez-Blanco, C. S.; Özkan, L. Modeling of reactive batch distillation processes for control. Computers & Chemical Engineering 2019, 121, 86–98. (10) Arimoto, S.; Kawamura, S.; Miyazaki, F. Bettering operation of robots by learning. Journal of Field Robotics 1984, 1, 123–140. (11) Wang, Y.; Gao, F.; Doyle, F. J. Survey on iterative learning control, repetitive control, and run-to-run control. Journal of Process Control 2009, 19, 1589–1600. (12) Lee, J. H.; Lee, K. S. Iterative learning control applied to batch processes: An overview. Control Engineering Practice 2007, 15, 1306–1318. (13) Adam, E. J.; González, A. H. Frontiers in Advanced Control Systems; InTech, 2012. (14) Shamma, J. S.; Athans, M. Analysis of gain scheduled control for nonlinear plants. IEEE Transactions on Automatic Control 1990, 35, 898–907. (15) Qin, S. J.; Badgwell, T. A. An overview of nonlinear model predictive control applications. Nonlinear model predictive control 2000, 369–392.

38

ACS Paragon Plus Environment

Page 38 of 42

Page 39 of 42 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

(16) Cueli, J.; Bordons, C. Iterative nonlinear control of a semibatch reactor. Stability analysis. Decision and Control, 2005 and 2005 European Control Conference. CDC-ECC’05. 44th IEEE Conference on. 2005; pp 2071–2076. (17) Cueli, J.; Bordons, C. Application of iterative nonlinear model predictive control to a batch pilot reactor. IFAC Proceedings Volumes 2005, 38, 63–68. (18) Cueli, J. R.; Bordons, C. Iterative nonlinear model predictive control. Stability, robustness and applications. Control Engineering Practice 2008, 16, 1023–1034. (19) Liu, X.; Kong, X. Nonlinear fuzzy model predictive iterative learning control for drumtype boiler–turbine system. Journal of Process Control 2013, 23, 1023–1040. (20) Shen, D.; Han, J.; Wang, Y. Convergence analysis of ILC input sequence for underdetermined linear systems. Science China Information Sciences 2017, 60, 099201. (21) Wang, Y.; Zhang, J.; Zeng, F.; Wang, N.; Chen, X.; Zhang, B.; Zhao, D.; Yang, W.; Cobelli, C. Learning can improve the blood glucose control performance for type 1 diabetes mellitus. Diabetes technology & therapeutics 2017, 19, 41–48. (22) Nagy, Z. K.; Braatz, R. D. Robust nonlinear model predictive control of batch processes. AIChE Journal 2003, 49, 1776–1786. (23) Lu, Y.; Arkun, Y. Quasi-min-max MPC algorithms for LPV systems. Automatica 2000, 36, 527–540. (24) Abbas, H. S.; Tóth, R.; Meskin, N.; Mohammadpour, J.; Hanema, J. A robust MPC for input-output LPV models. IEEE Transactions on Automatic Control 2016, 61, 4183–4188. (25) Hanema, J.; Tóth, R.; Lazar, M. Tube-based anticipative model predictive control for linear parameter-varying systems. Decision and Control (CDC), 2016 IEEE 55th Conference on. 2016; pp 1458–1463. 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

(26) Besselmann, T.; Lofberg, J.; Morari, M. Explicit MPC for LPV systems: Stability and optimality. IEEE Transactions on Automatic Control 2012, 57, 2322–2332. (27) Lofberg, J. Approximations of closed-loop minimax MPC. Decision and Control, 2003. Proceedings. 42nd IEEE Conference on. 2003; pp 1438–1442. (28) Casavola, A.; Famularo, D.; Franze, G. A feedback min-max MPC algorithm for LPV systems subject to bounded rates of change of parameters. IEEE Transactions on Automatic Control 2002, 47, 1147–1153. (29) Bemporad, A.; Borrelli, F.; Morari, M. Min-max control of constrained uncertain discrete-time linear systems. IEEE Transactions on automatic control 2003, 48, 1600– 1606. (30) Hovland, S.; Gravdahl, J. T.; Willcox, K. E. Explicit model predictive control for largescale systems via model reduction. Journal of guidance, control, and dynamics 2008, 31, 918–926. (31) Cisneros, P. S.; Werner, H. Parameter-dependent stability conditions for quasi-LPV Model Predictive Control. American Control Conference (ACC), 2017. 2017; pp 5032– 5037. (32) Cisneros, P. S.; Voss, S.; Werner, H. Efficient nonlinear model predictive control via quasi-lpv representation. Decision and Control (CDC), 2016 IEEE 55th Conference on. 2016; pp 3216–3221. (33) Lee, J. H.; Lee, K. S.; Kim, W. C. Model-based iterative learning control with a quadratic criterion for time-varying linear systems. Automatica 2000, 36, 641–657. (34) Moore, K. L.; Chen, Y.; Ahn, H.-S. Iterative learning control: A tutorial and big picture view. Decision and Control, 2006 45th IEEE Conference on. 2006; pp 2352–2357. (35) Oomen, T. Lecture notes. 40

ACS Paragon Plus Environment

Page 40 of 42

Page 41 of 42 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

(36) Bristow, D. A.; Tharayil, M.; Alleyne, A. G. A survey of iterative learning control. IEEE Control Systems 2006, 26, 96–114. (37) Xu, J.-X.; Wang, X.-W.; Heng, L. T. Analysis of continuous iterative learning control systems using current cycle feedback. American Control Conference, Proceedings of the 1995. 1995; pp 4221–4225. (38) Forgione, M.; Bombois, X.; Van den Hof, P. M. Data-driven model improvement for model-based control. Automatica 2015, 52, 118–124. (39) Hovakimyan, N.; Cao, C.; Kharisov, E.; Xargay, E.; Gregory, I. M. L1 Adaptive Control for Safety-Critical Systems. IEEE Control Systems 2011, 31, 54–104. (40) Amrhein, M.; Bhatt, N.; Srinivasan, B.; Bonvin, D. Extents of reaction and flow for homogeneous reaction systems with inlet and outlet streams. AIChE journal 2010, 56, 2873–2886. (41) Lehtonen, J.; Salmi, T.; Immonen, K.; Paatero, E.; Nyholm, P. Kinetic model for the homogeneously catalyzed polyesterification of dicarboxylic acids with diols. Industrial & engineering chemistry research 1996, 35, 3951–3963. (42) Marquez-Ruiz, A.; Ludlage, J.; Özkan, L. Computer Aided Chemical Engineering; Elsevier, 2018; Vol. 44; pp 577–582.

41

ACS Paragon Plus Environment

1 Industrial & EngineeringBATCH Chemistry Research Page 42 of 42 MPC

1 2 3 4 5 6 7 BATCH

FIRST PRINCIPLE MODEL H R1

+

O

R2

H

O

C

C

H O O H

R2

H

O

C

C

O

H

R1

FIRST PRINCIPLE MODEL

LEARNING REACTION AND SEPARATION

ACS Paragon Plus Environment INFORMATION TO BATCH 2

BATCH