Optimal feed rate policy for systems with two reactions

IFAC

Copyright © IFAC Dynamics and Control of Process Systems, Jejudo Island, Korea, 200 I

c:

0

C>

Publications www.elsevier.comllocatelifac

OPTIMAL FEED RATE POLICY FOR SYSTEMS WITH TWO REACTIONS B. Srinivasan t, O. Ubrich *, D. Bonvin t, F. Stoessel * t Institut d'Automatique, Ecole Poly technique Federale de

Lausanne, CH-l015 Lausanne, Switzerland ·Institut de Genie Chimique, Ecole Poly technique Federale de Lausanne, CH-l015 Lausanne, Switzerland Abstract: Optimization of batch processes operation is a natural choice for reducing production costs, improving product quality, and meeting safety and environmental regulations. The implementation of the optimal solution is straightforward if it is determined by the constraints of the optimization problem. This paper shows that such is the case for the majority of tw(}-reaction systems taking place in semi-batch reactors. Copyright © 2001lFAC Keywords: Dynamic optimization, Optimal control, Measurement-based optimization, Batch processes, Semi-batch reactors.

be implemented by tracking the active path constraints using on-line measurements (Visser et al., 2000) and the active terminal constraints using off-line measurements (Srinivasan et al., 2001a) . Whether or not the solution is determined by constraints depends on the compromises present. These compromises are classified in two categories: (a) the compromises intrinsic to the dynamic system, and (b) the compromises resulting from the formulation of the optimization problem. The first category corresponds to inputs having multiple conflicting influences on the states of the system. Since these influences act against each other, compromise values of the inputs need to be used for the sake of optimality. If there are no intrinsic compromises in the system, the optimal inputs are determined by the path constraints (input bounds and state constraints) . The compromises formulated in the optimization problem are achieved by proper sequencing of the various intervals and, more specifically, by the choice of the switching instants between them.

1. INTRODUCTION

Batch and semi-batch processes are of considerable importance in the batch chemical industry. A wide variety of specialty chemicals, pharmaceutical products, and polymers are manufactured in batch operations (Macchietto, 1998) . The operation of batch processes typically involves following recipes that have been developed in the laboratory. However, owing to differences in both equipment and scale, industrial production almost invariably necessitates modifications of these recipes in order to ensure productivity, quality, safety, and satisfaction of operational constraints (Wiederkehr, 1988) for which an optimization approach can be used. The operational decisions such as temperature and feed rate pr(}files are then determined from the solution to an optimization problem, where the objective is of economic nature and the various technical and operational constraints are considered explicitly. If the optimal solution is determined by the constraints of the problem (input bounds, state and terminal constraints), the optimal solution can

In this paper, conditions that guarantee the absence of compromise-seeking arcs in semi-batch

389

reactors with the feed rate used as manipulated input are presented. These conditions are based on the dynamic model only and not on the formulation of the optimization problem (terminal cost and constraints) . It will be shown that for the vast majority of semi-batch reaction systems with two reactions, the optimal feed rate is determined by the constraints of the problem.

Lagrange multipliers for the path constraints, and // 2 0 the Lagrange multipliers for the terminal constraints. The notation ab = ~~ is used.

2.3 Normal form of chemical reaction systems Consider a homogeneous, constant-density, semibatch chemical reaction system comprising S species and R independent reactions. Independent reactions are those which have both independent stoichiometries and independent kinetics (Srinivasan et al., 1998) . The system can be nonisothermal, in which case it it assumed that the reactor temperature is a manipulated varaible. The material balance equations and the continuity equation are described by:

The paper is organized as follows. Section 2 briefly reviews the formulation of the optimization problem and the reduction in dimensionality of chemical reaction systems. The main theoretical result is stated in Section 3, and its implication for two-reaction systems is discussed in Section 4. A catalogue of semi-batch reaction systems with two reactions is presented in Section 5. Three examples are provided in Section 6 to illustrate the theoretical developments, and conclusions are drawn in Section 7.

it = V K r n (n, V)

V=U

2.1 Optimization problem

A typical optimization problem that involves meeting certain specifications only at the end of the batch can be written as:

(1) :i; = F(x, u),

S(x, u)

~

0,

