Linear Quadratic Feasible Predictive Control

PII: S0005–1098(98)00106–X

Automatica, Vol. 34, No. 12, pp. 1583—1592, 1998  1998 Elsevier Science Ltd. All rights reserved Printed in Great Britain 0005-1098/98 $—see front matter

Brief Paper

Linear Quadratic Feasible Predictive Control* B. KOUVARITAKIS,- J. A. ROSSITER‡ and M. CANNONCharacterization of the class of stabilizing predictions is used to endow constrained predictive algorithms with the attributes of optimality and feasibility while retaining computational simplicity.

Key Words—Predictive control; stability; predictions; optimality; constraints; feasibility.

both predicted inputs and outputs thereby improving significantly upon problems of infeasibility due to the presence of constraints. Constraints present challenging problems (e.g. Sussman et al., 1994; Alvarez-Ramirez and Suarez, 1996; Rossiter et al., 1995; Kouvaritakis et al., 1996): for example in the case of open-loop unstable systems, depending on the initial conditions generated by earlier control moves, it is possible that all stabilizing future control trajectories require more control authority than constraints will allow, and this inevitably leads to closed-loop instability. Sussman et al. (1994), Alvarez-Ramirez and Suarez (1996), Rossiter et al. (1995) and Kouvaritakis et al. (1996), consider different aspects of the same problem. The first gives a characterization of the class of systems with bounded inputs which are globally stabilizable; the second shows that this class can be stabilized through the state feedback law that minimizes the usual LQR cost, J"#e(t)##j#u(t)#, where e(t) and u(t) denote the error and control input, providing that j is varied with time in a particular manner; the third establishes the necessary and sufficient conditions for the existence of feasible stabilizing solutions; whereas the fourth uses these necessary and sufficient conditions to ensure that feasible stabilizing solutions will exist at future times, providing that they exist at some initial time. The work of this paper bears some similarity to that of Alvarez and Suarez (1996) however, unlike Alvarez and Suarez (1996) which uses one degree of freedom and solves a nonlinear equation in order to derive sufficient feasibility conditions, here we use linear interpolation to define a Quadratic Program which guarantees feasibility as well as optimizes a prescribed performance cost over a finite number of degrees of freedom. The IHPC mechanism of generating predicted input IIRs is also exploited in Scokaert and Rawlings (1996) and Kouvaritakis et al. (1997), which combine LQ optimal control with constrained predictive control; this algorithm will be referred to as linear quadratic constrained predictive control (LQCPC). LQCPC has an attractive feature: at first sight optimality, in the presence of constraints, is not practicable because it leads to an infinite dimensional constrained optimization problem. LQCPC demonstrates that N, the number of degrees of freedom, need not be infinite; instead it is argued that if N is larger that some finite number N , optimality is guaranteed. However, N could   be arbitrarily large; the need for large Ns not only leads to implementation difficulties (i.e. increased computational demands), but also reduces the superiority of LQCPC over earlier algorithms; both SGPC and IHPC for example converge to optimality as N becomes arbitrarily large. Here we re-examine this issue starting from the class of stable predictions but define the degrees of freedom in a new way so that the unconstrained optimal is a member. Clearly this approach, referred to as LQIHPC, is related to LQCPC and we explore the connection. The relationship between the two is a simple one, but neither forms a subclass of the other; hence it is not possible to tell a priori which is most likely to be feasible or optimal. Given that in practice N, the number of degrees of freedom, has to be kept small, both LQCPC and LQIHPC can run into infeasibility especially in the presence of tight constraints. This is due to the

Abstract—Terminal constraints in predictive control provide a guarantee stability, but result in deadbeat predicted responses and can violate physical constraints. Two recent algorithms (Rossiter et al., 1996a; Scokaert and Rawlings, 1996) removed this restriction, but of these the former lacks the guarantee of l -optimality whereas the latter lacks the guarantee of feasibility.  Here we develop algorithms which overcome both these difficulties and illustrate their advantages by means of numerical examples.  1998 Elsevier Science Ltd. All rights reserved. 1. Introduction Generalized Predictive Control (GPC) (Clarke et al., 1987) has provided an effective means of controlling linear discrete-time time-invariant plant yet, in its original form it lacked a guarantee of stability. Later work (Clarke and Scattolini, 1991; Mosca and Zhang, 1992; Kouvaritakis et al., 1992; Rawlings and Muske 1993; Rossiter et al., 1996a; Rossiter et al., 1997; Scokaert and Rawlings, 1996; Kouvaritakis et al., 1997; Gilbert and Tan, 1991) provided a solution: for example Stable Generalized Predictive Control (SGPC) (Kouvaritakis et al., 1992) achieved stability through a pole/zero cancellation in the appropriate prediction equations. All relevant algorithms (Clarke and Scattolini, 1991; Mosca and Zhang, 1992; Kouvaritakis et al., 1992) however deploy deadbeat predictions and this often requires active control which in turn may come into conflict with hard physical limits thereby resulting in infeasibility and possibly instability. A recent state space approach (Rawlings and Muske, 1993), addressing a regulation problem, removed the need for cancellation of poles in the prediction equations and introduced instead the constraint that the state vector should be steered into a stable subspace of the state space. As a result it was possible to allow for infinite impulse responses (IIR) in terms of the output predictions. The predicted input behaviour however was still taken to be of a finite impulse response (FIR) type; IIR’s would lead to an infinite optimization problem, since in this approach the degrees of freedom are taken to be the future control moves. This difficulty was resolved by the infinite horizon predictive control (IHPC) (Gossner et al., 1997) which gave a characterization of the class of stable predictions and deployed the degrees of freedom available therein; IIR’s responses were deployed for *Received 26 April 1997; revised 1 November 1997, 9 March 1998, 12 May 1998; received in final form 2 June 1998. This paper was not presented at any IFAC meeting. This paper was recommended for publication in revised form by Associate Editor J. W. Polderman under the direction of Editor C. C. Hang. Corresponding author Professor B. Kouvaritakis. Tel. #44 1865 273000; fax #44 1865 273906; E-mail basil.kouvritakis@ eng.ox.ac.uk. -Oxford University, Department of Engineering Science, Parks Rd., Oxford OX1 3PJ, UK. ‡Loughborough University, Department of Mathematical Science, Loughborough, LE11 3TU, UK. 1583

