Simultaneous Model Predictive Control and Identification: Closed-Loop Properties S. Alper Eker and Michael Nikolaou1 Chemical Engineering Department University of Houston Houston, TX 77204-4792 Submitted to Automatica Abstract Model Predictive Control and Identification is an adaptive control technique which solves an on- line optimization problem to find process inputs for dual control problem. Its main goal is to bring robustness, to Model Predictive Control, increase the capability to handle uncertainties and time varying parameters in the processes. Theoretical properties, such as feasibility of the optimization problem, parameter convergence, are explored in this paper. In addition an algorithm combining recursive least squares objective with Sylvester Matrix constraint is proposed to obtain parameter estimates which do not exhibit pole zero cancellation in the identifier. Key Words-Adaptive Control, Persistence of Excitation, predictive Control, pole-zero cancellation 1 INTRODUCTION Simultaneous model predictive control and identification (MPCI), introduced by Genceli and Nikolaou (1996), is a class of adaptive control algorithms that employ on-line optimization at each time step, to generate process inputs which (a) satisfy ordinary operating constraints (e.g., lie within upper and lower bounds); (b) steer process outputs to their setpoints or within a desired range; and (c) satisfy excitation constraints compatible with the structure of the process model being identified. It is the last feature (c) of MPCI that places emphasis on model parameter identification over output regulation (the latter potentially relegated to satisfaction of output constraints) in order to develop a better model, which is at the heart of any model-based control algorithm. The focus of previous work on MPCI has been the formulation and numerical solution of the on-line optimization problem. In this work, we present a rigorous study of MPCI stability and parameter convergence properties, and suggest a new algorithm to avoid unstable pole-zero cancellation in the identification of over-parameterized models. 1 Author to whom all correspondence should be addressed. (713) 743 4309. [email protected] 1 The stability problem of adaptive control has challenged researchers for over three decades and there is a large body of related results, which we will not attempt to review. Rather, we will give a brief overview of results that are relevant to this work. Several early results on the stability of adaptive control have the common characteristic of employing a minimum-variance type of control law (Campi, 1996) resulting on cancellation of process zeros by controller poles, which requires processes to be minimum-phase. Using results of Sternby (1977) and Rootzen and Sternby (1984), Kumar (1990) established stability and optimality of a large class of adaptive regulators based on the recursive-least squares (RLS) algorithm, by exploiting the normal equations. The minimum-phase assumption is still needed in these results. Anderson and Johnson (1982) showed exponential parameter convergence of adaptive identification and control algorithms under deterministic or stochastic persistently exciting conditions on process inputs. The minimum-phase assumption was made in these results as well. Lozano and Zhao (1994) point out the difficulties of developing adaptive controllers for non-minimum-phase (NMP) processes: “The crucial problem in adaptive control of non-minimum phase processes arises when the process parameter estimates are such that they correspond to parameter values for which the process is not controllable.” Singularities such as unstable pole-zero cancellation in the control law create formidable problems in the analysis of NMP processes (Mosca, 1995, p. 334). For an adaptive pole-positioning algorithm, Anderson and Johnstone (1985) demonstrate such a problem that arises in the solution of the Diophantine equation for their control law. They point out that while it may be possible to give a with- probability-one argument that this problem will never happen, it seems impossible to completely rule out this possibility. Note that in this case even stable pole-zero cancellation is a problem. To completely prevent the possibility of unstable pole-zero cancellation in adaptive control schemes, research has taken two directions: (a) Use of excitation probing signals (Anderson and Johnstone 1985; Elliot et al., 1985; Goodwin and Teoh, 1985; Kreisselmeier and Smith, 1986; Girí et al., 1989; Polderman, 1989), and (b) Modification of the process parameter estimates away from values that create unstable pole-zero cancellation (Lozano and Zhao, 1994; Kreisselmeier, 1986; Ossman and Kamen, 1987; Middleton et al., 1988). In this paper our aim is to show closed-loop bounded-input-bounded-output and parameter convergence properties of the MPCI algorithm applied to stable processes. We show that, by the inherent structure of the algorithm, bounded and persistently exciting (PE) process inputs are guaranteed to exist when MPCI is started up with an appropriate procedure. PE process inputs are essential for the convergence of parameter estimates to their true values in most of the cases analyzed in this paper. No minimum-phase assumption for the identified process is needed. For the case where one cannot guarantee convergence to true parameter values even with PE inputs (over-parameterized Auto-Regressive-with-eXogenous-input (ARX) models), we propose a numerically tractable modification of the RLS identification algorithm, which prevents unstable pole-zero cancellation. The modified RLS algorithm may be used either within the MPCI framework or independently of it. The proposed algorithm relies on ensuring nonsingularity of the Sylvester matrix, an idea that has also appeared in Lozano and Zhao (1994) and Campi (1996). However, in contrast to these references, where nonsingularity is indirectly guaranteed through modification of least-squares parameter estimates, the proposed approach relies on direct enforcement of nonsingularity of the Sylvester matrix through addition of an 2 appropriate constraint to the RLS algorithm. The rest of this paper is organized as follows: The MPCI paradigm is briefly explained qualitatively and quantitatively in section 2. Then basic lemmas are given at the beginning of the main results in section 3, which guarantee the feasibility of the on-line optimization problem. Subsequently, the main results for deterministic and stochastic cases are given for both ARX and finite impulse response (FIR) stable process structures. For over parameterized ARX models, a modified RLS algorithm that prevents pole-zero cancellation is formulated and elaborated on in section 3.3. Finally, simulation examples illustrate the above results. 2 The MPCI Paradigm for Adaptive Control through On-line Optimization 2.1 MPCI qualitatively Focusing on the essential role of persistent excitation (PE) in parameter identification, Genceli and Nikolaou (1996) introduced the simultaneous MPC and identification (MPCI) paradigm for adaptive control. MPCI resorts to on-line optimization. An objective function over a moving horizon is minimized with respect to process inputs that satisfy (a) ordinary MPC constraints and (b) a constraint that data be persistently exciting with respect to the structure of the model being identified. The magnitude of the excitation may be set a priori or maximized by being included in the overall MPCI objective. The PE constraint (identification feature that distinguishes MPC from MPCI) may be turned on and off, as process identification needs arise. For example, the controlled process may be assessed to require MPCI because of changes in the equipment used in the process, drastic changes in raw materials, or poor controller performance detected by controller performance monitoring methods. Advantages of MPCI are that it does not require any explicit external dithering signals, PE, by construction of the algorithm, is trivially dependent on the feasibility of the on-line optimization, process excitation is maximized, constraints on process inputs are explicitly enforced, and constraints on process outputs are explicitly handled (Eker and Nikolaou, 1999). 2.2 MPCI quantitatively At time k, let the stable, linear, time-invariant model yt(| k) =f()t(T)qˆ k (1) be used to predict future output values yt(| k) of the linear time-invariant process yk()w( k=) f(k)Tq0 + (2) where w is white noise with zero mean. As an example, a finite-impulse-response (FIR) model has the form (cid:229) n yt(k|)b(u|)t(i|k=d)-t k ˆi+(k) (3) i=1 where: (cid:236) current or past output, yt(), measured at time t£ k y(t k)=(cid:237) (cid:238) future output at time t predicted at time k <t (cid:236) implementedinput, u(t), t<k u(t k)=(cid:237) (cid:238) current future potential input at time t predicted at time k <t 3 dt(|k )=estimate made at time k of disturbance introduced at time t bˆ =estimate made at time k of model coefficients (k) Based on the above, a variant of the MPCI on-line optimization can be formulated as follows. (cid:229) M maximize (vs - hm ) (4) i ukuk,kM,L, k ,+ -1 skk+k1M| kL, s + | i=1 (|ki+)k subject to uukik‡ui (1|)+1-‡,2M... = (5) maxmin D u ‡ D u(k +i- 1|k)‡ D u i =1,2...M (6) max min y +‡+m‡ y- k(|i )k y m max i min i, i =1, ,M (7) L m ‡ 0 i uk(|)(+|M+)=ikuki +k , i=0, ,M - 1 (8) L yk(|i +)k =f()k(i+) Tqkˆ i =1,2...M (9) P(1k)- f qˆqˆ()kf(k1k)=()(+1)qk - 1(+1f) TP k- (k)f -غ y(k) ˆ T - øß (10) ()k( ) k Pk(1P)-(k1)f f T - Pk()P(1k)=- - ()k( ) k (11) 1(+1f) TP k- f ()k( ) k (cid:229)s- 1 1s f(kj- i+ )f(kj- i+ )Tfs (|ki+)k If0 i =1,2,...M (12) j=0 where D u(k +i- 1|k)=ˆ u(k +i- 1|k)- u(k +i- 2|k) (13) [ ] Y(k)= y(k) y(k - 1) y(k - s+1) T (14) L f(-k j+)i =[uk(1ji)|-+(k)--+ukj- 1in|k L L ]T (15) v and h are weights on the PE lower bounds s and softening variables m on output bounds ki+k| i y and y , respectively; s refers to the length of the data set on which the PE condition is max min ˆ applied and used for parameter estimation; and eqns. (10) and (11) generate estimates q of the parameter q0 according to the RLS algorithm stemming from minimization of the quadratic cost function Jyk(q)(fq=(q)-+)q1(- q)(cid:229)(N q ) - T 2P 1 ˆ T - 1 ˆ (16) N 2 ()k(0)(0) 2 0 k=1 with respect to q , where P (=P(0) in equation (11)) is a “large” positive definite matrix and 0 ˆ q is the initial parameter estimate (Goodwin and Sin, 1984) (0) Loosely speaking, maximiziation of s , corresponding to the minimal eigenvalue of (|ki+)k the information matrix, ensures that maximal process information is generated (“small” variance and bias of parameter estimates). The information generated is known to be E-optimal (Ljung, 1987; Söderström and Stoica, 1989). In typical MPC fashion (Prett and Garcia, 1988), the above optimization problem is solved at time k, and the optimal u(k|k) is applied to the process. For more details on MPCI see 4 Eker and Nikolaou (1999). Throughout this paper we let s =M for notational simplicity (Figure 1). Note also that additional terms (e.g., related to set-point tracking, etc.) may be included in the objective of eqn. (4). 3 MAIN RESULTS Notice that only u(k|k) among all decision variables {u(),k,(k1u)k}M k+ - appears in L inequality (12) for the first window (i =1). The important implication of this observation is that PE is trivially guaranteed for the closed loop, regardless of the behavior of the true process, as long as the on-line optimization problem has a feasible solution with s >0. We show in the 1 following Lemmas that feasibility of the on-line optimization problem can be guaranteed through an appropriate start-up procedure. We then proceed to exploit the implications of these Lemmas in parameter convergence results for various cases. 3.1 Feasibility of MPCI ˆ Let the FIR model structure of eqn. (3) be used. At time k=l, a parameter estimate q of FIR (l) model (which is needed to make prediction of output and calculation of input via optimization) can be obtained as the result of the least-squares minimization Ø qø ØøØ y(1) uu(0)(1)n 1 ø - Œ 1œ L min ŒœŒ - œ Œ Mœ (17) ŒœŒ M M OM œ M Œ œ q q ŒœŒºßº y(l) ul(u1l)-()n 1 œß - Œ nœ L 14444244443º dß F (l) (l· n) Ø f Tø (1) Œ œ where F (l=) Œ f(2)Tœ . Similarly, when an ARX model structure is used, qˆ can be obtained as Œ œ (l) Œ M œ Œº f Tœß (l) the result of the least-squares minimization Ø q ø ØøØ y(1) - yu(0)(1)m 1 ø - Œ 1œ L min ŒœŒ - œ Œ M œ (18) ŒœŒ M M O œ M M Œ œ q q ŒœŒºßº y(l) -- yl(u1l)()m1 - œß Œ n+mœ L 14444244443º d ß F (l) ln· m(+ ) For unique solutions to the above minimization problems in eqns. (17) and (18) to exist, the following condition must be satisfied: F T(Fl)() l r0I (19) f l f Assuming that at time k=l the above condition is indeed satisfied and a unique result of the minimization in eqn. (17) or (18) has been obtained, it is natural to ask whether an input u(l|l) can be selected time k=l, such that the following inequality, containing u(l|l) as the only decision variable, can be satisfied: 5 F+T(Fl1)+(1) l r0 I (20) f l+1|l f Ø f Tø (2) Œ œ where F+(l1)= ˆ Œ f(3)Tœ and the input u(l|l) satisfies all input constraints. Satisfaction of the Œ œ Œ M œ Œº f(1l)+ Tœß ˆ above inequality (20) would ensure that a unique parameter estimate q(1l)+ would exist at the next time step l+1. We show below that such an input u(l|l) can always be found within the above context. In fact, we show that as long as we start MPCI from a feasible point (i.e., with a finite -length sequence of inputs that are persistently exciting), we will have a feasible solution for all subsequent time steps, i.e., the MPCI algorithm will be able to find a process input that satisfies the strong persistent excitation condition. We make that statement precise in the following Lemma 1 and Lemma 2. Lemma 1 - Feasibility of the MPCI on-line optimization problem for FIR model Assume that, a) The FIR structure of eqn. (3) is used. b) at time step k =l, the past process inputs u(1- n) to u(l1)- are such that the condition l (cid:229) f f T 0 (21) ()i() i f i=1 (cf. equation (19)) is satisfied. Then, assuming that all other constraints are feasible, there exists, at time k =l, a feasible value for the process input ul( l| ) such that l+1 (cid:229) f f T 0 (22) ()i() i f i=2 (cf. equation (20)) is satisfied. Remark Notice that, without loss of generality, the result is nothing but the feasibility of the PE constraint for the first window (i=1) in MPCI, equation (12), with s = M =l. Proof: See Appendix A. Corollary 1 Assume that at time step k =l, the process inputs u(1- n) to u(l1)- are such that the condition l (cid:229) f f T 0 (cf. equation (19)) is satisfied. Then, if the MPCI algorithm with FIR model ()i() i f i=1 structure is applied, there exists at any time k ‡ l a value for the process input uk( k| ) such that k+1 (cid:229) f f T 0 (cf. equation (20)). ()i() i f i=-k+l 2 Proof: Trivial, by repeated application of the above Lemma 1 for k =l+1, k =l+2, etc. 6 Lemma 2 - Feasibility of the MPCI on-line optimization problem for ARX model Assume that a) An LTI stable process is modeled with the following ARX structure A(z)yk-- 1(z)B()zu(k- =) ddel 1 + (23) with A(z)a--z1a1=z++ +- 1 n 1 L n Bz(bb-- 1z)b=+z - +1 m 0 1 K m where n,m and the time delay del are assumed to be known; the polynomials A(z- 1), B(z- 1) are coprime and d is a deterministic disturbance. b) at time step k =l, the past process inputs u(1- m) to u(l1)- are such that the condition l (cid:229) f f T 0 (24) ()i() i f i=1 (cf. equation (19)) is satisfied. Then, assuming that all other constraints are feasible, there exists, at time k =l, a feasible value for the process input ul( l| ) such that l+1 (cid:229) f f T 0 (25) ()i() i f i=2 (cf. equation (20)) is satisfied. Remark Notice again that, without loss of generality, the result is nothing but feasibility of the PE constraint for the first window (i=1) in MPCI, equation (12), with s = M =l. Proof: See Appendix B. Corollary 2 Assume that at time step k =l, the process inputs u(1- m) to u(l1)- are such that the condition l (cid:229) f f T 0 (cf. equation (19)) is satisfied. Then, if the MPCI algorithm with ARX model ()i() i f i=1 structure is applied, there exists at any time k ‡ l a value for the process input uk( k| ) such that k+1 (cid:229) f f T 0 (cf. equation (20)). ()i() i f i=-k+l 2 Proof: Trivial, by repeated application of the above Lemma 2 for k =l+1, k =l+2, etc. 3.2 Parameter Convergence for the Deterministic Case In this section, parameter convergence of the RLS algorithm, which is used for identification in MPCI, is studied with FIR and ARX structures for both known and unknown process order. Theorem 1 – Parameter convergence for FIR process of known order Assume that (a) An LTI stable process has the form yk()b( u=k-i(cid:229)n)+0d i i=1 7 where the order, n, of the system assumed to be known, and d is an unknown constant disturbance. (b) The MPCI adaptive control scheme, equations (4) – (15) , is used. (c) The MPCI algorithm starts in a feasible region such that there is enough excitation at the l starting point l = M , i.e., (cid:229) ff T r I>0. ()i( ) i f i=1 Then (a) The RLS parameter estimate error q%q =ˆ ˆq - 0, converges exponentially to zero. ()k( ) k (b) Predicted process output values converge to the real process output values. Proof: See Appendix C. Remarks (a) Notice that, by Lemma 1, feasible values of s >0, i =fi1¥ ,...,k are guaranteed to li+l+i1|+ exist. (b) It is important to note that the proof of Theorem 1, as well as the proofs of all subsequent theorems, do not assume that the MPCI algorithm reaches the global optimum of the on-line optimization. Rather, feasibility (guaranteed by Lemma 1 and Lemma 2) is the only requirement. Theorem 2 – Parameter convergence for ARX process of known order Assume that (a) An LTI stable process has the following ARX structure A(z)yk-- 1(z)B()zu(k- =) ddel 1 + (26) with A(z)a--z1a1=z++ +- 1 n (27) 1 L n Bz(bz-- 1b)=z 1- + m (28) 1 K m where n,m and the time delay del are known, and d is an unknown constant disturbance. (b) The MPCI adaptive control scheme for the above ARX structure is used, i.e. equations (4) – (14) along with the counterpart of equation (15) for the ARX structure, i.e., f(-k j+)i =--+[---+y--+(k--+1ji)|(-k)y(1kj)()in1|kukLji|kukj im|k L ]T(29) (c) The MPCI algorithm starts in a feasible region such that there is enough excitation at the l starting point l = M , i.e., (cid:229) ff T r I>0. ()i( ) i f i=1 Then (a) The parameter estimate error q%q =ˆ ˆq - 0 ,converges exponentially to zero ()k( ) k (b) Predicted output values converge to the real process output values. Proof: See Appendix C. Remarks The proofs of Theorem 1 and Theorem 2(Appendix C) start with arguments similar to those in Goodwin and Sin (1984, Equation 3.3.66, p. 61) and Anderson and Johnson (1982), to show 8 convergence (via a Lyapunov function argument) and exponential convergence of parameter estimates, respectively. However, the following important differences distinguish this work from the above references: (a) Instead of assuming that PE is a derivative property of the closed loop, resulting from external excitation and assumptions about the process and the controller, PE is directly enforced through the MPCI on-line optimization. (b) The process input u is constrained by the MPCI algorithm to be between upper and lower bounds. Therefore, the boundedness of u is trivial. Remarks on assumptions made in Theorem 1 and Theorem 2 (a) The values of n, m and del in equation (26) are known, and n is minimal, in order to prevent pole-zero cancellations and to provide convergence of model parameters to their actual system parameters. Without loss of generality, the time delay del is taken as 0 in our discussion. (b) The process is assumed to be stable. This is automatically ensured for FIR models. (c) The polynomials A(z- 1) and B(z- 1) for the process are assumed to be coprime. For FIR models coprimeness is trivially satisfied, since A()z- 11= . Coprimeness along with PE will prevent pole-zero cancellation during parameter estimation. (d) The assumption is often made in literature that the zeros of the polynomial B(z- 1) lie outside the unit disk (Anderson and Johnson, 1982), i.e. the process is minimum-phase. This assumption is to ensure that inputs u are bounded. However MPCI selects inputs satisfying input bounds, so the boundedness of u is trivially satisfied. Consequently, this assumption about the roots of B(z- 1) is not needed in our analysis. Theorem 3 – Parameter convergence for over-parameterized FIR model of FIR process Assume that (a) An LTI stable process has the form yk()b( u=k-i(cid:229)n)+0d with FIR structure, where an i i=1 upper bound of order, n, of the system is known, d is a deterministic disturbance which will be modeled as a constant in the model (b) The MPCI adaptive control scheme, equations (4)-(15) is used. (c) The MPCI algorithm starts in a feasible region such that there is enough excitation at the l starting point l = M , i.e., (cid:229) ff T r I>0. ()i( ) i f i=1 (d) The optimal values of s >0, produced by the MPCI algorithm for i=fi1¥ ,...,k , satisfy li+l+i1|+ the inequality s ‡ s > 0 (30) li+l+i1|+min Then (a) The parameter estimate error q%q =ˆ ˆq - 0 converges to zero, where the vector q0 contains ()k( ) k the true parameter vector augmented by zeros. (b) Predicted output values converge to the real process output values. Proof: See .Appendix D 9 Remark (a) Notice that because the LTI system model has higher order than the process, q0 contains the real process parameter vector augmented by a corresponding number of zeros, to satisfy equation (C-1) in Appendix C. (b) The proof of Theorem 3 uses a Lyapunov argument (equation (C-8) in Appendix C) to show eqn. (D-5) (cf. Lemma 3.3.6 in Goodwin and Sin (1984, p. 60)) which finally results in parameter convergence. Shown below is an extension of Theorem 3 for the case where an under-parameterized FIR model is used. Theorem 4 – Parameter convergence for under-parameterized FIR model of FIR process Assume that (a) The input-output data relationship is described by (cid:229)n¥ Process: yk()u(k=i) q - (31) i i=1 (cid:229) n Model: yk()u(k=i) q - (32) i i=1 where n<<n¥ . ¥ (b) The tail end of the kernel q decays exponentially as q <a e- i (33) i (c) The MPCI algorithm starts in a feasible region such that there is enough excitation at the l starting point l = M , i.e., (cid:229) ff T r I>0. ()i( ) i f i=1 (d) The optimal values of s >0, produced by the MPCI algorithm for i=fi1¥ ,...,k , satisfy li+l+i1|+ the inequality s ‡ s > 0 (34) li+l+i1|+min Then qˆ-<+q Ke+(e -+2- (1n)+2(2)1/2 en - 2n¥ ) (35) n· 1 3 L where the vector q contains the first n coefficients of the vector q and K <¥ is a constant. n· 1 3 Proof: See Appendix E. Remark (a) Notice that, by Lemma 1, feasible values of s >0, i =fi1¥ ,...,k are guaranteed to li+l+i1|+ exist. However an additional assumption (d) with, feasibility which is guaranteed by Lemma 1 and Lemma 2, is required in this case. Corollary 3 10