x(O) = Xo T(x(t f)

~

(3)

where J is the scalar performance index to be minimized, x the states with initial conditions Xo, u the scalar input, S the path constraints (which include state constraints and input bounds), T the terminal constraints, F the smooth vector field describing the dynamics of the system,

(

Ci;IT) cinl

n

IT

Cm

(7)

where the superscript + denotes the MoorePenrose pseudo inverse, takes the system (6) to the normal form (Srinivasan et al., 1998) :

~

Using Pontryagin 's Minimum Principle, the necessary conditions of optimality for (1 )-(3) can be obtained as:

/-LTS=O,

I _

-n-V I NTn

2.2 Necessary conditions

>.TFu+/-LTSu=O,

(6)

Let the rows of the matrix N E RSx(S-R-l) form the null space of [K Cin]T. Using singular value decomposition, N can be computed as the matrix of the right singular vectors corresponding to the zero singular values of [K Cin]T . Let I E R S be the null vector of [K N]T. The following transformation,

(2) 0

n(O) = no V(O) = Vo

where n is the S-dimensional vector of the number of moles, u the inlet volumetric flowrate, V the reactor volume, K the S x R stoichiometric matrix, rn the R-dimensional reaction rate vector, Cin the molar concentration of the inlet stream, no the initial number of moles, and Vo the initial volume. The terms V K rn and Cin U represent the effect of the reactions and the inlet stream on the number of moles, respectively. The molar concentrations are given by C = njV.

2. PRELIMINARIES

s.t .

+ Cin U

= V r,

(

= 0,

~

= 0, and V = u.

(8)

~ corresponds to the reaction variants, while ( and ~ are the reaction and flow invariants. The invariants can be removed and the reduced system

//TT=O , (4)

with

(9)

)..T=-)..TFx-/-LTSx , )..T(tf) = (
of dimension (R+ 1) considered for analysis. This model reduction, which uses no information on the kinetics except for the independence of reactions, eliminates the redundancy in n using information

where )..(t) =J 0 are the adjoints (Lagrange multipliers for the system equations), /-L(t) 2 0 the

390

present in the stoichiometric matrix. It is possible that the kinetics are such that the dimension of the space that is controllable from the input is less than (R + 1), thereby making further model reduction possible.

than 00 and yet the optimal solution be on the constraints.

3. SOLUTION ON THE BOUNDARY OF THE FEASIBLE REGION

4.1 Optimal feed rate for two-reaction systems

4. OPTIMAL SOLUTION FOR SYSTEMS WITH TWO REACTIONS

For two reactions the normal form (9) gives FT = [ V rl V r2 u] and FJ = [0 0 1] . Differentiation of Fu leads to:

The optimal input is obtained from the necessary conditions of optimality (4) . If the input does not appear explicitly in (4), the necessary conditions need to be differentiated with respect to time (Bryson, 1999) . If the input appears explicitly in one of the time derivatives AT Fu , the solution can be in the interior of the feasible region . Otherwise, the solution is necessarily on the boundary of the feasible region .

:t