1584

Brief Papers

fact that the implicit class of predictions allowed for by both algorithms is centred on the unconstrained optimal which can often be infeasible. Hence, in the absence of enough design freedom, the whole class of predictions can be infeasible. To overcome this problem, here we consider an alternative which ‘‘centres’’ predictions around a ‘‘mean level’’ control strategy and thus, under some weak assumptions, overcomes the feasibility difficulties (even for small N); for brevity the resulting algorithm will be referred to as LQFPC (LQ feasible predictive control). In addition to its feasibility attributes however, LQFPC uses a prediction class which contains the unconstrained optimal and thus converges to it as soon as this becomes feasible. In the case of loose constraints and/or a large number of degrees of freedom, LQFPC gives comparable performance to LQCPC and LQIHPC but clearly has the advantage that it retains feasibility even when the number of degrees of freedom is small and/or the constraints are tight. 2. Earlier work 2.1. Generalized predictive control. Consider the model b(z) y " *u ; "*u !*u "4*º R> a(z) R R     *u "u !u , R R R\ a(z)"*(z)a(z), *(z)"1!z\,

(1a) (1b)

x "Ax #B*u , e "Cx . I> I I I I

(1c)

b(z)*º(z)#P(z) ½(z)" a(z)



(2)

P(z)"[1 z\2z\L]P, P"H Du !H y @   ?   Du "[*u *u ]2, y "[y 2y ]2, 2   R\L> R\   R\L R where H is a Hankel matrix based on the coefficients of b(z); its @ structure as well as the structure of Toeplitz striped convolution matrices to be used in the sequel is

 

m  m H "  K $ mK L m  m C "  K $

2 2 $ 0 0 m  $

mK L \ m LK $ 2



mK L 0 ; $ 0



2.2. Guaranteeing stability via endpoint terminal equality constraints. GPC has no guarantee of stability and to remedy this algorithms were proposed which directly (Clarke and Scattolini, 1991) or indirectly (Kouvaritakis et al., 1992) force the predicted output to: (i) reach its target value of r after n steps; and W (ii) stay there subsequently. All future predicted error values beyond n are made zero, and since all future control increments W beyond n are zero, stability can be established by proving that S J is a monotonically decreasing function of time. The diffi%.! culty with this approach is that the predicted output response is deadbeat; this often leads to infeasibility. Below we review three recent approaches that tackle this problem. 2.3. Removing the restriction of error deadbeat responses. To avoid error deadbeat responses, an error state space model was deployed (Rawlings and Muske, 1993)

a(z)"a #a z\#2#z\L, (1d)   b(z)"b #b z\#2#b z\L> (1e)   L\ where a(z), b(z) are coprime; u , y denote the instantaneous input R R and output values, z is the z-transform variable, and zU when acting on u , y denotes the delay operator. For simplicity we R R assume that the problem is to track constant setpoints r and we only consider constraints on the control increments (other constraints can be treated similarly) and let *u "0. At time t,  equation (1a) can be used to relate the z-transforms of the output/input predictions:

½(z)"[1 z\z\2]Y, Y"[y y 2]2 R> R> Dº(z)"[1 z\ z\2]DU, *U"[*u *u 2]2 R R>

conditions and given by the inverse z-transform of P(z)/a(z). GPC [1] uses equation (4) to minimize J "#r1!Y ## %.!   j#DU # where 1 is a vector of 1’s and k'0. GPC imple  ments the first element of the optimal DU and repeats this   procedure at the next time instant.

2 0 2 0 ; $ $ 2 2 2 m m(z)"m #m z\#2#m z\LK . (3)   LK Inverse z-transformation of equation (2) gives the prediction equation Y "C DU #Y , (4)   E    where Y , DU are vectors comprising the first n , n ele    W S ments of Y and DU, with n , n being output and input predicW S tion horizons; all control increments beyond n are assumed to S be zero. [The subscript ‘‘fin’’ is used to emphasise the use of finite horizons.] C is the convolution matrix (as per equation (3) E based on the elements of model impulse response; correspondingly Y comprises the elements of the output due to the initial 

(5)

Consideration was given to a state regulation problem and, rather than cancel all the poles in the prediction equations, Rawlings and Muske (1993) invoked the terminal condition that, at the end of a given finite horizon n , the state vector should lie V in a subspace that is orthogonal to all the right eigenvectors of the state matrix A that are associated with unstable eigenvalues. Thus, the error response is no longer restricted to be a FIR, but is guaranteed to be stable. It is easy to show that under these circumstances even for infinite output horizons the performance cost is a computable quadratic function of the vector of future control increments:  L Q J" [e ]#j *u R>G R>G\ G GG "DU S DU #b2 DU #c , (6)             where S is positive definite, b a vector of conformal dimen    sion, and c a constant; initial conditions plus a Lyapunov   equation can be used to compute S , b , c .       2.4. Removing the restriction of error and input deadbeat response. The approach above still uses FIRs for future control sequences and may necessitate the use of active control inputs together with the concomitant feasibility difficulties. This problem is circumvented by IHPC (Rossiter et al., 1996a) which first characterizes the class of all stable predictions and then optimizes the performance index over the degrees of freedom available within the class. Thus subtracting equation (2) from the equation R(z)"r/(1!z), we have Q(z)!b(z)*º(z) E(z)"R(z)!½(z)" , a(z) Q(z)"a(z)r!P(z).

