- Email: [email protected]
Interpolating optimizing process control
J. Proc. Cont. Vol. 7, No. 2, pp. 129-138, 1997 © 1997 Elsevier Science Ltd Printed in Great Britain. All rights reserved 0959-1524/97 $17.00 + 0.00
•
~LSEVIEX
PIh S0959-1524(96)00023-6
Interpolating optimizing process control Bjarne A. Foss *t and S. Joe Qin* *Department of Engineering Cybernetics, Norwegian University of Science and Technology, Norway. email: [email protected] ~Department of Chemical Engineering, The University of Texas at Austin, USA. email: qin@che, utexas.edu Received 2 April 1996; revised 23 September 1996 A new model-based optimizing controller for a set of nonlinear systems is proposed. The nonlinear model set is based on a convex combination of two bounding linear models. An optimal control sequence is computed for each of the two bounding models. The proposed control algorithm is based on a convex combination of the two control sequences. A novel feature in these two optimizations is an added constraint related to the feasibility of the 'other' bounding model. The control algorithm can for example be used in model predictive control. We provide robust feasibility guarantees and an upper bound on the optimal criterion if the bounding models are linear FIR models. Further, simulation examples demonstrate significant feasibility improvements in the case where the bounding models are general linear state-space models• The proposed method guarantees robust feasibility for a l-step ahead prediction in the general case. This can be of interest in MPC applications. © 1997 Elsevier Science Ltd. All rights reserved
Keywords: model-based control, optimization, convexity, model predictive control, nonlinear systems, robustness, interpolation The combined use of dynamic models and optimization for process control offers a concept in which process knowledge can be linked to operational goals formulated by some optimization criterion. This concept has seen widespread use, particularly through the applications of model predictive control (MPC). MPC refers to a class of algorithms where an optimization problem is solved repetitively, at every new time-instant. Only the first part of the computed control sequence is applied to the system since a new optimal control sequence is computed and applied at the next time-step. Several reviews of MPC technology exist, see for example Lee ~, Rawlings et aL 2 and Qin and BadgwelP, the latter emphasizing industrial use of the technology. The interaction between control and optimization is discussed in an iliaminating way in Mayne 4. In this paper the divide between convex and non-convex problems is emphasized. In particular the problems arise by using nonlinear models since the optimization problem in these cases in general becomes non-convex, State-of-the-art optimizing control is based on linear dynamic models and linear constraints on the control inputs and system outputs. Nonlinear optimizing control has been studied by Rawlings et al. 2 and Oenceli and Nikalaou 5. Further, some approaches were
described by Bequette 6 in a somewhat earlier paper. These control strategies normally result in a nonconvex optimization problem. In this paper we explore optimization control based on nonlinear models. In particular, the goal is to derive an approach which can offer a smooth transition from linear optimizing control to nonlinear and robust optimizing control. By smooth transition we mean an approach which does not invoke the full mathematical 'machinery' of nonlinear optimizing control in the general case. This is particularly important from an industrial viewpoint as it simplifies the transition from the application of linear to nonlinear optimizing control. The smooth transition is accomplished by constraining the nonlinear optimization problem along three axes. First, the set of nonlinear models is limited. Second, constraints are added to the optimization problem to enhance feasibility of an optimal solution for a s e t of bounding linear models. Third, the control sequence is computed as a convex combination of control sequences computed on the basis of the bounding linear models. The remainder of this paper is structured as follows: in the next section we formulate the problem. Following this, the theoretical foundation for the proposed control algorithm is developed. In the following section three simulation examples are investigated to explore the theoretical findings, and to study the proposed method beyond the assumptions given by the theory. The pro-
*On sabbatical leave at The Department of Chemical Engineering,The University of Texas at Austin, USA 129
Interpolating optimizing process control: B.A. Foss and S.J. Oin
130
posed method is after this investigated in closed loop control. Finally, a comprehensive discussion is provided before the conclusions finalize the paper,
Problem definition Background
Assume an optimality criterion given by: (1)
¢ ( ~ , Z ) = ~.,l(x~+~,ui) ~:N
~ = arg min~,n¢(Zc, Z)
(6)
subject to the constraints Z ~ x; (5) using one of the models (2), (3) or (4). The optimization problem above is generally nonconvex when a nonlinear model within Y,,. is used. This can cause problems, especially in on-line control where computation time can be an issue. In MPC the minimization problem (6) is solved repetitively, at each time-step, with new initial conditions. Only the first control value of the zr sequence is actually applied to the process. Furthermore, it is typical to parameterize the control sequence as follows:
where zr = {u~ .... , u~ ..... u~}, M < N - 1
(7)
-- {.To,. ..,UN_I} , Ui E U cc_.R mu, I"[-~-U x . . . x U Z = {xT ..... xT}, x~ ~ X ~ R ' ~ , U, X - convex set
Z = X ×... × X
l • R m"×% --+ R ÷, l - convex function, i ~ I N I N = {0,..., N - 1}
The optimality criterion is defined on a time horizon from 0 to N. u~ denotes the control input to a system. The control input u~ is constant during the time span [i,i + 1). x~ denotes the system states at time i. Since U and X are convex sets, rI and Z are convex sets. Further, 0 is convex since l is convex, We define two linear state-space models (2~, Y~0)and a set of nonlinear models (Y~w)based on the interpolation of the two linear models. ~., 1: h l ( x i ÷ l , x ~ , u , ) = x~+l- f ( x ~ u ~ ) = B~u,. O, i ~ I,v =
xj+ x - A i x ~
(2)
0: h o ( x i . l , x i , u i ) = x i + x - fo(xi , ui) = x i ÷ I - Aoxi -
B0u~ = 0, i ~ I~.
(3)
Z,.---- {hw" hw(Xi+l, xi, ui) = w(xi, ui)hl(Xi*t,xi, ui) + (1 - w(x,, u,))ho(xi+ I, xi, ui) = 0 (4)
V i ~ Lv, w ~ W, W = { w ~ [0,1] Vx~ ~ X Vu~ ~ U}}
Cw(x~,u~)~
h, is a nonlinear function, constructed as convex cornbinations of h~ and h0. It should be noted that the continuous function w in general will depend on the states and the control inputs. We further observe that the interpolated model hw ~ ]~,,. is not generally a convex function since w may be any continuous function bounded by [0,1]. Finally, we define initial conditions, x0 - given
(5)
The problem we want to address is to minimize Equation (1) with respect to n based on the different constraints discussed above, hence we want to solve:
This means that the control input is constant during the last part of the control input sequence. The criterion function (1) does not cover all possible criteria; penalizing changes in the control input is for example not included. This type of change does not influence the results in this paper as long as the criterion function remains convex.
Interpolating control
We will in the following describe the proposed controller. First, we define two control sequences 7cI E 1-1 and Jr0 e FI and a set of interpolating controllers. H~ = {n~ " n~ = azr~ + (1 - a)n0,
V a ~ [0.1]}
(8)
The control sequence rc~ forms the basis for the controller. It should be noted that zr~ is feasible since it is based on interpolation within a convex set, i.e. FIo _C I-I. o~may in general vary from one time-instant to another. For later convenience we pursue with some definitions. Z~ and Z0 are the state sequences on the horizon { 1..... N} obtained by applying some control sequence n I --, Z ~N} I is the 1-I to Yl and 5".0, respectively. Z kl = {Z k~,state sequence obtained by applying the control sequence Jr~, I ~ { 1,0} to the system Zk, k ~ { 1, 0}. Further, we define the set o f states associated with the model set Y,, and the control sequence set I-I,. Z ,a ~ {Z,~~ " Z,~ = {Xw~l,...,Xwu},C~
Vh w ~ Ew,O~ ~ [0, 1]} (9)
Later we will in particular study the interpolated controller set I-I, where zq and zr0 are computed on the basis of the following two optimizing problems. zr~' = arg min=~ n ~(n',Z~ )
(10)
subject to the initial conditions (5), constraints Z~ ~ -K and Z~ given by Equation (2). o = arg m i n ~ n ~(zc, Z0 ) zc0
(ll)
Interpolating optimizing process control: B.A. Foss and S.J. Qin subject to the initial conditions (5), constraints Zo ~ X, and ~0 given by Equation (3). If we choose zq = zc[ and ~0 = gg the state sequences Zl~ and ~ are the optimal sequences for ~1 and ~0, respectively. Further, following the above notation 0~ is the criterion value obtained by applying the control sequence ~ , I s {1,0} to s y s t e m Y , , k s { 1 , 0 } . T h i s means that ~t and ~ are the optimal criterion values for El and ~;0Note that we assume that the state sequences discussed above are based on the initial value (5).
Analysis We will in this section analyze the use o f the interpolating controller on the model set Y~,~.. Guaranteeing feasibility There is no guarantee that the state constraints are satisfied on the horizon {I ..... N} when the control sequence n'~ e FI~ is applied to the system h,,. E ~,,,,. We will in the following include restrictions on the choice of zq and ~ , and on Y~,, so as to ensure that X~. C X. Before we do this we state the well-known connection between FIR-models and state space models. Proposition 1
Theorem 1 Assume two control sequences zc~ and zr0 such that ~ ~ X Vl e {0,1}, k s {0,1} for the initial conditions (5). Further, define a control input sequence set rI~ by Equation (8). If the linear models E1 and ~0 are FIR-models, and z~ ~ rI~, then Z~,. e x . ~ is defined in Equation (9). Proof." We first use Proposition 1. Since the bounding models (2) and (3) are FIR-models, we transform them into the equivalent input-output model structure (12). The model set Y~,,. contains nonlinear FIR-models and may, hence, also be transformed into a FIR-structure. Let z~ denote the j t h element of the predicted output vector y~+~ of the FIR-model (12) at some i e I N if we apply the control sequences rc~, l ~ { 1, 0} to the system Y~k, k ~ { I, 0f. The assumptions on Jr~ and zc0 guarantee that z~ ~ ~ where X = X~ × ... × Xmx" Let ~ denote the j t h element of output vector y,+~ of the FIR-model (12) by applying a control sequence zc~ lq~ to the system h,~ e Z,,. Let c~, co and Cw denote the j t h row vector of C~, C O and C,, w ~ W. respectively. We first show that c~T~" e X,, c~ ff~, = c~ ( a ~ + (1 - a)ff0~) = ac~ ff~i + -~ + (1 - a)z°~ (I - a)c~ ff0~ = m~
ke{l'O'W}'ieIN (12)
c ff~, = w(x,, ui)c I ff~, + (1 - w(x,, u,))c o ff~, is equivalent to a state-space model on the form (2), (3) or (4) if:
.
]
0O/ ~
=. A~
I Q~
...
I = .
C~L ~/ ,0
.
.
0 . " ...
.
.
(14)
We have above shown that c~ ~=: ~ Xj and c o ff=~ ~ X~. Since X~ is a convex set, the convex combination c ff=, must lie in X.
.
. " I
~ ! B~ =
(13)
c~ ~ , lies in ~ since it is a convex combination of two elements in the convex set Xr An identical argument may be formulated for co ff~,. We now choose some element cz e [0, 1] and show t h a t c u'~, e X ~ f o r a n y w ( x , , u , ) ~ [0,1].
An FIR-model of this form
)',÷~ = C ~ o u . ' + ' " + C ~ u , - c = C k U , ' 3', e Y ~ R my
131
k ~ {1, 0, W} '
The above result cannot be generalized to arbitrary A1 and A0 matrices. Further, it is critically dependent on the assumption that each of the control sequences n~ and ~ applied to either of the models ~t or ~0 obey the stateconstraintonthehorizon{1,...,N}. This theorem will form the basis for several corollaries. First, we use it in the case where ~c~ and ~ are based on the solution of some optimization problem. Corollary !
and y~ = (0 ..... 0, /)TX~. We note that C,,~, or C,, (w ~ W), and correspondingly A~ and B~ in general depend on the control inputs and states. If we limit the bounding linear models to be FIRmodels, and include feasibility constraints on the choice of zq and rc0 it is possible to guarantee feasibility in the sense that ~ C_ )2. This is shown below.
Given the following control sequences
~ = arg minion 0(z~,Z~)
(15)
subject to the initial conditions (5), constraints Z~ s x (Z~ is the state sequence obtained by applying ~c~ ~ FI to Y~t), and Z~ ~ X. ~ = arg min.~ n ¢(zc, Z0)
(16)
Interpolating optimizing process control: B.A. Foss and S.J. Qin
132
subject to the initial conditions (5), constraints ~ e 22 and ~ e .~ If the linear models ]~i and Y+0 are FIR-models, and zc~ ~ FI~, then Z~ ~ X.
Remark This corollary shows that the 1-step ahead prediction will satisfy the state constraints for any linear m o d e l s Y+I and Y+0 if zc~ and ~0 are computed using the additional constraints as in Equations (15) and (16). This implies that an MPC controller based on 7r~ and Equa-
Proof" Since the control sequences defined in Equations
tions (15) and (16) always will be feasible even though the constraints may be violated on the optimization horizon. The reason is that the MPC controller only
(15) and (16) satisfy the state constraints ~ e X V I e {0,1}, k e {0,1} in Theorem 1, the corollary follows directly from Theorem 1.
Remark This corollary states the important result that an arbitrary control sequence zr~ ~ I-I~ applied to an arbitrary model h,,. ~ Y+,,.will satisfy the state constraints on the horizon { 1.... , N} provided additional hard constraints are included in the optimization, i.e. including o O Zo ~ X for computing zrl, and Z1 ~ x for 7r0. If we limit our attention to 1-step ahead prediction, however, the result can be extended.
Given the control sequences tot, ~ and zc~ as defined in Equations (15), (16), and (8). The state sequence which arises by applying lr~ to h,. ~ Y+,,. is denoted ~ . = {~,.1..... ~tN}, cf. Equation (9). Then ~,.t e X. Corollary 2
Proof" Let z: denote element j of state vector x: which
k +1 arises by applying the control input u~0, l e { 1, 0} to the system z~k, k e { 1, 0}, and ~. denote elementj of state vector x ~t which arises by applying the control sequence zr~ ~ FI~ to the system h,,. ~ Y.,,. Let at, ao and a,. denote the j t h row vector of A1, Ao and A,,., and bl, b 0 and b,, denote the j t h row vector of B1, B0 and B,,.. Further, ul0 and u00 denote the control input of z~, and ~ at time 0. We need to show that z~ ~ ~ , where X =Xt × . . . ×
applies the first control input before recalculating the control input at the next time-step. Robust performance
Theorem 1 provides a possibility to obtain an upper bound on the following worst case scenario. The upper bound is defined by: ~ = suph,:z ,,.,,~ ~n~ ¢(zc,+,Y~)
(19)
This result is formulated by the following lemma. Corollary 3 Suppose the control sequence zc~ e FI~ is applied to the system h,,. e Y+,.with initial conditions (5). Y/ and z~0 are FIR-models. Further, suppose the optimality criterion function 1 is given by some separable norm function, i.e. l(X~+l, u,) = , , Ily~+llI?2 + Ilu~l[~, where y~ is the output of the FIR-model formulation (12). Then an upper bound for ~ in Equation (19) is given by:
(20)
~ <_ ~ Cu--~ + u~.
,~z~. where Cu,. = maxk.~,olHf~u,[[~ and ~
maxt~l,0~llu,I]~.
X,..
We first show that a 1 Xo + b 1 U~o ~ Xj.
Proof" We first compute an upper bound for
IlY~+II[~= IIC,~'oilI~= II(wCx + (1 - wf0))u+ell~~
a a Xo + b I U~o = a I (~Xo + (1 - a)Xo) + bl(aUto + (1 -
]ly~+lll~.
o:)U0o)
= a ( ~ 1 X 0 + btUl0 ) +
< [w]ICIUc~:HQ+ (1 - w)lIC0~+~,llo]2 < max~ u,0~[lf~u2t~
(1 - ot)(a~x o +
b 1 u0o) l + (1 -
cOz~
(17)
~" ~(~'1
~ =
a I x0 + b~ u~0 lies in X~ since it is a convex combination of two elements in the convex set X~. An identical argument may be formulated for aoX o + bou~o. We now choose some element ~ ~ [0,1] and show that a,.xo + b,. u~0 e Xj for any w(x~,u~) ~ [0,1].
maxk~ll,0~llC,[m71, + (1 -o0ff0i]]]~
< max~{,,01[o~lC,~l,[]~ + (1 < max,.:, i i.o~HC~u t~llb = Cu, A upper bound on
a,,xo + b,~uoo = w(x,,u,)(a:o + blu~0) + (1 - w(x,u,))(aoX o + bongo )
(18)
Since atx o + bt U~o e Xj and aox o + b o u=o ~ Xj and X] is a convex set, the convex combination a,,.Xo + b,u~o must lie in X~.
a)[IC~+0,ll~] = (21)
Itu,lt~ is easily computed.
Hu,IL~ = Ilu=~ll?~ = I[~ul, + (1 - c~)u0,ll~ < [c~lut+ll~+ (I - a)llu0:ll~]: < max]~ll,01lluz,ll~ = u:
(22)
Interpolating optimizing process control: B.A. Foss and S.J. Qin
133
By these expressions the corollary is proved. The important consequence of the above lemma is the fact that an easily computed upper bound on robust performance can be found. This can hence be used as a measure for robust performance. To compute the bound we use the output Yi instead of xi. This is no important limitation since y~ forms the output part of the states of a FIR-model. A weighting of y~ is, hence, the obvious choice for the criterion,
We will first study the feasibility of the proposed controller, Equation (8), by studying the output trajectoties. This is done in three steps.
Examples
2.
We will in this section further explore the findings in the above sections through three numerical examples. All the models, except for Y-0 in example 3, are typical models found in the process control problems. In the first example we will investigate the method on a nonlinear FIR-model. In examples 2 and 3 we will study the method beyond the assumptions given by the theory. In particular, we investigate situations where there is no theoretical guarantee for a feasible solution,
3.
1. The control sequences Jr~ = (0.90, 0.60, - 0.70, 0.33, 0.06, -0.11, 0.11) r and n o = (0.90, 0.60, -0.70, 0.13, 0.28, -0.17, 0.08) r are computed using Equations (15), (16), i.e. with the additional constraints in place. Six interpolating control sequences zr~ are generated by choosing a = {0.0, 0.2, 0.4, 0.6, 0.8,
1.0}. We choose six models within ~,,. by choosing constant interpolation weights w = {0.0, 0.2, 0.4, 0.6, 0.8, 1.0}. The six models are simulated using the six interpolating control sequences as inputs. This gives a total of 6 × 6 = 36 simulation runs. The system is simulated from 0 to 7, i.e. equal to the optimization horizon. The system outputs for each simulation run is checked against the constraints, i.e. [0.00, 1.01]. The result is shown in the matrix in Table 1. The 36 columns relate to each of the simulation runs, while the eight lines are related to time instants
{o,1,...,7',. E x a m p l e 1. F I R m o d e l s
We assume that a single-input single-output system is given by some model h,, ~ Y~,., cf. Equation (4), with bounding FIR-models given by: k ~ {1, 0} b~ = {-0.10, 0.30, 0.40, 0.30, 0.10} 50 = {0.20, 0.20, 0.15, 0.10, 0.05}
Y~+z = CkoU~+ ... ÷ C~4Ui-4 = CkU~,
(23)
Proposition 1 gives the equivalent state-space formulation. The two bounding models resemble a system with an inverse response and a steady state gain of 1.0, and a first order response and a steady state gain of 0.7, respectively, We use the well-known quadratic optimality criterion, add hard constraints on the system output and control input, and assign initial conditions. x 0( Jr, 2") = 0 . 5 ~ y~ + ru~_ 1 ~=~ N = 7 r = 0.01 U = [-1,1]
0 indicates no constraint violation, while 1 or 2 signifies a violation of the minimum or maximum allowable output, respectively. Table 1 shows that there is no violation of output constraints. This is in accordance with Theorem 1. Corollary 3 gives an upper bound on the optimality criterion, for this example we compute 0 = 1.93. The criterion values for the 36 simulation runs in Table 1 lie between 0.47 and 1.61. This indicates that the upper bound need not be very conservative. To study the impact of the additional constraint, i.e. Z0 ~ Y to compute 7r~ and Z ~ x to compute ~ro, we remove these constraints and perform the computations in items 1-3. This means computing lr~ and iv0 according to Equations (10) and (11). These results are shown in Table 2. There are in these cases many violations of the output constraints. In the latter case the control sequences are rc~ = (0.90, 0.60, -0.70, 0.33, 0.06, -0.12, 0.09) T and Jr0 = (-1.00, -0.34, 0.34, 0.17, 0.00, -0.17, 0.00) T. They are very different from each other. This is not the case in item 1 above, in that case the control sequences are identical except for the last control input. The reason for this are the added constraints which limit the possible choices of control inputs.
X = [-1, 11 × [-I, 11 x [-1, 11 x [1 - 1, 1] x [0.00, 1.01] x o = (1, I, 1,1, O) T
The state-space model has five states. The first four states are delayed control inputs while the last state is the system output, equivalent to y~ in Equation (23). Hence, the state constraints X limits the system output to [0, 1.01]. Note that the system output does not depend on the last element of x 0 since it does not have any auto-regressive terms.
E x a m p l e 2." Bounding second order models
We again assume that a model set Y~,, defines a singleinput single-output system. In this example the bounding linear models are given by
A,=
I'0.00 0.00 0.00
0.00 0.00 -0.79
0.00 / 1.00 / 1.78)
(1.001 /~ = / 0 . 0 0 [ ~0.01)
Interpolating optimizing process control: B.A. Foss and S.J. Qin
134
Table 1 Example 1 - Bounding FIR-models. Each column represents a simulation run using the control sequence zt~ on a model h,, e 5".. The
system outputs of each simulation run are checked against the constraints at each time-step. 0 indicates no output constraint violation. 1 or 2 signifies a violation of the minimum or maximum allowable output, respectively Time i=0 i=1 i=2 i=3 i=4 i=5 i=6 i--7
Model 1
w = 0.8
w = 0.6
w = 0.4
w = 0.2
Model 0
~o ..... =o0
=o, ..... ~ 00
=o ..... ~ o0
~o ..... ~ o0
~ ,0 ..... ~ o0
~ o ..... ~ 00
000000 000000 000000 000000 000000 000000 000000 000000
000000 000000 000000 000000 000000 000000 000000 000000
000000 000000 000000 000000 000000 000000 000000 000000
000000 000000 000000 000000 000000 000000 000000 000000
000000 000000 000000 000000 000000 000000 000000 000000
000000 000000 000000 000000 000000 000000 000000 000000
Table 2
Example 1 - Bounding FIR-models. ~ c h column represents a simulation run using the control sequence ~ on a model h~ ~ Z~. The ~stem outputs of each simulation run are checked against the constraints at each time-step. 0 indicates no output constraint violation. 1 or 2 signifies a violation of the minimum or m~imum allowable output, res~ctively Time
Model 1
w = 0.8
0
i=0 i=1 i=2 i=3 i=4 i=5 i=6 i=7
0
000000 022222 000000 000001 000111 000001 000000 000000
(0.00 0.00 A o = |0.00 0.00 ~,0.01 -0.855
w = 0.6 0
000000 000002 000000 000001 000011 000001 000000 000000
0.00) 1.00| 1.85)
The third state is the system output, hence Yi = (0, 0, 1)a x~. 21 is a second order model with two time constants o f a b o u t 6 and 15 time-steps, no time-delay and a steady state gain o f 1.0. 20 is a second order model with two time-constants o f a b o u t 10 and 20, a time delay o f 1 a n d a steady state gain o f 2.0. T o elaborate on the model set Y,. let w (x,, ui) ~ [0,1] be a constant. Hence, the model set Y~,,.consists o f only linear models. The steady state gain, time-delay and time-constants for this model set are shown in Figure 1. We observe that there is a close to linear change in all the characterizing variables. The optimality criterion, hard constraints on the systern o u t p u t and control inputs, and initial conditions are given below,
(24)
[0, 1.Ol].
0
000000 000000 000000 000000 000011 000001 000000 000000
Model 0
0
0
000000 000000 000000 000000 000011 000001 000000 000000
0
000000 000000 000000 000000 000000 000000 000000 000000
The control sequences r~ = (-7.48, -5.14, -3.28, - 1 . 8 7 , - 0 . 8 9 , - 0 . 2 8 ) r a n d ~ = (-7.49,--4.84, 2.80,-1.35, -0.43, 0.00) a- are c o m p u t e d using E q u a t i o n s (15), (16), i.e. with the additional constraints in place. We compute six interpolating control sequences zr~ by choosing o~ = {0.0, 0.2, 0.4, 0.6, 0.8, 1.0}, and choose six models within Y~,, by selecting constant interpolation weights w = {0.0, 0.2, 0.4, 0.6, 0.8, 1.0}. Again we simulate a total o f 6 x 6 = 36 runs a n d check the system outputs against the constraints. There is no violation at any time-step in any o f the 36 simulation runs. We also c o m p u t e the control sequences omitting the additional constraints, i.e. using Equations (10), (1 I). In this case we obtain identical control sequences as above. Hence, in this case the additional constraints m a k e no difference on the c o m p u t e d control sequences.
Example 3: Bounding second order models T o m a k e the p r o b l e m m o r e challenging we will in this example study a set o f systems where the b o u n d i n g models have different types o f behavior. Y'I is equal to ]~ in the above example, while Y0 is a second order model with oscillatory modes. The 'time-constant' of the oscillatory models is 29 time steps and the oscillation period is 5-6 time-steps. The steady state gain of Y~0 is 1.0 and the time-delay is 1. The b o u n d i n g models are described by:
~=(],,,,; The state constraints Z limits the system o u t p u t to
w = 0.2
0
000000 000000 000000 000001 000011 000001 000000 000000
(1.001 B0= |0.00| t0.00j
x ~(rc, Z) = 0 . 5 ~ y ~ +ru~_l ~=~ N = 6 r = 0.01 U = [-10,10] X -- (-oo, oo) × (-oo, oo) × [0.00, 1.01]
w = 0.4
0
(0.00 "
t,0.00
0.00
0.001
(,.001
0.00 1.00 /
B,--/0.00 !
-0.79
1.78)
~,0.01)
Interpolating optimizing process control: B.A. Foss and S.J. Qin a 1.8 1.6 1., 1.z
r,~ ~ , ~ , t l " o.lo ~ ~
i ~-.._~ 0.8 ~ 0.6 o.,
~m,,on,ta,,a'0.10
~d,~,y
~
0.2 oi* o'.s ole o17 o18 o'.~ wa~,~ngf~,,~on,* Figure 1 Example 2 - Bounding linear models. The gain, time-delay and time-constants of the models with constant w-valuesare depicted °o
o11
ola
o'.3
"] 0.00
0.00
0.62
-0.32
0=000 l
0.00
000 100! 0.70)
(
1 1.00
0--I000 / ~,0.00)
T o get a feel for the model set ]~,, we simulate a step response for six models with constant w-values, cf. Figure 2. The steady state gain does not change; the dynamics do, however, change a lot. The dynamics are particularly sensitive to w close to 1. The optimality criterion, constraints and initial conditions are equal to those in Equation (24). The output is again equal to y~ (0, 0, 1)Xx,. The control sequences z:] = (-0.61, 0.52, 0.00, 0.00, 0.00, -0.84) 7 and ~0 -- (-0.59, 0.50, 0.01, 0.00, 0.00, 0.00) ~ are computed using Equations (15), (16). Again we simulate a total of 6 × 6 = 36 runs and check the system outputs against the constraints. The result is shown in Table 3. We observe that in many of the late time-steps (4-6) there are violations of the minimum value o f the system output.
..
~.4
,o.~ 1.a 1
.~o., E
~0.6 o.4
~'*~
o.2 o
1'o '~'s
2'0
2's ~ 3's ~ ~s so r~,,,op, Figure 2 Example 3 - Bounding linear models. There is a step in u (u goes from 0 to 1) at time 1. The initial conditions are x = [0, 0, 0]v. Step responses are shown for the six models used in the computations in Example 3. i.e. w = {0.0, 0.2, 0.4, 0.6, 0.8, 1.0}
135
To study the impact of the additional constraints, we remove these constraints. This gives the control sequences ~ = (-7.48, -5.14, -3.28, -1.87, -0.89, -0.28) w and ~ = (-0.59, 0.50, 0.01, 0.00, 0.00, 0.00)L The results are shown in Table 4, the number of violations is substantially larger than in Table 3. The numbers behind these tables accentuate this difference. The minimum output constraint, i.e. 0, is only marginally violated in Table 3. Most of the violating outputs lie in the range [-0.10, -0.05] with a minimum value o f - 0 . 1 6 . In Table 4 most of the violating outputs lie in the range [-3, -2] with a minimum value of-5.92. To further illuminate the impact of the additional constraints we show the criterion values for the 36 simulation runs in the cases with and without the additional constraints in Tables 5 and 6, respectively. The (1, 1) and (6, 6) entries show the nominal cases for Yl and 5".0, i.e. ¢~] and ¢0. We observe, as should be expected, that the optimal criterion value in these nominal cases may increase due to the additional constraints. More interesting, however, is the fact that the added constraints seem to have a significant robustifying effect. The worst case criterion values show a dramatic improvement, 48.70 to 2.57. The improvements are especially pronounced for high c~-values, i.e. by comparing the lower left part of Tables 5 and 6.
Closed loop behavior As discussed earlier the proposed method can be used for closed loop control. In this section we will investigate this by extending the simulation examples and discuss issues related to feasibility, closed loop stability and robustness. Example 1 We extend Example 1 by using our method in MPC. The system, optimization criterion, constraints and initial values are unchanged. Figure 3 shows some typical results. In this simulation the system is given by Y~, i.e. w ~ 1, while the control input is computed by choosing cz = 0, cf. Equation (8), i.e. the control input is based on ~0. We observe that the output constraints are satisfied when the added constraints are used. This is, however, not the case if these added constraints are omitted. The above might lead to the conclusion that an initial feasible solution implies a feasible solution at any time. This cannot be guaranteed. This is seen by performing the same simulation as above but increasing the gain of the simulated system. We increase the gain of ~t by a factor of 5. Figure 4a shows a violation of the output constraints at certain time instances. The reason for this is that there exists some x s X where no feasible zr can be found, i.e. no n ~ H such that the predicted states satisfy the state constraints for both bounding linear models. It should be be noted that this shows a rather extreme situation as the first coefficient of the FIRmodel can vary between -0.50 and 0.20, cf. Equation
Interpolating optimizing process control: B.A. Foss and S.J. Qin
136
Table 3 Example 3 - Bounding second order models. Each column represents a simulation run using the control sequence zr~ on a model h,. ~w- The system outputs of each simulation run are checked against the constraints at each time-step. 0 indicates no output constraint violation. 1 or 2 signifies a violation of the minimum or maximum allowable output, respectively Time
Model 1 0
~o ..... ~0
i=0 i=1 i=2 i=3 i=4 i=5 i=6
000000 000000 000000 000000 000000 000000 000000
w = 0.8
w = 0.6 0
~o ..... ,~0
000000 000000 000000 000000 000000 000000 000000
0
w = 0.4 0
= , ..... ~0
000000 000000 000000 000000 000000 I 1111 I 11111 l
w = 0.2 0
~ ..... ~o
000000 000000 000000 000000 111111 I 11111 111111
Model 0 0
~o ..... =0
000000 000000 000000 000000 I 1111 I 11111 I 111111
0
~o ..... ~o
000000 000000 000000 000000 000000 000000 000000
Table 4 Example 3 - Bounding second order models. Each column represents a simulation run using the control sequence ~ on a model hw ~ . The system outputs of each simulation run are checked against the constraints at each time-step. 0 indicates no output constraint violation. 1 or 2 signifies a violation of the minimum or maximum allowable output, respectively Time
i=0 i=I i=2 i=3 i=4 i=5 i=6
Model 1 0 ~,-",~0 000000 000000 000000 000000 000000 000000 000000
w = 0.8 0 0 ~l,'",~0 000000 000000 110000 111100 I I 1110 111110 111111
w = 0.6 0 0 ~ i ..... ~0
w = 0.4 0 0 ~ l ..... ~0
w = 0.2 0 0 ~t,'"'~0
Model 0 0 ~"",~0
000000 000000 111100 111110 111110 111111 111111
000000 000000 111110 111110 111111 111111 I 11111
000000 000000 I 11110 111110 111111 111111 111111
000000 000000 111110 111110 111110 111110 111110
Table 5 Example 3 - Bounding second order models. The criterion values using the control sequence rc~ on model wht + (1 - w)he for different values of a and w. The extra constraints are included ¢ Model 1 w = 0.8 w = 0.6 w = 0.4 w = 0.2 Model 0
n°
a =0.8
a = 0.6
a = 0.4
~ = 0.2
zrOo
2.57 1.06 0.74 0.60 0.53 0.51
2.57 1.06 0.74 0.60 0.53 0.51
2.57 1.06 0.74 0.60 0.53 0.50
2.57 1.06 0.74 0.60 0.53 0.50
2.57 1.06 0.74 0.60 0.53 0.50
2.57 1.06 0.74 0.60 0.53 0.50
Table 6 Example 3 - Bounding second order models. The criterion values using the control sequence rt~ on model wh~ + (1 - w)ho for different values of a and w. The extra constraints are not included ¢ Model 1 w = 0.8 w = 0.6 w = 0.4 w = 0.2 Model 0
rt~ 1.67 19.28 38.52 45.94 48.08 48.70
(23) a n d r e m e m b e r t h a t t h e Comparing Figures 4a added constraints have though the constraints are
a =0.8
ce = 0.6
1.71 12.05 24.64 29.52 30.92 31.34
1.81 6.61 13.93 16.78 17.59 17.84
g a i n o f Y~ is m u l t i p l i e d b y 5. a n d 4b i n d i c a t e s t h a t t h e a robustifying effect even v i o l a t e d in b o t h c a s e s .
o: = 0.4
ct = 0.2
n Oo
2.00 2.97 6.37 7.71 8.08 8.20
2.25 1.12 1.98 2.32 2.39 2.42
2.57 1.06 0.74 0.60 0.53 0.50
satisfied w h e n t h e a d d e d c o n s t r a i n t s a r e i n v o k e d . T h i s is in a c c o r d a n c e w i t h C o r o l l a r y 2 since a f e a s i b l e s o l u t i o n c a n b e f o u n d a t all t h e s i m u l a t e d t i m e steps. S e c o n d , as in Example 1 the results indicate that the added constraints h a v e a r o b u s t i f y i n g effect o n t h e c l o s e d l o o p c o n t r o l l e r .
Example 2 T o f u r t h e r g a i n i n s i g h t w e u s e d o u r m e t h o d o n t h e syst e r n in E x a m p l e 2. T h e s y s t e m , o p t i m i z a t i o n c r i t e r i o n , c o n s t r a i n t s a n d i n i t i a l v a l u e s w e r e u n c h a n g e d . Figure 5 shows some typical results. In this simulation the system is g i v e n b y Y0, i.e. w ~ 0, w h i l e t h e c o n t r o l i n p u t is c o m p u t e d b y c h o o s i n g o~ = 1, cf. E q u a t i o n (8), i.e. t h e c o n t r o l i n p u t is b a s e d o n Y~. W e c a n m a k e t w o o b s e r v a t i o n s f r o m t h e results. First, t h e o u t p u t c o n s t r a i n t s a r e
Discussion T h e p r o p o s e d m e t h o d is a n a t t e m p t t o o f f e r a s m o o t h t r a n s i t i o n f r o m l i n e a r t o n o n l i n e a r o p t i m i z i n g c o n t r o l as for example nonlinear MPC. In our opinion a smooth t r a n s i t i o n is i m p o r t a n t t o p r o m o t e i n d u s t r i a l a p p l i c a tions of nonlinear model-based optimizing control. The f e a t u r e s t h a t u n d e r l i n e t h e s m o o t h t r a n s i t i o n are: (i)
Interpolating optimizing process control: B.A. Foss and S.J. Qin •
1.2
b
1.1~
a
o.a ~ ...............................
o.a[-
0.6
,I. I ! .I. i
0.6 ....................................... 0.4 ...........
: .........................
,i.
0.2 "'h
i .............
0.4 0.2
t
o ,.4.'
:-
~
....................................
i
o
I
;
-o.2 . . ,
.-o.2 " i
•- 0 . 4
.-0.4
.-o.e .-.f ..~ . . . . .
....................
.-0.6
....... I
, '
__
i
......... i .............
'l. . . . . . . . . . . i . . . . . . . . . . . . .
::
b
1.2
i
0.8
0.8
0.6
O.E
0.4
0.4
o.2
oa
: .............
:
:
y
137
.
.
: .............
: .......................................
c
! -0.8
i 10
-0.8 20
.30
I
-0.~ 10
lime
20
30
,: o.8..
.
2. 1.5 !i
.
' . ~
Y o.e . . . . . . . ~........................
0.5 •
O"4 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
0
0.2 ........... :........................ 0-. ~
-0"20
1'0
20
lime
30
-o.s -1
-1.5,
10
20
30
.me
Figure 4 Example 1 - MPC. Same simulation as in Figure 3, but the gain of the simulated system is increased by a factor of 5. In a the added constraints are invoked, while these are omitted in b. The figures show the system output and the control input (note the axis)
0
10
T h e c o n t r o l a l g o r i t h m uses w e l l - k n o w n o p t i m i z a U o n b a s e d on l i n e a r m o d e l s as its b u i l d i n g - b l o c k s . (ii) A controls sequence is c o m p u t e d b y s m o o t h i n t e r p o l a t i o n , a c o m m o n l y a p p l i e d engineering principle, T h e p r o p o s e d c o n t r o l l e r can be seen as an extension o f linear o p t i m i z i n g c o n t r o l f r o m b o t h a theoretical as well as a p r a c t i c a l perspective, First, we e l a b o r a t e on the t h e o r e t i c a l results. A key a n d novel feature o f the p r o p o s e d c o n t r o l l e r is the a d d e d c o n s t r a i n t s i n c l u d e d for the t w o minim i z a t i o n p r o b l e m s , i.e. c o n s t r a i n t z ~ ~ Y i n (15), a n d ~ X in E q u a t i o n (16). This m a k e s it possible to guarantee robust feasibility on some given o p t i m i z a t i o n h o r i z o n when the b o u n d i n g m o d e l s are linear F I R - m o d e l s p r o v i d e d t h a t there exists a feasible c o n t r o l i n p u t sequence for each o f the two b o u n d i n g linear models. It s h o u l d be n o t e d t h a t the
20
30
time
m o d e l set Y~,,.m a y have autoregressive terms even if the b o u n d i n g m o d e l s are F I R - m o d e l s , since the w e i g h t i n g f u n c t i o n m a y d e p e n d on xi. A n easily c o m p u t a b l e u p p e r b o u n d can be f o u n d for the o p t i m a l i t y criterion in the case o f n o n l i n e a r F I R - m o d e l s . This p r o v i d e s a tool to investigate r o b u s t p e r f o r m a n c e o f the p r o p o s e d c o n t r o l l e r scheme for this m o d e l class. C o r o l l a r y 2 shows r o b u s t feasibility for the 1-step a h e a d p r e d i c t i o n for a n y t y p e o f b o u n d i n g linear m o d e l s if zc~ a n d rc~ a r e c o m p u t e d using the a d d i t i o n a l c o n s t r a i n t s . This implies that an M P C cont r o l l e r b a s e d on t h e p r o p o s e d a p p r o a c h a l w a y s will be feasible if there exists a feasible c o n t r o l input sequence for each o f the t w o b o u n d i n g linear m o d els at each time-step. T h e r e a s o n for this is t h a t the o p t i m a l c o n t r o l sequence is r e c o m p u t e d at every time-step. It will be, however, neither o p t i m a l n o r necessarily stable if the c o n s t r a i n t s are v i o l a t e d on the o p t i m i z a t i o n h o r i z o n b e y o n d the first p r e d i c t i o n
•
•
30
Figure 5 Example 2 - MPC. The system is given by ,~, i.e. w -- 0, while the control input is computed by choosing a = 1 (8), i.e. the control input is based on ~ . In a the added constraints are invoked, while these are omitted in b. The figures show the system output
b •
20
time
Figure 3 Example 1 - MPC. The simulated system is given by 3~], i.e. w = 1, while the control input is computed by choosing a = 0 (8), i.e. the control input is based on 3"o. In a the added constraints are invoked, while these are omitted in b. The figures show the system output and the control input
,
--0.; 10
tklle
step. T h e theoretical results a p p l y w h e t h e r the u n d e r l y i n g system h,, e E,,. is linear or nonlinear. The richness o f the n o n l i n e a r m o d e l class covered by Y~,~a n d the relevance o f this m o d e l class for process c o n t r o l p r o b l e m s need f u r t h e r investigations.
S e c o n d , we discuss the s i m u l a t i o n results. •
•
E x a m p l e 1 with b o u n d i n g F I R - m o d e l s s u p p o r t the t h e o r e t i c a l findings. It also clearly d e m o n s t r a t e s , cf. Tables 1 a n d 2, the i m p o r t a n c e o f the a d d e d constraints. F u r t h e r , the c o m p u t e d u p p e r b o u n d on the o p t i m a l i t y criterion seems to be r e a s o n a b l e in the sense t h a t it is n o t overly conservative. E x a m p l e s 2 a n d 3 investigate the results b e y o n d the b o u n d s o f the t h e o r y . This is o f interest since practical use o f c o n t r o l l e r s u s u a l l y violates theoretical a s s u m p t i o n s at s o m e stage. In E x a m p l e 2 the b o u n d i n g m o d e l s are quite similar in the sense that
138
•
•
Interpolating optimizing process control: B.A. Foss and S . J , Qin
they both are second order damped systems, while the model characteristics are more divided in Exampie 3. The findings in Example 2 indicate that constraint handling need not be difficult when the bounding model are similar, in this example control c o m p u t a t i o n without the added constraints gave only feasible solutions. This picture changed significantly in Example 3, the positive effect of the added constraints was very pronounced here, even though there is no feasibility guarantee in this case. Corollary 1, i.e. the result on feasibility of the 1-step ahead prediction, is supported by the simulations, There were no constraint violation at time i = 1 with the added constraints, see Tables 1 a n d 3. This was not necessarily the case if these constraints were removed as can be seen in Table 2. The closed loop simulations substantiate the open loop findings. The added constraints do seem to improve robust stability and performance. It is important to note that the feasibility guarantee depends on the existence of a feasible control input sequence for each of the two bounding linear models. This is demonstrated in Figure 4a.
Finally, we elaborate on some additional issues, •
An open question on the proposed approach is the richness o f the model set Y.,,.. It might be necessary to include more than two bounding linear model, i.e. let the hw be a convex combination of several, say M, linear models. In this case the model set can comprise a large class o f nonlinear systems as shown in Johansen and Foss 7. The problem that arises is the complexity increase in our algorithm since the number o f added constraints to compute each of the
•
optimal control sequence increases. In addition, the control sequence zr~ will be a convex combination of M control sequences. There exists methods and tools for identifying mod-
•
els on a format like Equation (4), see e.g. Reference 8. An alternative approach to the one pursued in this paper is to compute a control sequence based on one nominal model with added constraints related to one or more bounding models.
•
This p a p e r does not provide a closed loop stability p r o o f for the proposed method. The findings, do, however, provide a basis for seeking such a result.
•
Nonlinear optimizing control provides a nonconvex optimization problem. By this, there is no guarantee
that the global o p t i m u m can be found. Hence, an alternative approach is to search for a suboptimal solution with some performance guarantee. This can be done by limiting the search to the set rI~ proposed in this paper, a e [0, 1] might be computed in one of three ways. (1) A gain-scheduling approach can be applied if it is possible to schedule a on some process variables. (2) a is assumed to be time-invariant, and it is computed using an off-line or on-line estimation scheme based on some performance measure. (3) a is time-varying and is computed using a recursive on-line algorithm. This m a y pave the way for a reliable adaptive M P C controller, since only one time-varying p a r a m e t e r is required. In practice, robustness of adaptive control is difficult to obtain when m a n y parameters are estimated on-line.
Conclusions A new model-based optimizing controller for a set of nonlinear systems is proposed. R o b u s t feasibility with respect to hard state and/or output constraints is guaranteed on the optimization horizon for a set of nonlinear FIR-models, while significant feasibility improvements are encountered in the case where the bounding models are general linear state-space models. Finally, the method guarantees robust feasibility for a 1-step ahead prediction in the general case. This can be of interest in MPC-applications.
References 1. Lee,J. H. Recent advances in model predictive control and other related areas. In Preprints CPC-V, Lake Tahoe, USA, 1996. 2. Rawlings, J. B., Meadows, E. S. and Muske, K. R. Nonlinear model predictive control: A tutorial and survey. In Preprints [FAC Symposium ADCHEM, Kyoto, Japan, 1994, pp. 203-214. 3. Qin, S. Joe and Badgwell, T. A. An overview of industrial model predictive control technology. In Preprints CPC-V, Lake Tahoe, USA, 1996. 4. Mayne,D. Q. Optimization in model based control. In Preprints IFAC Symposium DYCORD, Helsingor, Denmark, 1995. 5. Genceli, H. and Nikalaou, M. Design of robust constrained model-predictive controllers with volterra series. AIChEJournal, 1995, 41, 2098-2107. 6. Bequette, B. W. Nonlinear control of chemical processes: A review. Ind. Eng. Chem. Res., 1991, 30, 1391-1413. 7. Johansen, T. A. and Foss, B. A. Constructing NARMAX roodels using ARMAX models. Int. J. Control, 1993, 58, 1t25-1153. 8. Foss, B. A. and Johansen, T. A. Local modelling as a tool for semi-empirical or semi-mechanistic process modeling. In Neural Networks for Chemical Engineers, ed. A. Bulsari, Elsevier, Amsterdam, 1995, pp. 297-334.
•
~LSEVIEX
PIh S0959-1524(96)00023-6
Interpolating optimizing process control Bjarne A. Foss *t and S. Joe Qin* *Department of Engineering Cybernetics, Norwegian University of Science and Technology, Norway. email: [email protected] ~Department of Chemical Engineering, The University of Texas at Austin, USA. email: qin@che, utexas.edu Received 2 April 1996; revised 23 September 1996 A new model-based optimizing controller for a set of nonlinear systems is proposed. The nonlinear model set is based on a convex combination of two bounding linear models. An optimal control sequence is computed for each of the two bounding models. The proposed control algorithm is based on a convex combination of the two control sequences. A novel feature in these two optimizations is an added constraint related to the feasibility of the 'other' bounding model. The control algorithm can for example be used in model predictive control. We provide robust feasibility guarantees and an upper bound on the optimal criterion if the bounding models are linear FIR models. Further, simulation examples demonstrate significant feasibility improvements in the case where the bounding models are general linear state-space models• The proposed method guarantees robust feasibility for a l-step ahead prediction in the general case. This can be of interest in MPC applications. © 1997 Elsevier Science Ltd. All rights reserved
Keywords: model-based control, optimization, convexity, model predictive control, nonlinear systems, robustness, interpolation The combined use of dynamic models and optimization for process control offers a concept in which process knowledge can be linked to operational goals formulated by some optimization criterion. This concept has seen widespread use, particularly through the applications of model predictive control (MPC). MPC refers to a class of algorithms where an optimization problem is solved repetitively, at every new time-instant. Only the first part of the computed control sequence is applied to the system since a new optimal control sequence is computed and applied at the next time-step. Several reviews of MPC technology exist, see for example Lee ~, Rawlings et aL 2 and Qin and BadgwelP, the latter emphasizing industrial use of the technology. The interaction between control and optimization is discussed in an iliaminating way in Mayne 4. In this paper the divide between convex and non-convex problems is emphasized. In particular the problems arise by using nonlinear models since the optimization problem in these cases in general becomes non-convex, State-of-the-art optimizing control is based on linear dynamic models and linear constraints on the control inputs and system outputs. Nonlinear optimizing control has been studied by Rawlings et al. 2 and Oenceli and Nikalaou 5. Further, some approaches were
described by Bequette 6 in a somewhat earlier paper. These control strategies normally result in a nonconvex optimization problem. In this paper we explore optimization control based on nonlinear models. In particular, the goal is to derive an approach which can offer a smooth transition from linear optimizing control to nonlinear and robust optimizing control. By smooth transition we mean an approach which does not invoke the full mathematical 'machinery' of nonlinear optimizing control in the general case. This is particularly important from an industrial viewpoint as it simplifies the transition from the application of linear to nonlinear optimizing control. The smooth transition is accomplished by constraining the nonlinear optimization problem along three axes. First, the set of nonlinear models is limited. Second, constraints are added to the optimization problem to enhance feasibility of an optimal solution for a s e t of bounding linear models. Third, the control sequence is computed as a convex combination of control sequences computed on the basis of the bounding linear models. The remainder of this paper is structured as follows: in the next section we formulate the problem. Following this, the theoretical foundation for the proposed control algorithm is developed. In the following section three simulation examples are investigated to explore the theoretical findings, and to study the proposed method beyond the assumptions given by the theory. The pro-
*On sabbatical leave at The Department of Chemical Engineering,The University of Texas at Austin, USA 129
Interpolating optimizing process control: B.A. Foss and S.J. Oin
130
posed method is after this investigated in closed loop control. Finally, a comprehensive discussion is provided before the conclusions finalize the paper,
Problem definition Background
Assume an optimality criterion given by: (1)
¢ ( ~ , Z ) = ~.,l(x~+~,ui) ~:N
~ = arg min~,n¢(Zc, Z)
(6)
subject to the constraints Z ~ x; (5) using one of the models (2), (3) or (4). The optimization problem above is generally nonconvex when a nonlinear model within Y,,. is used. This can cause problems, especially in on-line control where computation time can be an issue. In MPC the minimization problem (6) is solved repetitively, at each time-step, with new initial conditions. Only the first control value of the zr sequence is actually applied to the process. Furthermore, it is typical to parameterize the control sequence as follows:
where zr = {u~ .... , u~ ..... u~}, M < N - 1
(7)
-- {.To,. ..,UN_I} , Ui E U cc_.R mu, I"[-~-U x . . . x U Z = {xT ..... xT}, x~ ~ X ~ R ' ~ , U, X - convex set
Z = X ×... × X
l • R m"×% --+ R ÷, l - convex function, i ~ I N I N = {0,..., N - 1}
The optimality criterion is defined on a time horizon from 0 to N. u~ denotes the control input to a system. The control input u~ is constant during the time span [i,i + 1). x~ denotes the system states at time i. Since U and X are convex sets, rI and Z are convex sets. Further, 0 is convex since l is convex, We define two linear state-space models (2~, Y~0)and a set of nonlinear models (Y~w)based on the interpolation of the two linear models. ~., 1: h l ( x i ÷ l , x ~ , u , ) = x~+l- f ( x ~ u ~ ) = B~u,. O, i ~ I,v =
xj+ x - A i x ~
(2)
0: h o ( x i . l , x i , u i ) = x i + x - fo(xi , ui) = x i ÷ I - Aoxi -
B0u~ = 0, i ~ I~.
(3)
Z,.---- {hw" hw(Xi+l, xi, ui) = w(xi, ui)hl(Xi*t,xi, ui) + (1 - w(x,, u,))ho(xi+ I, xi, ui) = 0 (4)
V i ~ Lv, w ~ W, W = { w ~ [0,1] Vx~ ~ X Vu~ ~ U}}
Cw(x~,u~)~
h, is a nonlinear function, constructed as convex cornbinations of h~ and h0. It should be noted that the continuous function w in general will depend on the states and the control inputs. We further observe that the interpolated model hw ~ ]~,,. is not generally a convex function since w may be any continuous function bounded by [0,1]. Finally, we define initial conditions, x0 - given
(5)
The problem we want to address is to minimize Equation (1) with respect to n based on the different constraints discussed above, hence we want to solve:
This means that the control input is constant during the last part of the control input sequence. The criterion function (1) does not cover all possible criteria; penalizing changes in the control input is for example not included. This type of change does not influence the results in this paper as long as the criterion function remains convex.
Interpolating control
We will in the following describe the proposed controller. First, we define two control sequences 7cI E 1-1 and Jr0 e FI and a set of interpolating controllers. H~ = {n~ " n~ = azr~ + (1 - a)n0,
V a ~ [0.1]}
(8)
The control sequence rc~ forms the basis for the controller. It should be noted that zr~ is feasible since it is based on interpolation within a convex set, i.e. FIo _C I-I. o~may in general vary from one time-instant to another. For later convenience we pursue with some definitions. Z~ and Z0 are the state sequences on the horizon { 1..... N} obtained by applying some control sequence n I --, Z ~N} I is the 1-I to Yl and 5".0, respectively. Z kl = {Z k~,state sequence obtained by applying the control sequence Jr~, I ~ { 1,0} to the system Zk, k ~ { 1, 0}. Further, we define the set o f states associated with the model set Y,, and the control sequence set I-I,. Z ,a ~ {Z,~~ " Z,~ = {Xw~l,...,Xwu},C~
Vh w ~ Ew,O~ ~ [0, 1]} (9)
Later we will in particular study the interpolated controller set I-I, where zq and zr0 are computed on the basis of the following two optimizing problems. zr~' = arg min=~ n ~(n',Z~ )
(10)
subject to the initial conditions (5), constraints Z~ ~ -K and Z~ given by Equation (2). o = arg m i n ~ n ~(zc, Z0 ) zc0
(ll)
Interpolating optimizing process control: B.A. Foss and S.J. Qin subject to the initial conditions (5), constraints Zo ~ X, and ~0 given by Equation (3). If we choose zq = zc[ and ~0 = gg the state sequences Zl~ and ~ are the optimal sequences for ~1 and ~0, respectively. Further, following the above notation 0~ is the criterion value obtained by applying the control sequence ~ , I s {1,0} to s y s t e m Y , , k s { 1 , 0 } . T h i s means that ~t and ~ are the optimal criterion values for El and ~;0Note that we assume that the state sequences discussed above are based on the initial value (5).
Analysis We will in this section analyze the use o f the interpolating controller on the model set Y~,~.. Guaranteeing feasibility There is no guarantee that the state constraints are satisfied on the horizon {I ..... N} when the control sequence n'~ e FI~ is applied to the system h,,. E ~,,,,. We will in the following include restrictions on the choice of zq and ~ , and on Y~,, so as to ensure that X~. C X. Before we do this we state the well-known connection between FIR-models and state space models. Proposition 1
Theorem 1 Assume two control sequences zc~ and zr0 such that ~ ~ X Vl e {0,1}, k s {0,1} for the initial conditions (5). Further, define a control input sequence set rI~ by Equation (8). If the linear models E1 and ~0 are FIR-models, and z~ ~ rI~, then Z~,. e x . ~ is defined in Equation (9). Proof." We first use Proposition 1. Since the bounding models (2) and (3) are FIR-models, we transform them into the equivalent input-output model structure (12). The model set Y~,,. contains nonlinear FIR-models and may, hence, also be transformed into a FIR-structure. Let z~ denote the j t h element of the predicted output vector y~+~ of the FIR-model (12) at some i e I N if we apply the control sequences rc~, l ~ { 1, 0} to the system Y~k, k ~ { I, 0f. The assumptions on Jr~ and zc0 guarantee that z~ ~ ~ where X = X~ × ... × Xmx" Let ~ denote the j t h element of output vector y,+~ of the FIR-model (12) by applying a control sequence zc~ lq~ to the system h,~ e Z,,. Let c~, co and Cw denote the j t h row vector of C~, C O and C,, w ~ W. respectively. We first show that c~T~" e X,, c~ ff~, = c~ ( a ~ + (1 - a)ff0~) = ac~ ff~i + -~ + (1 - a)z°~ (I - a)c~ ff0~ = m~
ke{l'O'W}'ieIN (12)
c ff~, = w(x,, ui)c I ff~, + (1 - w(x,, u,))c o ff~, is equivalent to a state-space model on the form (2), (3) or (4) if:
.
]
0O/ ~
=. A~
I Q~
...
I = .
C~L ~/ ,0
.
.
0 . " ...
.
.
(14)
We have above shown that c~ ~=: ~ Xj and c o ff=~ ~ X~. Since X~ is a convex set, the convex combination c ff=, must lie in X.
.
. " I
~ ! B~ =
(13)
c~ ~ , lies in ~ since it is a convex combination of two elements in the convex set Xr An identical argument may be formulated for co ff~,. We now choose some element cz e [0, 1] and show t h a t c u'~, e X ~ f o r a n y w ( x , , u , ) ~ [0,1].
An FIR-model of this form
)',÷~ = C ~ o u . ' + ' " + C ~ u , - c = C k U , ' 3', e Y ~ R my
131
k ~ {1, 0, W} '
The above result cannot be generalized to arbitrary A1 and A0 matrices. Further, it is critically dependent on the assumption that each of the control sequences n~ and ~ applied to either of the models ~t or ~0 obey the stateconstraintonthehorizon{1,...,N}. This theorem will form the basis for several corollaries. First, we use it in the case where ~c~ and ~ are based on the solution of some optimization problem. Corollary !
and y~ = (0 ..... 0, /)TX~. We note that C,,~, or C,, (w ~ W), and correspondingly A~ and B~ in general depend on the control inputs and states. If we limit the bounding linear models to be FIRmodels, and include feasibility constraints on the choice of zq and rc0 it is possible to guarantee feasibility in the sense that ~ C_ )2. This is shown below.
Given the following control sequences
~ = arg minion 0(z~,Z~)
(15)
subject to the initial conditions (5), constraints Z~ s x (Z~ is the state sequence obtained by applying ~c~ ~ FI to Y~t), and Z~ ~ X. ~ = arg min.~ n ¢(zc, Z0)
(16)
Interpolating optimizing process control: B.A. Foss and S.J. Qin
132
subject to the initial conditions (5), constraints ~ e 22 and ~ e .~ If the linear models ]~i and Y+0 are FIR-models, and zc~ ~ FI~, then Z~ ~ X.
Remark This corollary shows that the 1-step ahead prediction will satisfy the state constraints for any linear m o d e l s Y+I and Y+0 if zc~ and ~0 are computed using the additional constraints as in Equations (15) and (16). This implies that an MPC controller based on 7r~ and Equa-
Proof" Since the control sequences defined in Equations
tions (15) and (16) always will be feasible even though the constraints may be violated on the optimization horizon. The reason is that the MPC controller only
(15) and (16) satisfy the state constraints ~ e X V I e {0,1}, k e {0,1} in Theorem 1, the corollary follows directly from Theorem 1.
Remark This corollary states the important result that an arbitrary control sequence zr~ ~ I-I~ applied to an arbitrary model h,,. ~ Y+,,.will satisfy the state constraints on the horizon { 1.... , N} provided additional hard constraints are included in the optimization, i.e. including o O Zo ~ X for computing zrl, and Z1 ~ x for 7r0. If we limit our attention to 1-step ahead prediction, however, the result can be extended.
Given the control sequences tot, ~ and zc~ as defined in Equations (15), (16), and (8). The state sequence which arises by applying lr~ to h,. ~ Y+,,. is denoted ~ . = {~,.1..... ~tN}, cf. Equation (9). Then ~,.t e X. Corollary 2
Proof" Let z: denote element j of state vector x: which
k +1 arises by applying the control input u~0, l e { 1, 0} to the system z~k, k e { 1, 0}, and ~. denote elementj of state vector x ~t which arises by applying the control sequence zr~ ~ FI~ to the system h,,. ~ Y.,,. Let at, ao and a,. denote the j t h row vector of A1, Ao and A,,., and bl, b 0 and b,, denote the j t h row vector of B1, B0 and B,,.. Further, ul0 and u00 denote the control input of z~, and ~ at time 0. We need to show that z~ ~ ~ , where X =Xt × . . . ×
applies the first control input before recalculating the control input at the next time-step. Robust performance
Theorem 1 provides a possibility to obtain an upper bound on the following worst case scenario. The upper bound is defined by: ~ = suph,:z ,,.,,~ ~n~ ¢(zc,+,Y~)
(19)
This result is formulated by the following lemma. Corollary 3 Suppose the control sequence zc~ e FI~ is applied to the system h,,. e Y+,.with initial conditions (5). Y/ and z~0 are FIR-models. Further, suppose the optimality criterion function 1 is given by some separable norm function, i.e. l(X~+l, u,) = , , Ily~+llI?2 + Ilu~l[~, where y~ is the output of the FIR-model formulation (12). Then an upper bound for ~ in Equation (19) is given by:
(20)
~ <_ ~ Cu--~ + u~.
,~z~. where Cu,. = maxk.~,olHf~u,[[~ and ~
maxt~l,0~llu,I]~.
X,..
We first show that a 1 Xo + b 1 U~o ~ Xj.
Proof" We first compute an upper bound for
IlY~+II[~= IIC,~'oilI~= II(wCx + (1 - wf0))u+ell~~
a a Xo + b I U~o = a I (~Xo + (1 - a)Xo) + bl(aUto + (1 -
]ly~+lll~.
o:)U0o)
= a ( ~ 1 X 0 + btUl0 ) +
< [w]ICIUc~:HQ+ (1 - w)lIC0~+~,llo]2 < max~ u,0~[lf~u2t~
(1 - ot)(a~x o +
b 1 u0o) l + (1 -
cOz~
(17)
~" ~(~'1
~ =
a I x0 + b~ u~0 lies in X~ since it is a convex combination of two elements in the convex set X~. An identical argument may be formulated for aoX o + bou~o. We now choose some element ~ ~ [0,1] and show that a,.xo + b,. u~0 e Xj for any w(x~,u~) ~ [0,1].
maxk~ll,0~llC,[m71, + (1 -o0ff0i]]]~
< max~{,,01[o~lC,~l,[]~ + (1 < max,.:, i i.o~HC~u t~llb = Cu, A upper bound on
a,,xo + b,~uoo = w(x,,u,)(a:o + blu~0) + (1 - w(x,u,))(aoX o + bongo )
(18)
Since atx o + bt U~o e Xj and aox o + b o u=o ~ Xj and X] is a convex set, the convex combination a,,.Xo + b,u~o must lie in X~.
a)[IC~+0,ll~] = (21)
Itu,lt~ is easily computed.
Hu,IL~ = Ilu=~ll?~ = I[~ul, + (1 - c~)u0,ll~ < [c~lut+ll~+ (I - a)llu0:ll~]: < max]~ll,01lluz,ll~ = u:
(22)
Interpolating optimizing process control: B.A. Foss and S.J. Qin
133
By these expressions the corollary is proved. The important consequence of the above lemma is the fact that an easily computed upper bound on robust performance can be found. This can hence be used as a measure for robust performance. To compute the bound we use the output Yi instead of xi. This is no important limitation since y~ forms the output part of the states of a FIR-model. A weighting of y~ is, hence, the obvious choice for the criterion,
We will first study the feasibility of the proposed controller, Equation (8), by studying the output trajectoties. This is done in three steps.
Examples
2.
We will in this section further explore the findings in the above sections through three numerical examples. All the models, except for Y-0 in example 3, are typical models found in the process control problems. In the first example we will investigate the method on a nonlinear FIR-model. In examples 2 and 3 we will study the method beyond the assumptions given by the theory. In particular, we investigate situations where there is no theoretical guarantee for a feasible solution,
3.
1. The control sequences Jr~ = (0.90, 0.60, - 0.70, 0.33, 0.06, -0.11, 0.11) r and n o = (0.90, 0.60, -0.70, 0.13, 0.28, -0.17, 0.08) r are computed using Equations (15), (16), i.e. with the additional constraints in place. Six interpolating control sequences zr~ are generated by choosing a = {0.0, 0.2, 0.4, 0.6, 0.8,
1.0}. We choose six models within ~,,. by choosing constant interpolation weights w = {0.0, 0.2, 0.4, 0.6, 0.8, 1.0}. The six models are simulated using the six interpolating control sequences as inputs. This gives a total of 6 × 6 = 36 simulation runs. The system is simulated from 0 to 7, i.e. equal to the optimization horizon. The system outputs for each simulation run is checked against the constraints, i.e. [0.00, 1.01]. The result is shown in the matrix in Table 1. The 36 columns relate to each of the simulation runs, while the eight lines are related to time instants
{o,1,...,7',. E x a m p l e 1. F I R m o d e l s
We assume that a single-input single-output system is given by some model h,, ~ Y~,., cf. Equation (4), with bounding FIR-models given by: k ~ {1, 0} b~ = {-0.10, 0.30, 0.40, 0.30, 0.10} 50 = {0.20, 0.20, 0.15, 0.10, 0.05}
Y~+z = CkoU~+ ... ÷ C~4Ui-4 = CkU~,
(23)
Proposition 1 gives the equivalent state-space formulation. The two bounding models resemble a system with an inverse response and a steady state gain of 1.0, and a first order response and a steady state gain of 0.7, respectively, We use the well-known quadratic optimality criterion, add hard constraints on the system output and control input, and assign initial conditions. x 0( Jr, 2") = 0 . 5 ~ y~ + ru~_ 1 ~=~ N = 7 r = 0.01 U = [-1,1]
0 indicates no constraint violation, while 1 or 2 signifies a violation of the minimum or maximum allowable output, respectively. Table 1 shows that there is no violation of output constraints. This is in accordance with Theorem 1. Corollary 3 gives an upper bound on the optimality criterion, for this example we compute 0 = 1.93. The criterion values for the 36 simulation runs in Table 1 lie between 0.47 and 1.61. This indicates that the upper bound need not be very conservative. To study the impact of the additional constraint, i.e. Z0 ~ Y to compute 7r~ and Z ~ x to compute ~ro, we remove these constraints and perform the computations in items 1-3. This means computing lr~ and iv0 according to Equations (10) and (11). These results are shown in Table 2. There are in these cases many violations of the output constraints. In the latter case the control sequences are rc~ = (0.90, 0.60, -0.70, 0.33, 0.06, -0.12, 0.09) T and Jr0 = (-1.00, -0.34, 0.34, 0.17, 0.00, -0.17, 0.00) T. They are very different from each other. This is not the case in item 1 above, in that case the control sequences are identical except for the last control input. The reason for this are the added constraints which limit the possible choices of control inputs.
X = [-1, 11 × [-I, 11 x [-1, 11 x [1 - 1, 1] x [0.00, 1.01] x o = (1, I, 1,1, O) T
The state-space model has five states. The first four states are delayed control inputs while the last state is the system output, equivalent to y~ in Equation (23). Hence, the state constraints X limits the system output to [0, 1.01]. Note that the system output does not depend on the last element of x 0 since it does not have any auto-regressive terms.
E x a m p l e 2." Bounding second order models
We again assume that a model set Y~,, defines a singleinput single-output system. In this example the bounding linear models are given by
A,=
I'0.00 0.00 0.00
0.00 0.00 -0.79
0.00 / 1.00 / 1.78)
(1.001 /~ = / 0 . 0 0 [ ~0.01)
Interpolating optimizing process control: B.A. Foss and S.J. Qin
134
Table 1 Example 1 - Bounding FIR-models. Each column represents a simulation run using the control sequence zt~ on a model h,, e 5".. The
system outputs of each simulation run are checked against the constraints at each time-step. 0 indicates no output constraint violation. 1 or 2 signifies a violation of the minimum or maximum allowable output, respectively Time i=0 i=1 i=2 i=3 i=4 i=5 i=6 i--7
Model 1
w = 0.8
w = 0.6
w = 0.4
w = 0.2
Model 0
~o ..... =o0
=o, ..... ~ 00
=o ..... ~ o0
~o ..... ~ o0
~ ,0 ..... ~ o0
~ o ..... ~ 00
000000 000000 000000 000000 000000 000000 000000 000000
000000 000000 000000 000000 000000 000000 000000 000000
000000 000000 000000 000000 000000 000000 000000 000000
000000 000000 000000 000000 000000 000000 000000 000000
000000 000000 000000 000000 000000 000000 000000 000000
000000 000000 000000 000000 000000 000000 000000 000000
Table 2
Example 1 - Bounding FIR-models. ~ c h column represents a simulation run using the control sequence ~ on a model h~ ~ Z~. The ~stem outputs of each simulation run are checked against the constraints at each time-step. 0 indicates no output constraint violation. 1 or 2 signifies a violation of the minimum or m~imum allowable output, res~ctively Time
Model 1
w = 0.8
0
i=0 i=1 i=2 i=3 i=4 i=5 i=6 i=7
0
000000 022222 000000 000001 000111 000001 000000 000000
(0.00 0.00 A o = |0.00 0.00 ~,0.01 -0.855
w = 0.6 0
000000 000002 000000 000001 000011 000001 000000 000000
0.00) 1.00| 1.85)
The third state is the system output, hence Yi = (0, 0, 1)a x~. 21 is a second order model with two time constants o f a b o u t 6 and 15 time-steps, no time-delay and a steady state gain o f 1.0. 20 is a second order model with two time-constants o f a b o u t 10 and 20, a time delay o f 1 a n d a steady state gain o f 2.0. T o elaborate on the model set Y,. let w (x,, ui) ~ [0,1] be a constant. Hence, the model set Y~,,.consists o f only linear models. The steady state gain, time-delay and time-constants for this model set are shown in Figure 1. We observe that there is a close to linear change in all the characterizing variables. The optimality criterion, hard constraints on the systern o u t p u t and control inputs, and initial conditions are given below,
(24)
[0, 1.Ol].
0
000000 000000 000000 000000 000011 000001 000000 000000
Model 0
0
0
000000 000000 000000 000000 000011 000001 000000 000000
0
000000 000000 000000 000000 000000 000000 000000 000000
The control sequences r~ = (-7.48, -5.14, -3.28, - 1 . 8 7 , - 0 . 8 9 , - 0 . 2 8 ) r a n d ~ = (-7.49,--4.84, 2.80,-1.35, -0.43, 0.00) a- are c o m p u t e d using E q u a t i o n s (15), (16), i.e. with the additional constraints in place. We compute six interpolating control sequences zr~ by choosing o~ = {0.0, 0.2, 0.4, 0.6, 0.8, 1.0}, and choose six models within Y~,, by selecting constant interpolation weights w = {0.0, 0.2, 0.4, 0.6, 0.8, 1.0}. Again we simulate a total o f 6 x 6 = 36 runs a n d check the system outputs against the constraints. There is no violation at any time-step in any o f the 36 simulation runs. We also c o m p u t e the control sequences omitting the additional constraints, i.e. using Equations (10), (1 I). In this case we obtain identical control sequences as above. Hence, in this case the additional constraints m a k e no difference on the c o m p u t e d control sequences.
Example 3: Bounding second order models T o m a k e the p r o b l e m m o r e challenging we will in this example study a set o f systems where the b o u n d i n g models have different types o f behavior. Y'I is equal to ]~ in the above example, while Y0 is a second order model with oscillatory modes. The 'time-constant' of the oscillatory models is 29 time steps and the oscillation period is 5-6 time-steps. The steady state gain of Y~0 is 1.0 and the time-delay is 1. The b o u n d i n g models are described by:
~=(],,,,; The state constraints Z limits the system o u t p u t to
w = 0.2
0
000000 000000 000000 000001 000011 000001 000000 000000
(1.001 B0= |0.00| t0.00j
x ~(rc, Z) = 0 . 5 ~ y ~ +ru~_l ~=~ N = 6 r = 0.01 U = [-10,10] X -- (-oo, oo) × (-oo, oo) × [0.00, 1.01]
w = 0.4
0
(0.00 "
t,0.00
0.00
0.001
(,.001
0.00 1.00 /
B,--/0.00 !
-0.79
1.78)
~,0.01)
Interpolating optimizing process control: B.A. Foss and S.J. Qin a 1.8 1.6 1., 1.z
r,~ ~ , ~ , t l " o.lo ~ ~
i ~-.._~ 0.8 ~ 0.6 o.,
~m,,on,ta,,a'0.10
~d,~,y
~
0.2 oi* o'.s ole o17 o18 o'.~ wa~,~ngf~,,~on,* Figure 1 Example 2 - Bounding linear models. The gain, time-delay and time-constants of the models with constant w-valuesare depicted °o
o11
ola
o'.3
"] 0.00
0.00
0.62
-0.32
0=000 l
0.00
000 100! 0.70)
(
1 1.00
0--I000 / ~,0.00)
T o get a feel for the model set ]~,, we simulate a step response for six models with constant w-values, cf. Figure 2. The steady state gain does not change; the dynamics do, however, change a lot. The dynamics are particularly sensitive to w close to 1. The optimality criterion, constraints and initial conditions are equal to those in Equation (24). The output is again equal to y~ (0, 0, 1)Xx,. The control sequences z:] = (-0.61, 0.52, 0.00, 0.00, 0.00, -0.84) 7 and ~0 -- (-0.59, 0.50, 0.01, 0.00, 0.00, 0.00) ~ are computed using Equations (15), (16). Again we simulate a total of 6 × 6 = 36 runs and check the system outputs against the constraints. The result is shown in Table 3. We observe that in many of the late time-steps (4-6) there are violations of the minimum value o f the system output.
..
~.4
,o.~ 1.a 1
.~o., E
~0.6 o.4
~'*~
o.2 o
1'o '~'s
2'0
2's ~ 3's ~ ~s so r~,,,op, Figure 2 Example 3 - Bounding linear models. There is a step in u (u goes from 0 to 1) at time 1. The initial conditions are x = [0, 0, 0]v. Step responses are shown for the six models used in the computations in Example 3. i.e. w = {0.0, 0.2, 0.4, 0.6, 0.8, 1.0}
135
To study the impact of the additional constraints, we remove these constraints. This gives the control sequences ~ = (-7.48, -5.14, -3.28, -1.87, -0.89, -0.28) w and ~ = (-0.59, 0.50, 0.01, 0.00, 0.00, 0.00)L The results are shown in Table 4, the number of violations is substantially larger than in Table 3. The numbers behind these tables accentuate this difference. The minimum output constraint, i.e. 0, is only marginally violated in Table 3. Most of the violating outputs lie in the range [-0.10, -0.05] with a minimum value o f - 0 . 1 6 . In Table 4 most of the violating outputs lie in the range [-3, -2] with a minimum value of-5.92. To further illuminate the impact of the additional constraints we show the criterion values for the 36 simulation runs in the cases with and without the additional constraints in Tables 5 and 6, respectively. The (1, 1) and (6, 6) entries show the nominal cases for Yl and 5".0, i.e. ¢~] and ¢0. We observe, as should be expected, that the optimal criterion value in these nominal cases may increase due to the additional constraints. More interesting, however, is the fact that the added constraints seem to have a significant robustifying effect. The worst case criterion values show a dramatic improvement, 48.70 to 2.57. The improvements are especially pronounced for high c~-values, i.e. by comparing the lower left part of Tables 5 and 6.
Closed loop behavior As discussed earlier the proposed method can be used for closed loop control. In this section we will investigate this by extending the simulation examples and discuss issues related to feasibility, closed loop stability and robustness. Example 1 We extend Example 1 by using our method in MPC. The system, optimization criterion, constraints and initial values are unchanged. Figure 3 shows some typical results. In this simulation the system is given by Y~, i.e. w ~ 1, while the control input is computed by choosing cz = 0, cf. Equation (8), i.e. the control input is based on ~0. We observe that the output constraints are satisfied when the added constraints are used. This is, however, not the case if these added constraints are omitted. The above might lead to the conclusion that an initial feasible solution implies a feasible solution at any time. This cannot be guaranteed. This is seen by performing the same simulation as above but increasing the gain of the simulated system. We increase the gain of ~t by a factor of 5. Figure 4a shows a violation of the output constraints at certain time instances. The reason for this is that there exists some x s X where no feasible zr can be found, i.e. no n ~ H such that the predicted states satisfy the state constraints for both bounding linear models. It should be be noted that this shows a rather extreme situation as the first coefficient of the FIRmodel can vary between -0.50 and 0.20, cf. Equation
Interpolating optimizing process control: B.A. Foss and S.J. Qin
136
Table 3 Example 3 - Bounding second order models. Each column represents a simulation run using the control sequence zr~ on a model h,. ~w- The system outputs of each simulation run are checked against the constraints at each time-step. 0 indicates no output constraint violation. 1 or 2 signifies a violation of the minimum or maximum allowable output, respectively Time
Model 1 0
~o ..... ~0
i=0 i=1 i=2 i=3 i=4 i=5 i=6
000000 000000 000000 000000 000000 000000 000000
w = 0.8
w = 0.6 0
~o ..... ,~0
000000 000000 000000 000000 000000 000000 000000
0
w = 0.4 0
= , ..... ~0
000000 000000 000000 000000 000000 I 1111 I 11111 l
w = 0.2 0
~ ..... ~o
000000 000000 000000 000000 111111 I 11111 111111
Model 0 0
~o ..... =0
000000 000000 000000 000000 I 1111 I 11111 I 111111
0
~o ..... ~o
000000 000000 000000 000000 000000 000000 000000
Table 4 Example 3 - Bounding second order models. Each column represents a simulation run using the control sequence ~ on a model hw ~ . The system outputs of each simulation run are checked against the constraints at each time-step. 0 indicates no output constraint violation. 1 or 2 signifies a violation of the minimum or maximum allowable output, respectively Time
i=0 i=I i=2 i=3 i=4 i=5 i=6
Model 1 0 ~,-",~0 000000 000000 000000 000000 000000 000000 000000
w = 0.8 0 0 ~l,'",~0 000000 000000 110000 111100 I I 1110 111110 111111
w = 0.6 0 0 ~ i ..... ~0
w = 0.4 0 0 ~ l ..... ~0
w = 0.2 0 0 ~t,'"'~0
Model 0 0 ~"",~0
000000 000000 111100 111110 111110 111111 111111
000000 000000 111110 111110 111111 111111 I 11111
000000 000000 I 11110 111110 111111 111111 111111
000000 000000 111110 111110 111110 111110 111110
Table 5 Example 3 - Bounding second order models. The criterion values using the control sequence rc~ on model wht + (1 - w)he for different values of a and w. The extra constraints are included ¢ Model 1 w = 0.8 w = 0.6 w = 0.4 w = 0.2 Model 0
n°
a =0.8
a = 0.6
a = 0.4
~ = 0.2
zrOo
2.57 1.06 0.74 0.60 0.53 0.51
2.57 1.06 0.74 0.60 0.53 0.51
2.57 1.06 0.74 0.60 0.53 0.50
2.57 1.06 0.74 0.60 0.53 0.50
2.57 1.06 0.74 0.60 0.53 0.50
2.57 1.06 0.74 0.60 0.53 0.50
Table 6 Example 3 - Bounding second order models. The criterion values using the control sequence rt~ on model wh~ + (1 - w)ho for different values of a and w. The extra constraints are not included ¢ Model 1 w = 0.8 w = 0.6 w = 0.4 w = 0.2 Model 0
rt~ 1.67 19.28 38.52 45.94 48.08 48.70
(23) a n d r e m e m b e r t h a t t h e Comparing Figures 4a added constraints have though the constraints are
a =0.8
ce = 0.6
1.71 12.05 24.64 29.52 30.92 31.34
1.81 6.61 13.93 16.78 17.59 17.84
g a i n o f Y~ is m u l t i p l i e d b y 5. a n d 4b i n d i c a t e s t h a t t h e a robustifying effect even v i o l a t e d in b o t h c a s e s .
o: = 0.4
ct = 0.2
n Oo
2.00 2.97 6.37 7.71 8.08 8.20
2.25 1.12 1.98 2.32 2.39 2.42
2.57 1.06 0.74 0.60 0.53 0.50
satisfied w h e n t h e a d d e d c o n s t r a i n t s a r e i n v o k e d . T h i s is in a c c o r d a n c e w i t h C o r o l l a r y 2 since a f e a s i b l e s o l u t i o n c a n b e f o u n d a t all t h e s i m u l a t e d t i m e steps. S e c o n d , as in Example 1 the results indicate that the added constraints h a v e a r o b u s t i f y i n g effect o n t h e c l o s e d l o o p c o n t r o l l e r .
Example 2 T o f u r t h e r g a i n i n s i g h t w e u s e d o u r m e t h o d o n t h e syst e r n in E x a m p l e 2. T h e s y s t e m , o p t i m i z a t i o n c r i t e r i o n , c o n s t r a i n t s a n d i n i t i a l v a l u e s w e r e u n c h a n g e d . Figure 5 shows some typical results. In this simulation the system is g i v e n b y Y0, i.e. w ~ 0, w h i l e t h e c o n t r o l i n p u t is c o m p u t e d b y c h o o s i n g o~ = 1, cf. E q u a t i o n (8), i.e. t h e c o n t r o l i n p u t is b a s e d o n Y~. W e c a n m a k e t w o o b s e r v a t i o n s f r o m t h e results. First, t h e o u t p u t c o n s t r a i n t s a r e
Discussion T h e p r o p o s e d m e t h o d is a n a t t e m p t t o o f f e r a s m o o t h t r a n s i t i o n f r o m l i n e a r t o n o n l i n e a r o p t i m i z i n g c o n t r o l as for example nonlinear MPC. In our opinion a smooth t r a n s i t i o n is i m p o r t a n t t o p r o m o t e i n d u s t r i a l a p p l i c a tions of nonlinear model-based optimizing control. The f e a t u r e s t h a t u n d e r l i n e t h e s m o o t h t r a n s i t i o n are: (i)
Interpolating optimizing process control: B.A. Foss and S.J. Qin •
1.2
b
1.1~
a
o.a ~ ...............................
o.a[-
0.6
,I. I ! .I. i
0.6 ....................................... 0.4 ...........
: .........................
,i.
0.2 "'h
i .............
0.4 0.2
t
o ,.4.'
:-
~
....................................
i
o
I
;
-o.2 . . ,
.-o.2 " i
•- 0 . 4
.-0.4
.-o.e .-.f ..~ . . . . .
....................
.-0.6
....... I
, '
__
i
......... i .............
'l. . . . . . . . . . . i . . . . . . . . . . . . .
::
b
1.2
i
0.8
0.8
0.6
O.E
0.4
0.4
o.2
oa
: .............
:
:
y
137
.
.
: .............
: .......................................
c
! -0.8
i 10
-0.8 20
.30
I
-0.~ 10
lime
20
30
,: o.8..
.
2. 1.5 !i
.
' . ~
Y o.e . . . . . . . ~........................
0.5 •
O"4 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
0
0.2 ........... :........................ 0-. ~
-0"20
1'0
20
lime
30
-o.s -1
-1.5,
10
20
30
.me
Figure 4 Example 1 - MPC. Same simulation as in Figure 3, but the gain of the simulated system is increased by a factor of 5. In a the added constraints are invoked, while these are omitted in b. The figures show the system output and the control input (note the axis)
0
10
T h e c o n t r o l a l g o r i t h m uses w e l l - k n o w n o p t i m i z a U o n b a s e d on l i n e a r m o d e l s as its b u i l d i n g - b l o c k s . (ii) A controls sequence is c o m p u t e d b y s m o o t h i n t e r p o l a t i o n , a c o m m o n l y a p p l i e d engineering principle, T h e p r o p o s e d c o n t r o l l e r can be seen as an extension o f linear o p t i m i z i n g c o n t r o l f r o m b o t h a theoretical as well as a p r a c t i c a l perspective, First, we e l a b o r a t e on the t h e o r e t i c a l results. A key a n d novel feature o f the p r o p o s e d c o n t r o l l e r is the a d d e d c o n s t r a i n t s i n c l u d e d for the t w o minim i z a t i o n p r o b l e m s , i.e. c o n s t r a i n t z ~ ~ Y i n (15), a n d ~ X in E q u a t i o n (16). This m a k e s it possible to guarantee robust feasibility on some given o p t i m i z a t i o n h o r i z o n when the b o u n d i n g m o d e l s are linear F I R - m o d e l s p r o v i d e d t h a t there exists a feasible c o n t r o l i n p u t sequence for each o f the two b o u n d i n g linear models. It s h o u l d be n o t e d t h a t the
20
30
time
m o d e l set Y~,,.m a y have autoregressive terms even if the b o u n d i n g m o d e l s are F I R - m o d e l s , since the w e i g h t i n g f u n c t i o n m a y d e p e n d on xi. A n easily c o m p u t a b l e u p p e r b o u n d can be f o u n d for the o p t i m a l i t y criterion in the case o f n o n l i n e a r F I R - m o d e l s . This p r o v i d e s a tool to investigate r o b u s t p e r f o r m a n c e o f the p r o p o s e d c o n t r o l l e r scheme for this m o d e l class. C o r o l l a r y 2 shows r o b u s t feasibility for the 1-step a h e a d p r e d i c t i o n for a n y t y p e o f b o u n d i n g linear m o d e l s if zc~ a n d rc~ a r e c o m p u t e d using the a d d i t i o n a l c o n s t r a i n t s . This implies that an M P C cont r o l l e r b a s e d on t h e p r o p o s e d a p p r o a c h a l w a y s will be feasible if there exists a feasible c o n t r o l input sequence for each o f the t w o b o u n d i n g linear m o d els at each time-step. T h e r e a s o n for this is t h a t the o p t i m a l c o n t r o l sequence is r e c o m p u t e d at every time-step. It will be, however, neither o p t i m a l n o r necessarily stable if the c o n s t r a i n t s are v i o l a t e d on the o p t i m i z a t i o n h o r i z o n b e y o n d the first p r e d i c t i o n
•
•
30
Figure 5 Example 2 - MPC. The system is given by ,~, i.e. w -- 0, while the control input is computed by choosing a = 1 (8), i.e. the control input is based on ~ . In a the added constraints are invoked, while these are omitted in b. The figures show the system output
b •
20
time
Figure 3 Example 1 - MPC. The simulated system is given by 3~], i.e. w = 1, while the control input is computed by choosing a = 0 (8), i.e. the control input is based on 3"o. In a the added constraints are invoked, while these are omitted in b. The figures show the system output and the control input
,
--0.; 10
tklle
step. T h e theoretical results a p p l y w h e t h e r the u n d e r l y i n g system h,, e E,,. is linear or nonlinear. The richness o f the n o n l i n e a r m o d e l class covered by Y~,~a n d the relevance o f this m o d e l class for process c o n t r o l p r o b l e m s need f u r t h e r investigations.
S e c o n d , we discuss the s i m u l a t i o n results. •
•
E x a m p l e 1 with b o u n d i n g F I R - m o d e l s s u p p o r t the t h e o r e t i c a l findings. It also clearly d e m o n s t r a t e s , cf. Tables 1 a n d 2, the i m p o r t a n c e o f the a d d e d constraints. F u r t h e r , the c o m p u t e d u p p e r b o u n d on the o p t i m a l i t y criterion seems to be r e a s o n a b l e in the sense t h a t it is n o t overly conservative. E x a m p l e s 2 a n d 3 investigate the results b e y o n d the b o u n d s o f the t h e o r y . This is o f interest since practical use o f c o n t r o l l e r s u s u a l l y violates theoretical a s s u m p t i o n s at s o m e stage. In E x a m p l e 2 the b o u n d i n g m o d e l s are quite similar in the sense that
138
•
•
Interpolating optimizing process control: B.A. Foss and S . J , Qin
they both are second order damped systems, while the model characteristics are more divided in Exampie 3. The findings in Example 2 indicate that constraint handling need not be difficult when the bounding model are similar, in this example control c o m p u t a t i o n without the added constraints gave only feasible solutions. This picture changed significantly in Example 3, the positive effect of the added constraints was very pronounced here, even though there is no feasibility guarantee in this case. Corollary 1, i.e. the result on feasibility of the 1-step ahead prediction, is supported by the simulations, There were no constraint violation at time i = 1 with the added constraints, see Tables 1 a n d 3. This was not necessarily the case if these constraints were removed as can be seen in Table 2. The closed loop simulations substantiate the open loop findings. The added constraints do seem to improve robust stability and performance. It is important to note that the feasibility guarantee depends on the existence of a feasible control input sequence for each of the two bounding linear models. This is demonstrated in Figure 4a.
Finally, we elaborate on some additional issues, •
An open question on the proposed approach is the richness o f the model set Y.,,.. It might be necessary to include more than two bounding linear model, i.e. let the hw be a convex combination of several, say M, linear models. In this case the model set can comprise a large class o f nonlinear systems as shown in Johansen and Foss 7. The problem that arises is the complexity increase in our algorithm since the number o f added constraints to compute each of the
•
optimal control sequence increases. In addition, the control sequence zr~ will be a convex combination of M control sequences. There exists methods and tools for identifying mod-
•
els on a format like Equation (4), see e.g. Reference 8. An alternative approach to the one pursued in this paper is to compute a control sequence based on one nominal model with added constraints related to one or more bounding models.
•
This p a p e r does not provide a closed loop stability p r o o f for the proposed method. The findings, do, however, provide a basis for seeking such a result.
•
Nonlinear optimizing control provides a nonconvex optimization problem. By this, there is no guarantee
that the global o p t i m u m can be found. Hence, an alternative approach is to search for a suboptimal solution with some performance guarantee. This can be done by limiting the search to the set rI~ proposed in this paper, a e [0, 1] might be computed in one of three ways. (1) A gain-scheduling approach can be applied if it is possible to schedule a on some process variables. (2) a is assumed to be time-invariant, and it is computed using an off-line or on-line estimation scheme based on some performance measure. (3) a is time-varying and is computed using a recursive on-line algorithm. This m a y pave the way for a reliable adaptive M P C controller, since only one time-varying p a r a m e t e r is required. In practice, robustness of adaptive control is difficult to obtain when m a n y parameters are estimated on-line.
Conclusions A new model-based optimizing controller for a set of nonlinear systems is proposed. R o b u s t feasibility with respect to hard state and/or output constraints is guaranteed on the optimization horizon for a set of nonlinear FIR-models, while significant feasibility improvements are encountered in the case where the bounding models are general linear state-space models. Finally, the method guarantees robust feasibility for a 1-step ahead prediction in the general case. This can be of interest in MPC-applications.
References 1. Lee,J. H. Recent advances in model predictive control and other related areas. In Preprints CPC-V, Lake Tahoe, USA, 1996. 2. Rawlings, J. B., Meadows, E. S. and Muske, K. R. Nonlinear model predictive control: A tutorial and survey. In Preprints [FAC Symposium ADCHEM, Kyoto, Japan, 1994, pp. 203-214. 3. Qin, S. Joe and Badgwell, T. A. An overview of industrial model predictive control technology. In Preprints CPC-V, Lake Tahoe, USA, 1996. 4. Mayne,D. Q. Optimization in model based control. In Preprints IFAC Symposium DYCORD, Helsingor, Denmark, 1995. 5. Genceli, H. and Nikalaou, M. Design of robust constrained model-predictive controllers with volterra series. AIChEJournal, 1995, 41, 2098-2107. 6. Bequette, B. W. Nonlinear control of chemical processes: A review. Ind. Eng. Chem. Res., 1991, 30, 1391-1413. 7. Johansen, T. A. and Foss, B. A. Constructing NARMAX roodels using ARMAX models. Int. J. Control, 1993, 58, 1t25-1153. 8. Foss, B. A. and Johansen, T. A. Local modelling as a tool for semi-empirical or semi-mechanistic process modeling. In Neural Networks for Chemical Engineers, ed. A. Bulsari, Elsevier, Amsterdam, 1995, pp. 297-334.