where W is a vector field independent of u . Since u appears explicitly in the second differentiation , (7 = 2 provided the second term in 6.2 Fu is linearly independent of Fu and 6.Fu , i.e. ,

The above analysis can also be done independently of the adjoints by differentiating Fu along the trajectories of (2) using the following operator (Srinivasan et al., 2001 b):

8Fu 8F) 6.Fu = ( --F - -Fu 8x 8x

8Fu . + -u 8u

= 8(V rd 8 2(V r2) _ 8(V r2) 8 2 (V rd #- 0

c-

(10)

The above operator is also studied in the systems literature using tools of Lie algebra.

8V

8V2

8V

8V2

If C = 0, further differentiations of Fu would required. However, these can be avoided since the rank of M = [Fu 6.Fu 6.2Fu] (defined as p) provides information on whether or not there is a possibility of u appearing explicitly and independently. Three subcases are considered:

Let (7 be the minimum number of time differentiations of Fu required for the input u to appear explicitly with 8f:l.;uF" being independent of 6. j Fu , j = 0, 1, . .. , (7 - 1. (7 is referred to as the order of singularity in the literature (Palanki et al. , 1993) .

(1) p = 3 (everywhere in the operating domain) : Since full rank is already attained, u cannot appear independently of Fu, 6.Fu , and 6. 2 Fu in further differentiations and so (7 = 00 . (2) p = 2 (everywhere in the operating domain) : Since 6. 2 Fu is a linear combination of Fu and 6.Fu , the rank of M cannot be augmented by further differentiations of Fu . This indicates a redundancy in the states and, thus further model reduction is possible. As with the earlier case, since maximum rank is already attained, u cannot appear explicitly and independently in further differentiations of Fu and so (7 = 00 . (3) p = 2 on a manifold and p = 3 elsewhere: In the manifold where M loses rank, it is possible that u appears explicitly and independently with (7 > 2. Though such a case is theoretically possible, it was not encountered in the examples considered in this paper.

Theorem 1. If (7 = 00 , the optimal solution of (1)(3) is on the boundary of the admissible region . This theorem is a restatement and generalization of the results presented in (Palanki et al., 1993; Benthack, 1997) where the reader is referred to for the proof. (7 = 00 occurs if the system is feedback linearizable or differentially flat 1 . Intuitively, the system can be transformed into a chain of integrators, implying that the input cannot have multiple conflicting effects on the states of the system. So, the solution is determined by the constraints of the optimization problem (1)-(3), whatever the cost function and the constraints might be. On the other hand , if (7 < 00, depending on the cost function and the constraints, the solution may have intervals where the input is not determined by the constraints. However, it is interesting to note that the condition (7 = 00 is sufficient, but not necessary, for the optimal solution to be determined by problem constraints, i.e., (7 can be less

4.2 Influence of variations in temperature Since it was assumed that the temperature is a manipulated variable, the above results are valid

I For single input systems, the notions of feedback Iinearizability and differential flatness are equivalent.

391

#

Kinetic model

Description

1

Reversible reaction

(7

p

A+B~C

00

2

A+B~ 2A

00

2

manipulated input, then a < 00. The competition results in a compromise and, thus, the optimal feeding policy will not necessarily be determined by the constraints. The following notations are used: • A is the reactant in the reactor at time t = o. • B is the reactant to be fed. • C and D are products. • lA and IB are impurities present with A and B, respectively.

k2

2

A utocatalyzed reversible reaction

k2

3

Production of isomers

A+B~C A+B~D

00

2

4

Consecutive reactions

A+B~C C~D

00

3

5

Autocatalyzed consecutive reactions

A+B~C C+D~2D

00

3

6

Consecutive parallel reactions

A+B~C A+C~D

00

3

7

Consecutive parallel reactions

A+B~C C+B~D

00

2

8

Reactions with impurity

A+IB ~ D

00

3

9

Reactions with impurity

00

2

lA

10

Reactions with impurity

C+I A ~ D

00

3

11

Reactions with impurity

C+IB~D

00

3

12

Reactions with decomposition

A+B~C aA~D

00

3

13

Reactions with decomposition

A+B~C bB~D

2

3

14

Parallel

A+B~C aA+B~ D

2

3

reactions

(a i- 1)

15

A+B~C

A+B~C +B ~ D A+B~C

A+B~C

A+B~C 2 A+bB ~ D (b i- 1) Table 1. AnalYSIS of reactlOn systems by manipulating the feed rate of B: a is the order of singularity and p the dimension of the reachable state space. Parallel

reactions

Table 1 summarizes the results found for some generic chemical reaction schemes. Most of them give rise to a = 00, for which the optimal feeding policy will necessarily be determined by the constraints. In Cases 1-3, there is only one independent reaction and no intrinsic compromise. If B and I B do not intervene in the second reaction (Cases 4, 5, 6, 10, 12), the absence of compromise is obvious. Interestingly, if B reacts with a product (Case 7) or with an impurity present in the reactor (Case 9), the two reactions are dependent (p = 2) and there is no compromise possible between them. More interestingly, if IB reacts with a reactant (Case 8) or a product (Case 11), though the two reactions are independent, there is no possiblity for a compromise. This can be explained intuitively by the fact that the dynamics of Band IB have "simillar" nonlinearities. However, compromise solutions are possible if B decomposes (Case 13) or reacts with A in more than one way (Cases 14-15) since, in such a case, the problem of selectivity arises . 6. EXAMPLES

3

Instead of providing a detailed study of all cases, three examples are considered corresponding to: (1) a = 00, p = 3, (2) a = 00, p = 2, and (3) a = 2, p= 3. 6.1 Consecutive reactions: Case

as long as the temperature is an independent variable, i.e., the influence of the feed rate on the temperature is either compensated or negligible. When the temperature is considered as an input variable, the optimum temperature profile, in most cases, corresponds to a compromise and is in the interior of the feasible region. Indeed, since T appears nonlinearly in F, FT contains T explicitly, leading to aT = O.

Reaction system: A Kinetics:

+B

4

~ C !:.; D

Transformations : n ..... ~ :

6 = nAO -

~ ..... n :

5. A CATALOG OF SYSTEMS WITH TWO REACTIONS

C = 0;

{

nB

=

ne

= ~1 -

det(M ) =

6 = nD 6, ~l + CBi" (V - Vo)

nA,

n A =n AO nBO -

6,

nD

kfk2CB

V4

=

6

2

.n

nA (nB - cBin V)

M loses rank for det(M) = 0 <=} nA = 0 or = CBin V . However, nA = 0 occurs only at the end of the reaction. Likewise, nB = CBin V implies CB = CBin, which is not possible when nA =I O.

In this section, several semi-batch reaction systems with two reactions are considered and checked whether or not a = 00. If the two reactions compete from the point of view of the

nB

392

So, p = 3 throughout the chemical reaction, and a = 00 . Thus, for any optimization problem, the solution is determined by the constraints.

Here C = 0 and det(M) = 0, i.e, p = 2 implying a reduandancy in the state representation. It can be checked that n::\(k 2 /k d ((k2 - kt}ne - k1 n A) remains constant throughout the reaction. Thus, 6 can be eliminated from the state space and computed algebraically using:

Objective: Maximize production of C. Constraints: Input bounds, constraint on the volume in the reactor . max J = ne(t,)

c c .,,2=.,,1

u(t},tj

s.t . System dynamics V(t,) ~ Vmax

a=

Optimal Solution:

max J

CB

5

in

UTnin U Tnax

VTnax

l/h

nA o

0.01 l/h 1 I

0

nBo Vo

= nc{t,)

S.t. System dynamics

Umin ~ U ~ u max ,

V(t,) ~ Vmax

nD(t,) ~ nD ,max = 0.05 mol

Optimal solution:

For the numerical values given in Table 2, the optimal cost is J = 0040 mol. Note that there is a compromise related to the choice of the terminal time but no compromise value for u(t) which is always on either one of the input bounds.

l/h g/ l

(((kl-k2} / kd (k2 / kd ) n AO nA -nA •

00.

u(t},t,

It

0.5 l/ (molh)

1

Objective: Maximize production of C. Constraints: Input bounds, constraint on the volume in the reactor, constraint on the maximum number of moles of D .

• As seen in Figure 1, the input is initially at the upper bound, U = u max , in order to increase nB as quickly as possible and maximize the rate of reaction . • Once the volume reaches Vmax , the input is set to Umin = O. • The number of moles of C decreases after a certain time due to the second reaction. So, to maximize ne, the terminal time t, is determined such that ~ I = O.

0.02

2 -

So, there is only one independent reaction and

Umin ~ U ~ Umax ,

kl k2

+ -k1k k

1 g/ l g/ l 0.9 I 0

Table 2. Model parameters, operating constraints and initial conditions

• As seen in Figure 2, the input is initially at the upper bound, u = u max . • For the numerical values given in Table 2, the constraint nD(t,) ~ nD,max is more restrictive than V(t,) ~ Vmax . So, the input switches to Umin = 0 so as to be able to satisfy nD(t,) = nD ,max . • The number of moles of C decreases after a certain time due to the second reaction. So, to maximize ne, the terminal time t, is determined such that ~ t I = O.

I

....

" ' r - - - - - - - '- , . . -- -- ---,

The optimal cost is J = 0.33 mol. Here, the compromises intrinsic to the optimization problem are reflected in the choice of the switching and final times. Despite the fact that B appears in both reactions, there is no compromise value for u(t), and the solution is on the input bounds.

0 .009

0 .000 0.007

0.000 u(~l

0.005 0.00<

0.003

O.O'r - - - - - - -'''''''' y------,

0002

0 .000 0.00 1

1._10lh] \_'4.M,,,) °o~----:----....:----:,':-o---'--~"

0.008 0 .007

0 .000 u lllh )

Fig. l. Optimal feedrate for Example 1

0 .005

0.004

6.2 Consecutive parallel reactions: Case 7

0.003

. system: A R eactzon Kinetics:

+B

kl -->

C, B

0.1)02

k2 D . + C -->

0.001 t• • •. 15(h

\_13.01 [1'1]

~~~-~-~.~~.-~'~O--'-~I2~~ TII'M[h]

Fig. 2. Optimal feedrate for Example 2

Transformations: = nAO - nA, 6 = nD nA = nAO - 6 n : nB = nBO - ~1 - 6 + CB'n (V - VO) { ne = 6 - 6, nD = 6

n

--> ~ :

~

-->

6

6.3 Reactions with decomposition: Case 13 Reaction system: A

393

+ B !::. C,

B ~ D.

Kinetics:

the solution will not have had the arc u comp , i.e. , = Urnax till ts = 10 hand U = Urnin for the rest of the batch with t f = 20 h. U

Transformations: n--;~:6=nAo-nA,

~ --; n :

= nAO - ~l = nBO - 6 - 6 + CB,n (V = ~l' nD = 6

nA nB

{

7. CONCLUSIONS

6=nD

ne

This paper investigates whether the optimal feed rate policy will be on the constraints of the optimization problem for semi-batch reaction systems. The majority of reaction schemes for two reactions is shown to possess this desirable property.

- VO)

c=

2 k 1k2nACBm(CBin V - nB) Vs The value of er is 2 since C i O. The compromiseseeking feed rate can take values inside the feasible region and is computed from det(M) = 0:

Since for most two-reaction systems, the optimal feed rate policy is determined by the constraints of the optimization problem, a measurement-based scheme for tracking the active constraints is an interesting alternative to numerical model-based optimization, especially in the presence of uncertainty (Bonvin et al., 2001).

Objective: Maximize production of C. Constraints: Input bounds, constraint on the volume in the reactor , constraint on the maximum number of moles of D , constraint on the final time. max J = ne (t f)

8. REFERENCES Benthack, C. (1997). Feedback-Based Optimization of a Class of Constrained Nonlinear Systems: Application to a Biofilter. Phd thesis 1717. Swiss Federal Institute of Technology. Lausanne, Switzerland. Bonvin , D., B. Srinivasan and D. Ruppen (2001) . Dynamic optimization in the batch chemical industry. In: Chemical Process Control - 6. Tucson, AZ . Bryson, A. E. (1999). Dynamic Optimization. Addison-Wesley, Menlo Park, California. Macchietto, S. (1998). Batch process engineering revisited: Adding new spice to old recipes. In: IFAC DYCOPS-S. Corfu, Greece. Palanki, S. , C. Kravaris and H. Y. Wang (1993). Synthesis of state feedback laws for end-point optimization in batch processes. Chem . Engng. Sci. 48(1) , 135-152. Srinivasan , B., C. J Primus, D. Bonvin and N. L. Ricker (2001 a). Run-to-run optimization via generalized constraint control. Control Eng. Practice In Press. Srinivasan, B. , M. Amrhein and D. Bonvin (1998). Reaction and flow variants/ invariants in chemical reaction systems with inlet and outlet streams. AIChE 1. 44(8) , 1858-1867. Srinivasan, B .. S. Palanki and D. Bonvin (2001b). A tutorial on the optimization of batch processes: I. Characterization of the optimal solution. Comp. Chem. Eng. Submitted. Visser , E .. B. Srinivasan . S. Palanki and D. Bonvin (2000). A feedback-based implementation scheme for batch process optimization . J. Process Contr. 10, 399-410. \Viederkehr , H. (1988 ). Examples of process improvements in the fine chemicals industry. Comp . Chem. Eng. 43. 1783-179l.

(13)

u (t).t{

S.t. System dynamics

SuS umax , V(tf) S nD(tf) S nD ,max = 0.02 mol Umin

tf

S

Vmax

tf ,max = 20h

Optimal solution: • As seen in Figure 3, the input is initially at the upper bound, U = u max . • The compromise between the production of C and D is reached through u comp . • Since C is not consumed by another reaction , the optimal value of t f is t f ,max· For the numerical values given in Table 2, the optimal cost is J = 0.41 mol. '.....' 0 0 1- - . - - - - - - - - - - 1

oooo l

I

OOOO ~

OOO' r o ()()8 ~ u[""1 . ooos·

t... \ .g

~hl

,,

o 003r 0002 "

000 ' ·

I ~~-~~-~,-c---,~,--~~ T_ i";

Fig. 3. Optimal feedrate for Example 3

Remark: Even when there is the possibility of a compromise-seeking arc. its presence or absence in the optimal solution depends on the formulat ion of the optimization problem. For example. if the constraint on nD (t f) had not been there in (13) ,

394