(7a) (7b)

Then denoting by a (z), b (z) the strictly Hurwitz factors of a(z) \ \ and b(z), and by a (z)b (z) the unstable factors, we conclude > > from the above that the predicted error will be stable iff: Q(z)!b(z)*º(z)"a (z) (z) > Q(z)!a (z) (z) > 0 *º(z)" , b(z)

(8a) (8b)

where a (z)"*(z)a (z), and (z) is a polynomial or stable > > transfer function. Equation (8b) on the other hand, implies that the future control increment sequence is stable if and only if Q(z)!a (z) (z)"b (z)t(z) > > 0 a (z) (z)#b (z)t(z)"Q(z) > >

(9a) (9b)

Brief Papers with t(z) is a polynomial or stable transfer function. If + (z), t(z), is a particular solution pair to equation (9b) then the whole class of solution pairs is + (z)#b (z)F(z), t(z)!a (z)F(z),, > > where F(z) is a stable transfer function, or equivalently the infinite-order polynomial corresponding to its impulse response; in practice F(z) is taken to be a polynomial of finite order. The general solution pair above with equations (7a) and (8b) give the class of stable predictions as *º(z)"G (z)#G (z)F(z);  

(10a)

E(z)"H (z)#H (z)F(z),  

(10b)

t(z) G (z)" ,  b (z) \

(10c)

a (z) G (z)"! > ,  b (z) \

(10d)

(z) H (z)" ,  a (z) \

(10e)

b (z) H (z)" > .  a (z) \

(10f)

With F(z) a polynomial and F the vector of its coefficients it can be shown that the l -norm of the predicted error and control  increment can be written as quadratic functions of F #E#"F2S F# b2F#c    #DU#"F2S F#b2F#c   





N

  J" [r!y ]#j *u "F2SF#b2F#c, R>G R>G\ G GG (11) S"S #jS , b"b #jb , c"c #jc .      

Hence it is possible to perform the minimization of J explicitly when there are no constraints or via QP otherwise; to invoke the constraints note that equation (10) implies that DU"g #C F,  E

(12)

where g is the vector of the elements of the impulse response of  G (z) and C  is the convolution matrix associated with g , the  E  vector of elements of the impulse response of G (z). For general  initial conditions the QP optimization above may be infeasible but as mentioned earlier throughout this paper we make the necessary/sensible assumption that: Assumption 2.1. At start time there exists at least one F (which could be infinite dimensional) for which the stabilizing prediction solution of equation (12) is feasible. If f denotes the first coefficient of the optimal F(z), and g,   g the first elements of g , g then the current optimal control    increment is given as *u "g f #g R   

(13)

which can be implemented and the procedure repeated at the next time instant. 2.5. Introducing ¸Q optimal state feedback control law. For no input constraints, the predicted law that minimizes the cost J of equation (11) is *u "!K x , k50, R>I R>I

1585

LQCPC (LQ constrained predictive control) is  min J, J" [e ]#j[*u ] subject to R>G R>G\ DU   G (i) equation (14) ∀k5N, (ii) equation (1b) ∀t,



where DU "[*u *u 2*u ]2.   R R> R>,\ The subscript ‘‘fin’’ distinguishes the vector above from that of all future control increments. LQCPC has the property that there exists an finite integer N such that for all N'N the minimization of equation (15)   yields the optimal solution to the infinite-dimensional optimization:  min J, J" [e ]#j[*u ] subject to (1b) ∀t. R>G R>G\ DU G (16) Thus, the seemingly infinite-dimensional problem turns out to be finite dimensional and equivalent to that of equation (15) providing of course that N'N .  This is an appealing feature and suggests that LQCPC cannot be improved upon. However N , though finite, can be arbitrar ily large. Thus often N must be chosen to be large, and then there is little difference between the results of LQCPC and IHPC: if F(z) were to be taken to be a polynomial of arbitrarily large order, then IHPC would utilize almost the entire allowable class of stable predictions. Given the convergent nature of the solution to equation (16) it follows that deviations from the ideal optimal solution due to taking F(z) to be of finite, albeit large, order can be made arbitrarily small. Besides, often the use of a large N is impracticable due to the excessive computational load demanded by QP. For this reason in the current paper we explore alternatives to LQCPC and first of all, we approach the problem from the characterization of all stable predictions. 3. ¸Q infinite horizon predictive control (¸QIHPC) 3.1. ¹he algorithm. The choice of polynomial F(z) in IHPC implies the need of high orders in order to retrieve the optimal solution, even in the unconstrained case. An obvious way to remove this difficulty is to compute at each t the unconstrained optimal F(z), say F (z), and use this as the ‘‘centre’’ of the class  of stable control increment prediction, namely replace equation (10) by *º(z)"G (z)#G (z)F(z)   E(z)"H (z)#H (z)F(z)   a (z) t(z) > G (z)" ! F (z),  b (z) b (z)  \ \ a (z) G (z)"! > ,  b (z) \

(z) b (z) H (z)" # > F (z),  a (z) a (z)  \ \ b (z) H (z)" > ,  a (z) \ where F(z) is a finite-order polynomial.

(17a) (17b) (17c) (17d) (17e) (17f )

¹heorem 3.1. The stable F (z) which minimizes the cost of  equation (11) subject to equation (10) is n (z) F (z)"   "(z)

(14)

where K is an optimal LQ state feedback controller. In the presence of constraints however such a predicted control law may be infeasible, and can only be invoked after N time instants, providing that N is large enough so that the predicted control increments from t to t#N!1 can be deployed to steer the x sufficiently close to the origin. The implied algorithm R>, (Scokaert and Rawlings, 1996), which will be referred to as

(15)

(18a)

n (z)"!b (z)S (z)#ja (z)S (z),  \  \  "(z)"C (z)"b(z)bC (z)#ja(z)aC (z),

(18b)

C

(18d)

C

" (z)S (z)#a (z)¹ (z)"z\b (z) (z),  \  "C (z)S (z)#b (z)¹ (z)"aC (z)t(z),  \ 

(18c)

(18e)

1586

Brief Papers

where "(z) is strictly Hurwitz, and (.)C denotes reversal of the order of the coefficients. Proof: Invoking Parseval’s theorem to equation (11) we have  J" [e #j*u ] R>G R>G G 1 dz ["E(z)"#j"*º(z)"] " j2n z X

(z) b(z) 1 a (z) "(z) \ " # j2n t(z) a(z) X j !j b (z) "(z) \

 

3.2. Relationship between ¸QCPC and ¸QIHPC. To accomplish this comparison, here we derive a transfer function formulation for LQCPC and to do this we assume that the control law of equation (14) is applied not only for k5N but for all nonnegative k, and introduce the degrees of freedom as perturbations c [Rossiter et al., 1997; Kouvaritakis et al., 1997], to the I control increments thus generated. Then providing c "0 for R>I k5N, the state space equations:

 F(z)

dz , z (19a)

"(z) F(z)" F(z), (19b) a (z)b (z) \ \ here "(z) is as defined in the statement of the theorem. Then using an appropriate unitary transformation (for details see Kouvaritakis et al., 1997) it is easy to show that:

 



z\bI(z) (z) aI(z)t(z)  !j #F(z) "I(z)a (z) "I(z)b (z) X \ \ j a (z) (z)#b (z)t(z) dz dz > > ; # j2n "(z) z z X ¹ (z) 1 ¹ (z)  !j  " KI(z) j2n KI(z) X S (z)  S (z) #  !j  #F(z) a (z) b (z) \ \ dz j a (z) (z)#b (z)t(z)  dz > > ; # (20) z j2n "(z) z X where S (z), S (z), ¹ (z) and ¹ (z) are as given in equation (18).     The first term (in square brackets) is antistable, whereas the second is stable. The optimal F(z) is defined by setting the second term equal to zero; solving for F(z) from equation (19b) returns the F (z) of the theorem.  1 J" j2n

      

 



tive performance between the three algorithms is dependent on the model, constraints and initial conditions. Intuitively however, it is obvious that LQIHPC bears a resemblance to LQCPC, and in the following section we explore this through the use of the state space model of equation (5).



Algorithm 3.1 (LQIHPC). Step 1. Compute the stable factor "(z) of b(z)bI(z)# ja(z)aI(z) and set t"0. Step 2. Compute the minimal order solution pairs + (z), t(z),, +S (z), ¹ (z),, and +S (z), ¹ (z), of the Bezout identities     (9b), (18c—e). Step 3. Apply QP to minimize over F the cost of equation (11) subject to equation (17, 1b). Step 4. Introduce the first element of F, f , together with the  first elements g, g of the impulse response of the G (z), G (z) of     equation (17), into equation (13) to compute the optimal current control increment Step 5. Increment t by one and return to Step 2. Remark 3.1. A simple alternative which avoids the explicit computation of F (z) is presented later in Theorem 3.2.  Clearly, this algorithm has the important attributes of IHPC, namely guaranteed stability and convergence, but in addition in the unconstrained case it recovers the optimal solution, despite the fact that F(z) is chosen to be a polynomial of finite order. It is reasonable to conjecture therefore that even in the constrained case it will give a better approximation to the ideal optimal solution (i.e. that of the minimization problem of equation (16)) than does IHPC. Unfortunately, given the nonlinear nature of the optimization problem, it is not possible to draw analytical comparisons between any of LQCPC, IHPC, LQIHPC (unless of course N for LQCPC is chosen such that N'N , in which  case LQCPC will retrieve the ideal optimal solution). As mentioned earlier, due to the heavy computation demands of QP, N has to be chosen to be small, and in such a case the compara-

x "(A!BK)x #Bc , e "Cx R>I> R>I R>I R>I R>I c "0 ∀k5N (21) R>I generate identical predictions to those of LQCPC. Thus we can state the following result. ¸emma 3.1. At a time instant t, the solution pair + (z), t(z), of equation (9) and the transfer function F (z) of Theorem 3.1,  satisfy the conditions

(z) b (z) # > F (z)"zC(zI![A!BK])\[A!BK]xt , a (z) a (z) MNR \ \ t(z) a (z) ! > F (z)"!K(zI![A!BK])\xt (22) b (z) b (z) MNR \ \ Proof. When the perturbation sequence c is taken to be R>I identically zero (for all nonnegative k), the state space model of equation (21) generates the optimal error predictions for the given initial equation xt which therefore must have the same z-transform as that given by equation (17) when F(z) is taken to be identically zero. Equation (22) follows from setting the ztransform of the predicted sequences e and Kx , ∀k50 of R>I R>I equation (21) (for c "0 ∀k50) equal to E(z), and *º(z) as R>I given by equation (17) (when F(z)"0). ¹heorem 3.2. The prediction equation for LQIHPC and those for LQCPC are given correspondingly as k (z) b (z) E(z)"  # > F(z), "(z) a (z) \ m (z) a (z) *º(z)"  ! > F(z), b (z) "(z) \ k (z) b(z) E(z)"  # c(z), "(z) "(z)

(23b)

m (z) a(z) *º(z)"  ! c(z), "(z) "(z)

(23d)

(23a)

(23c)

where the polynomial k (z), m (z) are defined by the Bezout   identity: a(z)k (z)#b(z)m (z)""(z)Q(z). (24)   Proof. For F(z)"0 and c(z)"0 the prediction equations for LQIHPC and LQCPC coincide with the unconstrained optimal predictions of equation (22) both of which have as denominator the determinant of (zI![A!BK]), the optimal closed-loop pole polynomial "(z). Let the respective numerator polynomials for E(z), *º(z) be k (z) and m (z). Then introducing these   expressions for the z-transforms of the optimal error and control increment predictions into equation (7) we derive equation (24). Equatiosn (23a) and (b) then follow directly from equation (17), whereas equation (23c) follows from the z-transforms of the e sequence of equation (21) for c O0, k"0, 1, 2 , N!1: R>I R>I E(z)"zC(zI![A!BK])\[A!BK]x R #C(zI![A!BK])\Bc(z).

(25)

The first term is clearly k (z)/"(z), whereas the numerator of  the second is the zero polynomial, b(z), of the open-loop system; state feedback does not affect the system zeros. Equation (23d)

Brief Papers can be obtained by substituting equation (23c) into equation (7) in conjunction with equation (24). Remark 3.2. Equations (17) and (23) clearly produce the same class of stable predictions, however, from a computational viewpoint, equation (23) is more efficient: it requires the solution of one Bezout identity, that given in equation (24), as opposed to the two Bezout identities required by equation (17) and given in equation (18). This change can be incorporated trivially in Step 2 of Algorithm 3.1. Equation (23) affords useful insight into the relationship between LQCPC and LQIHPC. It shows for example that for F(z) and c(z) both polynomials of arbitrarily large order, the two approaches become identical. It also shows that for the same number of degrees of freedom, when we can set F(z)"c(z), the only difference between the control increment predictions in LQCPC and LQIHPC is the filtering of the perturbation term; in the case of LQCPC c(z) is filtered by G "[!a(z)/"(z)]  whereas in LQIHPC they are filtered by G "[!a (z)/b (z)].  > \ In general, for jO0 the poles of G (z) will be faster than those of  G (z) and so the LQCPC filtering will be more aggressive. As  a consequence, it is reasonable to conjecture that LQCPC is more likely to run into infeasibility than LQIHPC. Given the nonlinear nature of the constrained optimization, it is not possible to provide an analysis that will support this conjecture, but Monte-Carlo simulations (discussed in the examples section) provide some evidence. One fact however is beyond doubt: both LQCPC and LQIHPC can run into infeasibility whenever the constraints are tight and/or the number of degrees of freedom is small, and what is worse, for some examples this infeasibility problem is not overcome unless one is prepared to use an arbitrarily large number of degrees of freedom (illustrations of this are included in the Examples section). It is noted that if infeasibility is the result of a stringent setpoint change, then it is possible to replace the new demanded setpoint by a slack variable (Rossiter et al., 1996b; Bemporad and Mosca, 1994) which can be steered (as fast as the constraints will allow) towards the desired setpoint. However, infeasibility may also arise due to disturbances, in which case one must anticipate the worst case disturbance and invoke an artificially tight set of control constraints to ensure that feasibility at t will be preserved at the next sampling instant (Gossner et al., 1997). Such a strategy of course decreases the control authority that is available for the purposes of performance optimization. In this perspective it would be advantageous to devise algorithms which are not susceptible to problems of feasibility (in the disturbance-free case in any event). It is the purpose of the following section to propose one such algorithm.

1587

matrix associated with the elements of the impulse response of the transfer function G (z) of equation (10). Unfortunately, the  computation of such an F(z) is impracticable (because it may be necessary to use an F(z) of arbitrarily high order). To avoid this problem in the sequel we make the following (weak) assumption. Assumption 4.1. If at any time equation (26) admits a feasible finite-dimensional F then the minimization of #g !C F#  E  over infinite-dimensional F yields a feasible minimizer F :   #g #C F # 4*º . (27)     E    In contrast to equation (26) the optimal solution to equation (27) is known explicitly as stated below. ¹heorem 4.1. The solution of the optimization problem in the left-hand side of equation (27) and corresponding control-increment prediction equation are given by R(z) F (z)" ,   *(z)aI(z) > ¹(z) *º (z)" .   aI(z) > where R(z) and ¹(z) are defined by the bezout identity

(28a) (28b)

b (z)¹(z)#a (z)R(z)"aI(z)t(z) (29) \ > > Hence, if DU denotes the vector of elements of the impulse   response of *º (z), Assumption 4.1 implies   ¹(z) l "#DU # " Z\ 4*º (30)         aI(z)  > where ZU is the inverse z-transform operator that produces the vector of elements of the impulse response of a transfer function.

  

Proof. Using the l optimization procedure we invoke Par seval’s theorem to equation (10) to get

 



 1 t(z) a (z) dz #DU#" *u " ! > F(z) R>G j2n b (z) b (z) z X \ \ G 1 aI(z)t(z) aI(z) dz 1 > " ! > F(z) " j2n a (z)b (z) b (z) z j2n X > \ X \ aI(z)t(z) aI(z)  dz ! > ! > F(z) a (z)b (z) b (z) z > \ \ 1 ¹(z) R(z) aI(z)  dz " # # > F(z) (31) j2n a (z) b (z) b (z) z X > \ \ which is minimized for the F(z) of equation (28a) since that sets the square bracket equal to zero; the decomposition into terms involving ¹(z), R(z) used above is possible via the bezout identity (29). Substitution of F (z) into equation (10a) and use of the   bezout identity yields (28b).



   











4. ¸Q predictive control with improved feasibility The weakness of IHPC is that for F(z) of finite degree, the class of predictions of equation (10a) may not include the unconstrained optimal. Hence, the instantaneous optimal cost J can be unduly conservative even when the constraints are not overly active. A large predicted cost J does not necessarily imply a large actual running cost; after all at the next time instant further degrees of freedom will be available to reduce the running cost. However, the instantaneous cost provides an upper bound for the running cost, and nonconservative instantaneous costs are desirable. Therefore, it is desirable that the class of allowable predictions should include the unconstrained optimum. Both LQCPC and LQIHPC were designed to meet this objective, but on account of it are likely to run into problems of infeasibility. To avoid this situation a further obvious objective to aim for is that the class includes a stabilizing predicted control increment trajectory which is feasible; by Assumption 2.1 it is known that at start time at least one such exists, and that therefore equation (16) admits an ideal optimal solution. The existence of such a solution implies that an F(z) of finite but large order exists (at start time) such that:

Remark 4.1. F(z) is not in l , on account of the integrator  introduced through *(z), it nevertheless yields *º which is in   l and therefore provides an admissible solution. However to  avoid technical difficulties, and in particular to derive an F(z) which yields error predictions in l , henceforth we shall replace  the F (z) of equation (28a) and the corresponding prediction   equation for the control increments by

min #g #C F# 4*º (26)  E     $ where F is the vector of coefficients of F(z), g the vector of the  elements of the impulse response of the transfer function G (z) of  equation (10) and C  denotes the Toeplitz striped convolution E

where e'0 is chosen to be sufficiently small so as not to have an appreciable effect on l (equation (30); equation (32b) is derived   by substituting equation (32a) into equation (10) and invoking the identity of equation (29). Note that an alternative choice for F (z) (not used here) is that of equation (18) for j large.  

R(z) F (z)" ,   D(z)aI(z) > ¹(z) [D(z)!*(z)]a (z)R(z) > *º (z)" # ,   aI (z) b (z)D(z)aI (z) > \ > D(z)"1!(1!e)z\,

(32a) (32b) (32c)

1588

Brief Papers

F (z) will be referred as the ‘‘mean level ’’ solution and will be   used in place of the F (z) of Section 3, equation (17). However,  our preoccupation here is not only with feasibility: we also require our class of predicted control increments to include the unconstrained optimal. The perturbation on F (z) that re  trieves F (z) is given by  n (z) B$ , dF(z)"F (z)!F (z)"    "(z)D(z)aI (z) > n (z)"n (z)D(z)aI (z)!"(z)R(z), B$  >

(33)

where n (z) is as given in equation (18). Thus, we have the  following obvious result: Corollary 4.2. Denote by l the degree of n (z), and let p (z) be B$  a factor of n (z), of degree n ; furthermore let c(z) be a polyB$  nomial of degree n . Then providing that n 5l!n the class of A  A predictions of equation (10) for p (z)c(z)  F(z)"F (z)# (34)   "(z)D(z)aI (z) > contains both the unconstrained optimal and the mean level solutions. Proof. The choice of c(z)"0 recovers the mean level solution, whereas since n 5l!n , c(z) can be chosen to be n (z)/p (z) A  B$  thereby recovering the unconstrained optimal solution. Remark 4.2. In practice the factorization of n (z) is not necesB$ sary. Instead with n "N, p (z) can be taken equal to n (z) for A  B$ N(l, and 1 otherwise without affecting the corollary; the implied prediction class contains both the unconstrained optimal and mean level solutions. The general form of the F(z) used above is F(z)"cF (z)#(1!c)F (z)#'(z), where c is    a constant and '(z) denotes a stable perturbation. The implied particular choice of '(z) is by no means unique, but it is a convenient one and produces good results (see Section 5). Remark 4.3. Corollary 4.2/Assumption 4.1 suggest an LQ predictive algorithm which is always feasible. Algorithm 4.1 (LQFPC - LQ feasible predictive control). Step 1. Compute the stable factor "(z) of b(z)bI (z)#ja(z)aI(z) and set t"0; select the number of degrees of freedom, n , to be used in the LQ optimization. A Step 2. Compute the minimal order solution pairs + (z), t(z),, +S (z), ¹ (z),, and +S (z), ¹ (z), of the Bezout identities     (9b), (18c,d) and hence from equations (18) and (33) compute n (z) and n (z); if l, the degree of n (z) is less than n set  B$ B$ A p (z)"1, otherwise select a suitable factor p (z) of n (z).   B$ Step 3. Apply QP to minimize over F (as defined by equation (34)) the cost equation (11) subject to equations (10) and (1b) Step 4. Introduce the first element f of F together with the  first elements g , g of the impulse response of the G (z), G (z) of     equation (10), into equation (13) to compute the optimal current control increment Step 5. Increment t by one and return to Step 2. Clearly Algorithm 4.1 has all the attributes of stability and convergence of IHPC, LQCPC and LQIHPC, but the guarantee of these attributes does not require any strong (restrictive) assumptions (over and above Assumption 2.1 which is clearly needed for any predictive control algorithm): we merely assume that the weak condition of Assumption 4.1 is satisfied. This

implies some considerable advantages which are illustrated in the next section. 5. Numerical examples 5.1. ¹esting optimality and comparing costs. The advantage of LQCPC is that if N is large enough, and in particular larger than N , then it yields optimal results. However, it was claimed that  N can often be prohibitively large and that if N is chosen too  small and/or the constraints are tight, LQCPC can result in infeasibility. In this section simulations are performed to illustrate this point. Models are generated with four random poles, one of which is unstable, and three random zeros, one of which is placed outside the unit circle; thus all the models considered are unstable and non-minimum phase. In addition, random initial conditions are ensured by driving the models (prior to simulation) for a period of 10 samples by a random input; subsequently, in this subsection, the various algorithms discussed of the paper are applied at the first time instant only. For each run, the infinity norms l ,  l of the predicted control increments for the unconstrained   optimal solution and the mean level solution are computed. This provides a mechanism for defining tight or loose constraints through *º "k l #(1!k) l . (35)       Values of k close to 1 result in stringent rate constraints, whereas values close to 0 should present little challenge to all the algorithms described in this paper. Throughout Section 5, for convenience we assume a common value for the j parameter of the predictive control cost J, namely we set j"1, and unless otherwise stated D(z)"1!0.95z\. Unsurprisingly for N"2 (small number of degrees of freedom) and k"0.9 LQCPC, LQIHPC and IHPC run into infeasibility often, yet the costs (including that for LPFPC) averaged over all cases where no algorithms run into infeasibility difficulties were similar as shown in Table 1. What is surprising is that at times LQCPC failed even though the other algorithms achieved a cost comparable to the unconstrained as illustrated by the results of trials listed in Table 2. As is apparent from the tables above, it is not possible to tell a priori which of LQCPC, LQIHPC or IHPC are likely to fail, or which of LQFPC, LQIHPC or IHPC will give the smallest cost when LQCPC fails. However when LQCPC did not fail, it gave the best cost, though the difference was not significant (less than 3%). The incidence of infeasibility drops considerably with a significant relaxation of the constraints or a more modest constraints Table 1. Incidence of infeasibility (N"2, k"0.9, 200 trials)

Failure rate (%) Average cost

LQFPC

LQCPC

LQIHPC

0 1.5430

80 1.3374

76 1.6072

Trials"200; N"2; k"0.3 Average cost Trials"200; N"4; k"0.5 Average cost

0% 1.0988 0% 1.2335

66 1.7380

Table 2. Comparison of costs (N"2, k"0.9) Unconst. opt

LQFPC

LQCPC

LQIHPC

IHPC

27.0413 10.8410 1.1089

27.7796 16.0966 3.7983

Failed Failed Failed

Failed Failed 2.4016

Failed 21.2061 2.0631

Table 3. Infeasibility incidence for greater N and/or smaller k LQFPC

IHPC

LQCPC 5% 1.0479 18% 1.2281

LQIHPC 7% 1.0825 25% 1.2359

IHPC 23% 1.262 26% 1.3051

Brief Papers Table 4. Infeasibility and comparison of costs for N"10; k"0.9

Failure rates (%) Average costs

LQFPC

LQCPC

LQIHPC

0 1.9871

20 1.9758

22 1.9996

IHPC 26 2.0154

relaxation combined with an increase of N as shown in Table 3. However, the observations made above still hold true: (i) LQCPC failed even when LQFPC produced costs close to the unconstrained optimal; e.g. for a particular run the

1589

unconstrained cost was 5.8292, LQFPC produced a cost equal to 6.1198, yet LQCPC failed (ii) It is not known a priori which of LQCPC, LQIHPC and IHPC will fail (iii) When LQCPC did not fail it yielded the lowest cost but the difference was not significant (e.g. less than 0.5% for N"4, k"0.5). Stringent constraints result in a large failure rate for LQCPC, LPIHPC and IHPC even if N is chosen to be large (though of course the incidence rates are not as high as those of Table 1) as shown by the results of 200 runs summarized in Table 4. In this case none of LQCPC, LQIHPC and IHPC have any advantages to offer: their incidence of infeasibility is high, and

Fig. 1. Infeasibility of LQCPC and comparison of cost with LQFPC; D(z)"1!0.95z\, k"0.9, *º "0.1291.   

Fig. 2. Infeasibility of LQCPC and comparison of cost with LQFPC; D(z)"1!0.99z\, k"0.95, *º "0.0452.   

1590

Brief Papers

even in the case of no failure, they produce average costs which are to within 0.6% of the average cost produced by LQFPC.

a(z)"1!1.7361z\!0.74z\#0.3324z\ #0.05761z\,

5.2. ¹esting convergence to optimality and closed-loop responses. When N is small, LQCPC can produce very good (near opti mal) results. Tight constraints however often cause N to be come impractically large, and in cases like this the advantage of the optimality associated with LQCPC is no longer apparent. To test this here we consider a fixed randomly generated model and initial conditions, but varying N. The results plotted in Figure 1 illustrate a typical situation; they correspond to a model for which:

b(z)"0.3211!2.1571z\#1.6084z\!0.316z\, u "10\[0.6517,!0.5578, 0.5956,!0.3559,   !0.086, 0.7151,!0.735, 0.5964, 0.8643], *º "0.1291, k"0.9, (36)    where u denotes the random input sequence used to generate   the random initial conditions. The value of J"!1 for LQCPC

Fig. 3. Closed loop LQFPC input and output responses for the model of equation (36).

Fig. 4. Infeasibility of LQCPC and comparison of cost with LQFPC for model of equation (37).

Brief Papers when N"2,3, 2 , 12 shows that for those N the algorithm was infeasible. Note that the LQFPC cost initially is large, but for N'7 is within 5% of the constrained optimum value of 5.4618. By contrast LQCPC remains infeasible up to N"12, assumes the value of 5.4745 at N"13 and reaches the optimum for N"15. Thus, by the time LQCPC becomes feasible the difference in costs is insignificant, whereas LQFPC homes on values close to the optimal for values of N for which LQCPC is infeasible. The plots of Fig. 2 are less typical (due to the use of tight constraints) but illustrate that N can be arbitrarily large. The  same model and initial conditions are used, but now D(z)"1!0.99z\, k"0.95, Dº "0.0452. LQCPC re  

1591

mains infeasible for all N(43; still larger values arise out of tighter limits which occur as D(z) tends to 1!z\ and k tends to 1. Returning to the data of equation (36) we next consider closed-loop responses. LQFPC does not run into feasibility problems for any value of N, and indeed even for small values of N gives satisfactory responses as the ones shown in Fig. 3 for N"4. In contrast to this LQCPC remains unstable (due to infeasibility) for N(12 and produces almost identical responses to those of Fig. 3 for N"12; the total running cost here is 5.462 and is within 0.7% of the corresponding value obtained by LQFPC for N"4!.

Fig. 5(a). Closed loop LQCPC input and output responses for the model of equation (37). (b) Closed loop LQFPC input and output responses for the model of equation (37).

1592

Brief Papers

Infeasibility need not always cause such instability problems as is demonstrated by the example below for which a(z)"1#0.9026z\!1.7924z\#0.113z\#0.1261z\, b(z)"0.247!0.8957z\#0.5027z\!0.0178z\, u "10\[!0.9309, !0.8931, 0.0594, 0.3423, !0.9846,   !0.2332,!0.8663, !0.165, 0.3735] *º "0.1667, k"0.9, (37)    For this data LQCPC remains infeasible for all values N(12 (see Fig. 4), yet the closed-loop responses for even small values of N, (e.g. for N"2 see Fig. 5) are stable (albeit erratic for 25 samples). The plots were obtained using MATLAB which, for infeasible QP problems minimizes constraint violations by considering a subset of active constraints and returning a suitable least squares solution. The corresponding predicted control trajectory is of course still infeasible, but the first predicted control move may or may not exceed the constraint limits thereby, respectively, obviating or necessitating the use of ‘‘clipping’’ prior to implementation. LQFPC on the other hand, does not run into infeasibility and produces better responses (for the same number of degrees of freedom), these are shown in Fig. 5b. 6. Conclusions In summary, despite their appealing properties, IHPC, LQCPC, and LQIHPC often run into infeasibility difficulties for small numbers of degrees of freedom, and therefore result in instability or poor closed-loop responses. In instances like this it is of course possible to invoke setpoint conditioning techniques to overcome some of the relevant problems associated with large setpoint changes; infeasibility however would still remain a problem with respect to unknown (but bounded disturbances) and would necessitate artificial tightening of the input constraints which in turn would result in unduly conservative results. The algorithm proposed in this paper offers an alternative to setpoint conditioning and overcomes the infeasibility problem at the expense of a modest (often insignificant) increase in cost which nevertheless is close to optimal except under very stringent conditions. Acknowledgements—The authors thank the Engineering and Physical Sciences Research Council for financial support. References Alvarez-Ramirez, J. and R. Suarez (1996). Global stabilization of discrete-time linear systems with bounded inputs, Int. J. Adaptive Control and Signal Process., 10, 409—416.

Bemporad, A. and E. Mosca (1994). Constrained fulfilment of feedback control via prediction management reference, Proc. 3rd IEEE CCA, Glasgow, pp. 1909—1914. Clarke, D. W., C. Mohtadi and P. S. Tuffs (1987). Generalized predictive control, Parts 1 and 2. Automatica, 23, 137—160. Clarke, D. W. and R. Scattolini (1991). Constrained receding horizon predictive control., Proc. IEE, Pt.-D, 138 (4), 347—354 Gilbert, E. G. and K. T. Tan (1991). Linear systems with state and control constraints: The theory and practice of maximal admissible sets. IEEE ¹rans. AC, 36(9), 1008—1020. Gossner, J. R., B. Kouvaritakis and J. A. Rossiter (1997). Stable generalized predictive control with constraints and bounded disturbances. Automatica, 33(4), 551—568. Kouvaritakis, B., J. A. Rossiter and A. O. T. Chang (1992). Stable generalized predictive control: an algorithm with guaranteed stability. Proc. IEE, Pt-D, 139(4), 349—362. Kouvaritakis, B., J. R. Gossner and J. A. Rossiter (1996). A priori stability conditions for an arbitrary number of unstable poles. Automatica, 32(10), 1441—1446 Kouvaritakis, B., J. A. Rossiter, M. Cannon and C. K. J. Wong (1997). Constrained stable predictive control: an l -optimal  approach. 5th IEEE MED conf., Cyprus. Mosca, E. and J. Zhang (1992). Stable redesign of predictive control. Automatica, 28(6), 1229—1233. Rawlings, J. B. and K. R. Muske (1993). The stability of constrained receding horizon control. IEEE ¹rans. AC, 38(10), 1512—1516 Rossiter, J. A., B. Kouvaritakis and J. R. Gossner (1995). Feasibility and stability results for constrained stable generalized predictive control. Automatica, 31(6), 863—877. Rossiter, J. A., J. R. Gossner and B. Kouvaritakis (1996a). Infinite horizon stable predictive control. IEEE ¹rans. AC, 41(10), 1522—1527 Rossiter, J. A., B. Kouvaritakis and J. R. Gossner (1996b). Guaranteeing feasibility in constrained stable generalized predictive control. Proc. IEE, Pt-D, 143(5), 463—469. Rossiter, J. A., M. J. Rice, B. Kouvaritakis (1997). A robust statespace approach to stable predictive control strategies, ACC ’97, Albuquerque. Scokaert, P. O. M., and J. B. Rawlings (1996). Infinite horizon linear quadratic control with constraints. IFAC Congress ’96 Sussman, H. J., E. D. Sontag and Y. Yang (1994). A general result on the stabilization of linear systems using bounded control. IEEE ¹rans. AC, 39(12), 2411—2424.